Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
32 changes: 23 additions & 9 deletions config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90
Original file line number Diff line number Diff line change
Expand Up @@ -112,8 +112,10 @@ module MOM_surface_forcing_gfdl
!! salinity to a specified value.
logical :: restore_temp !< If true, the coupled MOM driver adds a term to restore sea
!! surface temperature to a specified value.
real :: Flux_const_salt !< Piston velocity for surface salt restoring [Z T-1 ~> m s-1]
real :: Flux_const_temp !< Piston velocity for surface temp restoring [Z T-1 ~> m s-1]
real :: Flux_const_salt !< Piston velocity for surface salinity restoring [Z T-1 ~> m s-1]
real :: Flux_const_temp !< Piston velocity for surface temperature restoring [Z T-1 ~> m s-1]
real :: rho_restore !< The density that is used to convert piston velocities into salt
!! or heat fluxes with salinity or temperature restoring [R ~> kg m-3]
logical :: trestore_SPEAR_ECDA !< If true, modify restoring data wrt local SSS
real :: SPEAR_dTf_dS !< The derivative of the freezing temperature with
!! salinity [C S-1 ~> degC ppt-1].
Expand Down Expand Up @@ -268,7 +270,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
isr = is-isd+1 ; ier = ie-isd+1 ; jsr = js-jsd+1 ; jer = je-jsd+1

kg_m2_s_conversion = US%kg_m2s_to_RZ_T
if (CS%restore_temp) rhoXcp = CS%Rho0 * fluxes%C_p
if (CS%restore_temp) rhoXcp = CS%rho_restore * fluxes%C_p
open_ocn_mask(:,:) = 1.0
fluxes%vPrecGlobalAdj = 0.0
fluxes%vPrecGlobalScl = 0.0
Expand All @@ -281,7 +283,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
! flux type has been used.
if (fluxes%dt_buoy_accum < 0) then
call allocate_forcing_type(G, fluxes, water=.true., heat=.true., ustar=.true., press=.true., &
fix_accum_bug=CS%fix_ustar_gustless_bug)
fix_accum_bug=CS%fix_ustar_gustless_bug, tau_mag=.true.)

call safe_alloc_ptr(fluxes%sw_vis_dir,isd,ied,jsd,jed)
call safe_alloc_ptr(fluxes%sw_vis_dif,isd,ied,jsd,jed)
Expand Down Expand Up @@ -363,7 +365,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
do j=js,je ; do i=is,ie
delta_sss = data_restore(i,j) - sfc_state%SSS(i,j)
delta_sss = sign(1.0,delta_sss) * min(abs(delta_sss), CS%max_delta_srestore)
fluxes%salt_flux(i,j) = 1.e-3*US%S_to_ppt*G%mask2dT(i,j) * (CS%Rho0*CS%Flux_const_salt)* &
fluxes%salt_flux(i,j) = 1.e-3*US%S_to_ppt*G%mask2dT(i,j) * (CS%rho_restore*CS%Flux_const_salt)* &
(CS%basin_mask(i,j)*open_ocn_mask(i,j)*CS%srestore_mask(i,j)) * delta_sss ! R Z T-1 ~> kg Salt m-2 s-1
enddo ; enddo
if (CS%adjust_net_srestore_to_zero) then
Expand All @@ -386,7 +388,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
delta_sss = sfc_state%SSS(i,j) - data_restore(i,j)
delta_sss = sign(1.0,delta_sss) * min(abs(delta_sss), CS%max_delta_srestore)
fluxes%vprec(i,j) = (CS%basin_mask(i,j)*open_ocn_mask(i,j)*CS%srestore_mask(i,j))* &
(CS%Rho0*CS%Flux_const_salt) * &
(CS%rho_restore*CS%Flux_const_salt) * &
delta_sss / (0.5*(sfc_state%SSS(i,j) + data_restore(i,j)))
endif
enddo ; enddo
Expand Down Expand Up @@ -717,7 +719,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
! mechanical forcing type has been used.
if (.not.forces%initialized) then
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., &
press=.true.)
press=.true., tau_mag=.true.)

