Mercurial > op > rk_sakamoti

comparison rkmod.f90 @ 3:2347cff808b0

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
chage a parameter of rkmod::calc
author "uncorrelated zombie" <uncorrelated@yahoo.co.jp>
date Sun, 21 Jul 2024 04:11:29 +0900
parents 206e90e91191
children
comparison
equal deleted inserted replaced
2:206e90e91191 3:2347cff808b0
14 rhside(1)=a*(wk(2)-wk(1)) 14 rhside(1)=a*(wk(2)-wk(1))
15 rhside(2)=b*wk(1)-wk(2)-wk(1)*wk(3) 15 rhside(2)=b*wk(1)-wk(2)-wk(1)*wk(3)
16 rhside(3)=wk(1)*wk(2)-c*wk(3) 16 rhside(3)=wk(1)*wk(2)-c*wk(3)
17 end function 17 end function
18 !Rから呼び出すサブルーチン 18 !Rから呼び出すサブルーチン
19 subroutine calc(x, max_t, r) 19 subroutine calc(x, max_i, r)
20 implicit none 20 implicit none
21 real(kind=dpkind),dimension(3),intent(in) :: x !状態変数 21 real(kind=dpkind),dimension(3),intent(in) :: x !状態変数
22 real(kind=dpkind),intent(in) :: max_t 22 integer, intent(in) :: max_i
23 real(kind=dpkind),allocatable,dimension(:,:),intent(out) :: r 23 real(kind=dpkind) :: max_t
24 real(kind=dpkind),dimension(max_i, 4),intent(out) :: r
24 real(kind=dpkind) :: zero = 0 25 real(kind=dpkind) :: zero = 0
25 integer i, max_i 26 integer i
26 !微分方程式の右辺を計算する関数、解法、時間刻みの指定 27 !微分方程式の右辺を計算する関数、解法、時間刻みの指定
27 call rktype%setup(rhside,"rk4",1d-5) 28 call rktype%setup(rhside,"rk4",1d-5)
28 max_i = int(max_t / 0.01) 29 max_t = max_i * 0.01
29 allocate(r(0:max_i, 4))
30 !時刻max_t秒まで、0.01秒毎に配列に記録。 30 !時刻max_t秒まで、0.01秒毎に配列に記録。
31 !ただし、内部的には刻み幅1e-5秒としている。 31 !ただし、内部的には刻み幅1e-5秒としている。
32 r(0, :) = [zero, x] ! 開始点 32 r(1, :) = [zero, x] ! 開始点
33 do i = 1, max_i 33 do i = 2, max_i
34 call rktype%solver(rktype%time+1e-2,x) !積分実行 34 call rktype%solver(rktype%time+1e-2,x) !積分実行
35 r(i, :) = [rktype%time, x] 35 r(i, :) = [rktype%time, x]
36 end do 36 end do
37 end subroutine 37 end subroutine
38 end module 38 end module