company_banner

«Задачка-то сошлась с ответом!»



    Я очень часто по работе слышу вопрос, задаваемый из, большей частью, академической среды, ввиду огромного количества выполняемых вычислений именно там: «Почему наши результаты разные от запуска к запуску одного и того же приложения? Мы же ничего не меняем в нем». Стоит отметить, разговор про это уже был, но лишь частично отвечающий на вопрос. Попробую рассказать про эту проблему ещё чуть-чуть.

    Почему мы можем видеть следующую картину на одной и той же системе:


    Или, скажем, следующую, на двух разных машинах:


    Основная причина, из-за которой всё это происходит, состоит в том, что закон ассоциативности, который в математике всегда работает, в вычислениях на машинах с ограниченной разрядной сеткой невыполним, то есть (a+b)+c не равно a+(b+c).
    По этой теме много интересного можно найти в статье What Every Computer Scientist Should Know About Floating-Point Arithmetic, автор David Goldberg. Я же приведу простой наглядный пример, когда в вычислениях встречаются очень маленькие и относительно большие числа:


    Вроде бы пока всё ожидаемо, если мы считаем как математики. Теперь предположим, что мы работаем с числами с одинарной точностью, и сложим маленькое число с большим:


    С точки зрения машинной арифметики, мы получили правильный результат, согласно стандарту IEEE 754. И следующий результат тоже верен:


    Именно отсюда и начинают вылезать все проблемы. Ведь компилятор вполне может оптимизировать наш код так, что порядок сложений в выражениях поменяется, да и ещё много чего другого. И здесь помогут опции –fp-model, про которые я уже подробно рассказывал.

    Вот только не всё зависит от компилятора. Помните, оригинальный вопрос то звучал так, что вроде бы ничего не менялось в коде, а значит его не «пересобирали», но результаты всё равно разные.

    Так вот, на сцену выходят другие факторы. Для последовательного приложения, равно как и для параллельного, одним из важных аспектов для получения воспроизводимых результатов является выравнивание данных. Кстати, на производительность оно так же имеет огромное влияние. И если данные у нас в коде никто не выравнивал, то в циклах, которые векторизованы, разное количество итераций будут распознаны как пролог, ядро цикла и эпилог, а, соответственно, и порядок вычислений будет изменён. Для того, чтобы предотвратить подобные проблемы, можно пойти двумя путями:

    1) Сохранить оптимизации, но при этом явно выравнивать данные, например, используя _mm_malloc() и _mm_free(). В Фортране есть директива !DIR$ ATTRIBUTES ALIGN:64 :: arrayname, а так же компиляторная опция -align array64byte
    2) Отказаться от векторизации редукций и ряда других оптимизаций, используя опцию –fp-model precise, что повлечет за собой потерю в производительности.

    Отлично! Допустим, что мы скомпилировали наш код с –fp-model precise, но результаты до сих не отличаются постоянством. Что же теперь не так? На самом деле, есть ещё ряд причин.

    Во-первых, если наше приложение параллельно, то это число потоков. Скажем, редукции в OpenMP никоим образом не затрагиваются опциями компилятора. Таким образом, изменяя число потоков, а на разных системах оно вполне может быть разным по умолчанию, мы неявно изменяем порядок сложений в редукциях (равно как и их количество). Больше того, никто и не говорит о том, что промежуточные суммы в OpenMP будут суммироваться в одинаковой последовательности от запуска к запуску, соответственно, даже при одинаковом количестве потоков мы можем получать разный порядок суммирования, а значит, и разные результаты. Кстати, это свойственно и icc, и gcc. При этом результаты вполне могут быть более точными, чем в последовательной версии редукции.

    Нам нужно контролировать это самим, и в компиляторе Intel, начиная с версии 13.0, есть переменная окружения, которая позволяет выполнять детерминированные редукции. Это KMP_DETERMINISTIC_REDUCTION. Выставляем её в yes, и можем быть спокойными, но только при условии, что мы используем статическое распределение работы между потоками (static scheduling), и число потоков неизменно.
    Если же мы используем Intel Threading Building Blocks, то и там мы можем контролировать редукции, через вызов функции parallel_deterministic_reduce(). Кстати, в Intel TBB нет зависимости от числа потоков, и результат может отличаться от последовательной версии редукции.

    Хорошо. Мы добавили и это, но опять не помогает.
    Тогда есть подозрение, что используются библиотеки, которые знать не знают про ваше стремление к «постоянности». Хочу сначала обратить внимание на вызовы математических функций. Опция –fp-model позволяет добиться консистентности вызовов функций, но для получения идентичных результатов на разных архитектурах (в том числе и не Интеловских), рекомендуется использовать ключ –fimf-arch-consistency.

    Более «мощный» пример — использование библиотеки Intel Math Kernel Library, которая по умолчанию не ориентирована на воспроизводимость результатов, но этого можно добиться, используя дополнительные настройки, и при соблюдении ряда ограничений.

    Например, при одинаковом количестве потоков, статическом планировщике (OMP_SCHEDULE=static, выставлен по умолчанию), одной и той же ОС и архитектуре, можно получать воспроизводимые результаты. Эта возможность называется в MKL Conditional Numerical Reproducibility. Я думаю понятно, почему она условная. Кстати, раньше так же было требование к выравниванию данных, но в последней версии инженеры чего-то «нашаманили», и об этом условии можно больше не заботиться при работе с MKL. Так вот при соблюдении этих требований, выставляя переменную окружения MKL_CBWR_BRANCH (например, в значение «COMPATIBLE»), или вызывая функцию mkl_cbwr_set(MKL_CBWR_COMPATIBLE) в коде, мы получим долгожданную «стабильность».

    Таким образом, работая с большим количеством вычислений мы должны понимать, что добиться «постоянных»результатов возможно, но лишь с соблюдением ряда условий, о которых мы говорили. Почти во всех инструментах сегодня существует способ контролировать этот вопрос, нужно быть внимательным и выставлять нужные настройки и опции. Ну и ещё – часто причина банальна, и кроется в обычных ошибках. Иногда интересно видеть код, вычисляющий серьёзные научные формулы, используя, скажем, неинициализированные переменные. А потом мы удивляемся, чего это очередная ракета полетела не туда. Будьте внимательны при написании кода и используйте правильные средства для его создания и отладки.

    Кстати, такой вопрос. Почему же в примере на системах с разными Xeon'ами (картинка наверху), мы получали разные результаты? Я про это явно так и не написал, было бы интересно услышать ваши идеи.
    Intel
    Company

    Comments 43

      0
      набор и реализация инструкций самого МП тоже может вносить свои коррективы. А так же различные помехи извне. Это довольно хорошо может быть заметно на МК.
        0
        В ксеонах вероятно есть аппаратные оптимизации float-point операций, которые реализованы по-разному в двух моделях. Или числа больше разрядной сетки в одном случае округляются до единиц, в другом до нулей
          +1
          Варианты интересные, но пока не то.
            0
            А формула расчетов есть? Может я и не внимателен, но я не нашел откуда результат взялся тот, что на скринах…
              0
              Нет, формулы нет. Есть просто факт — программка выдаёт такие результаты, что бы могло быть причиной. Кстати, ниже уже написано что.
              0
              На мой взгляд именно округление (точность) и сыграло ту роковую роль в расчетах. Причиной тому, что расчеты разбиваются на части и здесь много факторов, апаратные и программные алгоритмы оптимизации (компилятор, процессор, ОС).

              Допустим есть формула 1/3 + 1/3 + 1/3 = 1

              Но если считать так (1/3 + 1/3) + 1/3 != 1 или так 1/3 + (1/3 + 1/3) != 1

              (1/3 считаем в первую очередь, как в математике, не брал в дужки, что б не засорять формулу)

              Более сложную формулу может разбивать на части и одна из частей считаться быстрее.

              К примеру есть формула:
              1/3 + 2/3 + 3/3 + 4/3 + 5/3 = 5

              Допустим апаратный или программный алгоритм разбил задачу на части так:
              (1/3 + 2/3) и (3/3 + 4/3) (5/3 еще не затронуло)

              В первом варианте первая часть посчиталась быстрее и вышло так ((1/3 + 2/3) + 5/3) + (3/3 + 4/3)
              В другом варианте вторая часть посчиталась быстрее и вышло так ((3/3 + 4/3) + 5/3) + (1/3 + 2/3)

              Вот и вышло — оба варианта имеют разный ответ.
            0
            о, интересный момент, а почему не писать код так, чтобы избежать подобных казусов?
            наверняка, есть стандартные приемы подавляющие подобные ошибки?
              +1
              Есть, разумеется.
              Существует целый раздел вычислительной математики, называется «методы решения некорректно поставленных задач». Для вхождения нужно осилить средних размеров томик, но если на пальцах — сглаживаем, округляем и по возможности сортируем. Но самое главное — сами заранее устанавливаем предел машинной точности (эпсилон), скажем. не в 17 десятичных знаков double, а например, вполне достаточные почти во всех реальных приложениях 10 знаков. И тогда все эти проблемы уходят сами собой.

              Подобные вопросы ожидаемы от программистов, но как-то очень странно, что математики не знают математики.
                0
                «Алгеброики» не знают арифметики
                0
                Они есть, другой вопрос, почему, когда пишется код, ими пренебрегают. Ещё один — не все хорошо знают возможности средств разработки, в частности, компилятора, библиотек, и именно на это я делал фокус. Ну и, как следует из поста, не следует ожидать, что ваше приложение всегда будет выдавать одни и те же результаты на разных конфигурациях, но минимизировать требования можно.
                0
                Самое главное скажите, как воспроизвести проблему?
                  0
                  Смотря что является причиной проблемы. Вот, скажем, на картинках у меня показано, как. Просто запускаем на системах с разными процессорами. У разных приложений — разные причины. Скажем, невыравненные данные могут стать причиной редко воспроизводимой проблемы.
                  0
                  Спасибо Вам огромное за статью. Я как раз вчера столкнулся с такой проблемой, правда на куде.
                  Вы случайно не знаете аналогичные опции компилятора для куды?
                    0
                    Спасибо за интерес к данной теме. К сожалению, с компилятором для куды помочь не могу, но причины проблем везде одинаковы.
                    +6
                    А потом такие заголовки появляются…
                      0
                      Кстати, в Intel TBB нет зависимости от числа потоков, и результат может отличаться от последовательной версии редукции.

                      Мне кажется, что в этой фразе выпущено «не». Должно читаться:
                      Кстати, в Intel TBB нет зависимости от числа потоков, и результат не может отличаться от последовательной версии редукции.

                      Или я неправильно понял параграф про TBB?
                        0
                        Нет, всё правильно было написано. В TBB нет зависимости от числа потоков (partitioner делит работу одним и тем же способом всегда), но вот только результат последовательной версии и параллельной, возможно, будут не совпадать.
                          0
                          Тогда логичнее будет вместо «и» поставить противопоставление «но».
                        0
                        > Почему же в примере на системах с разными Xeon'ами, мы получали разные результаты?
                        > в циклах, которые векторизованы, разное количество итераций будут распознаны как пролог, ядро цикла и эпилог, а, соответственно, и порядок вычислений будет изменён

                        Разная разрядность SIMD векторов (SSE vs AVX)?
                          0
                          Или же FMA?
                            0
                            Поддержка FMA появилась начиная с Haswell'а, так что тут версия неверна.
                            0
                            Да, суть верна. Осталось только понять, почему один и тот же бинарник выполнился на одной машине с набором инструкций SSE, а на другой — с AVX, но я думаю ответ очевиден.
                              0
                              > но я думаю ответ очевиден.
                              Спешу вас удивить, что не для всех :).

                              > почему один и тот же бинарник выполнился на одной машине с набором инструкций SSE, а на другой — с AVX
                              Я так думаю, что если бы CPU без поддержки AVX увидел бы регистр по имени ymm (или же AVX-инструкцию), то выдал бы что то вроде Illegal instruction. Поэтому очевидно (мне), что приложение скомпилировано с опцией -msseX (например для gcc). А в AVX-регистр SSE-код попал, так как xmm является младшей (младшими 128 битами) частью ymm. Я думаю так… Но как это могло повлиять на результат мне непонятно. Ведь в обоих случаях разрядность регистра в котором происходит вычисление одинакова и равна 128 битам.
                                +1
                                Да, правда… при выполнении неподдерживаемой инструкции позитива было бы мало, но её не будет.
                                Причина в том, что приложение было скомпилировано с опцией -axAVX, которая говорит компилятору о том, чтобы создавать не только специфичные инструкции для AVX, но и дефолтную (baseline) ветку (code path) для SSE2. Вообще, их можно перечислять в этой опции, например -axSSE4.2,AVX, причем дефолтная ветка создается всегда (можно контролировать минимальный набор поддерживаемых инструкций). Тем самым гарантируется запуск приложения практически на любом процессоре.
                                При запуске мы знаем cpuid и выполняется нужный code path.
                                Потому и результаты разные — в одном случае выполнялся код, сгенерированный c SSE, в другом — с AVX. Размеры векторов разные, а значит и суммирование поменялось в редукциях. Вот как-то так.
                                  0
                                  Немаловажное дополнение — с опцией -ax, доступной в компиляторе Intel.
                                    0
                                    Большое спасибо за развернутый ответ!

                                    > -ax, доступной в компиляторе Intel
                                    А для gcc, я так понимаю для таких целей прийдется вручную писать разные реализации и использовать Function Multiversioning. Или в gcc тоже есть такая волшебная опция компиляции?
                            0
                            А что, кто-то считает в float-ах/double-ах в случае, если нужно сохранять точность? Есть же специально предназначенные для этого типы данных.
                              +1
                              О, я вас уверяю, считают и ещё как. Кстати, что бы вы предложили для Fortran'а, какие типы?
                                0
                                Ошибся веткой… Вопрос адресован eugenius_nsk
                                0
                                > Есть же специально предназначенные для этого типы данных.
                                Подскажите пожалуйста, что это за типы? Речь идет о целых типах и вычислениях с фиксированной точкой?
                                +1
                                У нас такая проблема возникла при переносе программы на CUDA (там вообще был винегрет из OpenMP, MPI и CUDA). Удалось добиться побитового совпадения результатов сравнительно малой кровью — ключевым решением была реализация подмножества стандартной библиотеки. В результате скорость вычисления упала на полтора порядка, но этого было достаточно, чтобы пройти тесты и убедиться, что мы не накосячили с Кудой.
                                  0
                                  Забавно, когда щупал OpenCL, в котором негодяи из Кроноса вырезали сишный тип complex, то оказалось, если такой тип навелосипедить через структуру и определить все основные операции самому, то результат будет отличаться от сишного. Очевидно это было потому, что в нормальном complex.h перед действиями double переводится в long double для сохранения точности. А т.к. в OpenCL такого типа нет, то пришлось бы велосипедить и его.
                                    0
                                    Очень интересно. Сам сейчас бюсь с такой проблемой на CUDA, слава Богу пока без MPI и OpenMP.
                                    До сих пор не удалось добиться стабильности.
                                    Узнать, каким образом Вы решили эту задачу очень интересно, но я не понял, что Вы имеете ввиду под реализацией подмножества стандартной библиотеки.
                                    Вы что, функции из math.h вручную для куды переписали что ли?
                                      +1
                                      Да, примерно так и есть. Была взята сановская библа (сейчас попытался было ее найти, но из нее тащат код все кому не лень, поэтому находятся производные проекты) и обернута в класс, реализующий операции над double (класс нужен, чтобы пользовательский исходный код не переписывать).

                                      Собственно, не обязательно всю библиотеку переносить, достаточно лишь те функции (sqrt, log, ...), которые используются.

                                      Еще был какой-то нюанс, связанный с тем, что в Куде есть операция типа A=B*C+D (забыл как называется), не имеющая эквивалента на CPU, и компилятор может довольно агрессивно менять порядок вычислений, чтобы вместо отдельных умножений и сложений использовать эту операцию. Поэтому нужно заменить обычную арифметику на специальные CUDA-вызовы вроде __dadd_rn().
                                        0
                                        > в Куде есть операция типа A=B*C+D
                                        > не имеющая эквивалента на CPU

                                        FMA4?
                                          0
                                          Да, оно. Только, во-первых, FMA на разных ядрах может отличаться реализацией (или вообще отсутствовать, как в нашем случае), во-вторых, его использование сильно зависит от причуд компилятора — даже не знаю, можно ли в принципе заставить nvcc генерировать код с FMA для CPU.
                                          0
                                          Спасибо!
                                      0
                                      Может ли instruction reordering на современных процессорах влиять на наблюдаемый результат fp-вычислений? Борются ли с этим как-нибудь на суперскалярных архитектурах (например, введением зависимостей по данным)?
                                        +3
                                        В вычислительной математике (вернее, в области физического моделирования, за другие не скажу — не знаю) есть две точности — «машинная» точность (то количество знаков после запятой, которое в выполненном вычислении гарантированно соответствует истине) и «физическая точность». Цель программиста-моделиста, очевидно, в том, чтобы вторая нигде не превзошла первую.

                                        Вот эта, вторая точность, является результатом хитрых оценок рассматриваемой модели и как правило составляет 3-4 знака, не более. Бывают, конечно, и другие случаи, но в физике всё равно экспериментальный результат ценнее аналитического, а для того, чтобы поставить эксперимент (то есть «знать, куда целиться»), обычно хватает и такой точности.

                                        Поэтому для любой (или почти любой) практической цели просто хватает числа double.

                                        А если у вас в формуле возникает ситуация, когда вы складываете 10^(-20) с единицей, значит у вас худая модель. И очень хорошо, что программа в этой ситуации явно врёт — меньше шансов, что вы дурацкий результат примете за правду.
                                          0
                                          А по моему это отличные новости для крипто-алгоритмов.
                                          Возможно уже кто то реализовал?
                                          Дополнительная случайность не повторяющаяся на других процессорах, или даже в последующем запуске на том же процессоре никому не повредит.
                                          Надо собирать особенности и делать на них библиотеку «Дополнительный генератор случайности» (копирайт за имя мой :) )
                                            0
                                            Любой генератор случайных чисел должен быть «предсказуемо случайным». Иными словами, он должен основываться на величинах, функция распределения которых известна (желательно, равномерна). В противном случае ценность его низка. Это интересно описано в книге «Искусство Программирования», например. Там автор рассказывал, как пытался придумать генератор случайности сам, пользуясь подобным подходом, и что из этого вышло.

                                            А еще там написано, что большая часть алгоритмов генерации случайных чисел (причем хороших) на сегодня основана на формуле типа: x[k+1] = (x[k] * a + b) % c, где a, b и c специальным образом подобраны. Как ни странно, но самое случайное, что удалось найти математикам на сегодня — взятие остатка от деления.

                                            Тот же алгоритм (вернее просто остаток от деления, без умножения) используется при вычислении контрольной суммы CRC
                                            0
                                            Эээ я же не предлагаю заменть то что есть, дополнить.
                                            Все что я читал про генерацию случайных чисел и что осталось в голове, это то что без физических приборов, чисто логикой сложно генерить случайные числа. Тут сложно сказать что совсем без физики (некая привязка к процессорному железу есть) но вполне интересное дополнение, например тот же остаток от деления + еше такую псевдослучайную величину. На тему гарантий так их вообше нет, на то и есть основные сертифицированные алгоритмы. Я бы рассматривал эту особенность как виагру, искали одно помогло в другом, не излечивает, но текушую задачу помогает выполнить качественнее…
                                              0
                                              > Эээ я же не предлагаю заменть то что есть, дополнить.
                                              Т.е. использовать CPU в качестве источника энтропии?

                                            Only users with full accounts can post comments. Log in, please.