call safe_alloc_ptr(forces%p_surf,isd,ied,jsd,jed)
call safe_alloc_ptr(forces%p_surf_full,isd,ied,jsd,jed)
Expand Down Expand Up @@ -1276,6 +1278,8 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, wind_stagger)
! Local variables
real :: utide ! The RMS tidal velocity [Z T-1 ~> m s-1].
real :: Flux_const_dflt ! A default piston velocity for restoring surface properties [m day-1]
real :: rho_TKE_tidal ! The constant bottom density used to translate tidal amplitudes into the
! tidal bottom TKE input used with INT_TIDE_DISSIPATION [R ~> kg m-3]
logical :: new_sim ! False if this simulation was started from a restart file
! or other equivalent files.
logical :: iceberg_flux_diags ! If true, diagnostics of fluxes from icebergs are available.
Expand Down Expand Up @@ -1501,6 +1505,11 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, wind_stagger)
"The derivative of the freezing temperature with salinity.", &
units="deg C PSU-1", default=-0.054, scale=US%degC_to_C*US%S_to_ppt, &
do_not_log=.not.CS%trestore_SPEAR_ECDA)
call get_param(param_file, mdl, "RESTORE_FLUX_RHO", CS%rho_restore, &
"The density that is used to convert piston velocities into salt or heat "//&
"fluxes with RESTORE_SALINITY or RESTORE_TEMPERATURE.", &
units="kg m-3", default=CS%Rho0*US%R_to_kg_m3, scale=US%kg_m3_to_R, &
do_not_log=.not.(CS%restore_temp.or.CS%restore_salt))

! Optionally read tidal amplitude from input file [Z T-1 ~> m s-1] on model grid.
! Otherwise use default tidal amplitude for bottom frictionally-generated
Expand All @@ -1525,6 +1534,11 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, wind_stagger)
"The constant tidal amplitude used with INT_TIDE_DISSIPATION.", &
units="m s-1", default=0.0, scale=US%m_to_Z*US%T_to_s)
endif
call get_param(param_file, mdl, "TKE_TIDAL_RHO", rho_TKE_tidal, &
"The constant bottom density used to translate tidal amplitudes into the tidal "//&
"bottom TKE input used with INT_TIDE_DISSIPATION.", &
units="kg m-3", default=CS%Rho0*US%R_to_kg_m3, scale=US%kg_m3_to_R, &
do_not_log=.not.(CS%read_TIDEAMP.or.(CS%utide>0.0)))

call safe_alloc_ptr(CS%TKE_tidal,isd,ied,jsd,jed)
call safe_alloc_ptr(CS%ustar_tidal,isd,ied,jsd,jed)
Expand All @@ -1537,13 +1551,13 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, wind_stagger)
rescale=US%m_to_Z*US%T_to_s)
do j=jsd, jed; do i=isd, ied
utide = CS%TKE_tidal(i,j)
CS%TKE_tidal(i,j) = G%mask2dT(i,j)*CS%Rho0*CS%cd_tides*(utide*utide*utide)
CS%TKE_tidal(i,j) = G%mask2dT(i,j)*rho_TKE_tidal*CS%cd_tides*(utide*utide*utide)
CS%ustar_tidal(i,j) = sqrt(CS%cd_tides)*utide
enddo ; enddo
else
do j=jsd,jed; do i=isd,ied
utide = CS%utide
CS%TKE_tidal(i,j) = CS%Rho0*CS%cd_tides*(utide*utide*utide)
CS%TKE_tidal(i,j) = rho_TKE_tidal*CS%cd_tides*(utide*utide*utide)
CS%ustar_tidal(i,j) = sqrt(CS%cd_tides)*utide
enddo ; enddo
endif
Expand Down
4 changes: 2 additions & 2 deletions config_src/drivers/mct_cap/mom_surface_forcing_mct.F90
Original file line number Diff line number Diff line change
Expand Up @@ -276,7 +276,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
! flux type has been used.
if (fluxes%dt_buoy_accum < 0) then
call allocate_forcing_type(G, fluxes, water=.true., heat=.true., ustar=.true., &
press=.true., fix_accum_bug=CS%fix_ustar_gustless_bug)
press=.true., fix_accum_bug=CS%fix_ustar_gustless_bug, tau_mag=.true.)

call safe_alloc_ptr(fluxes%sw_vis_dir,isd,ied,jsd,jed)
call safe_alloc_ptr(fluxes%sw_vis_dif,isd,ied,jsd,jed)
Expand Down Expand Up @@ -649,7 +649,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
! mechanical forcing type has been used.
if (.not.forces%initialized) then

