Сейчас я буду рассказывать и показывать на примерах, почему физики тоже должны использовать инструменты статического анализа кода. Мне бы хотелось, чтобы этим инструментом был PVS-Studio. Но, конечно, любой другой инструмент тоже будет полезен. Анализатор кода сократит время на отладку приложений и уменьшит головные боли от тупых ошибок. Лучше побольше думать о физике и поменьше об ошибках в программах на языке Си++.
Введение
Я продолжаю серию заметок, посвященных анализу кода, используемого в научных целях. Предыдущие заметки:
В этот раз на очереди проект Geant4. Приведу описание проекта, взятое из Wikipedia:
Geant4 (англ. GEometry ANd Tracking — геометрия и трекинг) — программа для моделирования прохождения элементарных частиц через вещество с использованием методов Монте-Карло. Разработана в CERN на объектно-ориентированном языке программирования С++. Является дальнейшим развитием предыдущих версий GEANT, существенно переработанным и дополненным.
Как заявлено на официальном сайте проекта, «области применения включают в себя физику высоких энергий и исследование ядерных реакций, медицину, ускорители частиц, и космические физические исследования». ПО используется во многих исследовательских проектах по всему миру, в том числе и в России.
Сайт проекта: http://geant4.org. Проект имеет средний размер: 76 мегабайт исходного кода. Для сравнения:
- проект VirtualDub — 13 мегабайт;
- проект Apache HTTP Server — 26 мегабайт;
- проект Chromium (вместе с дополнительными библиотеками) — 710 мегабайт.
Для анализа был использован анализатор кода PVS-Studio. Так как проект Geant4 не маленький, шанс найти в нём интересное был очень велик. Не найти ошибок можно только в маленьких проектах (см. о нелинейной плотности ошибок). Бывают, конечно, и крупные проекты, где PVS-Studio ничего не находит. Жаль, но это только исключение.
Сразу хочу попросить прощения, если случайно написал какую-то глупость, касающуюся физики. Но обратите внимание, что я смог найти ошибки в программном обеспечении, не понимая суть партонов, да и вообще почти всего касающегося ядерных реакций!
Примечание. В статье перечислено далеко не всё, что удалось заметить. Более полный список предупреждений, на которые я обратил внимание: geant4.txt.
Давайте же посмотрим, что удалось найти интересного в коде Geant4.
Copy-Paste и мюоны
Про мюоны я знаю только то, что это элементарная частица. Зато я гораздо лучше понимаю, что такое Copy-Paste. Вот красивый пример ошибки, когда строка кода была размножена. Скопированные строки подверглись изменениям, но не везде. Результат:
void G4QMessenger::SetNewValue(G4UIcommand* aComm, G4String aS) { if(photoDir) { if (aComm==theSynchR) thePhoto->SetSynchRadOnOff(aS); else if(aComm==minGamSR) thePhoto->SetMinGammaSR(.... else if(aComm==theGamN) thePhoto->SetGammaNuclearOnOff(.... else if(aComm==theMuoN) thePhoto->SetElPosNuclearOnOff(.... else if(aComm==theMuoN) thePhoto->SetMuonNuclearOnOff(aS); else if(aComm==theMuoN) thePhoto->SetTauNuclearOnOff(aS); else if(aComm==biasPhotoN)thePhoto->SetPhotoNucBias(.... } .... }Предупреждение PVS-Studio: V517 The use of ‘if (A) {…} else if (A) {…}’ pattern was detected. There is a probability of logical error presence. Check lines: 195, 196. G4phys_builders g4qmessenger.cc 195
Обратите внимание на три раза повторяющуюся проверку (aComm==theMuoN).
Распад барионов
Тяжело изучать радиоактивный распад или пытаться выявить распад протона. Особенно это сложно делать, когда программы содержат ошибки.
void G4QEnvironment::DecayBaryon(G4QHadron* qH) { .... else if(qM<mSzPi) // Only Lambda+PiM is possible { fQPDG=lQPDG; // Baryon is Lambda fMass=mLamb; sQPDG=pimQPDG; // Meson is Pi- sMass=mPi; } else if(qM<mSzPi) // Both Lambda+PiM & Sigma0+PiM are possible { if(G4UniformRand()<.6) { .... }Предупреждение PVS-Studio: V517 The use of ‘if (A) {…} else if (A) {…}’ pattern was detected. There is a probability of logical error presence. Check lines: 8373, 8380. G4hadronic_body_ci g4qenvironment.cc 8373
Два раза используется одно и то же условие (qM<mSzPi). Это значит, что мы никогда не попадем в ветку, к которой относится комментарий "// Both Lambda+PiM & Sigma0+PiM are possible".
Посчитаем количество партонов
Партон — точечноподобная составляющая адронов, проявляющаяся в экспериментах по глубоко неупругому рассеянию адронов на лептонах и других адронах.
Жаль, но с подсчетом их количества есть проблема:
G4bool G4CollisionMesonBaryonElastic:: IsInCharge(const G4KineticTrack& trk1, const G4KineticTrack& trk2) const { G4bool result = false; G4ParticleDefinition * p1 = trk1.GetDefinition(); G4ParticleDefinition * p2 = trk2.GetDefinition(); if( (GetNumberOfPartons(p1) != 2 || GetNumberOfPartons(p2) != 3) ||(GetNumberOfPartons(p1) != 3 || GetNumberOfPartons(p2) != 2) ) { result = false; } .... }Предупреждение PVS-Studio: V547 Expression is always true. Probably the ‘&&’ operator should be used here. G4had_im_r_matrix g4collisionmesonbaryonelastic.cc 53
Сразу суть ошибки может быть незаметна. Поэтому я упрощу выражение:
A = GetNumberOfPartons(p1); B = GetNumberOfPartons(p2); if ( (A != 2 || B != 3) || (A != 3 || B != 2) ) Скобки для простоты тоже можно отбросить. Получаем:
if ( A != 2 || B != 3 || A != 3 || B != 2 )Это условие всегда истинно. Переменная ‘A’ всегда не равна 2 или не равна 3. Аналогично и с переменной ‘B’. Видимо, где-то произошла путаница. Скорее всего, здесь должен присутствовать оператор ‘&&’.
Кулоновская блокада и выход за границы массива
Кулоновская блокада — блокирование прохождения электронов через квантовую точку, включённую между двумя туннельными контактами, обусловленное отталкиванием электронов в контактах от электрона на точке, а также дополнительным кулоновским потенциальным барьером, который создаёт электрон, усевшийся на точке.
С функцией SetCoulombEffects() что-то не так. В массив указателей ‘sig’ помещаются адреса несуществующих подмассивов. Работа с ‘sig’ приведёт к неопределённому поведению программы. Хорошо, если программа завершится в аварийном режиме. Гораздо хуже, если программа продолжит работу и будет хаотично читать и писать в память, занятую другими массивами и переменными.
enum { NENERGY=22, NANGLE=180 }; class G4LEpp : public G4HadronicInteraction { .... G4float * sig[NANGLE]; static G4float SigCoul[NENERGY][NANGLE]; .... }; G4LEpp::SetCoulombEffects(G4int State) { if (State) { for(G4int i=0; i<NANGLE; i++) { sig[i] = SigCoul[i]; } elab = ElabCoul; } .... }Предупреждение PVS-Studio: V557 Array overrun is possible. The value of ‘i’ index could reach 179. g4lepp.cc 62
Массив ‘sig’ содержит 180 указателей. Эти указатели должны ссылаться на разные строки двумерного массива ‘SigCoul’. Но в массиве ‘SigCoul’ только 22 строки. Получается, что большинство указателей в массиве ‘sig’ будут указывать непонятно куда.
Затруднительно сказать, в чем именно состоит ошибка. Возможно, неправильно объявлен массив ‘SigCoul’. Например, надо поменять местами количество строк и столбцов:
SigCoul[NENERGY][NANGLE] —>> SigCoul[NANGLE][NENERGY]
Опечатка и Twist Surfaces
void G4VTwistSurface::GetBoundaryLimit(G4int areacode, G4double limit[]) const { .... if (areacode & sC0Min1Max) { limit[0] = fAxisMin[0]; limit[1] = fAxisMin[1]; } else if (areacode & sC0Max1Min) { limit[0] = fAxisMax[0]; limit[1] = fAxisMin[1]; } else if (areacode & sC0Max1Max) { limit[0] = fAxisMax[0]; limit[1] = fAxisMax[1]; } else if (areacode & sC0Min1Max) { limit[0] = fAxisMin[0]; limit[1] = fAxisMax[1]; } .... }Предупреждение PVS-Studio: V517 The use of ‘if (A) {…} else if (A) {…}’ pattern was detected. There is a probability of logical error presence. Check lines: 793, 802. G4specsolids g4vtwistsurface.cc 793
Существует 4 переменные:
- sC0Min1Max
- sC0Max1Min
- sC0Min1Min
- sC0Max1Max
При работе с границами используются только 3 из них. При этом проверка (areacode & sC0Min1Max) выполняется 2 раза. Если присмотреться, то видно, что после первой проверки выбираются минимумы: fAxisMin[0], fAxisMin[1]. Скорее всего, первая проверка должна была быть такой:
if (areacode & sC0Min1Min) { limit[0] = fAxisMin[0]; limit[1] = fAxisMin[1]; }
Рентгеновское излучение и неправильный IF
Рентгеновское излучение — электромагнитные волны, энергия фотонов которых лежит на шкале электромагнитных волн между ультрафиолетовым излучением и гамма-излучением.
Как я понимаю, класс G4ForwardXrayTR имеет отношение к X-Ray. Ошибка закралась в сравнение индексов.
G4VParticleChange* G4ForwardXrayTR::PostStepDoIt(....) { .... if (iMat == jMat || ( (fMatIndex1 >= 0 && fMatIndex1 >= 0) && ( iMat != fMatIndex1 && iMat != fMatIndex2 ) && ( jMat != fMatIndex1 && jMat != fMatIndex2 ) ) .... }Предупреждение PVS-Studio: V501 There are identical sub-expressions to the left and to the right of the ‘&&’ operator: fMatIndex1 >= 0 && fMatIndex1 >= 0 G4xrays g4forwardxraytr.cc 620
Два раза проверяется индекс с именем ‘fMatIndex1’. А про ‘fMatIndex2’ забыли. Наверное, должно быть написано так:
(fMatIndex1 >= 0 && fMatIndex2 >= 0)
Опечатка и нейтроны
G4double G4MesonAbsorption:: GetTimeToAbsorption( const G4KineticTrack& trk1, const G4KineticTrack& trk2) { .... if(( trk1.GetDefinition() == G4Neutron::Neutron() || trk1.GetDefinition() == G4Neutron::Neutron() ) && sqrtS>1.91*GeV && pi*distance>maxChargedCrossSection) return time; .... }Предупреждение PVS-Studio: V501 There are identical sub-expressions ‘trk1.GetDefinition() == G4Neutron::Neutron()’ to the left and to the right of the ‘||’ operator. G4had_im_r_matrix g4mesonabsorption.cc 285
Да, я не знаю, что делает эта функция. Но, как мне кажется, на вход она получает две траектории движения элементарных частиц. Функция должна обработать особым образом ситуацию, если хотя бы одна частица являются нейтронами. Но, на самом деле, проверятся только первая частица.
Видимо, на самом деле, планировали написать:
trk1.GetDefinition() == G4Neutron::Neutron() || trk2.GetDefinition() == G4Neutron::Neutron() Аналогичную опечатку можно обнаружить здесь: g4scatterer.cc 138
Добавление партона
Существует функция InsertParton(), которая вставляет в контейнер партон. Можно указать, после какого партона нужно вставить новый элемент. Если место вставки не указано, то, наверное, можно вставить куда угодно. Но как раз этот случай реализован некорректно.
typedef std::vector<G4Parton *> G4PartonVector; inline void G4ExcitedString::InsertParton( G4Parton *aParton, const G4Parton * addafter) { G4PartonVector::iterator insert_index; .... if ( addafter != NULL ) { insert_index=std::find(thePartons.begin(), thePartons.end(), addafter); .... } thePartons.insert(insert_index+1, aParton); }Предупреждение PVS-Studio: V614 Potentially uninitialized iterator ‘insert_index’ used. g4excitedstring.hh 193
Если указатель ‘addafter’ равен нулю, то итератор «insert_index» остаётся неинициализированным. В результате, вставка нового элемента может привести к непредсказуемым последствиям.
Обработка не всех нуклонов
Нуклоны — общее название для протонов и нейтронов.
Имеется функция packNucleons(), которая обрабатывает не все элементы. Дело в том, что цикл завершается сразу после первой итерации. В конце тела цикла имеется оператор ‘break’, но при этом нигде нет оператора ‘continue’.
void G4QMDGroundStateNucleus::packNucleons() { .... while ( nmTry < maxTrial ) { nmTry++; G4int i = 0; for ( i = 1 ; i < GetMassNumber() ; i++ ) { .... } if ( i == GetMassNumber() ) { areTheseMsOK = true; } break; } .... }Предупреждение PVS-Studio: V612 An unconditional ‘break’ within a loop. g4qmdgroundstatenucleus.cc 274
Мне кажется, оператор ‘break’ в конце просто лишний и написан случайно.
Lund string model и опечатка в индексе
In particle physics, the Lund string model is a phenomenological model of hadronization.
Когда приходится работать с отдельными элементами массивов, очень легко допустить опечатку. Именно это и произошло в конструкторе класса G4LundStringFragmentation. В приведенном ниже коде, ошибка сразу видна. Два раза присваиваем значения одной и той же ячейке. Однако, это очень большая функция, где инициализируется множество элементов массива. И заметить ошибку, рассматривая функцию, крайне сложная задача. Этот тот случай, где статический анализ кода просто незаменим.
G4LundStringFragmentation::G4LundStringFragmentation() { .... BaryonWeight[0][1][2][2]=pspin_barion*0.5; .... BaryonWeight[0][1][2][2]=(1.-pspin_barion); .... }Предупреждение PVS-Studio: V519 The ‘BaryonWeight[0][1][2][2]’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 205, 208. g4lundstringfragmentation.cc 208
Примечание. Хочу отметить, что в коде мне встретилось множество мест, где переменной два раза подряд присваиваются различные значения. Многие места вполне безобидны. Например, в начале, переменным присваивается 0, а затем уже настоящее значение. Однако, многие из таких мест могут содержать ошибку. Поэтому я рекомендую авторам Geant4 внимательно отнестись к диагностическим предупреждениям V519. Сам я просмотрел эти предупреждения очень поверхностно.
Кстати, не очень понятен подход, когда вначале переменная инициализируется значением по умолчанию, а затем уже нужным значением. Зачем это нужно? Проще сразу объявлять переменную там, где она нужна, и инициализировать нужным числом.
Некоторые другие предупреждения V519
Подозрительно смотрится оператор копирования в классе G4KineticTrack:
const G4KineticTrack& G4KineticTrack::operator=( const G4KineticTrack& right) { .... the4Momentum = right.the4Momentum; the4Momentum = right.GetTrackingMomentum(); .... }Предупреждение PVS-Studio: V519 The ‘the4Momentum’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 451, 452. g4kinetictrack.cc 452
Кстати, на конструкторы вообще было много предупреждений V519. Быть может, в таком коде есть смысл? Например, для целей отладки. Непонятно. Так что думаю, стоит привести ещё пару примеров подозрительного кода:
void G4IonisParamMat::ComputeDensityEffect() { .... fX0density = 0.326*fCdensity-2.5 ; fX1density = 5.0 ; fMdensity = 3. ; while((icase > 0)&&(fCdensity < ClimiG[icase])) icase-- ; fX0density = X0valG[icase]; fX1density = X1valG[icase]; .... }Предупреждения PVS-Studio: V519 The ‘fX0density’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 245, 247. g4ionisparammat.cc 247
V519 The ‘fX1density’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 245, 247. g4ionisparammat.cc 247
void G4AdjointPhotoElectricModel::SampleSecondaries(....) { .... pre_step_AdjointCS = totAdjointCS; post_step_AdjointCS = AdjointCrossSection(aCouple, electronEnergy,IsScatProjToProjCase); post_step_AdjointCS = totAdjointCS; .... }Предупреждение PVS-Studio: V519 The ‘post_step_AdjointCS’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 76, 77. g4adjointphotoelectricmodel.cc 77
И последний, замеченный мной подозрительный фрагмент. Обратите внимание на член ‘erecrem’.
void G4Incl::processEventIncl(G4InclInput *input) { .... varntp->mzini = izrem; varntp->exini = esrem; varntp->pxrem = pxrem; varntp->pyrem = pyrem; varntp->pzrem = pzrem; varntp->mcorem = mcorem; varntp->erecrem = pcorem; varntp->erecrem = erecrem; .... }Предупреждение PVS-Studio: V519 The ‘varntp->erecrem’ variable is assigned values twice successively. Perhaps this is a mistake. Check lines: 896, 897. g4incl.cc 897
Счет элементов массива с 1
Мне кажется, здесь кто-то на мгновение забыл, что элементы массивов в Си++ нумеруются с нуля. В результате, может произойти выход за границу массива. Плюс отсутствует сравнение со значением 1.4, расположенным в самом начале массива.
void G4HEInelastic::MediumEnergyClusterProduction(....) { .... G4double alem[] = {1.40, 2.30, 2.70, 3.00, 3.40, 4.60, 7.00}; .... for (j = 1; j < 8; j++) { if (alekw < alem[j]) { jmax = j; break; } } .... }Предупреждение PVS-Studio: V557 Array overrun is possible. The value of ‘j’ index could reach 7. g4heinelastic.cc 4682
Физика и undefined behavior
Си++ жестокий язык. Чуть-чуть расслабился — и тут же отстрелил себе протоном ногу. Причем, сразу это можно и не заметить. Вот пример неправильного реализованного оператора сложения:
template <typename T> GMocrenDataPrimitive<T> & GMocrenDataPrimitive<T>::operator + (const GMocrenDataPrimitive<T> & _right) { GMocrenDataPrimitive<T> rprim; .... return rprim; }Предупреждение PVS-Studio: V558 Function returns the reference to temporary local object: rprim. G4GMocren g4gmocrenio.cc 131
Функция возвращает ссылку на локальный объект. Попытка работать с этой ссылкой приведёт к неопределенному поведению.
Пора закругляться
К сожалению, экскурсию в мир физики пора заканчивать. Дело в том, что это всё-таки статья, а не многостраничный отчёт. Жаль, что я не успел рассказать про многие другие ошибки. Более подробно познакомиться с подозрительным кодом, который я успел заметить, просматривая диагностические предупреждения PVS-Studio, можно в этом файле — geant4.txt.
Только прошу не делать вывод, что это все ошибки, который смог найти PVS-Studio. Я смотрел отчёт достаточно поверхностно и мог много чего пропустить. Прошу разработчиков проекта проверить код с помощью PVS-Studio самостоятельно. Напишите нам, и мы выделим для этого ключ на некоторое время.
Ну и конечно, как всегда повторюсь, что эффективнее использовать статический анализ кода регулярно, а не от случая к случаю. Чтобы осознать важность регулярного использования, очень прошу почитать здесь и здесь.
Как я уже сказал, в файле приведено намного больше подозрительных фрагментов кода, чем в этой статье. Я думаю все ситуации достаточно понятны. В крайнем случае, хочу напомнить, что для каждого кода ошибки в документации есть подробное описание с примерами.
Единственно, хочу кратко пояснить диагностики V636 и V624. Иногда они могут свидетельствовать о неточности в вычислениях. Мне кажется, эти диагностики важны, когда речь идёт о вычислительных программах.
Пример диагностики V636:
G4double G4XAqmTotal::CrossSection( const G4KineticTrack& trk1, const G4KineticTrack& trk2) const { .... G4int sTrk1 = ....; G4int qTrk1 = ....; G4double sRatio1 = 0.; if (qTrk1 != 0) sRatio1 = sTrk1 / qTrk1; .... }Предупреждение PVS-Studio: V636 The ‘sTrk1 / qTrk1’ expression was implicitly casted from ‘int’ type to ‘double’ type. Consider utilizing an explicit type cast to avoid the loss of a fractional part. An example: double A = (double)(X) / Y;. g4xaqmtotal.cc 103
Результат деления «double X = 3/2» равен 1, а не 1.5 как может показаться из-за невнимательности. В начале, происходит целочисленное деление, а уже затем результат превращается в тип ‘double’. Рассматривая незнакомый код, трудно сказать, нужно ли действительно использовать целочисленное деление или это ошибка. Думаю, такие места в Geant4 заслуживают внимания.
Примечание. Если требуется именно целочисленное деление, такие места стоит сопровождать комментарием.
Пример диагностики V624:
dSigPodT = HadrTot*HadrTot*(1+HadrReIm*HadrReIm)* (....)/16/3.1416*2.568;Предупреждение PVS-Studio: V624 The constant 3.1416 is being utilized. The resulting value could be inaccurate. Consider using the M_PI constant from <math.h>. g4elastichadrnucleushe.cc 750
Непонятно, почему массово используются жесткие константы для обозначения числа Пи, Пи/2 и так далее. Готов поверить, что точности вполне хватает. Но всё равно, это не объяснение. Такие константы — лишний повод для опечаток и неточностей. Есть масса предопределенных констант, таких как M_PI, M_PI_4, M_LN2. PVS-Studio выдаёт соответствующие рекомендации по замене жестких констант.
Чуть не забыл. В файле geant4.txt я также привёл сообщения, относящиеся к микрооптимизациям. Например, по поводу инкремента итераторов:
class G4PhysicsTable : public std::vector<G4PhysicsVector*> { .... }; typedef G4PhysicsTable::iterator G4PhysicsTableIterator; inline void G4PhysicsTable::insertAt (....) { G4PhysicsTableIterator itr=begin(); for (size_t i=0; i<idx; ++i) { itr++; } .... }Предупреждение PVS-Studio: V803 Decreased performance. In case ‘itr’ is iterator it’s more effective to use prefix form of increment. Replace iterator++ with ++iterator. g4physicstable.icc 83
Подобнее, почему этот стоит изменить, описано в статье: Есть ли практический смысл использовать для итераторов префиксный оператор инкремента ++it, вместо постфиксного it++.
Заключение
Смиритесь, что ошибки и опечатки допускают все. И да, даже вы, уважаемые читатели. Это неизбежно. Используя инструменты статического анализа кода, можно исправить множество ошибок на самых ранних этапах. Это сократит время, необходимое на поиск дефектов и позволит большее внимание уделить непосредственно решаемой задачи.
Дополнительные ссылки
- Андрей Карпов. Мифы о статическом анализе. Миф второй — профессиональные разработчики не допускают глупых ошибок.
- Андрей Карпов. Ответы на вопросы, которые часто задают после прочтения наших статей.
- Новости о языке Си++, интересные статьи и информация о проверенных нами проектах: @Code_Analysis.
- Первое знакомство с анализатором: PVS-Studio для Visual C++.
ссылка на оригинал статьи http://habrahabr.ru/company/pvs-studio/blog/200938/
Добавить комментарий