Skip to content

Commit

Permalink
Introduced common subroutine for adaptive integrators
Browse files Browse the repository at this point in the history
  • Loading branch information
jcandy committed Sep 17, 2024
1 parent 0152e30 commit f8d6711
Show file tree
Hide file tree
Showing 4 changed files with 43 additions and 57 deletions.
33 changes: 31 additions & 2 deletions cgyro/src/cgyro_step.f90
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ module cgyro_step
real, parameter :: eps = 2e-12

integer :: itrk,istep,conv

! adaptive timesteps
real :: delta_t_gk
real :: delta_t_last,delta_t_last_step
Expand All @@ -34,6 +34,35 @@ module cgyro_step

real, dimension(2) :: error_sum,error_x

logical :: is_first

end module cgyro_step


subroutine get_delta_t(c0)

use cgyro_globals
use cgyro_step

implicit none

real, intent(in) :: c0

delta_t_gk = max(delta_t_last,c0*deltah2)

if (delta_t_last_step <= small) delta_t_last_step = delta_t_last

if (delta_t_last_step < 0.1*delta_t_last) then
delta_t_gk = delta_t_last+delta_t_last_step
else
if (delta_t_last_step/itrk < 0.1*delta_t_last) then
delta_t_gk = delta_t_last+delta_t_last_step/itrk
endif
endif

delta_t_gk = min(delta_t,delta_t_gk)
total_local_error = var_error

end subroutine get_delta_t



22 changes: 4 additions & 18 deletions cgyro/src/cgyro_step_gk_bs5.f90
Original file line number Diff line number Diff line change
Expand Up @@ -67,8 +67,6 @@ subroutine cgyro_step_gk_bs5
real, parameter :: e5 = 243.d0/12800.d0
real, parameter :: e6 = -1.d0/95.d0

logical :: is_first

tol = error_tol

itrk = 0
Expand Down Expand Up @@ -195,7 +193,9 @@ subroutine cgyro_step_gk_bs5
!------------------

call timer_lib_in('str')
! using a multiplication by 0 in one element is still efffienct, since the matrix element was read for the 1st equation
! using a multiplication by 0 in one element is still efficient, since
! the matrix element was read for the 2nd equation

call cgyro_vel_solution_werror(5, h_x, &
h0_x, &
deltah2*b1, rhs(:,:,:,1), &
Expand Down Expand Up @@ -253,20 +253,6 @@ subroutine cgyro_step_gk_bs5

! caller expects the fields to be updated on exit
call cgyro_field_c(.TRUE.)
call get_delta_t(4.0/5.0)

delta_t_gk = max(delta_t_last,4.0/5.0*deltah2)

if (delta_t_last_step <= small) delta_t_last_step = delta_t_last

if (delta_t_last_step < 0.1*delta_t_gk) then
delta_t_gk = delta_t_last+delta_t_last_step
else
if (delta_t_last_step/itrk < 0.1*delta_t_gk) then
delta_t_gk = delta_t_gk + delta_t_last_step/itrk
endif
endif

delta_t_gk = min(delta_t,delta_t_gk)
total_local_error = var_error

end subroutine cgyro_step_gk_bs5
23 changes: 4 additions & 19 deletions cgyro/src/cgyro_step_gk_ck.f90
Original file line number Diff line number Diff line change
Expand Up @@ -24,8 +24,6 @@ subroutine cgyro_step_gk_ck
! Apar -> field(2)
! Bpar -> field(3)

logical :: is_first

tol = error_tol

! total iteration counter
Expand Down Expand Up @@ -153,7 +151,9 @@ subroutine cgyro_step_gk_ck
!-------------------

call timer_lib_in('str')
! using a multiplication by 0 in one element is still efffienct, since the matrix element was read for the 2nd equation
! using a multiplication by 0 in one element is still efficient, since
! the matrix element was read for the 2nd equation

call cgyro_vel_solution_werror(4, h_x, &
h0_x, &
deltah2*( 37.d0/ 378.d0), rhs(:,:,:,1), &
Expand Down Expand Up @@ -221,21 +221,6 @@ subroutine cgyro_step_gk_ck

! caller expects the fields to be updated on exit
call cgyro_field_c(.TRUE.)

delta_t_gk = max(delta_t_last,deltah2)

if (delta_t_last_step <= small) delta_t_last_step = delta_t_last

if (delta_t_last_step < 0.1*delta_t_gk) then
delta_t_gk = delta_t_last+delta_t_last_step
else
if (delta_t_last_step/itrk < 0.1*delta_t_gk) then
delta_t_gk = delta_t_gk+delta_t_last_step/itrk
endif
endif

delta_t_gk = min(delta_t,delta_t_gk)

total_local_error = var_error
call get_delta_t(1.0)

end subroutine cgyro_step_gk_ck
22 changes: 4 additions & 18 deletions cgyro/src/cgyro_step_gk_v76.f90
Original file line number Diff line number Diff line change
Expand Up @@ -89,8 +89,6 @@ subroutine cgyro_step_gk_v76
real, parameter :: b9h = 0.d0
real, parameter :: b10h = 0.20295184663356282227670547938d-1

logical :: is_first

tol = error_tol

itrk = 0
Expand Down Expand Up @@ -254,7 +252,9 @@ subroutine cgyro_step_gk_v76
!-------------------

call timer_lib_in('str')
! using a multiplication by 0 in one element is still efffienct, since the matrix element was read for the 2nd equation
! using a multiplication by 0 in one element is still efficient, since
! the matrix element was read for the 2nd equation

call cgyro_vel_solution_werror(7, h_x, &
h0_x, &
deltah2*b1, rhs(:,:,:,1), &
Expand Down Expand Up @@ -316,20 +316,6 @@ subroutine cgyro_step_gk_v76

! caller expects the fields to be updated on exit
call cgyro_field_c(.TRUE.)

delta_t_gk = max(delta_t_last,6.0/7.0*deltah2)
call get_delta_t(6.0/7.0)

if (delta_t_last_step <= small) delta_t_last_step = delta_t_last

if (delta_t_last_step < 0.1*delta_t_last) then
delta_t_gk = delta_t_last+delta_t_last_step
else
if (delta_t_last_step/itrk < 0.1*delta_t_last) then
delta_t_gk = delta_t_last+delta_t_last_step/itrk
endif
endif

delta_t_gk = min(delta_t,delta_t_gk)
total_local_error = var_error

end subroutine cgyro_step_gk_v76

0 comments on commit f8d6711

Please sign in to comment.