step-rk4.f90

四阶RK方法计算单波方程

!!!!!!!!!!!!!!!!!!!!!!!!!
!!!    阶跃单波方程   !!!
!!!!!!!!!!!!!!!!!!!!!!!!!
program step
implicit none

!!! 输入网格点数
integer,parameter::N=301
!!! 输入迭代步数
integer,parameter::M=300
real k1(N),k2(N),k3(N),k4(N),u(N),u1(N)
integer i,j
real ot,ox,c
c=1.0
ox=0.01
ot=ox/c/2

!!!初始化
do i=1,N
if (i<=101)then
u(i)=1
else
u(i)=0
end if
end do

do i=1,M

  do j=2,N
    k1(j)=-c*(u(j)-u(j-1))/ox
  end do
  k1(1)=2*k1(2)-k1(3)

  do j=2,N
    k2=-(c*(u(j)-u(j-1))/ox+c*(k1(j)-k1(j-1))/ox*ot/2)
  end do
  k2(1)=2*k2(2)-k2(3)

  do j=2,N
    k3=-(c*(u(j)-u(j-1))/ox+c*(k2(j)-k2(j-1))/ox*ot/2)
  end do
  k3(1)=2*k3(2)-k3(3)

  do j=2,N
    k4=-(c*(u(j)-u(j-1))/ox+c*(k3(j)-k3(j-1))/ox*ot)
  end do
  k4(1)=2*k4(2)-k4(3)

  do j=1,N
    u1(j)=u(j)+(k1(j)+2*k2(j)+2*k3(j)+k4(j))*ot/6
  end do

  do j=1,N
    u(j)=u1(j)
  end do

write(*,*) i

end do

do i=1,N
write(*,*) i,u(i)
open(unit=10,file='out.txt')
write(10,*) u(i)
end do

stop
end