call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true.)
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true., tau_mag=.true.)

call safe_alloc_ptr(forces%p_surf,isd,ied,jsd,jed)
call safe_alloc_ptr(forces%p_surf_full,isd,ied,jsd,jed)
Expand Down
4 changes: 2 additions & 2 deletions config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90
Original file line number Diff line number Diff line change
Expand Up @@ -308,7 +308,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
if (fluxes%dt_buoy_accum < 0) then
call allocate_forcing_type(G, fluxes, water=.true., heat=.true., ustar=.true., &
press=.true., fix_accum_bug=CS%fix_ustar_gustless_bug, &
cfc=CS%use_CFC, hevap=CS%enthalpy_cpl)
cfc=CS%use_CFC, hevap=CS%enthalpy_cpl, tau_mag=.true.)

call safe_alloc_ptr(fluxes%sw_vis_dir,isd,ied,jsd,jed)
call safe_alloc_ptr(fluxes%sw_vis_dif,isd,ied,jsd,jed)
Expand Down Expand Up @@ -716,7 +716,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
! mechanical forcing type has been used.
if (.not.forces%initialized) then

call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true.)
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true., tau_mag=.true.)

call safe_alloc_ptr(forces%p_surf,isd,ied,jsd,jed)
call safe_alloc_ptr(forces%p_surf_full,isd,ied,jsd,jed)
Expand Down
16 changes: 11 additions & 5 deletions config_src/drivers/solo_driver/MESO_surface_forcing.F90
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ module MESO_surface_forcing
use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe
use MOM_file_parser, only : get_param, log_version, param_file_type
use MOM_forcing_type, only : forcing, mech_forcing
use MOM_forcing_type, only : allocate_forcing_type, allocate_mech_forcing
use MOM_forcing_type, only : allocate_forcing_type
use MOM_grid, only : ocean_grid_type
use MOM_io, only : file_exists, MOM_read_data, slasher
use MOM_time_manager, only : time_type, operator(+), operator(/)
Expand All @@ -30,6 +30,8 @@ module MESO_surface_forcing
real :: Rho0 !< The density used in the Boussinesq approximation [R ~> kg m-3].
real :: G_Earth !< The gravitational acceleration [L2 Z-1 T-2 ~> m s-2].
real :: Flux_const !< The restoring rate at the surface [Z T-1 ~> m s-1].
real :: rho_restore !< The density that is used to convert piston velocities into salt
!! or heat fluxes with salinity or temperature restoring [R ~> kg m-3]
real :: gust_const !< A constant unresolved background gustiness
!! that contributes to ustar [R L Z T-2 ~> Pa]
real, dimension(:,:), pointer :: &
Expand Down Expand Up @@ -166,14 +168,14 @@ subroutine MESO_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
! call MOM_error(FATAL, "MESO_buoyancy_surface_forcing: " // &
! "Temperature and salinity restoring used without modification." )

rhoXcp = CS%Rho0 * fluxes%C_p
rhoXcp = CS%rho_restore * fluxes%C_p
do j=js,je ; do i=is,ie
! Set Temp_restore and Salin_restore to the temperature (in degC) and
! salinity (in ppt or PSU) that are being restored toward.
if (G%mask2dT(i,j) > 0.0) then
fluxes%heat_added(i,j) = G%mask2dT(i,j) * &
((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const)
fluxes%vprec(i,j) = - (CS%Rho0 * CS%Flux_const) * &
fluxes%vprec(i,j) = - (CS%rho_restore * CS%Flux_const) * &
(CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / &
(0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)))
else
Expand All @@ -188,7 +190,7 @@ subroutine MESO_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
"Buoyancy restoring used without modification." )

! The -1 is because density has the opposite sign to buoyancy.
buoy_rest_const = -1.0 * (CS%G_Earth * CS%Flux_const) / CS%Rho0
buoy_rest_const = -1.0 * (CS%G_Earth * CS%Flux_const) / CS%rho_restore
do j=js,je ; do i=is,ie
! Set density_restore to an expression for the surface potential
! density [R ~> kg m-3] that is being restored toward.
Expand Down Expand Up @@ -272,7 +274,11 @@ subroutine MESO_surface_forcing_init(Time, G, US, param_file, diag, CS)
"variable NET_SOL.", fail_if_missing=.true.)
call get_param(param_file, mdl, "INPUTDIR", CS%inputdir, default=".")
CS%inputdir = slasher(CS%inputdir)

