로렌츠 방정식

program mmm

real a(20000), b(20000), c(20000) 

!----------------------
!세 배열 a,b,c => x,y,z 좌표
!--------------------------

real rk1, rk2, rk3, rk4

real ak1, ak2, ak3, ak4

real ck1, ck2, ck3, ck4,h

integer i

a(1)=2

b(1)=1

c(1)=0.03

h=0.02


do i=1,15000,1

rk1=10.0*(b(i)-a(i))
rk2=10.0*(b(i)-a(i)-rk1*h/2.0)
rk3=10.0*(b(i)-a(i)-rk2*h/2.0)
rk4=10.0*(b(i)-a(i)-h*rk3)

ak1=a(i)*(28.0-c(i))-b(i)
ak2=a(i)*(28.0-c(i))-b(i)-h*ak1/2.0
ak3=a(i)*(28.0-c(i))-b(i)-h*ak2/2.0
ak4=a(i)*(28.0-c(i))-b(i)-h*ak3

ck1=a(i)*b(i)-8.0/3.0*c(i)
ck2=a(i)*b(i)-8.0/3.0*(c(i)+h*ck1/2.0)
ck3=a(i)*b(i)-8.0/3.0*(c(i)+h*ck2/2.0)
ck4=a(i)*b(i)-8.0/3.0*(c(i)+h*ck3)

a(i+1)=a(i)+(h/6.0)*(rk1+2.0*(rk2+rk3)+rk4)

b(i+1)=b(i)+(h/6.0)*(ak1+2.0*(ak2+ak3)+ak4)

c(i+1)=c(i)+(h/6.0)*(ck1+2.0*(ck2+ck3)+ck4)





end do

open (unit=5,file='lorenz.txt',status='replace')


do i=1,10000

write(5,*)a(i),b(i) ,c(i)

end do

close(5)








end program mmm