Skip to content

Commit

Permalink
+Added the INTERFACE_IC_QUANTA runtime parameter
Browse files Browse the repository at this point in the history 
  Added a new parameter, INTERFACE_IC_QUANTA, that is used in the seamount and
dumbbell test case to ensure that the initial interface heights can be easily
represented in output files and are not sensitive to order of arithmetic changes
in the setup.  By default all answers are bitwise identical, but the
MOM_parameter_doc.short files change for the seamount and dumbbell test cases.
  • Loading branch information
@Hallberg-NOAA
Hallberg-NOAA committed Nov 6, 2018
1 parent 6ccb6b6 commit 0874ac4
Show file tree
Hide file tree
Showing 2 changed files with 28 additions and 13 deletions.
20 changes: 14 additions & 6 deletions src/user/dumbbell_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ module dumbbell_initialization
use MOM_domains, only : sum_across_PEs
use MOM_dyn_horgrid, only : dyn_horgrid_type
use MOM_error_handler, only : MOM_mesg, MOM_error, FATAL, is_root_pe
use MOM_file_parser, only : get_param, param_file_type
use MOM_file_parser, only : get_param, log_version, param_file_type
use MOM_get_input, only : directories
use MOM_grid, only : ocean_grid_type
use MOM_sponge, only : set_up_sponge_field, initialize_sponge, sponge_CS
Expand Down Expand Up @@ -86,12 +86,14 @@ subroutine dumbbell_initialize_thickness ( h, G, GV, US, param_file, just_read_p
! negative because it is positive upward.
real :: eta1D(SZK_(G)+1)! Interface height relative to the sea surface
! positive upward, in depth units (Z).
integer :: i, j, k, is, ie, js, je, nz
real :: x
real :: delta_h
real :: min_thickness, S_surf, S_range, S_ref, S_light, S_dense
real :: min_thickness ! The minimum layer thicknesses, in Z.
real :: S_surf, S_range, S_ref, S_light, S_dense ! Various salinities, in ppt.
real :: eta_IC_quanta ! The granularity of quantization of intial interface heights, in Z-1.
! This include declares and sets the variable "version".
# include "version_variable.h"
character(len=20) :: verticalCoordinate
logical :: just_read ! If true, just read parameters but set nothing.
integer :: i, j, k, is, ie, js, je, nz

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

Expand All @@ -100,6 +102,7 @@ subroutine dumbbell_initialize_thickness ( h, G, GV, US, param_file, just_read_p
if (.not.just_read) &
call MOM_mesg("MOM_initialization.F90, initialize_thickness_uniform: setting thickness")

if (.not.just_read) call log_version(param_file, mdl, version, "")
call get_param(param_file, mdl,"MIN_THICKNESS", min_thickness, &
'Minimum thickness for layer',&
units='m', default=1.0e-3, do_not_log=just_read, scale=US%m_to_Z)
Expand All @@ -125,6 +128,10 @@ subroutine dumbbell_initialize_thickness ( h, G, GV, US, param_file, just_read_p
call get_param(param_file, mdl, "S_REF", S_ref, default=35.0, do_not_log=.true.)
call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", S_light, default = S_Ref, do_not_log=.true.)
call get_param(param_file, mdl, "TS_RANGE_S_DENSE", S_dense, default = S_Ref, do_not_log=.true.)
call get_param(param_file, mdl, "INTERFACE_IC_QUANTA", eta_IC_quanta, &
"The granularity of initial interface height values \n"//&
"per meter, to avoid sensivity to order-of-arithmetic changes.", &
default=2048.0, units="m-1", scale=US%Z_to_m, do_not_log=just_read)
if (just_read) return ! All run-time parameters have been read, so return.

do K=1,nz+1
Expand All @@ -137,7 +144,8 @@ subroutine dumbbell_initialize_thickness ( h, G, GV, US, param_file, just_read_p
e0(K) = - G%max_depth * ( ( S_light - S_surf ) + ( S_dense - S_light ) * &
( (real(K)-1.5) / real(nz-1) ) ) / S_range
! Force round numbers ... the above expression has irrational factors ...
e0(K) = nint(2048.*US%Z_to_m*e0(K)) / (2048.*US%Z_to_m)
if (eta_IC_quanta > 0.0) &
e0(K) = nint(eta_IC_quanta*e0(K)) / eta_IC_quanta
e0(K) = min(real(1-K)*GV%Angstrom_Z, e0(K)) ! Bound by surface
e0(K) = max(-G%max_depth, e0(K)) ! Bound by bottom
enddo
Expand Down
21 changes: 14 additions & 7 deletions src/user/seamount_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -85,11 +85,13 @@ subroutine seamount_initialize_thickness ( h, G, GV, US, param_file, just_read_p
logical, optional, intent(in) :: just_read_params !< If present and true, this call will
!! only read parameters without changing h.

real :: e0(SZK_(G)+1) ! The resting interface heights, in m, usually !
! negative because it is positive upward. !
real :: eta1D(SZK_(G)+1)! Interface height relative to the sea surface !
! positive upward, in m. !
real :: min_thickness, S_surf, S_range, S_ref, S_light, S_dense
real :: e0(SZK_(G)+1) ! The resting interface heights, in depth units (Z), usually
! negative because it is positive upward.
real :: eta1D(SZK_(G)+1)! Interface height relative to the sea surface
! positive upward, in depth units (Z).
real :: min_thickness ! The minimum layer thicknesses, in Z.
real :: S_surf, S_range, S_ref, S_light, S_dense ! Various salinities, in ppt.
real :: eta_IC_quanta ! The granularity of quantization of intial interface heights, in Z-1.
character(len=20) :: verticalCoordinate
logical :: just_read ! If true, just read parameters but set nothing.
integer :: i, j, k, is, ie, js, je, nz
Expand All @@ -102,7 +104,7 @@ subroutine seamount_initialize_thickness ( h, G, GV, US, param_file, just_read_p
call MOM_mesg("MOM_initialization.F90, initialize_thickness_uniform: setting thickness")

call get_param(param_file, mdl,"MIN_THICKNESS",min_thickness, &
'Minimum thickness for layer',&
'Minimum thickness for layer', &
units='m', default=1.0e-3, do_not_log=just_read, scale=US%m_to_Z)
call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE",verticalCoordinate, &
default=DEFAULT_COORDINATE_MODE, do_not_log=just_read)
Expand All @@ -126,6 +128,10 @@ subroutine seamount_initialize_thickness ( h, G, GV, US, param_file, just_read_p
call get_param(param_file, mdl, "S_REF", S_ref, default=35.0, do_not_log=.true.)
call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", S_light, default = S_Ref, do_not_log=.true.)
call get_param(param_file, mdl, "TS_RANGE_S_DENSE", S_dense, default = S_Ref, do_not_log=.true.)
call get_param(param_file, mdl, "INTERFACE_IC_QUANTA", eta_IC_quanta, &
"The granularity of initial interface height values \n"//&
"per meter, to avoid sensivity to order-of-arithmetic changes.", &
default=2048.0, units="m-1", scale=US%Z_to_m, do_not_log=just_read)
if (just_read) return ! All run-time parameters have been read, so return.

do K=1,nz+1
Expand All @@ -138,7 +144,8 @@ subroutine seamount_initialize_thickness ( h, G, GV, US, param_file, just_read_p
e0(K) = - G%max_depth * ( ( S_light - S_surf ) + ( S_dense - S_light ) * &
( (real(K)-1.5) / real(nz-1) ) ) / S_range
! Force round numbers ... the above expression has irrational factors ...
e0(K) = nint(2048.*US%Z_to_m*e0(K))/(2048.*US%Z_to_m)
if (eta_IC_quanta > 0.0) &
e0(K) = nint(eta_IC_quanta*e0(K)) / eta_IC_quanta
e0(K) = min(real(1-K)*GV%Angstrom_Z, e0(K)) ! Bound by surface
e0(K) = max(-G%max_depth, e0(K)) ! Bound by bottom
enddo
Expand Down

0 comments on commit 0874ac4

Please sign in to comment.