【1階微分 ルンゲクッタ法にて】


csvファイル出力版はこちら
****** main ******
      external rk1,f1
      
      write(*,*)'h='
      read(*,*)h
      
      write(*,*)'t0='
      read(*,*)t0
      
      write(*,*)'y0='
      read(*,*)y0
      
      write(*,*)'tend='
      read(*,*)tend
      
      call rk1(f1,h,t0,y0,tend)
      stop
      end
      
****** subroutine program*******
      subroutine rk1(f1,h,t0,y0,tend)
      real*4 k1,k2,k3,k4
      t=t0
      y=y0
      
   20 continue
      if(t.gt.tend)go to 30
      k1= f1(y,t)
      k2= f1(y+h*k1/2,t+h/2)
      k3= f1(y+h*k2/2,t+h/2)
      k4= f1(y+h*k3,t+h)
      y = y+(h/6)*(k1+2*k2+2*k3+k4)
      t=t+h
      goto 20
      
   30 write(*,*)' y  =',y
      return
      end
      
****** function *******
      function f1(y,t)
      f1 = 1/y
      return
      end