FORTRAN90実行時エラーについて

はじめまして。

FORTRANを用いて微分方程式をルンゲクッタ法で解こうとしています。
しかし、実行時にエラーが発生してしまい、解消できずに困っています。
ソースコードを添付しますので、どのようにしたら解消されるかお分かりの方、教えてください。
よろしくお願いします。

◆症状
計算を実行すると、
f ( 0.430) = NaN
と表示される。

◆症状を解消するために試したこと
do while文の初回の計算のときに、0での除算が含まれることが原因と推測し、
fやgに0.1などの0以外の値を代入して計算を実行したが、解消されず。

◆ソースコード
implicit none

!==================
!DEFINE
!==================
real ita,theta,beta
real f,g,h
real k1,k2,k3,k4
real l1,l2,l3,l4
real m1,m2,m3,m4

!===================
!INITIALIZE
!===================
ita=0.0
theta=0.1 !step
beta=2.0
f=0.0
g=0.0
h=1.6

!===================
!SOLVE
!===================
do while(ita<5.0)

k1 = theta * g
l1 = theta * h
m1 = theta * (-(fh)-((1-(gg))*beta))

k2 = theta * (g+(l1/2))
l2 = theta * (h+(m1/2))
m2 = theta * (-(f+k1/2)(h+m1/2)-((1-(g+l1/2)(g+l1/2))*beta))

k3 = theta * (g+(l2/2))
l3 = theta * (h+(m2/2))
m3 = theta * (-(f+k2/2)(h+m2/2)-((1-(g+l2/2)(g+l2/2))*beta))

k4 = theta * (g+l3)
l4 = theta * (h+m3)
m4 = theta * (-(f+k3)(h+m3)-((1-(g+l3)(g+l3))*beta))

f = f + (k1+2k2+2k3+k4)/6
g = g + (l1+2l2+2l3+l4)/6
h = h + (m1+2m2+2m3+m4)/6

ita=ita+theta

open(1,file='hw.dat',status='replace')
write (1,*) ita,g
close(1)

end do

end

mie

2015/07/22 02:20

＞ do while文の初回の計算のときに、0での除算が含まれることが原因と推測し、除算は /2 と/6 しかないですね。０除算にはならないし、イテレーションもたかだか50回くらいなのでこの程度の計算ではオーバーフローにもならないと思います。ソースコードは添付のもので間違いないですか？＞ fやgに0.1などの0以外の値を代入して計算を実行したが、解消されず。ｆやｇは累積を計算しているので、初期値を意味なく変えるのは正しくありません。