makris
06-29-2007, 08:06 AM
Hi all,
I am trying to solve the famous Lorenz 3x3 system of ODEs using rk4.f from NR.
I studied an example from the NR example book and then I tried to code my program.
Here is the problem
C................................................. .....
C
C NUMERICAL SOLUTION OF THE FOLLOWING 3x3 SYSTEM
C OF ODEs
C LORENZ
C
C DX
C ---- = 10 * ( Y - X )
C DT
C
C DY
C ---- = X * (28 -Z) - Y
C DT
C
C DZ
C ---- = X * Y - (8/3) * Z
C DT
C
C WITH INITIAL CONDITIONS
C X = 3
C 0
C
C Y = 15
C 0
C
C Z = 1
C 0
C................................................. .....
But when I run the following program the solution is not satisfcatory. Could you please give me a hint of what is goint wrong?
INTEGER N
PARAMETER(N=3)
INTEGER i,j
REAL h,x,y(N),dydx(N),yout(N)
EXTERNAL derivs
x=0.0
y(1)=3
y(2)=15
y(3)=1
call derivs(x,y,dydx)
WRITE(*,*) ' time X Y Z '
WRITE(*,*) '================================================= ===='
do 11 i=1,100
h=0.0001*i
call rk4(y,dydx,N,x,h,yout,derivs)
write(*,*) h , yout(1), yout(2), yout(3)
11 continue
END
SUBROUTINE derivs(x,y,dydx)
REAL x,y(*),dydx(*)
dydx(1)= 10.0 * ( y(2) - y(1) )
dydx(2)= y(1) *( 28.0 - y(3) ) - y(2)
dydx(3)= y(1) * y(2) - (2.666666666*y(3))
return
END
Thanks
Makris
I am trying to solve the famous Lorenz 3x3 system of ODEs using rk4.f from NR.
I studied an example from the NR example book and then I tried to code my program.
Here is the problem
C................................................. .....
C
C NUMERICAL SOLUTION OF THE FOLLOWING 3x3 SYSTEM
C OF ODEs
C LORENZ
C
C DX
C ---- = 10 * ( Y - X )
C DT
C
C DY
C ---- = X * (28 -Z) - Y
C DT
C
C DZ
C ---- = X * Y - (8/3) * Z
C DT
C
C WITH INITIAL CONDITIONS
C X = 3
C 0
C
C Y = 15
C 0
C
C Z = 1
C 0
C................................................. .....
But when I run the following program the solution is not satisfcatory. Could you please give me a hint of what is goint wrong?
INTEGER N
PARAMETER(N=3)
INTEGER i,j
REAL h,x,y(N),dydx(N),yout(N)
EXTERNAL derivs
x=0.0
y(1)=3
y(2)=15
y(3)=1
call derivs(x,y,dydx)
WRITE(*,*) ' time X Y Z '
WRITE(*,*) '================================================= ===='
do 11 i=1,100
h=0.0001*i
call rk4(y,dydx,N,x,h,yout,derivs)
write(*,*) h , yout(1), yout(2), yout(3)
11 continue
END
SUBROUTINE derivs(x,y,dydx)
REAL x,y(*),dydx(*)
dydx(1)= 10.0 * ( y(2) - y(1) )
dydx(2)= y(1) *( 28.0 - y(3) ) - y(2)
dydx(3)= y(1) * y(2) - (2.666666666*y(3))
return
END
Thanks
Makris