Remove RKCK (time_stepper == 4)

Author: Diego Vaca

src/common/m_checker_common.fpp

@@ -62,7 +62,7 @@ contains
         !! Called by s_check_inputs_common for simulation and post-processing
     subroutine s_check_inputs_time_stepping
         if (cfl_dt) then
-            @:PROHIBIT((cfl_target < 0 .or. cfl_target > 1._wp) .and. .not. rkck_adap_dt)
+            @:PROHIBIT(cfl_target < 0 .or. cfl_target > 1._wp)
             @:PROHIBIT(t_stop <= 0)
             @:PROHIBIT(t_save <= 0)
             @:PROHIBIT(t_save > t_stop)

src/common/m_constants.fpp

@@ -49,35 +49,6 @@ module m_constants
     integer, parameter :: mapCells = 3 !< Number of cells around the bubble where the smoothening function will have effect
     real(wp), parameter :: R_uni = 8314._wp ! Universal gas constant - J/kmol/K
 
-    ! RKCK constants
-    integer, parameter :: num_ts_rkck = 6 !< Number of time-stages in the RKCK stepper
-    ! RKCK 4th/5th time stepper coefficients based on Cash J. and Karp A. (1990)
-    real(wp), parameter :: rkck_c1 = 0._wp, rkck_c2 = 0.2_wp, rkck_c3 = 0.3_wp, rkck_c4 = 0.6_wp, &     ! c1 c2 c3 c4 c5 c6
-                           rkck_c5 = 1._wp, rkck_c6 = 0.875_wp
-    real(wp), dimension(6), parameter :: rkck_coef1 = (/0.2_wp, 0._wp, 0._wp, 0._wp, 0._wp, 0._wp/)     ! a21
-    real(wp), dimension(6), parameter :: rkck_coef2 = (/3._wp/40._wp, 9._wp/40._wp, 0._wp, 0._wp, &     ! a31 a32
-                                                        0._wp, 0._wp/)
-    real(wp), dimension(6), parameter :: rkck_coef3 = (/0.3_wp, -0.9_wp, 1.2_wp, 0._wp, 0._wp, 0._wp/)  ! a41 a42 a43
-    real(wp), dimension(6), parameter :: rkck_coef4 = (/-11._wp/54._wp, 2.5_wp, -70._wp/27._wp, &       ! a51 a52 a53 a54
-                                                        35._wp/27._wp, 0._wp, 0._wp/)
-    real(wp), dimension(6), parameter :: rkck_coef5 = (/1631._wp/55296._wp, 175._wp/512._wp, &          ! a61 a62 a63 a64 a65
-                                                        575._wp/13824._wp, 44275._wp/110592._wp, &
-                                                        253._wp/4096._wp, 0._wp/)
-    real(wp), dimension(6), parameter :: rkck_coef6 = (/37._wp/378._wp, 0._wp, 250._wp/621._wp, &       ! b1 b2 b3 b4 b5 b6
-                                                        125._wp/594._wp, 0._wp, 512._wp/1771._wp/)
-    real(wp), dimension(6), parameter :: rkck_coefE = (/37._wp/378._wp - 2825._wp/27648._wp, 0._wp, &   ! er1 er2 er3 er4 er5 er6 (4th/5th error)
-                                                        250._wp/621._wp - 18575._wp/48384._wp, &
-                                                        125._wp/594._wp - 13525._wp/55296._wp, &
-                                                        -277._wp/14336._wp, 512._wp/1771._wp - 0.25_wp/)
-    ! Adaptive rkck constants
-    real(wp), parameter :: verysmall_dt = 1e-14_wp !< Very small dt, stop stepper
-    real(wp), parameter :: SAFETY = 0.9_wp !< Next dt will be maximum 0.9*dt if truncation error is above tolerance.
-    real(wp), parameter :: RNDDEC = 1e8_wp !< Round calculated dt (16th digit) to avoid the inclusion of random decimals
-    real(wp), parameter :: PSHRNK = -0.25_wp !< Factor to reduce dt when truncation error above tolerance
-    real(wp), parameter :: SHRNKDT = 0.5_wp !< Factor to reduce dt due to negative bubble radius
-    real(wp), parameter :: ERRCON = 1.89e-4_wp !< Limit to slightly increase dt when truncation error is between ERRCON and 1
-    real(wp), parameter :: PGROW = -0.2_wp !< Factor to increase dt when truncation error is between ERRCON and 1
-
     ! Relativity
     integer, parameter :: relativity_cons_to_prim_max_iter = 100
 

src/post_process/m_global_parameters.fpp