call get_param(param_file, mdl, "RESTORE_FLUX_RHO", CS%rho_restore, &
"The density that is used to convert piston velocities into salt or heat "//&
"fluxes with RESTORE_SALINITY or RESTORE_TEMPERATURE.", &
units="kg m-3", default=CS%Rho0*US%R_to_kg_m3, scale=US%kg_m3_to_R, &
do_not_log=(CS%Flux_const==0.0).or.(.not.CS%restorebuoy))
endif

end subroutine MESO_surface_forcing_init
Expand Down
38 changes: 23 additions & 15 deletions config_src/drivers/solo_driver/MOM_surface_forcing.F90
Original file line number Diff line number Diff line change
Expand Up @@ -79,9 +79,11 @@ module MOM_surface_forcing

real :: Rho0 !< Boussinesq reference density [R ~> kg m-3]
real :: G_Earth !< gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
real :: Flux_const !< piston velocity for surface restoring [Z T-1 ~> m s-1]
real :: Flux_const_T !< piston velocity for surface temperature restoring [Z T-1 ~> m s-1]
real :: Flux_const_S !< piston velocity for surface salinity restoring [Z T-1 ~> m s-1]
real :: Flux_const = 0.0 !< piston velocity for surface restoring [Z T-1 ~> m s-1]
real :: Flux_const_T = 0.0 !< piston velocity for surface temperature restoring [Z T-1 ~> m s-1]
real :: Flux_const_S = 0.0 !< piston velocity for surface salinity restoring [Z T-1 ~> m s-1]
real :: rho_restore !< The density that is used to convert piston velocities into salt
!! or heat fluxes with salinity or temperature restoring [R ~> kg m-3]
real :: latent_heat_fusion !< latent heat of fusion times [Q ~> J kg-1]
real :: latent_heat_vapor !< latent heat of vaporization [Q ~> J kg-1]
real :: tau_x0 !< Constant zonal wind stress used in the WIND_CONFIG="const"
Expand Down Expand Up @@ -250,9 +252,9 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US

if (CS%first_call_set_forcing) then
! Allocate memory for the mechanical and thermodynamic forcing fields.
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true.)
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true., tau_mag=.true.)

call allocate_forcing_type(G, fluxes, ustar=.true., fix_accum_bug=CS%fix_ustar_gustless_bug)
call allocate_forcing_type(G, fluxes, ustar=.true., fix_accum_bug=CS%fix_ustar_gustless_bug, tau_mag=.true.)
if (trim(CS%buoy_config) /= "NONE") then
if ( CS%use_temperature ) then
call allocate_forcing_type(G, fluxes, water=.true., heat=.true., press=.true.)
Expand Down Expand Up @@ -837,7 +839,7 @@ subroutine wind_forcing_by_data_override(sfc_state, forces, day, G, US, CS)
call callTree_enter("wind_forcing_by_data_override, MOM_surface_forcing.F90")

if (.not.CS%dataOverrideIsInitialized) then
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true.)
call allocate_mech_forcing(G, forces, stress=.true., ustar=.true., press=.true., tau_mag=.true.)
call data_override_init(G%Domain)
CS%dataOverrideIsInitialized = .True.
endif
Expand Down Expand Up @@ -953,7 +955,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS)

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec

if (CS%use_temperature) rhoXcp = CS%Rho0 * fluxes%C_p
if (CS%use_temperature) rhoXcp = CS%rho_restore * fluxes%C_p

