Фантом
Буду очень признательна если Вы напишете как это сделать, а то я в хелпе как-то не нашла...

FCM
Специально для Вас провела ряд экспериментов. Буду писать сначала кусочек кода, потом результат:

*******************
rho3 = rho
print *, rho 3

   1.0240488621215573    
*******************
rho3 = rho
print *, rho 3, rho

 -1.15385900643022870E+137 -1.15385900643022870E+137
*******************
rho3 = rho
print *, rho1, rho 3, rho

   1.0000000000000000        1.0107194919764675        1.0107194919764675 
*******************
rho3 = rho
print *, rho3, rho1, rho 3, rho

  1.0107194919764675        1.0000000000000000        1.0107194919764675        1.0107194919764675 
*******************
rho3 = rho
print *, rho3, rho1

 -1.15385900643022870E+137   1.0000000000000000    
*******************
rho3 = rho
print *, rho1, rho3

1.0000000000000000   -1.15385900643022870E+137     
*******************
rho3 = rho
print *, rho3, rho3

   1.0240488621215573        1.0240488621215573 
*******************

Т.е. результат не зависит от повторений, от порядка переменных, а только от их наличия. 

Ну как, нравится?   

Нравится. Но для полного счастья нужна подробная информация о rho - возможно проблемы именно в нем.

Это сообщение отредактировал(а) FCM - 5.10.2011, 17:51

   real*8 function V(q)
    implicit real*8 (a-h,o-z)
    common /zn/Znucl
    common /screen/i_scr,iks,nwf_a,ka_a
    parameter (kMaxL = 9)
    parameter (NR1_g = 208)
    parameter (NR2_g = 2*(NR1_g+kMaxL))
    parameter (acl = 137.03599911d0)
    parameter (pi = 3.14159265358979323846)

    common /type_of_splines/type_of_splines
    common /wfgrid/break(NR1_g) /wfconst/nwfunc,na,kmax
    common /gaus/xx(64,13),cc(64,13) /Nst/Nstor(13)
    common /rn/rnucl
    common /inucl/inucl
    dimension r(0:114)
    dimension DI_B(0:10,0:10)

    do iii = 1,1                                             (это цикл включен для проверки в нескольких точках)
    q = 1.d-12 *10.d0**iii
    r_n = rnucl
    p = q
c    if (inucl.eq.1) then
c---------------------------c
c           Shell           c
c---------------------------c
          pr = q * r_n
          rho = dsin(pr) / pr
          rho1 = rho
c    elseif (inucl.eq.2) then
c---------------------------c
c           Sphere          c
c---------------------------c
    pr = q*r_n*dsqrt(5.d0/3.d0)
    if (pr.le.0.00001) then
    V1 = 1./3. - pr**2/30. + pr**4/840. - pr**6/45360.
    else
    V1 =  (dsin(pr)-pr*dcos(pr))/(pr)**3
     endif
    rho = V1*3.
    rho2 = rho
c    elseif (inucl.eq.3) then
c---------------------------c
c           Fermi           c
c---------------------------c
        ru_fermi = 52917.721d0 / 137.03599911d0
        af = 2.30d0 / 4d0 / dlog(3d0) / ru_fermi
        cf = dsqrt ( 5d0/3d0 * rnucl**2 - 7d0/3d0 * (pi*af)**2 )
    pac = pi * ac
    S_3_0 = S_k(-cf/af,3)
        S_5_0 = S_k(-cf/af,5)
    qf = 1d0 + pac**2 - 6d0 * ac**3 * S_3_0
          pc = q * cf
          pa = q * af
          spc = dsin(pc)
          cpc = dcos(pc)
          rho = 3.d0 * ( spc - pc * cpc ) / pc**3
          if (pc.lt.1d-5) then
            rho = 1.d0 - 0.1d0 * pc**2
          endif
          do m = 1, 1000
          tmp = 6d0 * ac**2 * (-1)**m / pc / ( m**2 + pa**2 )
     &        * ( ( m**2 - pa**2 ) / ( m**2 + pa**2 ) * spc + pc * cpc
     &            + m * pa / ( m**2 + pa**2 ) * dexp ( - m / ac ) )
            rho = rho - tmp 
          enddo
          rho = rho / qf
     rho3 = rho
     print *, rho3, rho3