@@ -314,7 +314,7 @@ module m_global_parameters
 
     !> @name Lagrangian bubbles
     !> @{
-    logical :: bubbles_lagrange, rkck_adap_dt
+    logical :: bubbles_lagrange
     !> @}
 
     real(wp) :: Bx0 !< Constant magnetic field in the x-direction (1D)
@@ -443,7 +443,6 @@ contains
 
         ! Lagrangian bubbles modeling
         bubbles_lagrange = .false.
-        rkck_adap_dt = .false.
 
         ! IBM
         num_ibs = dflt_int

src/post_process/m_mpi_proxy.fpp

@@ -170,8 +170,7 @@ contains
             & 'file_per_process', 'relax', 'cf_wrt',                           &
             & 'adv_n', 'ib', 'cfl_adap_dt', 'cfl_const_dt', 'cfl_dt',          &
             & 'surface_tension', 'hyperelasticity', 'bubbles_lagrange',        &
-            & 'rkck_adap_dt', 'output_partial_domain', 'relativity',           &
-            & 'cont_damage' ]
+            & 'output_partial_domain', 'relativity', 'cont_damage' ]
             call MPI_BCAST(${VAR}$, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr)
         #:endfor
 

src/post_process/m_start_up.f90

@@ -83,7 +83,7 @@ subroutine s_read_input_file
             polydisperse, poly_sigma, file_per_process, relax, &
             relax_model, cf_wrt, sigma, adv_n, ib, num_ibs, &
             cfl_adap_dt, cfl_const_dt, t_save, t_stop, n_start, &
-            cfl_target, surface_tension, bubbles_lagrange, rkck_adap_dt, &
+            cfl_target, surface_tension, bubbles_lagrange, &
             sim_data, hyperelasticity, Bx0, relativity, cont_damage
 
         ! Inquiring the status of the post_process.inp file
@@ -113,7 +113,7 @@ subroutine s_read_input_file
 
             nGlobal = (m_glb + 1)*(n_glb + 1)*(p_glb + 1)
 
-            if (cfl_adap_dt .or. cfl_const_dt .or. rkck_adap_dt) cfl_dt = .true.
+            if (cfl_adap_dt .or. cfl_const_dt) cfl_dt = .true.
 
         else
             call s_mpi_abort('File post_process.inp is missing. Exiting.')

src/pre_process/m_global_parameters.fpp

@@ -261,10 +261,6 @@ module m_global_parameters
     integer :: chemxb, chemxe
     !> @}
 
-    !> @ lagrangian solver parameters
-    logical :: rkck_adap_dt
-    !> @}
-
     integer, allocatable, dimension(:, :, :) :: logic_grid
 
     type(pres_field) :: pb
@@ -526,9 +522,6 @@ contains
             fluid_pp(i)%G = 0._wp
         end do
 
-        ! Lagrangian solver
-        rkck_adap_dt = .false.
-
         Bx0 = dflt_real
 
     end subroutine s_assign_default_values_to_user_inputs

src/pre_process/m_start_up.fpp

@@ -143,8 +143,8 @@ contains
             palpha_eps, ptgalpha_eps, ib, num_ibs, patch_ib, &
             sigma, adv_n, cfl_adap_dt, cfl_const_dt, n_start, &
             n_start_old, surface_tension, hyperelasticity, pre_stress, &
-            rkck_adap_dt, elliptic_smoothing, elliptic_smoothing_iters, &
-            viscous, bubbles_lagrange, Bx0, relativity, cont_damage
+            elliptic_smoothing, elliptic_smoothing_iters, bubbles_lagrange, &
+            viscous, Bx0, relativity, cont_damage
 
         ! Inquiring the status of the pre_process.inp file
         file_loc = 'pre_process.inp'
@@ -171,7 +171,7 @@ contains
 
             nGlobal = (m_glb + 1)*(n_glb + 1)*(p_glb + 1)
 
-            if (cfl_adap_dt .or. cfl_const_dt .or. rkck_adap_dt) cfl_dt = .true.
+            if (cfl_adap_dt .or. cfl_const_dt) cfl_dt = .true.
 
         else
             call s_mpi_abort('File pre_process.inp is missing. Exiting.')

src/simulation/m_acoustic_src.fpp

@@ -158,7 +158,6 @@ contains
         integer, parameter :: mass_label = 1, mom_label = 2
 
         sim_time = t_step*dt
-        if (time_stepper == 4) sim_time = rkck_time_tmp ! Probably create a time_stepper == 5 for the rkck stepper
 
         !$acc parallel loop collapse(3) gang vector default(present)
         do l = 0, p