! Read the buoyancy forcing file
call get_time(day, seconds, days)
Expand Down Expand Up @@ -1152,7 +1154,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS)
if (G%mask2dT(i,j) > 0.0) then
fluxes%heat_added(i,j) = G%mask2dT(i,j) * &
((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const_T)
fluxes%vprec(i,j) = - (CS%Rho0*CS%Flux_const_S) * &
fluxes%vprec(i,j) = - (CS%rho_restore*CS%Flux_const_S) * &
(CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / &
(0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)))
else
Expand All @@ -1164,7 +1166,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS)
do j=js,je ; do i=is,ie
if (G%mask2dT(i,j) > 0.0) then
fluxes%buoy(i,j) = (CS%Dens_Restore(i,j) - sfc_state%sfc_density(i,j)) * &
(CS%G_Earth * CS%Flux_const / CS%Rho0)
(CS%G_Earth * CS%Flux_const / CS%rho_restore)
else
fluxes%buoy(i,j) = 0.0
endif
Expand Down Expand Up @@ -1220,7 +1222,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US
is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec
isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed

if (CS%use_temperature) rhoXcp = CS%Rho0 * fluxes%C_p
if (CS%use_temperature) rhoXcp = CS%rho_restore * fluxes%C_p

if (.not.CS%dataOverrideIsInitialized) then
call data_override_init(G%Domain)
Expand Down Expand Up @@ -1258,7 +1260,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US
if (G%mask2dT(i,j) > 0.0) then
fluxes%heat_added(i,j) = G%mask2dT(i,j) * &
((CS%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoXcp * CS%Flux_const_T)
fluxes%vprec(i,j) = - (CS%Rho0*CS%Flux_const_S) * &
fluxes%vprec(i,j) = - (CS%rho_restore*CS%Flux_const_S) * &
(CS%S_Restore(i,j) - sfc_state%SSS(i,j)) / &
(0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)))
else
Expand All @@ -1270,7 +1272,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US
do j=js,je ; do i=is,ie
if (G%mask2dT(i,j) > 0.0) then
fluxes%buoy(i,j) = (CS%Dens_Restore(i,j) - sfc_state%sfc_density(i,j)) * &
(CS%G_Earth * CS%Flux_const / CS%Rho0)
(CS%G_Earth * CS%Flux_const / CS%rho_restore)
else
fluxes%buoy(i,j) = 0.0
endif
Expand Down Expand Up @@ -1457,8 +1459,8 @@ subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, US, CS)
S_restore = CS%S_south + (CS%S_north-CS%S_south)*y
if (G%mask2dT(i,j) > 0.0) then
fluxes%heat_added(i,j) = G%mask2dT(i,j) * &
((T_Restore - sfc_state%SST(i,j)) * ((CS%Rho0 * fluxes%C_p) * CS%Flux_const))
fluxes%vprec(i,j) = - (CS%Rho0*CS%Flux_const) * &
((T_Restore - sfc_state%SST(i,j)) * ((CS%rho_restore * fluxes%C_p) * CS%Flux_const))
fluxes%vprec(i,j) = - (CS%rho_restore*CS%Flux_const) * &
(S_Restore - sfc_state%SSS(i,j)) / &
(0.5*(sfc_state%SSS(i,j) + S_Restore))
else
Expand All @@ -1472,7 +1474,7 @@ subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, US, CS)
!do j=js,je ; do i=is,ie
! if (G%mask2dT(i,j) > 0.0) then
! fluxes%buoy(i,j) = (CS%Dens_Restore(i,j) - sfc_state%sfc_density(i,j)) * &
! (CS%G_Earth * CS%Flux_const / CS%Rho0)
! (CS%G_Earth * CS%Flux_const / CS%rho_restore)
! else
! fluxes%buoy(i,j) = 0.0
! endif
Expand Down Expand Up @@ -1874,6 +1876,12 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, tracer_flow_C
"at the southern end of the domain toward which to "//&
"to restore.", units="PSU", default=35.0, scale=US%ppt_to_S)
endif
call get_param(param_file, mdl, "RESTORE_FLUX_RHO", CS%rho_restore, &
"The density that is used to convert piston velocities into salt or heat "//&
"fluxes with RESTORE_SALINITY or RESTORE_TEMPERATURE.", &
units="kg m-3", default=CS%Rho0*US%R_to_kg_m3, scale=US%kg_m3_to_R, &
do_not_log=(((CS%Flux_const==0.0).and.(CS%Flux_const_T==0.0).and.(CS%Flux_const_S==0.0))&
.or.(.not.CS%restorebuoy)))
endif
call get_param(param_file, mdl, "G_EARTH", CS%G_Earth, &
"The gravitational acceleration of the Earth.", &
Expand Down
Loading