Comments 58
Небольшое замечание по моему определению SSA формы (спасибо читателям, которые это заметили): на самом деле, оно, возможно, слишком уж строгое. Существуют реализации IR на основе SSA, в которых CFG в явном виде не существует, по крайней мере до какого-то момента. Одним из них является IR в Graal VM. Там всё равно есть финоды, и на фундаментальном уровне всё очень похоже, но мне проще объяснять основные принципы на примере LLVM IR, который имеет явный CFG.
В дальнейшем мы будем предполагать наличие явного CFG, а про более сложные вариации, где зависимости по управлению представляются иначе, можете поискать информацию самостоятельно.
Вторая часть про доминирование тут.
Спасибо за статью! Давно хотел разобраться в SSA хотя бы немного, но всё руки не доходили. Понял (надеюсь, что правильно) наконец как работают Ф-функции. На самом деле, оказалось даже проще, чем мне казалось. В этом примере, где Ф-ноды зависят от других Ф-нод, сразу понял что дело в том, что значения из предыдущего блока берутся и нет никакой разницы в порядке вычислений внутри блока.
Спасибо за статью! Понятно о сложных вещах - это верный путь.
Спасибо за статью, появился ряд вопросов, если не трудно, то интересно прочесть ответы, в связи с тем что SSA это как-бы "RISC с бесконечным числом регистров" и с 3-адресной арифметикой:
Умеют ли бэк-энды оптимизорать код из трех-адресной, скажем в одно или двух адресную или даже стековую арифметику (Java Jit в частности - стековый)?
Как могут быть реализованы виртуальные таблицы методов, на отличной от типовых ветвлений архитектурах, например "табличной" (внутренняя регистровая таблица переходов)
Реализация switch конструкций, где как-бы "данные" (регистры) становятся точками переходов, как?
Возможно вопросы нубские, прошу прощения. Компиляторами занимался уж сильно давно (до 1990-го)..
Умеют ли бэк-энды оптимизорать код из трех-адресной, скажем в одно или двух адресную или даже стековую арифметику (Java Jit в частности - стековый)?
Ну, строго говоря, LLVM-инструкции не всегда трёхадресные, там иногда вообще нет def'а, или же аргумент всего один. Но если коротко - да, умеют, это задача instruction selector'а и register allocator'а. Я не шибко большой эксперт конкретно в этой теме, но если есть такой запрос, то чего-нибудь написать смогу.
Java JIT в Фалконе переводит джава-байткод в LLVM IR и оптимизирует уже его, таким образом стековость там исчезает уже на этапе парсинга. Потом, на этапе кодогенерации, это превращается в обычный asm (например, x86). С2, насколько я знаю, занимается примерно тем же самым, только там Sea of nodes вместо LLVM IR.
Как могут быть реализованы виртуальные таблицы методов, на отличной от типовых ветвлений архитектурах, например "табличной" (внутренняя регистровая таблица переходов)
Я не уверен, что понял вопрос. Если речь о том, как моделируется виртуальный call, то для этого есть специальная инструкция, которая просто зовёт другую фукнцию. Она стоит где-нибудь в середине блока и не рассматривается как часть CFG данной функции. У неё просто есть сайд-эффект, говорящий о том, что вот из этой инструкции мы, например, можем и не вернуться обратно (а уйти в бесконечный цикл или вообще жёстко прервать программу). Реализована же она с помощью обычной таблицы виртуальных вызовов.
Реализация switch конструкций, где как-бы "данные" (регистры) становятся точками переходов, как?
switch аналогичен инструкции br, просто у br аргумент типа i1 и всего два константных варианта перехода (по true и false). Switch в LLVM разрешает аргумент любого целочисленного типа и большее количество переходов (но все они также идут по константным значениям). Подробно можно почитать здесь:https://llvm.org/docs/LangRef.html#switch-instruction
Только что смотрел и анализировал ассемблерный листинг большой программы на Си, скомпилированный свежим gcc с опцией -O3. Сначала было трудно, почти невозможно, понять в нём логику управления потоком исполнения команд. Потом попробовал -O2, -O1, -O0 и стал немного разбираться, как всё это устроено - стало полегче... А сейчас, после статьи, многое стало вообще понятно. Приятно чувствовать себя более разбирающимся в предмете своего интереса. Хотя я ассемблерный текст смотрел просто из любопытства - хотел посмотреть как ассемблерная вставка в сишной inline-функции вставилась...
Спасибо за статью, на моём не особенно глубоком уровне понимания оптимизации - зашло на ура!
Может быть на основе вашего опыта, вы можете дать какие-то общие рекомендации по стилю (архитектуре?) написания программ на Си?
Спасибо за отзыв, я рад, если моя статья правда помогает людям. :)
Что касается стиля написания программ на тех или иных языках - думаю, на этот счёт написано немало книг, и всё об одном и том же. Можете почитать "Совершенный код" и иже с ним, плюс книгу банды четырёх про паттерны, ну и следуйте тем конвенциям, которые есть в вашем проекте. Если нет никаких - то заведите какие-нибудь (лучше всего, опять же, слизать с какого-нибудь качественного опенсорсного проекта :) ).
Спасибо за рекомендации, но я уже пользуюсь теми, что предназначены для лучшего понимания человеками. Я имел в виду стиль, способствующий глубокой автоматической оптимизации. Или таковой не существует?
Более конкретно - что-то из области:
"не создавайте много промежуточных переменных, а лучше лепите всё в одно выражение или делайте оператор в операторе (по возможности)"; или "не ходите много по меткам, а лучше избавьтесь от них", "вместо констант используйте #define, а лучше перечисления" (ну или наоборот). Чтобы это способствовало оптимизации, так сказать.
А, понял. Пару советов дать могу, но учтите, что компиляторы разные и ведут себя по-разному. То, о чём я говорю, основано на опыте с компиляторами на основе LLVM:
Чего НЕ обязательно или не нужно делать:
Особо возиться с константами. Это очень базовый класс оптимизаций, и компиляторы это давно умеют, по крайней мере в С++. enum не enum - тоже по большому счёту разницы нет, это всё равно константы. Единственное - могут быть сложности, если константа передаётся через какой-нибудь extern в другой модуль, но это какие-то очень вычурнутые случаи. Пишите как удобно.
Лепить выражения в одну строку вместо разбиения на разные времянки согласно логике. Компилятору всё равно, на уровне IR программа будет, скорее всего, одинаковая. Вы всегда можете в этом убедиться при помощи чудо-сайта https://godbolt.org/. Так что пишите так, как лучше читается человеку.
Применять "хакерские секретики" типа замены "i++" на "++i" в for-циклах или замены умножения/делений на 2 на сдвиги. Наверное, лет 20 назад были компиляторы, которые этого не умели, но современный clang, да и другие современные компиляторы, это делают сами. Если это не помогает читаемости кода, заниматься этим бессмысленно.
Чего НУЖНО делать:
Избегать goto. Особенно в/из циклов. goto может привести к появлению т.н. non-loop cycles в CFG, которые являются циклами (с точки зрения теории графов), но не распознаются как таковые компилятором. На них не применяются цикловые оптимизации.
В сигнатурах функций для параметров-ссылок и указателей писать const везде, где возможно. В случае, если функция не заинлайнится, отсутствие const не позволит компилятору доказать, что никто по этой ссылке ничего не пишет. Это может поломать оптимизации.
Циклы с положительным шагом в целом оптимизируются лучше, чем с отрицательным. Это связано с некоторыми особенностями реализации анализов. В теории разницы нет, но на практике она ещё как есть. Я довольно много времени потратил, чтобы поддержать отрицательные шаги в разных цикловых оптимизациях, но эта работа далека от завершения. Поэтому, если всё равно как бежать по циклу, бегите с положительным шагом.
Всякую сложную дебажную логику (особенно с принтами) лучше уводить под #ifdef, чтобы в релизе её не было. Наличие сложного и развесистого кода с сайд-эффектами может поломать оптимизации много где.
Если сомневаетесь, может ли ваш компилятор сделать то или это, можете проверить в https://godbolt.org/ и зафайлить соотв. issue. Иногда на них кто-то смотрит (но это не точно). :)
- В теории разница тоже есть, только она связана с кешированием памяти и предсказанием переходов в процессоре.
Судя по тому, что видел я, префетчер на х86 одинаково хорошо работает как с итерированием как вперёд, так и назад. Так что проблем с кэшами быть не должно. Что касается предсказателя переходов, то он всегда считает, что backedge - likely taken, независимо от условия цикла. Так что разницы именно по этим причинам я бы не ожидал.
По направлению чтения/записи в память сам лично не проверял, спорить не буду, просто отложилось в памяти из какого-то руководства по оптимизации. По предсказанию переходов — это к тому, что даже один goto вперёд внутри цикла (например по if) заметно влияет на производительность.
Проблема любых книжек в том, что они устаревают уже на момент публикации. :)
А это касается только книжек или статей тоже?)
Ну тут всё зависит от того, насколько длинный процесс издания. У книг это год или годы, статья с шансами опубликуется за несколько месяцев (если это научный журнал) или сразу (если это сайт типа хабра), и с шансами ещё не устареет. :) К тому же, академики довольно часто пишут о вещах, которые ещё никто даже не начал реализовывать.
Вот это здОрово! Как раз то, что доктор прописал нужно!
Небольшое замечание по моему определению SSA формы (спасибо читателям, которые это заметили): на самом деле, оно, возможно, слишком уж строгое. Существуют реализации IR на основе SSA, в которых CFG в явном виде не существует, по крайней мере до какого-то момента. Одним из них является IR в Graal VM. Там всё равно есть финоды, и на фундаментальном уровне всё очень похоже, но мне проще объяснять основные принципы на примере LLVM IR, который имеет явный CFG.
В дальнейшем мы будем предполагать наличие явного CFG, а про более сложные вариации, где зависимости по управлению представляются иначе, можете поискать информацию самостоятельно.
Напишите, пожалуйста побольше о том, что дает IR для оптимизаций потом.
Огромное спасибо за статью! С нетерпением жду продолжения серии!
Очень интересно! Давно интересуюсь компиляторами, в студенчестве даже написал парочку простеньких, но до оптимизаций там дело не дошло. С нетерпением жду продолжения. Спасибо.
Вопрос-оффтопик про ошибку оптимизации
К теме статьи не относится, но вдруг Вы сможете подсказать? У меня есть программа на фортране, которую я собираю компилятором Интел. При выключенной оптимизации она работает корректно, а при включении ключей /MP /O2 /QaxSSE3 /QxSSE3 /Qparallel совершенно обычное присвоение с преобразованием типов
R8=I8
(64-битное целое преобразуется в 64-битное число с плавающей точкой) ИНОГДА (очень редко! - один раз за десятки миллионов вызовов) я получаю R8=Nan в случайных местах. Пришлось даже сделать подпорку - вместо простого присвоения я теперь вызываю функцию, которая проверяет результат и в случае R8=Nan начинает добавлять-вычитать единичку из I8, пока не найдет такое значение, которое конвертируется корректно. Увы, но эту подпорку я могу вставить только в собственный код, а при вызове встроенных функций баг все равно происходит...
И второй интересный нюанс: баг возникает, только если где-то в моей программе есть запрос системного времени, например:
call getdat(i4Year, i4Month, i4Day)
При этом не важно - какая именно функция используется для запроса времени/даты (в фортране их несколько). Если программа вообще не интересуется текущим временем/датой, то баг не возникает.
Ну и третий нюанс = баг возникает только на некоторых процессорах. На других программа с оптимизацией работает корректно. (Список тех и других лежит в файле FORTRAN_URAND_BUG.doc).
Впервые я обнаружил этот баг (Nan появлялся в генераторе случайных чисел) чуть больше года назад. Тогда же я задал этот вопрос на Хабре, и получил целую кучу советов от разных людей. Но раскрутить эту проблему до конца (найти причину бага и как его устранить) мы так и не смогли. Кроме того, объем собранной по багу информации намного превысил те 10000 символов, которые отведены для вопроса. Поэтому подробности пришлось вынести в отдельный файл FORTRAN_URAND_BUG.doc.
P.S. Сборка программы:
Ключи компиляции
<Tool Name="VFFortranCompilerTool" SuppressStartupBanner="true" DebugInformationFormat="debugEnabled" MultiProcessorCompilation="true" GenAlternateCodePaths="codeForSSE3" UseProcessorExtensions="codeExclusivelySSE3" Parallelization="true" BufferedIO="true" EnableEnhancedInstructionSet="codeArchSSE3" FixedFormLineLength="fixedLength132" ErrorLimit="60" DebugParameter="debugParameterAll" WarnDeclarations="true" WarnUnusedVariables="true" WarnIgnoreLOC="true" WarnTruncateSource="true" WarnInterfaces="true" ByteRECL="true" InitLocalVarToNAN="true" LocalSavedScalarsZero="true" ModulePath="$(INTDIR)/" ObjectFile="$(INTDIR)/" PdbFile="$(OutDir)\vc90.pdb" Traceback="true" BoundsCheck="true" UninitializedVariablesCheck="true" RuntimeLibrary="rtMultiThreadedDebug" GlobalOptimizations="false"/>
<Tool Name="VFLinkerTool"
OutputFile="$(OUTDIR)/my.exe" LinkIncremental="linkIncrementalNo" SuppressStartupBanner="true" IgnoreDefaultLibraryNames="LIBCMT.LIB" GenerateDebugInformation="true" ProgramDatabaseFile="$(OUTDIR)/Lib_Test.pdb" SubSystem="subSystemConsole" />
Настройки Floating-Point
P.P.S. В любом случае, спасибо за статью! Я про компиляторы почти ничего не знаю (общаюсь с ними только как пользователь). Но для расширения кругозора всегда стараюсь читать подобные материалы, даже если понимаю не все. В этот раз понял довольно много - было очень интересно (и полезно ;-)
Был у меня похожий случай на Борландовском компиляторе c++. программа падала в странных местах с ошибками (то ли NAN, то ли Division by zero, точно не помню). Но помню, что связано было с плавающей точкой и падала она в том месте, где таких ошибок просто не могло произойти. Дебаг не помогал. Добавление отладочной выдачи внезапно "вылечивало" баг и он пропадал.
В итоге, как я понял, дело было в том, что компилятор, видимо в силу оптимизаций, не всегда очищал регистр флагов процессора. И получалось так, что какая-то операция пыталась работать с плавающей точкой (в моём случае это было чтение из неинициализированной переменной типа double). Иногда там оказывался мусор и процессор взводил флаг некорректности числа. Но т.к. Операция простого чтения значения не предполагала, что что-то могло пойти не так (не совершались никакие арифметические операции), то далее компилятор не вставлял инструкции проверки этого флага и не сбрасывал его. А вот потом следовала безобидная операция (например, присвоение нормального значения в тип double). После этой операции, компилятор по всем своим правилам лез в регистр флагов, чтобы убедиться, что всё прошло хорошо и внезапно обнаруживал там флаг некорректной работы с Double! И на этом месте радостно вываливался в исключение. И естественно при этом показывал пальцем на безобидную инструкцию.
Добавление отладочной выдачи внезапно "вылечивало" баг и он пропадал.
Да, у меня очень похоже.
Если вставить оператор print сразу после оператора присвоения Real_8=Integer_8, то Nan-ы не появляются. И некоторые другие операторы, использующие значение Real_8, тоже купируют баг, если их поставить сразу же после присвоения. При этом он может в одном месте исчезнуть, а в другом появиться...
Но если бы это были просто флаги, то откуда связь с запросом системного времени? Если время не запрашивать, то бага не будет...
А еще у меня проблема осложняется тем, что обычно у меня все ряды данных с пропусками, и Nan - это совершенно легальное значение при расчетах. Поэтому обнаружить их "самозарождение" я могу только по косвенным признакам уже после окончания вычислений. Например, если в исходном сигнале Nan-ов не было, а после прибавления к каждому значению ряда некоторой константы они вдруг появились. Поэтому заметить этот баг мне было непросто. Сперва я это явление обнаружил при вызове ГСЧ (когда Nan-ов в принципе быть не должно). Потом обнаружил то же самое при вычислении БПФ (для такого расчета я все пропуски в сигнале заранее заполняю). А уже потом начал целенаправленно тестить, и тут оказалось, что и некоторые другие операции небезгрешны...
Могу предположить, что где-то в недрах каких-то функций взводится какой-то флаг (я не большой специалист в низком уровне, просто знаю, что такое есть), но не анализируется и не сбрасывается. В нормальном состоянии компилятор предполагает, что этот флаг обычно сброшен и когда включены оптимизации, выкидывает лишние команды сброса этого флага (которые могли быть натыканы просто на всякий случай перед каждой работой с плавающей точкой).
В вашем случае, я бы постарался локализовать проблему и сделать максимально короткий тестовый пример, который бы повторял баг. Ну или если в вашем языке программирования можно делать ассемблерные вставки, то можно попытаться напрямую читать состояние регистра флагов (или сбрасывать его, чтобы понять его влияние).
Спасибо за совет, идея понятна. Но я пока не понимаю, как прочитать регистр флагов. Готовой функции для этого в моем компиляторе вроде бы нет (я не нашел). А asm-вставки я последний раз делал в DOS-фортране... умеет ли это компилятор от Интел, я пока не знаю. Быстрым поиском в справке я ничего не нашел. Возможно, что не умеет... а скомпилировать нужную функцию отдельно и прицепить на этапе линковки мне пока сложновато.
Nan — это совершенно легальное значение при расчетахИспользовать NaN в качестве маркера отсутствия данных — плохая идея. Потому что:
1) его битовое представление разное для разных типов float,
2) процессор работает с ним очень медленно, прям очень,
3) проверку на NaN часто делают через сравнение переменной с самой собой, а оптимизирующий компилятор может счесть такое сравнение избыточным и выкинуть его.
Использовать NaN в качестве маркера отсутствия данных — плохая идея
А что, есть что-то лучше?
Изначально у меня вместо Nan была предопределенная собственная константа. При этом после каждого вычисления приходилось проверять - не совпадает ли результат с этой константой, и корректировать его, если вдруг совпадает. В конечном счете я пришел к выводу, что удобнее (и быстрее) использовать стандартную константу IEEE_QUITE_NAN, которая самопроизвольно появляться вроде бы не должна. Скорость работы программы в результате не изменилась (мои тесты даже небольшое ускорение показали), а логика стала проще: значение проверяется на Nan не при вычислении, а только перед использованием. Это уменьшает число проверок и алгоритм легче читается. А перед использованием мне значение в любом случае надо на Nan проверять, т.к. пропуск - это вполне легальные данные, и их наличие - это не аварийная, а вполне штатная ситуация. Например, при подсчете среднего по N измерениям я суммирую только не-пропуски, а потом либо делю сумму на количество не-пропусков M, и возвращаю число, либо возвращаю Nan, если процент пропусков (N-M)/N превышает некоторый известный предел P (где P является одним из обязательных параметров функции). И аналогично при любых других вычислениях.
Ну и проверку на NaN я делаю не вручную, а вызовом функции isNun(). Так что компилятор ее заведомо не выкидывает ;-) Это же фортран, он на такие ситуации специально заточен ;-)
Да, есть хорошие практики программирования, а есть плохие. Полагаться на потенциальный UB — совершенно точно плохая практика. Ну как бы сам повод нашего общения совершенно непрозрачно на это намекает. И заточенность Фортрана на вычисления никак в борьбе с багом вам не помогла.
Я сам тоже иногда NaN использую, но только когда могу полностью контролировать формат и процесс. QUITE в IEEE_QUITE_NAN означает, что исключения выкидываться не будут. Но для FPU QUITE обеспечивается флагами на глобальном уровне (о чём я уже говорил), а вовсе не самой константой.
Хорошая практика — это сделать новый тип из double и bool, в котором можно заодно и выход за пределы контролировать, и точность вычислений, и точность преобразований, и поведение при нештатных ситуациях явно прописать. Если у вас временные ряды — туда же можно и штамп времени засунуть, и вообще не факт, что в этом случае имеет смысл фиксировать пропущенное значение.
Хорошая практика — это сделать новый тип из double и bool, в котором можно заодно и выход за пределы контролировать, и точность вычислений, и точность преобразований, и поведение при нештатных ситуациях явно прописать. Если у вас временные ряды — туда же можно и штамп времени засунуть, и вообще не факт, что в этом случае имеет смысл фиксировать пропущенное значение.
В идеальном сферическом мире такая практика, наверное, неплоха. Только вот какая производительность будет у этого типа? Если что, у меня гигабайтные ряды - обычное дело. При том, что АЦП обычно 24- или даже 16-битные, поэтому данные при обработке хранятся, как Real*4 (а в базе вообще иногда ужимаются в Real*2). Надстройка из bool либо увеличит потребный размер памяти вдвое, либо катастрофически замедлит скорость из-за проблем с выравниванием. Я уж не говорю про штамп даты, который в моем случае потребует минимум 8 байт на каждое 4-байтное значение (пакет обрабатывает сигналы с частотой дискретизации от 1МГц до нескольких лет, и с длиной до миллионов лет. Хотя конечно, сигналы с высокой частотой всегда намного короче ;-).
А точность вычислений я и так контролирую - она у меня в конфиге прописывается и затем везде проверяется, чтобы случайно деления на ноль или переполнения не случилось - там, где это потенциально возможно. К примеру, в конфиге есть параметр Zеro, который говорит, что все числа, меньшие +-Zеro, надо считать нулем. Поэтому деление на такое число у меня вернет Nun, и т.д. Сейчас эти "ручные" проверки, наверно, можно было бы как-то оптимизировать, но у меня они были сделаны еще в DOS фортране, и с тех пор рефакторятся по-минимуму (работает - не трогай). Хотя всякие функции типа nearest или huge я, конечно, тоже цепляю, по мере того, как они в языке (компиляторе) появляются...
В общем, в моем случае кодировать пропущенные наблюдения особым числом - оптимально. Ряды данных в этом случае выглядят просто, как наборы значений, а вся остальная информация хранится в паспорте ряда. Помимо прочего, это позволяет при совместной обработке пересчитывать (синхронизировать) ряды с разными шкалами времени (совершенно типичная у нас ситуация) простым циклом суммирования/интерполяции. Что гораздо эффективнее, чем, например, в матлабе, где для этого как раз используются календарные вектора, и каждой операции с несинхронными временными рядами предшествует не всегда тривиальная операция с этими векторами. Впрочем, я матлабом интересовался последний раз несколько лет назад - может, с тех пор там что-то более интересное прикрутили...
Так если у вас выход с АЦП 16/24-битный — то зачем преобразовывать его в 32/64-битный float? Откуда берутся пропущенные значения, если выход с АЦП всегда имеет конкретное целочисленное значение? Если данные пишутся через равные промежутки времени — то конечно, привязывать индивидуальный штамп времени к каждому значению смысла нет. Если нет — то вовсе не обязательно писать его вместе с датой, дата же не будет меняться через каждую секунду (у меня, например, это количество миллисекунд от начала месяца или года в UINT32).
Так если у вас выход с АЦП 16/24-битный — то зачем преобразовывать его в 32/64-битный float?
Тут все просто. Выход с АЦП идет в условных единицах, а работать-то хочется в градусах, миливольтах и пр. Это только кажется ерундой, а на самом деле эрогономически очень важно.Особенно если это разведочный анализ, и в обработке всегда участвует человек, а не просто какие-то формальные алгоритмы. Для этого данные надо пересчитывать по различным формулам, иногда не очень простым. А еще мне при обработке всегда нужен однородный сигнал, вне зависимости от того, что там где-то какой-то шунт поменяли. Это все делается на этапе первичной обработки еще до загрузки в базу. Поэтому в базе по-любому будут real-значения. Только, конечно, не 64-битные (при 24-битном АЦП это бессмысленно), а 32 или 16 бит. (Для краткости я их называю real*4 и real*2).
Откуда берутся пропущенные значения, если выход с АЦП всегда имеет конкретное целочисленное значение?
Во-вторых, АЦП работает не всегда. Сплошь и рядом бывают технические перерывы и сбои, а также всякое техобслуживание. Ну и просто что-то может выйти из строя. Или прибор зашкаливает, так как в сигнале есть сильный дрейф (выбрали слишком высокую чувствительность при настройке), а вернуть его в рабочий диапазон некому. Это же не лаборатория, где все под контролем, а какая-то полевая точка, которая обслуживается в лучшем случае раз в неделю, а в худшем - раз в несколько месяцев.
Во-вторых, в экспериментальных данных достаточно регулярно появляется всякий брак. Иногда это помеха, иногда аппаратурный сбой. Я твердо придерживаюсь той точки зрения, что подобные "аномалии" надо заменять Nan-ами, а не какими-то искусственными (вымышленными) значениями. А затем обрабатывать с учетом этого факта.
дата же не будет меняться через каждую секунду (у меня, например, это количество миллисекунд от начала месяца или года в UINT32
Да, это возможный вариант... только вот в году 3*Е+7 сек, или 3Е+10мс. Поэтому даже для мс-ряда придется дополнительно хранить номер года. Который в моем случае может достигать очень больших величин (палеоклимат и др). А еще и шаг опроса может быть задан не в мс, а в мкс. В общем, 4 байт никак не хватает, если хочется сделать универсальную хронологию. Я вот недавно столкнулся с угрозой переполнения 8-битного счетчика при синхронизации скважности (там сначала умножение, а потом деление - пришлось специальные проверки прикручивать и в сложных случаях в два этапа считать).
Тут все просто. Выход с АЦП идет в условных единицах, а работать-то хочется в градусах, миливольтах и пр.Так не обязательно работать с данными именно в том виде, как они хранятся у вас на магнитной ленте или перфокарте. Хорошей практикой считается не смешивать логики. Хранение данных — это одно, обработка (например удаление) шума — другое, анализ (статистический/Фурье/etc) третье.
У вас же должно прилагаться к этим данным хоть какое-то описание. Туда можно добавить масштабирующий фактор и единицы измерения, которые будут учитываться при конвертации.
Во-вторых, АЦП работает не всегдаНо ведь можно же просто не писать данные, когда их нет. А когда появляются, продолжать писать в другом файлике.
Я твердо придерживаюсь той точки зрения, что подобные «аномалии» надо заменять Nan-ами, а не какими-то искусственными (вымышленными) значениями.А можно их вообще не заменять. Шум, выходы за пределы и прочее — это тоже информация, которую можно измерять и анализировать.
Да, это возможный вариант… только вот в году 3*Е+7 сек, или 3Е+10мс.Для этого давно придумали структуру под названием «chunk», когда данные хранятся не одним непрерывным куском, а разбиваются на одинаковые небольших размеров. Перемотка в видео- и аудио-файлах достижима именно таким подходом.
А затем обрабатывать с учетом этого факта.По вашей ссылке чужих не пускают.
Так не обязательно работать с данными именно в том виде, как они хранятся у вас на магнитной ленте или перфокарте. Хорошей практикой считается не смешивать логики. Хранение данных — это одно, обработка (например удаление) шума — другое, анализ (статистический/Фурье/etc) третье.
Разумеется, хранение и обработка - это разные вещи. Но отсюда вовсе не следует, что хранить данные надо строго в том виде, как они получены с АЦП. Многие операции по первичной обработке (включая, например, пересчет из единиц АЦП в килограммы или выбраковка аппаратурных дефектов) должны выполняться ВСЕГДА, при любой обработке. Поэтому нет абсолютно никакого смысла откладывать такие процедуры на какой-то более поздний этап анализа. Их надо делать один раз при получении данных и их записи в базу. Во всяком случае, для наших наблюдений это именно так. Да и в других похожих случаях это
скорее система, чем исключение
Именно из-за того, что некоторые операции надо перед обработкой выполнять всегда, при режимных геофизических наблюдениях первичные данные уровня "0" (выход АЦП) чаще всего хранятся лишь "для архива". Иногда бывает до 3-4 уровней первичной обработки. Для примера могу привести службу гравиметрического мониторинга на самых современных сверхпроводящих гравиметрах (см. www.gfz-potsdam.de). Она объединяет десятки станций по всему миру, там наблюдениями занимаются сотни людей, первичной обработкой - тоже десятки (сотни?) людей, а потом все это выкладывается на один общий сайт. Так вот, хотя там на этом сайте выложены данные разных уровней, но на практике пользователям обычно нужен только верхний уровень либо предпоследний. А на некоторых станциях некоторые промежуточные уровни вообще не выкладываются.
У вас же должно прилагаться к этим данным хоть какое-то описание.
Разумеется. Именно для этого и нужна база данных - чтобы каждый сигнал в ней был документирован и паспортизирован так, как удобно пользователю сигнала. Именно поэтому в эту базу заносятся не первичные данные, а уже обработанные и поправленные. Благодаря чему описание сигнала в БД - это не многометровые таблицы с названиями каналов и коэффициентов, включая температурные поправки датчика и т.д., а что-то вроде "Относительная влажность воздуха в %, датчик N13". И в этих данных уже нету ни отрицательных значений, ни 146%, ни просто единичных выбросов, когда влажность, которая неделю держалась на уровне 90-94%, на одну секунду вдруг падает до нуля... а потом опять всю неделю 90-94%.
В общем, тут можно провести прямую аналогию между "публичным интерфейсом" (готовый сигнал) и ""приватной реализацией" (те самые АЦП, коэффициенты и чистка заведомо дефектных значений ;-)
Туда можно добавить масштабирующий фактор и единицы измерения, которые будут учитываться при конвертации.
Вы, вероятно, очень далеки от практических режимных наблюдений "в полях", поэтому, кажется, смотрите на вещи, как теоретик. Увы, но на практике мой сферический конь очень часто оказывается с одной стороны кубическим, с другой - плоским, а шея и голова у него и вовсе выполнены новомодно по Клейну...
Поясню на примере наших наблюдений в Баксане (см. рис.5, полный текст статьи доступен вот здесь). Первичные данные в этом случае - это ana-файл длиной 1 час, с заголовком вида "год-мес-день-час-регистратор_1388-канал_13". Хранить в базе миллион таких файлов, конечно можно. Только вот "базовая единица" при любой обработке у нас - это не часовой фрагмент, а временной ряд длиной в годы. А главное, при любой обработке мне (или коллеге) всегда нужны градусы температуры или нанорадианы наклона. А вовсе не условные значения от 0 до 32767 с многоуровневыми таблицами формул для их пересчета в разумные единицы. А еще в систему наблюдений иногда добавляются/убираются датчики. И номера каналов при этом, соответственно, меняются. Вы предлагаете дополнить основную базу таблицей переключений номеров каналов? А еще у наклономера есть такая особенность - наклон постоянно дрейфует (это типичный фликкер-шум). Поэтому раз в полгода специально обученный человек садится на электричку и едет в самый конец тоннеля, где стоит наклономер, чтобы выставить ноль. После чего другой специально обученный человек сшивает старый сигнал с новым, чтобы обеспечить однородность ряда. Причем это не совсем тривиальная операция, т.к. для такой сшивки надо сначала выкинуть переходный процесс (когда наблюдатель приближается к наклономеру, то он, естественно, это чувствует и наклон "убегает"), потом вычесть прилив и тренд, потом определить фактическую амплитуду скачка уровня, и только потом можно ввести поправку в исходный сигнал. Причем это делается не автоматически, а под контролем "эксперта", так как наблюдения экспериментальные, и почти каждый случай имеет свои особенности. В том числе это относится и к будущим случаям, которые еще не наступили. особенно интересно бывает, когда наблюдения только начинаются, и мы про эти "особенности" еще вообще ничего не знаем, т.к. у нас все без исключения подобные случаи - будущие..
И это я еще не начал говорить про аппаратурный брак - начиная от импульсных помех и кончая всякими "невозможными" глюками, которых в данных в принципе быть не может и не должно, если, конечно, производитель аппаратуры ничего не наврал в техдоках... С которыми, тоже, кстати, обычно приходится разбираться творчески. Поскольку данные на первой странице паспорта датчика противоречат данным на третьей странице того же паспорта, а формулы в обновленной версии паспорта (на тот же датчик!) противоречат обоим страницам первого паспорта. Зато если подставить в эти формулы все приведенные коэффициенты, включая температурные поправки датчика и т.д., то относительная влажность воздуха может запросто получиться 104%... И это не какая-то поделка от сомалийского жир-мол-треста, а небезызвестная Honeywell, чьи датчики аналогичных серий стоят и в нефтянке, и на самолетах... В общем, этот мир "железяк" совсем не так идеален, как это может показаться с позиций сферически чистого кодера в вакууме.
Так вот, если я правильно понял, Вы предлагаете делать все вышеперечисленное каждый раз, когда мне зачем-то понадобился ряд, или просто захотелось на него посмотреть? Нет уж, спасибо, я лучше как-нибудь все это выполню один раз при первичной обработке ;-) Но чтобы в базе потом была именно "относительная влажность" в процентах - очищенная от дефектов настолько, насколько это вообще возможно... См. рисунок 6 из той же статьи.
Впрочем, базу с первичными данными (обработка уровня 1 - сигнал, сшитый в один большой ряд из множества часовых фрагментов) мы при этом тоже храним. Но если честно, за три года наблюдений мы ей еще ни разу не пользовались...
Вы, вероятно, очень далеки от практических режимных наблюдений «в полях», поэтому, кажется, смотрите на вещи, как теоретикВ полях — да. Потому что у меня производство, всё под боком. Ну и специфика — а) производство непрерывное, б) опасное, в) сигналов на один агрегат тысячи и десятки тысяч, которые г) в идеале должны храниться вечно.
Так вот, если я правильно понял, Вы предлагаете...Не, вы неправильно поняли, я ничего не предлагаю. Я просто поделился своим опытом по схожей задаче. Насколько он вам или кому-то ещё релевантен/интересен вопрос уже другой.
Но ведь можно же просто не писать данные, когда их нет. А когда появляются, продолжать писать в другом файлике.
Можно (если честно, то они с АЦП именно так и пишутся ;-), но обрабатывать потом будет очень неудобно. В соседнем комментарии подробнее написал. А если на каждый пропуск новый файл начинать - то и вообще катастрофа случится. Т.к. этих пропусков в ВЧ сигнале всегда миллионы, а каждый файл - это не только данные, но и еще пустой фрагмент в конце сектора...
А можно их вообще не заменять. Шум, выходы за пределы и прочее — это тоже информация, которую можно измерять и анализировать.
Ну, тут есть шум и есть шум. Если у меня динамический диапазона прибора Е-5 радиан, то выброс с амплитудой Е-4 - это заведомый брак системы регистрации. На большинстве датчиков есть ограничения на скорость изменения сигнала, вытекающие из физики, и т.д. Поэтому такие заведомые дефекты лучше устранить на этапе первички. Ну а когда генезис проблемы не ясен - то конечно, такой шум идет в основную БД для последующего анализа.
Для этого давно придумали структуру под названием «chunk», когда данные хранятся не одним непрерывным куском, а разбиваются на одинаковые небольших размеров.
В разных случаях оптимальны разные подходы. Нас пока плоская база вполне устраивает. Нет никакого смысла (да и возможности) кратно усложнять алгоритмы ради совсем неочевидного выигрыша. Например, мне поступила на вход базы новая порция данных. Там могут быть фрагменты произвольной длины и с произвольными датами. Если у меня плоский ряд (столбец), я их просто накладываю в нужное место по определенным правилам (в конфиге есть всякие настройки, связанные с наложением пропусков на данные и т.д.). А так пришлось бы переконфигурировать все эти куски. Которые вдобавок в нашем случае пришлось бы делать либо разных размеров, либо очень маленькими. Ну и от пропусков внутри куска (в том числе и в его начале) это все равно не спасает, т.к. куски с совсем уж произвольными границами обрабатывать неудобно.
А еще не знаю я заранее, какую выборку потребует юзер. Например, первичный сигнал может быть 10-герцовый, а юзер захочет выборку с шагом 2.5 сек, причем первый дискрет должен начинаться в 22.2 секунды и т.д. В общем, плоские ряды (в одном большом файле) гораздо удобнее - и работает намного быстрее (все накопители идеально умеют кэширование вперед), и алгоритмы проще гораздо.
Перемотка в видео- и аудио-файлах достижима именно таким подходом.
А вот тут я не понял немного. Кто мешает перемотать в произвольную точку один большой файл? Или речь про магнитные ленты с последовательным доступом?
По вашей ссылке чужих не пускают.
Ой, извините, пожалуйста... Вроде бы там для просмотра текста регистрация не нужна была... Но сейчас сайт http://elibrary.ru и правда не открывается вовсе - глюк какой-то. Я-то туда по ночам обычно хожу, а днем, похоже, у них что-то там не справляется :-(
Вообще, если интересен наш подход к работе с пропущенными наблюдениями, то упомянутую статью можно скачать с Я-диска. Чуть больше технических подробностей про конкретные алгоритмы можно найти вот здесь (файл WinABD_Help.chm, разделы А4 и 3I1).
Но общая идея обработки сигналов с пропущенными наблюдениями там очень простая:
1) Чтобы аккуратно заполнить пропуски, нам обязательно нужна хорошая модель сигнала. Которая правильно подсказывает - какое именно значение мы бы могли получить в данный момент, если бы измерение прошло штатно.
2) Но при разведочном анализе у нас обычно еще нет этой модели сигнала. Наоборот, ее построение обычно как раз и является одной из целей такого анализа.
3) Если же мы будем заполнять пропуски на основе какой-то произвольно выбранной априорной модели, то велика вероятность, что потом при анализе мы эту априорную модель и получим (подтвердим). Особенно если пропусков было много... Попросту говоря, получим желаемое вместо действительного. Зачем тогда вообще измерения проводить?
Ну и некоторые практические следствия отсюда:
1) Пропуски надо заполнять только и исключительно тогда, когда без этого ну совсем никак. Перед БПФ, например
2) В остальных случаях надо свести алгоритм к расчету матожиданий (ну или других статистик, настолько же терпимых к пропускам). Так, для стационарного процесса можно выкинуть некоторое число измерений (которые мы хотели провести... но вот не смогли), и все равно получить несмещенную оценку. Да, дисперсия вырастет - ну так бесплатный сыр только в мышеловке, увы...
3) Аналогично при расчете ВКФ для стационарных процессов X и Y. Если мы выкинем некоторую часть пар {X(t),Y(t+dt)}, то все равно можем получить несмещенную оценку по оставшимся парам.
4) И так далее. То есть, во всех случаях, где мы сумели, мы заменяем классический алгоритм (работающий с эквидистантным сигналом) на его суррогатный аналог, где часть данных просто исключена из расчетов. По мере уменьшения числа пропусков, алгоритм сходится к классике.
А теперь некоторые контрапункты к этой идее.
1) А зачем вообще так мучиться-то? Есть же алгоритмы для работы с неэквидистантными (по t) временными рядами. Чем они не устраивают?
Ответ: алгоритмы-то есть, только вот их на порядок меньше, чем для эквидистантных сигналов. И расчетов там кратно больше, а теории почти нет. Так как "общий случай" для неэквидистантных измерений - это за гранью добра и зла. В общем, не получается обойтись только ими. Совсем. Нужно что-то выдумывать (компромиссы искать).
2) Математическая теория предлагает неисчислимое количество строгих методов с доказанной оптимальностью для статистической обработки данных. Вы что, самые умные - вместо этого кладезя мудрости предлагаете какие-то странные самодельные подходы использовать? А как же воспетая классиками строгость математическая?!
Ответ: мы не против математики, только вот в жизни экспериментатора-практика этой строгости нет, и не может быть почти никогда. Так как кроме самого метода, для нее необходимо точное соответствие данных той модели, в рамках которой построен наш строгий математический метод. А такое соответствие бывает чуть чаще, чем никогда. (И даже если оно действительно есть, доказать его все равно не получится, скорее всего).
Но если данные не удовлетворяют условиям применимости критерия/метода, то о какой строгости вообще может идти речь? В лучшем случае это будет иллюзия строгости, причем с превосходными шансами сесть в лужу со значимостью. Так что Вы можете верить в сказки... но если Вам нужны практически полезные результаты, то лучше все-таки повзрослеть...
3) Хорошо, пусть у вас все сводится к матожиданиям. Но на практике все ряды в геофизике - нестационарные. Эргодичность присутствует только в мультяшках для первоклашек и к взрослой жизни отношения не имеет. Как в вашем случае вообще можно что-то говорить про оценку матожиданий?!
Ответ: а все верно - никак нельзя. Просто альтернатив-то нету. Данные таковы, как они есть - это факт. Вот и приходится какие-то трюки выдумывать (например, декомпозировать сигнал на стационарную и нестационарную часть, хотя на выходе все равно получается суррогат). Правильного решения нет в принципе. Как сказал классик, все модели неверны, но некоторые из них полезны. Вот и у нас то же самое: анализ экспериментальных временных рядов в геофизике - это не о математической строгости. Совсем. Попросту говоря, это не наука, а скорее искусство. Вот так и живем...
А вот тут я не понял немного. Кто мешает перемотать в произвольную точку один большой файл?Это достижимо только если если замеры равноотстоящие во времени и не используется сжатие данных. Тогда индекс можно посчитать просто умножением. А если нет — придётся сканировать с самого начала (таблица индексов конечно ускорит этот процесс, но не решает проблему, с которой мы начали — как минимизировать кол-во байт под индивидуальную метку времени для каждого замера).
Можно пропустить данные через дельта-кодек, а потом сжать. А рядом положить файл "засечек", описывающий смещения в сжатых данных, чтобы быстро делать seek и разжимать только нужные блоки.
Если конвейер обработки данных фиксирован и заранее известен, то можно очень эффективно используя SIMD прогонять через него больше потоки. Если возможности зафиксировать его полностью нет - например, оператор в интерактивном режиме вводит аналитический запрос, то можно разбить его на стадии, а самые тугие циклы прогнать через JIT.
Сжатие данных экономит пропускную способность шин и позволяет записывать их в удобном человеку формате. Выигрыш в скорости часто оказывается больше, чем затраты ресурсов CPU на работу распаковщика.
Похожие симптомы могут быть при наличии неопределённого поведения (UB) в программе. Современные компиляторы стараются по полной использовать его и часто на его основе делают оптимизации в местах, казалось бы, не связанных с той инструкцией, которая непосредственно вызвала UB. Раскапывать такие баги бывает непросто (шаг влево, шаг вправо, добавили принт или дебаг - и эффект пропадает). Пожалуй, это заслуживает отдельного очерка (если доберусь) :)
Похожие симптомы могут быть при наличии неопределённого поведения (UB) в программе.
Я немного читал про UB в программах на Си, но в фортране такого вроде бы быть не должно? Кроме того, что язык к этому подталкивает (у меня включен максимальный уровень диагностики, включена проверка выхода индекса за пределы массива, и т.д.), я еще и сам в принципе не использую никаких трюков, которые могут хоть как-то к UB приводить - стараюсь чтобы все было однозначно. Например, не использую указатели, избавляюсь от устаревших конструкций вроде оператора common, счетчики циклов у меня 64-битные (чтобы переполнения не было), для real-операций в подозрительных местах стоят проверки диапазона данных (хотя по идее там по умолчанию должно получаться не UB, а inf, которое я сразу же конвертирую в Nun) и т.д. Тем более, в тестовой проге...
К сожалению, про UB в Си на Хабре пишут достаточно регулярно, а вот про UB в фортране я вообще ни одной статьи не нашел :-(
Потому что в настоящее время на Фортране пишут два с половиной человека. Формат float сам по себе практически UB и от языка программирования не зависит, равно как и другие виды UB.
ООП же совсем не случайно придумали. Как, например, в процедурном стиле защититься от складывания метров и километров, если и там и там тип float? Правильно — никак. А в ООП можно отдельно задать типы Meter и KiloMeter и их нельзя будет сложить до тех пор, пока явно не перегрузить операции сложения и приведения к одному типу.
Формат float сам по себе практически UB и от языка программирования не зависит, равно как и другие виды UB.
Справедливости ради, в фортране есть довольно много встроенных функций, которые предназначены как раз для контроля точности вычислений и предотвращения float UB. Конечно, это не панацея, но если я правильно понимаю, в других языках их аналоги появились заметно позже. А чтобы выйти за границы массива, надо специально постараться, так как, к примеру, при его передаче в функцию в качестве параметра туда одновременно передается размер и форма. Включая и динамические массивы. Аналогично со строками. Ну и не надо использовать указатели там, где они не нужны (а в фортране они почти нигде не нужны). Кстати, чтобы попасть на UB в фортране, мало все это проигнорировать - надо еще и выключить соответствующие проверки в компиляторе. Которые, если я правильно понимаю, реализованы достаточно эффективно - так как их выключение в типичных случаях на производительность почти не влияет (а значит и незачем выключать).
Конечно, выстрелить себе в ногу можно на любом языке. Но сделать это неумышленно в фортране достаточно сложно. Если, конечно, не писать в стиле F-77, со всеми этими COMMОN, EQUIVALENCE, ENTRY и прочими GOTO.
Потому что в настоящее время на Фортране пишут два с половиной человека.
Я бы даже сказал, что Вы сильно завысили эту цифру, если речь идет о типичных коммерческих приложениях вроде бухгалтерии или склада. Но, тем не менее, существуют также отдельные узкие ниши, где таких людей около половины от всех программистов. Научные вычисления - это как раз одна из таких ниш.
ООП же совсем не случайно придумали
Кстати, да! Очень уместное замечание, так как начиная с фортрана 90, а тем более - с фортрана-2003, основные элементы ООП вполне себе поддерживаются языком. Конечно, я лишь очень поверхностно наслышан про ООП, так как почти не знаю других языков, кроме фортрана. Но если судить про ООП по всяческим википедиям, то в фортране практически все это есть, причем пользоваться этим легко и просто. Производные типы, объединение данных и работающих с ними методов в модулях (фактически это классы), публичные интерфейсы и приватное внутреннее устройство, конструкторы (чуть хуже с деструкторами), наследование, перегрузка функций и операторов - все это есть из коробки "на блюдечке". Даже такие динозавры, как я, все вышеперечисленное (кроме деструкторов) используют постоянно. А конструкции на основе FORALL и WHERE позволяют на 90% избавиться от явно прописанных циклов при обработке многомерных массивов, причем без потери эффективности (наоборот, обещается, что компилятор якобы сам все это распараллелит - и есть подозрения, что это у него даже иногда получается ;-). В сумме все это дает возможность непрофессионалам писать достаточно эффективные программы для научных расчетов. Что при использовании, например, С++ было бы невозможно (там для этого требуется совершенно другой уровень программистской квалификации).
Чаще всего floating point реализован по IEEE 754 и не зависит от языка. Но к такому может приводить undefined behavior вообще в другой части кода, а не в какой-нибудь FP-инструкции. Я на фортране не пишу, но гугл говорит, что UB в нём есть.
Чаще всего floating point реализован по IEEE 754
Да, в интел-фортране именно так (во всяком случае, по умолчанию)
и не зависит от языка.
От языка, может, и не зависит, а вот от компилятора - очень даже. В том же интел-фортране есть тьма опций для более тонкой донастройки floating point. Беда в том, что перебором проверить их все - нереально, а где именно может скрываться причина моего бага, я понять не могу...
Я на фортране не пишу, но гугл говорит, что UB в нём есть.
Ну, разумеется, есть ;-)))
Я, конечно, читал, что в некоторых современных языках
UB невозможно.
Однако мне очень сложно представить компилятор, который бы действительно давал 100%-ную гарантию отсутствия UB для не совсем тривиального (не кастрированного) языка (впрочем, это может быть от ограниченности моего кругозора - я новыми языками интересуюсь только сугубо теоретически и о многом даже не подозреваю наверняка).
Но вот сомневаться в наличии UB в языке, декларирующем совместимость с кодом 60-летней давности, это было бы
очень странно
Если нужно, я сам могу на интел-фортране программу с UB написать и такие опции компилятора подобрать, чтобы он ее пропустил ;-)
Вопрос скорее в другом - что я стараюсь опасные конструкции не использовать, и стелю соломку везде, где умею. Но, конечно, такой шанс всегда остается. Даже больше того, у меня кое-где конкретные подозрения есть. Только вот проверить их я пока не сумел :-(
Нормально вы так заморочились) А причина бага понятна, если посмотреть на ассемблерный код — это использование FPU несмотря на то, что явно указано использовать SSE. Это может быть только по двум причинам:
1) в библиотеке рандома явно указан расширенный тип с плавающей точкой (80-битный) и у компилятора нет другого выбора, или
2) разработчики компилятора слегка схалтурили.
В чём проблема использования FPU:
1) он имеет стековую организацию регистров,
2) на регистры с плавающей точкой маппятся регистры MMX и их нельзя использовать одновременно,
3) он спроектирован для однопоточного использования и имеет глобальные статусы, определяющие точность вычисления, сторону округления и прочее.
Собственно сам код FPU из вашего листинга просто умножает число с плавающей точкой на целое число, в каком формате сохраняет, без опкодов непонятно. Если в целом — то возможно исключение, если результат "не влезет". При этом в коде присутствует загрузка в стек числа, которое затем затирается — тоже непонятно, зачем и на оптимизирующую компиляцию похоже мало.
Правильное решение — это:
1) Отключите /Qparallel. Если вам нужна многопоточность — её лучше делать вручную точно зная, какие ресурсы разделяемы, какие нет, и как их правильно синхронизировать;
2) напишите свой генератор случайных чисел;
3) перейдите на другой, более современный язык программирования. Чем больше популярен язык, тем меньше багов у компилятора.
Правильное решение — это:
1) Отключите /Qparallel. Если вам нужна многопоточность — её лучше делать вручную точно зная, какие ресурсы разделяемы, какие нет, и как их правильно синхронизировать;
Пробовал. Скорость резко падает практически на всех машинах, где удалось провести тесты, а баг не исчезает (проверено на одной из машин) :-((
2) напишите свой генератор случайных чисел;
Уже давно написал обертку к нему (просто вызываю ГСЧ повторно, если он вернул Nun). Но к сожалению, баг появляется не только при вызове ГСЧ. Просто там я его впервые заметил, т.к. ГСЧ вообще не должен Nan продуцировать. У меня ведь Nun - это вполне штатное значение данных. В экспериментальных рядах мониторинга они сплошь и рядом встречаются, особенно после чистки аппаратурного брака. Поэтому эти спонтанные "десятимиллионники" не так просто заметить, если специально не заморачиваться. Но главное - они появляются в непредсказуемом месте, в том числе и при вызове библиотечных функций. А к каждой функции обертку не сделаешь...
3) перейдите на другой, более современный язык программирования. Чем больше популярен язык, тем меньше багов у компилятора.
Хороший совет, только вот нереалистичный. У меня все задачи расчетные, там фортрану нет равных по простоте языка, при том, что по производительности он в топе. Одни массивные операторы и сечения чего стоят. Но главное, у меня несколько сотен тысяч строк кода собственного legacy, и весь этот код постоянно работает и используется. Новые методы просто встраиваются в этот пакет. Для многих коллег это основной рабочий инструмент, хотя и весьма узконишевый. Переписывать его на другом языке - задача неподъемная и бессмысленная. Да и где гарантия, что новые компиляторы будут менее подвержены багам? Скорее наоборот - чем более бурно язык развивается, тем больше шансов, что там где-то есть глюки. Я за разными языками не особо слежу, но вот только недавно здесь на Хабре ругали один из современных языков за многочисленные ошибки, которые годами не устраняются. Или можно вспомнить про другой очень популярный современный язык, стандарт которого при выпуске очередной версии поменялся настолько, что старые программы проще запускать в среде 15-летней давности, чем переписать под новый стандарт (впрочем, этот язык для моих задач по-любому не актуален).
Впрочем, Вы правы, что мне стоило бы на наиболее свежую версию компилятора перейти. Возможно, там этот баг поправлен. Но тут уже не все от меня зависит...
Если в программе в произвольный момент времени и в произвольном месте может появиться ничем не обусловленный NaN — это не баг, а эпик фейл. Я бы не доверился ни такой программе, ни компилятору (если проблема действительно в нём, а не в программе), ни программисту, который решает такие проблемы костылями типа while(isNaN(x)) x=random(); (ничего личного). Если компилятор приобретён легально — то обращаться нужно к самим Intel. Если нет — обращаться надо к рутрекеру за другими, не обязательно более новыми версиями.
Я бы не доверился ни такой программе, ни компилятору ... ни программисту
Да, ситуация скользкая. Однако у меня ведь не бухгалтерия, где исходные данные имеют абсолютную точность, и надо каждую копейку без погрешностей посчитать. У меня экспериментальный сигнал, который всегда измеряется с приличной ошибкой. Иногда там из 16-ти бит только 4-5 точных, а все остальное - шум. Шум же не белый, а фликкерный... Поэтому сигнал постоянно гуляет и выбрать оптимальный диапазон очень непросто.
Кроме того, у меня на каждый такой эпик-фейл-Nan есть примерно миллион "своих" Nan-ов. Которые всегда обрабатываются, так что на результат вычислений их наличие почти никак не влияет (конечно, пока "хороших значений заметно больше, чем Nan-ов - но это условие, к счастью, почти всегда выполнено). Фактически из-за этого "эпик-фейла" я просто теряю еще одно значение данных из 10 миллионов в дополнение к миллиону уже потерянных раньше.
Для любого разумного алгоритма такая потеря не должна быть критичной. К примеру, при вычислении матожидания потеря даже 10% случайно выбранных значений де-факто почти не скажется на дисперсии результата. Тем более, если потеряны малые доли процента.
Вообще, принцип робастности настоятельно требует, чтобы при потере 0.0001% значений данных результат совсем не менялся. Если вдруг это не так, то для анализа экспериментальных сигналов такой алгоритм (который не слишком устойчив к малому изменению входных данных) категорически непригоден. Мы ведь обычно ищем устойчивые закономерности на фоне шумов, а они не должны опираться на единичные артефакты...
Возможно, Вы имели в виду ситуацию, когда я использую Nan, как обычное число при расчетах? Понятно, что даже единичный гигантский выброс может катастрофически исказить результат. Но как я уже написал, это абсолютно исключено, т.к. проверка на Nan выполняется абсолютно при любых вычислениях. Иначе у меня бы вообще ничего не работало.
Если компилятор приобретён легально — то обращаться нужно к самим Intel. Если нет — обращаться надо к рутрекеру за другими, не обязательно более новыми версиями.
Если честно, то тут у меня тоже все в "серой" зоне. Компилятор был куплен легально, но очень давно, на организацию, но после истечения срока лицензии она официально не продлевалась. То есть, изначально это был "белый" продукт, но теперь уже нет. Поэтому я, с одной стороны, не хочу от него отказываться и переходить на что-то совсем пиратское, но с другой - пока не могу убедить организацию купить за немаленькие деньги актуальную версию. Хотя и надеюсь, что этот вопрос все же решится за какое-то разумное время.
Если честно, то тут у меня тоже все в «серой» зоне. Компилятор был куплен легально, но очень давно, на организацию, но после истечения срока лицензии она официально не продлевалась.Я ещё обратил внимание, что код 32-битный, заметно старый компилятор. Проблемы могут быть ещё и от того, что он не учитывает особенности современных процессоров с технологиями типа Hyper-Threading.
И на Интеле же свет клином не сошёлся, есть же ещё как минимум GNU Fortran. А экстремальная оптимизация в любом случае достигается лишь ручным программированием на ассемблере.
Я ещё обратил внимание, что код 32-битный, заметно старый компилятор.
Да, компилятор действительно старый (надеюсь организация сможет это исправить).
Проблемы могут быть ещё и от того, что он не учитывает особенности современных процессоров с технологиями типа Hyper-Threading.
Да, тут Вы в точку попали, это было одним из основных подозрений сначала. Правда, после тестов примерно на десятке компов оказалось, что баг возникает (или не возникает) на разных машинах почти вне связи с их возрастом. Но все равно эта версия пока одна из приоритетных. Только вот компетенций, чтобы копнуть ее вглубь, у нас ни у кого нету :-(
И на Интеле же свет клином не сошёлся, есть же ещё как минимум GNU Fortran. А экстремальная оптимизация в любом случае достигается лишь ручным программированием на ассемблере.
Он, к сожалению, заметно медленнее. Мы делали специальные тесты (именно под наши типичные вычисления) при выборе компилятора. Правда, все это было очень давно, поэтому не исключено, что современный GNU Fortran будет не хуже старого Intel Fortran-а. Но я пока что надеюсь на покупку более современного компилятора Intel - , как для скорости, так и чтобы программы не переделывать. (Например, у меня для мнемоники имена констант начинаются с $, а это уже расширение стандарта языка, и др.).
Отключите /Qparallel. Пробовал.Из этого, кстати, программа вовсе не обязана стать однопоточной, особенно если остаётся возможность распараллеливать задачи вручную. Реальное количество потоков в программе можно увидеть даже через стандартный диспетчер задач.
В статье не хватает описания одной интересной особенности.
Операции с памятью пришиваются как-то с боку к SSA-форме. Операции load не производят значения, но имеют сторонний эффект и посему не могут быть так просто выкинуты. Операции store же зависят не только от входных значений, но и от состояния внешнего по отношению к SSA-форме мира - памяти.
Вторая часть про доминирование тут.
Поговорим об оптимизирующих компиляторах. Сказ первый: SSA-форма