Как устранить ошибку сегментации при вызове подпрограммы в подпрограмме в FORTRAN 77?
Я борюсь с проблемой ошибки сегментации в FORTRAN 77 уже несколько месяцев, и я действительно не могу понять, где я ошибаюсь. Я новичок в Фортране.
Я пишу программу (main) в котором я вызываю подпрограмму (source). Эта подпрограмма выполняет метод Рунге-Кутты 4-го порядка (как предложено в " Численных рецептах" на Фортране 77 издательствами Press, Teukolsky, Vetterling, Flannery) для решения дифференциального уравнения и продвижения определенных значений входных данных от момента времени к следующему.
Чтобы выполнить метод Рунге-Кутты, в source другая подпрограмма deriv называется, который численно вычисляет производную в дифференциальном уравнении. Подпрограмма deriv должен быть вызван четыре раза для метода Рунге-Кутты, и я должен выполнить его для значений данных в обоих направлениях x и y: так, глобально, deriv называется восемь раз в течение source,
Значения для продвижения через Рунге-Кутта собраны в векторы максимального размера mxnn=200000 (параметр включен в param.fi), но в ходе программы только часть из них (размером nodes) будет заполнен и так продвинут. Эти векторы передаются source и их компоненты передаются deriv для метода Рунге-Кутта.
Это выдержка из main (это действительно большая программа, поэтому я вырезал то, что мне показалось ненужным):
PROGRAM main
INCLUDE 'param.fi'
INTEGER nodes
DOUBLEPRECISION z0
DOUBLEPRECISION deltaxnp1(mxnn),zxnp1(mxnn),u0xnp1(mxnn)
DOUBLEPRECISION deltaynp1(mxnn),zynp1(mxnn),u0ynp1(mxnn)
DOUBLEPRECISION deltaxn(mxnn),zxn(mxnn),u0xn(mxnn)
DOUBLEPRECISION deltayn(mxnn),zyn(mxnn),u0yn(mxnn)
...
...
CALL source(nodes,
+ z0,deltaxn,u0xnp1,u0xn,zxn,zxnp1,deltaxnp1,
+ deltayn,u0ynp1,u0yn,zyn,zynp1,deltaynp1,
+ [...])
...
...
END program main
Это подпрограмма source, который я написал в другом модуле:
SUBROUTINE source(nodes,
+ z0,deltaxin,u0xnew,u0xin,zx,zxout,deltaxout,
+ deltayin,u0ynew,u0yin,zy,zyout,deltayout,
+ [...])
INTEGER nodes, i
DOUBLEPRECISION z0
DOUBLEPRECISION deltaxin(*),u0xnew(*),u0xin(*)
DOUBLEPRECISION zx(*), zxout(*),deltaxout(*)
DOUBLEPRECISION deltayin(*),u0ynew(*)
DOUBLEPRECISION u0yin(*),zy(*),zyout(*),deltayout(*)
DOUBLEPRECISION duxdt(mxnn),duydt(mxnn)
*** Variables needed for Runge-Kutta ***
DOUBLEPRECISION dtmezzo,d6,th
DOUBLEPRECISION dzini(mxnn),dzt(mxnn),dzm(mxnn)
DOUBLEPRECISION zt(mxnn)
DO i=1,nodes
...
dtmezzo=dt*0.5D0
d6=dt/6.D0
th = time + dtmezzo
*** Runge-Kutta for values in x direction
*** INPUT: zx(i)
*** OUTPUT: zxout(i)
CALL deriv(zx(i),u0xnew(i),z0,duxdt(i),dzini(i))
zt(i)=zx(i)+dtmezzo*dzini(i)
CALL deriv(zt(i),u0xnew(i),z0,duxdt(i),dzt(i))
zt(i)=zx(i)+dtmezzo*dzt(i)
CALL deriv(zt(i),u0xnew(i),z0,duxdt(i),dzm(i))
zt(i)=zx(i)+dt*dzm(i)
dzm(i)=dzt(i)+dzm(i)
CALL deriv(zt(i),u0xnew(i),z0,duxdt(i),dzt(i))
zxout(i)=zx(i)+t6*(dzini(i)+dzt(i)+2.*dzm(i))
*** Runge-Kutta for values in y direction
*** INPUT: zy(i)
*** OUTPUT: zyout(i)
CALL deriv(zy(i),u0ynew(i),z0,duydt(i),dzini(i)) !!!!!
zt(i)=zy(i)+dtmezzo*dzyini(i)
CALL deriv(zt(i),u0ynew(i),z0,duydt(i),dzt(i))
zt(i)=zy(i)+dtmezzo*dzyt(i)
CALL deriv(zt(i),u0ynew(i),z0,duydt(i),dzm(i))
zt(i)=zy(i)+dt*dzm(i)
dzm(i)=dzt(i)+dzm(i)
CALL deriv(zt(i),u0ynew(i),z0,duydt(i),dzt(i))
zyout(i)=zy(i)+t6*(dzini(i)+dzt(i)+2.*dzm(i))
...
ENDDO
RETURN
END
А вот и подпрограмма deriv:
SUBROUTINE deriv(Z,u0,z0,dudt,der)
DOUBLEPRECISION Z,u0,z0,dudt,der
der = ( (((0.4D0**2)*u0)/z0) - Z*(dexp(Z)-Z-1)
+ * (dudt/u0) ) / ( dexp(Z)*(Z-1)+1 )
END
При компиляции программы (с GNUfortran) это не доставляет мне проблем: но после запуска программы Segmentation fault (core dumped) появляется ошибка и программа прерывается. После небольшого анализа кажется, что проблема лежит, когда я звоню deriv в первый раз для значений в направлении у (т. е. пятый вызов, который я пометил с !!!!! в коде).
Это странно для меня, так как, если я пытаюсь скомпилировать и запустить программу без блока направления y (решение Runge-Kutta только для направления x), то все идет хорошо. Форма deriv одинаково для обоих направлений: меняются только входы и выходы.
Я проверил размеры массива и ничего не кажется неправильным (по крайней мере, мне), и я боюсь, что в моем скрипте есть концептуальная ошибка. Я прошел через различные ответы здесь, в Переполнении Стека, но безрезультатно.
Я думал об отладке с dbg и -g флаг, но программа сделана через make-файл, и я не мог найти, как использовать флаг с ним.