これをFortranプログラムで書くと、第１段は
do j=1,jm
rho_h(j)=0.5d0*(rho(j)+rho(j+1))
* -0.5d0*dt/dx*(px(j+1)-px(j))
px_h(j)=0.5d0*(px(j)+px(j+1))
* -0.5d0*dt/dx*(p(j+1)-p(j)+px(j+1)*v(j+1)-px(j)*v(j))
e_h(j)=0.5d0*(e(j)+e(j+1))
* -0.5d0*dt/dx*((e(j+1)+p(j+1))*v(j+1)
* -(e(j)+p(j))*v(j))
end do
ここでhのついた量は
に関するものを表す。
また変数名は、普通に使われているものにしてあるが、pxはx方向の
運動量でpx=
vをあらわす。
これから、
に必要な他の量は次のようにして求められる。
do j=1,jm
v(j)=px_h(j)/rho_h(j)
p(j)=(gamma-1)*(e_h(j)-0.5d0*px_h(j)**2/rho_h(j))
end do
これらの量は、もう
の値は必要ないので、第１段のhつきの変数ではなく、
同じ変数に代入している。
第２段は、
do j=2,jm
rho(j)=rho(j)-dt/dx*(px_h(j)-px_h(j-1))
px(j)=px(j)-dt/dx*(p(j)-p(j-1)+px_h(j)*v(j)-px_h(j-1)*v(j-1))
e(j)=e(j)-dt/dx*((e_h(j)+p(j))*v(j)-(e_h(j-1)+p(j-1))*v(j-1))
end do
の様にして計算される。
これから、
に必要な他の量は次のようにして求められる。
do j=2,jm
v(j)=px(j)/rho(j)
p(j)=(gamma-1)*(e(j)-0.5d0*px(j)**2/rho(j))
end do
これを必要な回数繰り返すことによって時間経過が計算される。