Попался вот такой кусок кода для вычисления тривиального числа Pi
n = 1000000000
62 h = 1.0 / (double) n;
63 sum = 0.0;
64
65 for (i = 1; i <= n; i++)
66 {
67 x = h * ((double)i - 0.5);
68 sum +=(4.0 / (1.0 + x*x));
69 }
70
71 mypi = h * sum;
изначально все строилось для MPI , правда захотелось перенести на OpenMP
казалось бы все просто:
53 #pragma omp parallel for reduction(+: sum)
54 for (int i = 1; i <= n; i++)
55 {
56
57 x = h * ((double)i - 0.5);
58 sum +=(4.0 / (1.0 + x*x));
59 }
(гусары - молчать , пока что ... ).
(компилируем как gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=sse)
Практически все (окромя точности, но об этом позже) устраивает и почти красиво, пример Thread Profiler for Linux:
tprofile_cl ./pim
Building project
Instrumenting
16% pim ( API Imports ):..
Running: /mpi/pi/pi/pim
n = 1000000000
Num threads before = 1
myid: 1 i: 0: 0 c:0
mypi: 3.1415926535898211
Application finished
Intel® Thread Profiler 3.1 Summary Report
application: /mpi/pi/pi/pim
collection: Wed Jan 7 17:52:15 2009
runtime: 3.43215s
# of processors: 4
# of threads: 4
# of waits: 3
wait frequency: 0.874088
average concurrency: 3.96682
Concurrency:
0 [..........] 0% 0
1 [..........] 0.913% 0.0301119
2 [..........] 0.0051% 0.000168224
3 [..........] 0.349% 0.0114959
4 [##########] 98.7% 3.25637
Паралельность просто фантастическая!
однако сюрпризы только начинаються если присмотреться подольше и поближе :) ,
вот пример нескольких последовательных, и не только, запусков:
[maxd@fc9i pi]$ time ./pim
n = 1000000000
mypi: 3.1415926535898211
real 0m3.345s
user 0m13.125s
sys 0m0.032s
[maxd@fc9i pi]$ time ./pim
n = 1000000000
mypi: 3.1415926535898211
real 0m55.433s
user 3m38.935s
sys 0m0.033s
В итоге - иногда время работы ~3 сек и очень часто - ~55 сек - крута, 18х slowdown. Что то мы не подкрутили!
посмотрим поглубже и видим что то типа
400710: f2 0f 2a c0 cvtsi2sd %eax,%xmm0
400714: 83 c0 01 add $0x1,%eax
400717: 66 0f 28 ea movapd %xmm2,%xmm5
40071b: 39 c2 cmp %eax,%edx
40071d: f2 0f 5c c4 subsd %xmm4,%xmm0
400721: f2 41 0f 59 44 24 08 mulsd 0x8(%r12),%xmm0
400728: f2 41 0f 11 44 24 10 movsd %xmm0,0x10(%r12)
40072f: f2 0f 59 c0 mulsd %xmm0,%xmm0
400733: f2 0f 58 c3 addsd %xmm3,%xmm0
400737: f2 0f 5e e8 divsd %xmm0,%xmm5
40073b: f2 0f 58 cd addsd %xmm5,%xmm1
40073f: 7d cf jge 400710 <main.omp_fn.0+0x70>
(btw AT&T синтаксис вовсе не дурен)
запустив VTune/PTU увидим впечатляющий CPI (циклов на инструкцию) ~ 45 :) при желаемом -как можно ближе к 0.25 но меньше 1 уже хорошо . Кстати - саму ошибку уже видно с разных сторон. Для разнообразия запустим Thread Checker for Linux и уже наглядно видим что то типа
tcheck_cl ./pim
Intel® Thread Checker 3.1 command line instrumentation driver (26185)
Copyright (c) 2007 Intel Corporation. All rights reserved.
Building project
Instrumenting
16% pim ( All Functions ):..
Running: /mpi/pi/pi/pim
n = 1000000000
Num threads before = 1
myid: 1 i: 0: 0 c:0
mypi: 3.1415926535898211
Application finished
....
(форматирование несколько "поехало") - кстати Thread Checker оказался самым быстрым способом борьбы.
Тут даже придумано красивое название - False Sharing: несколько потоков мешают друг другу работая с одним ресурсом. Как другой вариан объяснения - читайте blog недавно были интересные блоги насчет красивых и иногда жизненых примеров для общих проблем.
На четырех ядерной машинке 18х ухудшения - интересно.
Как вылечить ? Например так :
55 #pragma omp parallel for private(x) reduction(+: sum)
56 for (i = 1; i <= n; i++)
57 {
58 x = h * ((double)i - 0.5);
59 sum +=(4.0 / (1.0 + x*x));
60 }
добавленое - выделено.
Кстати - исходный код был С где declaration для переменной x было вне цикла, как можно было сделать для С++ и тогда проблемы не было бы вообще.
в итоге наш цикл превращаеться в более простую конструкцию типа:
400710: f2 0f 2a c0 cvtsi2sd %eax,%xmm0
400714: 83 c0 01 add $0x1,%eax
400717: 66 0f 28 f3 movapd %xmm3,%xmm6
40071b: 39 c2 cmp %eax,%edx
40071d: f2 0f 5c c5 subsd %xmm5,%xmm0
400721: f2 0f 59 c2 mulsd %xmm2,%xmm0
400725: f2 0f 59 c0 mulsd %xmm0,%xmm0
400729: f2 0f 58 c4 addsd %xmm4,%xmm0
40072d: f2 0f 5e f0 divsd %xmm0,%xmm6
400731: f2 0f 58 ce addsd %xmm6,%xmm1
400735: 7d d9 jge 400710 <main.omp_fn.0+0x70>
которая выполняеться и быстрее и стабильнее.
Насчет точности:
Используемый алгоритм подвержен проблеме округления в частности при использовании OpenMP.
например для использования SSE можно попробывать использовать
gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=sse
для FPU:
gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=387
результаты будут разные:
sse: 3.1415926535898211
387: 3.1415926535897936
org: 3.14159265358979323846
где org - просто Pi взятое из какой либо википедии, кстати FPU более точнее но об этом я уже как то писал.
Время работы FPU кстати неплохое
sse: real 0m3.280s
fpu: real 0m3.894s
Присутствующая векторизация в gcc для этого цикла не сработала и пришлось бы менять код.
Intel Compiler сделал вещи шустрее (используя packed arguments instead of scalar ~ векторизация , loop unrolling etc ) и с уже приведенным C кодом :
icc ./pim.c -O3 -xSSSE3 -openmp -vec-report5 -o ./pim.icc
и результат неплох
real 0m1.857s
user 0m6.528s
sys 0m0.033s
За сим прощаюсь с успехом поздравив Вас в удачной оптимизации :)
n = 1000000000
62 h = 1.0 / (double) n;
63 sum = 0.0;
64
65 for (i = 1; i <= n; i++)
66 {
67 x = h * ((double)i - 0.5);
68 sum +=(4.0 / (1.0 + x*x));
69 }
70
71 mypi = h * sum;
изначально все строилось для MPI , правда захотелось перенести на OpenMP
казалось бы все просто:
53 #pragma omp parallel for reduction(+: sum)
54 for (int i = 1; i <= n; i++)
55 {
56
57 x = h * ((double)i - 0.5);
58 sum +=(4.0 / (1.0 + x*x));
59 }
(гусары - молчать , пока что ... ).
(компилируем как gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=sse)
Практически все (окромя точности, но об этом позже) устраивает и почти красиво, пример Thread Profiler for Linux:
tprofile_cl ./pim
Building project
Instrumenting
16% pim ( API Imports ):..
Running: /mpi/pi/pi/pim
n = 1000000000
Num threads before = 1
myid: 1 i: 0: 0 c:0
mypi: 3.1415926535898211
Application finished
Intel® Thread Profiler 3.1 Summary Report
application: /mpi/pi/pi/pim
collection: Wed Jan 7 17:52:15 2009
runtime: 3.43215s
# of processors: 4
# of threads: 4
# of waits: 3
wait frequency: 0.874088
average concurrency: 3.96682
Concurrency:
0 [..........] 0% 0
1 [..........] 0.913% 0.0301119
2 [..........] 0.0051% 0.000168224
3 [..........] 0.349% 0.0114959
4 [##########] 98.7% 3.25637
Паралельность просто фантастическая!
однако сюрпризы только начинаються если присмотреться подольше и поближе :) ,
вот пример нескольких последовательных, и не только, запусков:
[maxd@fc9i pi]$ time ./pim
n = 1000000000
mypi: 3.1415926535898211
real 0m3.345s
user 0m13.125s
sys 0m0.032s
[maxd@fc9i pi]$ time ./pim
n = 1000000000
mypi: 3.1415926535898211
real 0m55.433s
user 3m38.935s
sys 0m0.033s
В итоге - иногда время работы ~3 сек и очень часто - ~55 сек - крута, 18х slowdown. Что то мы не подкрутили!
посмотрим поглубже и видим что то типа
400710: f2 0f 2a c0 cvtsi2sd %eax,%xmm0
400714: 83 c0 01 add $0x1,%eax
400717: 66 0f 28 ea movapd %xmm2,%xmm5
40071b: 39 c2 cmp %eax,%edx
40071d: f2 0f 5c c4 subsd %xmm4,%xmm0
400721: f2 41 0f 59 44 24 08 mulsd 0x8(%r12),%xmm0
400728: f2 41 0f 11 44 24 10 movsd %xmm0,0x10(%r12)
40072f: f2 0f 59 c0 mulsd %xmm0,%xmm0
400733: f2 0f 58 c3 addsd %xmm3,%xmm0
400737: f2 0f 5e e8 divsd %xmm0,%xmm5
40073b: f2 0f 58 cd addsd %xmm5,%xmm1
40073f: 7d cf jge 400710 <main.omp_fn.0+0x70>
(btw AT&T синтаксис вовсе не дурен)
запустив VTune/PTU увидим впечатляющий CPI (циклов на инструкцию) ~ 45 :) при желаемом -как можно ближе к 0.25 но меньше 1 уже хорошо . Кстати - саму ошибку уже видно с разных сторон. Для разнообразия запустим Thread Checker for Linux и уже наглядно видим что то типа
tcheck_cl ./pim
Intel® Thread Checker 3.1 command line instrumentation driver (26185)
Copyright (c) 2007 Intel Corporation. All rights reserved.
Building project
Instrumenting
16% pim ( All Functions ):..
Running: /mpi/pi/pi/pim
n = 1000000000
Num threads before = 1
myid: 1 i: 0: 0 c:0
mypi: 3.1415926535898211
Application finished
_______________________________________________________________________________
|ID|Short De|Seve|C|Contex|Description |1st Acc|2nd Acc|
| |scriptio|rity|o|t[Best| |ess[Bes|ess[Bes|
| |n |Name|u|] | |t] |t] |
| | | |n| | | | |
| | | |t| | | | |
_______________________________________________________________________________
....
|3 |Write ->|Erro|2|[pim, |Memory write at [pim, 0x728] |[pim, |[pim, |
| |Write da|r |1|0x6a0]|conflicts with a prior memory write|0x728] |0x728] |
| |ta-race | |4| |at [pim, 0x728] (output dependence)| | |
| | | |8| | | | |
| | | |7| | | | |
| | | |6| | | | |
| | | |1| | | | |
| | | |1| | | | |....
(форматирование несколько "поехало") - кстати Thread Checker оказался самым быстрым способом борьбы.
Тут даже придумано красивое название - False Sharing: несколько потоков мешают друг другу работая с одним ресурсом. Как другой вариан объяснения - читайте blog недавно были интересные блоги насчет красивых и иногда жизненых примеров для общих проблем.
На четырех ядерной машинке 18х ухудшения - интересно.
Как вылечить ? Например так :
55 #pragma omp parallel for private(x) reduction(+: sum)
56 for (i = 1; i <= n; i++)
57 {
58 x = h * ((double)i - 0.5);
59 sum +=(4.0 / (1.0 + x*x));
60 }
добавленое - выделено.
Кстати - исходный код был С где declaration для переменной x было вне цикла, как можно было сделать для С++ и тогда проблемы не было бы вообще.
в итоге наш цикл превращаеться в более простую конструкцию типа:
400710: f2 0f 2a c0 cvtsi2sd %eax,%xmm0
400714: 83 c0 01 add $0x1,%eax
400717: 66 0f 28 f3 movapd %xmm3,%xmm6
40071b: 39 c2 cmp %eax,%edx
40071d: f2 0f 5c c5 subsd %xmm5,%xmm0
400721: f2 0f 59 c2 mulsd %xmm2,%xmm0
400725: f2 0f 59 c0 mulsd %xmm0,%xmm0
400729: f2 0f 58 c4 addsd %xmm4,%xmm0
40072d: f2 0f 5e f0 divsd %xmm0,%xmm6
400731: f2 0f 58 ce addsd %xmm6,%xmm1
400735: 7d d9 jge 400710 <main.omp_fn.0+0x70>
которая выполняеться и быстрее и стабильнее.
Насчет точности:
Используемый алгоритм подвержен проблеме округления в частности при использовании OpenMP.
например для использования SSE можно попробывать использовать
gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=sse
для FPU:
gcc ./pim.c -o ./pim -O3 -fopenmp -mfpmath=387
результаты будут разные:
sse: 3.1415926535898211
387: 3.1415926535897936
org: 3.14159265358979323846
где org - просто Pi взятое из какой либо википедии, кстати FPU более точнее но об этом я уже как то писал.
Время работы FPU кстати неплохое
sse: real 0m3.280s
fpu: real 0m3.894s
Присутствующая векторизация в gcc для этого цикла не сработала и пришлось бы менять код.
Intel Compiler сделал вещи шустрее (используя packed arguments instead of scalar ~ векторизация , loop unrolling etc ) и с уже приведенным C кодом :
icc ./pim.c -O3 -xSSSE3 -openmp -vec-report5 -o ./pim.icc
и результат неплох
real 0m1.857s
user 0m6.528s
sys 0m0.033s
За сим прощаюсь с успехом поздравив Вас в удачной оптимизации :)
Комментарии
Ох уж мне эти линуксоиды :)
тег "code" лечит проблемы форматирования.
Кстати, Максим меня несколько опередил - я как раз подумывал о следующем этапе "объяснений на пальцах"... Хотел провести небольшой конкурс, где нужно было бы привести пример кода с той или иной ошибкой, и, соответственно, варианты решения. Впрочем, почему "хотел"? :) - все еще хочу.
спасибочки
а то вот пробывал "man code" - не работает :D
(юмор с уклоном в сторону Линукса, man - отличная справка в Линуксе)
а вот на чем я как то попался, правда простой последовательный код:
for(unsigned char uc = 0; uc <= 255; ++uc)
{
// do something
.....
}
ой сколько же я ждал завершения программы, пока не понял что (uc <= 255) выполняется всегда
2 Ильнар: это ж классика :)
Что-то я не понимаю ни откуда возникает проблема, ни почему именно так она решается.
Видно, что в ассемблерном коде при умножении на шаг h возникает чтение из памяти (по-видимому, чтение h) и затем запись в память (чего?). Почему h не оптимизируется в регистр, непонятно. Что и куда пишется, непонятно абсолютно. Причём тут False Sharing и где он возникает? x - автоматическая переменная на стэке потока, и в разных потоках локальные "копии" x вряд ли разделяют одну кэш линию. И уж если Thread Сhecker ругается на состязание (data race), тут точно не false, а true sharing - но на какую переменную?
Похоже, ситуация снова из разряда "компилятор не любит параллельный код", но какая-то совсем уж удивительная, с использованием автоматической переменной на стэке, которая по всем правилам (включая стандарт OpenMP) должна быть локальной для потока, и оптимизироваться в регистр.
Не понимаю :)
рад что тема вызвала такой хороший вопрос - пока еще оставлю на - "подумать глубже" :)
только скромно добавлю - data race с переменной h не может быть - ибо эта переменная не меняеть по всему циклу.
Плюс могу добавить еще один хинт - каково должно быть значение x после выхода из цикла (это глобальная переменная и ее видимость после цикла я думаю очевидна )
(перечитывая позже увидел опечатку внутри for () { double x = } - это как бы решение для С++ причем лишеное проблемы , естественно поменял , спасибо за хинт ! )
Максим, ваши записи в блоге для меня - одни из самых интересных. Жаль что вы редко пишете. Но это наверно потому, что качество и количество имеют связь :-)
спасибо, ksili, за оценку , рад что нравиться.
Писать чаще надо, пока что стремлюсь писать регулярно :)
уже есть идеи и к следующей теме, надо ее только еще хорошенько проработать :)
> (перечитывая позже увидел опечатку внутри for () { double x = } - это как бы решение для С++ причем лишеное проблемы , естественно поменял , спасибо за хинт ! )
А, ну вот теперь почти всё встало на свои места :) Теперь понятен и код ассемблера (запись в глобальную переменную, разделяемую между потоками), и диагностика чекера (ну да, состязание на эту самую запись в x), и откуда ложное разделение данных (false sharing - для глобальных переменных x и h), и что приватизация переменной решает проблему.
Что по-прежнему не совсем понятно, так это что помешало компилятору в начальном варианте поместить h в регистр. Подозрение на то, что запись в x может изменить h? Так вроде очевидно из адресов и размера переменных, что этого не происходит. Ну да ладно, это уже мелочи.
регистров хватает , но gcc категорически не хочет юзать xmm6, xmm7
что интерестно объявление const помогает вообще колосально:
например const double h = 1.0/n;
400718: f2 0f 2a c8 cvtsi2sd %eax,%xmm1
40071c: 83 c0 01 add $0x1,%eax
40071f: 66 0f 28 f3 movapd %xmm3,%xmm6
400723: 39 c2 cmp %eax,%edx
400725: f2 0f 5c cd subsd %xmm5,%xmm1
400729: f2 0f 59 4c 24 08 mulsd 0x8(%rsp),%xmm1
40072f: 66 0f 28 c1 movapd %xmm1,%xmm0
400733: f2 0f 59 c1 mulsd %xmm1,%xmm0
400737: f2 0f 58 c4 addsd %xmm4,%xmm0
40073b: f2 0f 5e f0 divsd %xmm0,%xmm6
40073f: f2 0f 58 d6 addsd %xmm6,%xmm2
400743: 7d d3 jge 400718
но больше 3х констант на цикл всеравно видеть не хочет
Using built-in specs.
Target: x86_64-redhat-linux
Configured with: ../configure --prefix=/usr --mandir=/usr/share/man --infodir=/usr/share/info --with-bugurl=http://bugzilla.redhat.com/bugzilla --enable-bootstrap --enable-shared --enable-threads=posix --enable-checking=release --with-system-zlib --enable-__cxa_atexit --disable-libunwind-exceptions --enable-languages=c,c++,objc,obj-c++,java,fortran,ada --enable-java-awt=gtk --disable-dssi --enable-plugin --with-java-home=/usr/lib/jvm/java-1.5.0-gcj-1.5.0.0/jre --enable-libgcj-multifile --enable-java-maintainer-mode --with-ecj-jar=/usr/share/java/eclipse-ecj.jar --disable-libjava-multilib --with-cpu=generic --build=x86_64-redhat-linux
Thread model: posix
gcc version 4.3.0 20080428 (Red Hat 4.3.0-8) (GCC)
Страницы