c     print *, rho1, rho3, rho

        enddo
    stop

c    else
c    stop 'wrong nuclear model'
c     rho = 1.
c    endif
c    
    V = - Znucl/acl/q**2 * rho

    
 444    return
    end
c_____________________________________

Мне тут в чужом твитте довелось увидет дюже прелестный issue в одном баг-трекере.  

Говнокод настолько адский, что разбираться тут никто не будет. Видимо, даже компилятор не способен разобрать, что вы от него хотите. Единственное решение подсказано выше.  

IMPLICIT типизация, COMMON блоки, бросающиеся в глаза случаи некорректной арифметики с плавающей точкой (здесь parameter (pi = 3.14159265358979323846) точно должно быть d0 на конце), приправленные кодом из 70-х (фиксированный формат, великолепные названия переменных).

Хотя, даже интересно.  Поковырять что-ли... gfortran, говорите...

Фуф.  Жесть. Вот чем хорошо модульное программирование, о котором, Вы, боюсь ничего не знаете, так тем, что модули (хорошо написанные, безусловно) можно тестировать изолированно. А тут... Ну выложили Вы этот кусок кода? И? Толку то?

Итак, я поубирал все COMMON блоки, повставлял фиктивных значений, для переменных, заданных этими блоками, покоментил чутка.


У меня никакого "бага" при выводе естественно нет. Чтобы и в каком порядке я не выводил rho3 у меня всегда 2.2325875331982603. Просто "для галочки". 

Ещё раз - Вы несёте откровенную ахинею. Таких багов нет и быть не может.

Симптоматически можно предположить порчу памяти. НО! Если это "чистый" Fortran 77, то там, вроде никак не запортить. Однако исключать возможность использования всяких расширений (Cray Pointers) нельзя. Кроме того тут может быть ведь и Тандем. Не тот, который давече определился, а другой: сишечка + фортран. Но сразу после присвоения... Как-то с трудом верится.

Но код воняет. Не просто пахнет, воняет. Даже в этих 20 строчках кроме уже указанных недостатков, есть совершенно эпический цикл со счётчиком iii, где iii - число с плавающей точкой и оператор STOP, который не последний в теле функции. За такие вещи надо людей к стенке ставить...

Баг где-то есть. И не один, и не два. Их тонны в таком коде.

Nutsy, попробуй, как уже отметил kemiisto, везде, где встречаются floating point- литералы (т.е. буквально заданные вещественные константы) добавить суффикс d0.
На будущее - никогда не пользуйся неявной типизацией - всегда ставь в начале любой программной единицы инструкцию  IMPLICIT NONE. При неявной типизации невозможно на этапе компиляции отследить опечатки в именах объектов данных в программе.


Это сообщение отредактировал(а) FCM - 6.10.2011, 11:32

FCM
Спасибо за советы. Отказываться от неявной типизации, конечно, очень тяжко, но смысл этого действия вполне ясен. Про d0 тоже ясно, боюсь, тут мне нет оправданий - не ставила исключительно из-за лени. Может, будут еще какие-нибудь советы? Выслушаю внимательно и благодарно.

kemiisto
В Ваших сообщениях слишком много эмоций и слишком мало информации, чего можно ожидать скорее от человека, желающего самоутвердиться на форуме, чем от умного опытного програмиста, желающего помочь другим.  Я верю, что Вы очень умный и это получилось случайно как результат шока - поэтому не могли бы Вы написать поконкретнее?

Чем так плох код? (Не надо отвечать всем - это неинформативно)

Как названия переменных влияют на качество кода? Если есть какая-то связь, буду рада узнать - и тогда все исправлю. А то мне всегда казалось, что компьютеру все равно как что назвать, и это скорее для человека, а я свои названия понимаю хорошо - они все очень логичны (для меня). 

Чем можно заменить коммон блоки? Явная передача этих величин как параметров функции или подпрограммы невозможна.

Почему Вы решили что iii  - с плавающей точкой? Нам с компилятором казалось что оно вполне целое, где ошибка?

И что такого ужасного в операторе stop? Мне казалось, что для данной задачи (остановить исполнение кода на конкретном месте) он и придуман. Я не права? Есть альтернативы? В чем отличие альтернативных операторов от столь нелюбимого Вами? 

Я совершенно не сомневалась, что этот отдельный кусок, запущенный на другой машине не выдаст никакого аномального результата, потому что как бы неэлегантно не был написан код, он тем не менее не ошибочный. И проблема, к сожалению, глубже. Тем не менее я благодарна Вам за Ваш труд.

Это сообщение отредактировал(а) Nutsy - 6.10.2011, 14:37

Вы так говорите, как будто в желании самоутвердится есть что-то плохое.  А эмоции... Так они вообще и отличают нас от мерзких машин!


Я писал уже.


Это надо смотреть, зачем они используются. Единого рецепта тут нет.


  Что значит "невозможна"? Почему?


Потому что забыл латинский алфавит.  


Такой задачи возникать не должно. Поэтому и оператор STOP не нужен.


Я уже писал по этому поводу. Если Вас всё в коде устраивает, тогда зачем Вы здесь?

Сначала сделай с d0 и посмотри будут ли плавать соответствующие результаты.
Там еще в одном месте if (pr.le.0.00001) then
Вещественные числа обычно сравнивают только на строгое больше или меньше. 

обоже...

от них зависит читабельность. с именми навроде rho, rho1, V, V1 получается writeonly
ну да, главное чтобы компилер да ты всё понимали 
http://habrahabr.ru/blogs/htranslations/111348 пункт 2

FCM
Сделала в этом куске все с d0 и исправила сравнение как ты советовал. Не помогло. Теперь ползаю по всему коду, правлю d0 везде. Вернусь не скоро 

bems
Я вообще-то пыталась вежливо намекнуть одному из советующих, что его советы не по делу. Кажется, не удалось... А что касается названий - большинство из переменных определенные физические величины, и им всегда сообветствует одна и та же буква (в данном случае - rho).  Так что назвать ее по-другому точно было бы усложнением себе жизни 

Ссылка понравилась, спасибо  Пункт 2 представить сложно - потому как у нас все программы исключительно для индивидуального пользования - но я постараюсь. Должно помочь 

Еще вопрос для всех - а как же тогда правильно писать сложные выражения? Ну вот допустим, все что начинается с d - real*8, все что начинается с i - integer. Что из ниженаписанного правильно, а что нет?

dsum = d1 + 1.d0 + dble(i)
dsum = d2**i1
dsum = d2**d3

i1 = int(d1)
i2 = i1 + int (d1 + d2/d3 + d4**d5)

Надо ли всегда ставить dble  и int, или иногда в этом нет необходимости?

P.S. А может кто-то порекомендует ХОРОШУЮ книжку почитать?

Не всегда. В ряде случаев работают неявные преобразования. Но нужно понимать, что "иногда" (когда речь идет "вещественных" числах разной точности) они дают не то, чтобы хотелось.
Например, 


 
Насчет книжки. Если интересует именно сохранение точности в вычислениях, можно посмотреть соответсвующий параграф ("Точность вычислений") в <Немнюгин, Стесик "Современный Фортран"2004>. Возможно и в других книгах или интернет-источниках этот вопрос обсуждается.


Это сообщение отредактировал(а) FCM - 11.10.2011, 11:14

FCM
Спасибо, смысл поняла. 
За рекомендацию книжки спасибо. Вы не поверите, но один из авторов и был моим преподавателем по Фортрану ;) 


Проблема все еще со мной, но общее направление понятно. Вычищать все мелкие некорректности, проверять типы переменных и т.п. Встретила человека, который не удивился на мой рассказ, а сказал, что у него было то же самое. В чем была проблема тоже, к сожалению, не знает. 

Спасибо всем, кто старался мне помочь!