MOM6
Neverland_initialization.F90
Go to the documentation of this file.
1 !> Initialization for the "Neverland" configuration
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
10 use mom_get_input, only : directories
11 use mom_grid, only : ocean_grid_type
17 
18 use random_numbers_mod, only: initializerandomnumberstream, getrandomnumbers, randomnumberstream
19 
20 implicit none ; private
21 
22 #include <MOM_memory.h>
23 
26 
27 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
28 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
29 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
30 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
31 
32 contains
33 
34 !> This subroutine sets up the Neverland test case topography.
35 subroutine neverland_initialize_topography(D, G, param_file, max_depth)
36  type(dyn_horgrid_type), intent(in) :: g !< The dynamic horizontal grid type
37  real, dimension(G%isd:G%ied,G%jsd:G%jed), &
38  intent(out) :: d !< Ocean bottom depth in the units of depth_max
39  type(param_file_type), intent(in) :: param_file !< Parameter file structure
40  real, intent(in) :: max_depth !< Maximum ocean depth in arbitrary units
41 
42  ! Local variables
43  real :: pi ! 3.1415926... calculated as 4*atan(1)
44  real :: d0 ! A constant to make the maximum !
45  ! basin depth MAXIMUM_DEPTH. !
46  real :: x, y
47  ! This include declares and sets the variable "version".
48 # include "version_variable.h"
49  character(len=40) :: mdl = "Neverland_initialize_topography" ! This subroutine's name.
50  integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
51  real :: nl_roughness_amp, nl_top_amp
52  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
53  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
54 
55  call mom_mesg(" Neverland_initialization.F90, Neverland_initialize_topography: setting topography", 5)
56 
57  call log_version(param_file, mdl, version, "")
58  call get_param(param_file, mdl, "NL_ROUGHNESS_AMP", nl_roughness_amp, &
59  "Amplitude of wavy signal in bathymetry.", default=0.05)
60  call get_param(param_file, mdl, "NL_CONTINENT_AMP", nl_top_amp, &
61  "Scale factor for topography - 0.0 for no continents.", default=1.0)
62 
63  pi = 4.0*atan(1.0)
64 
65 ! Calculate the depth of the bottom.
66  do j=js,je ; do i=is,ie
67  x = (g%geoLonT(i,j)-g%west_lon) / g%len_lon
68  y =( g%geoLatT(i,j)-g%south_lat) / g%len_lat
69 ! This sets topography that has a reentrant channel to the south.
70  d(i,j) = 1.0 - 1.1 * spike(y-1,0.12) - 1.1 * spike(y,0.12) - & !< The great northern wall and Antarctica
71  nl_top_amp*( &
72  (1.2 * spike(x,0.2) + 1.2 * spike(x-1.0,0.2)) * spike(min(0.0,y-0.3),0.2) & !< South America
73  + 1.2 * spike(x-0.5,0.2) * spike(min(0.0,y-0.55),0.2) & !< Africa
74  + 1.2 * (spike(x,0.12) + spike(x-1,0.12)) * spike(max(0.0,y-0.06),0.12) & !< Antarctic Peninsula
75  + 0.1 * (cosbell(x,0.1) + cosbell(x-1,0.1)) & !< Drake Passage ridge
76  + 0.5 * cosbell(x-0.16,0.05) * (cosbell(y-0.18,0.13)**0.4) & !< Scotia Arc East
77  + 0.4 * (cosbell(x-0.09,0.08)**0.4) * cosbell(y-0.26,0.05) & !< Scotia Arc North
78  + 0.4 * (cosbell(x-0.08,0.08)**0.4) * cosbell(y-0.1,0.05)) & !< Scotia Arc South
79  - nl_roughness_amp * cos(14*pi*x) * sin(14*pi*y) & !< roughness
80  - nl_roughness_amp * cos(20*pi*x) * cos(20*pi*y) !< roughness
81  if (d(i,j) < 0.0) d(i,j) = 0.0
82  d(i,j) = d(i,j) * max_depth
83  enddo ; enddo
84 
86 ! -----------------------------------------------------------------------------
87 
88 !> Returns the value of a cosine-bell function evaluated at x/L
89 real function cosbell(x, L)
90  real , intent(in) :: x !< non-dimensional position
91  real , intent(in) :: l !< non-dimensional width
92  real :: pi !< 3.1415926... calculated as 4*atan(1)
93 
94  pi = 4.0*atan(1.0)
95  cosbell = 0.5 * (1 + cos(pi*min(abs(x/l),1.0)))
96 end function cosbell
97 
98 !> Returns the value of a sin-spike function evaluated at x/L
99 real function spike(x, L)
100 
101  real , intent(in) :: x !< non-dimensional position
102  real , intent(in) :: l !< non-dimensional width
103  real :: pi !< 3.1415926... calculated as 4*atan(1)
104 
105  pi = 4.0*atan(1.0)
106  spike = (1 - sin(pi*min(abs(x/l),0.5)))
107 end function spike
108 
109 !> This subroutine initializes layer thicknesses for the Neverland test case,
110 !! by finding the depths of interfaces in a specified latitude-dependent
111 !! temperature profile with an exponentially decaying thermocline on top of a
112 !! linear stratification.
113 subroutine neverland_initialize_thickness(h, G, GV, US, param_file, eqn_of_state, P_ref)
114  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
115  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
116  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
117  real, intent(out), dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h !< The thickness that is being
118  !! initialized [H ~> m or kg m-2].
119  type(param_file_type), intent(in) :: param_file !< A structure indicating the open
120  !! file to parse for model
121  !! parameter values.
122  type(eos_type), pointer :: eqn_of_state !< integer that selects the
123  !! equation of state.
124  real, intent(in) :: p_ref !< The coordinate-density
125  !! reference pressure [Pa].
126  ! Local variables
127  real :: e0(szk_(g)+1) ! The resting interface heights, in depth units [Z ~> m],
128  ! usually negative because it is positive upward.
129  real, dimension(SZK_(G)) :: h_profile ! Vector of initial thickness profile [Z ~> m].
130  real :: e_interface ! Current interface position [Z ~> m].
131  real :: x,y,r1,r2 ! x,y and radial coordinates for computation of initial pert.
132  real :: pert_amp ! Amplitude of perturbations measured in Angstrom_H
133  real :: h_noise ! Amplitude of noise to scale h by
134  real :: noise ! Noise
135  type(randomnumberstream) :: rns ! Random numbers for stochastic tidal parameterization
136  character(len=40) :: mdl = "Neverland_initialize_thickness" ! This subroutine's name.
137  integer :: i, j, k, k1, is, ie, js, je, nz, itt
138 
139  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
140 
141  call mom_mesg(" Neverland_initialization.F90, Neverland_initialize_thickness: setting thickness", 5)
142  call get_param(param_file, mdl, "INIT_THICKNESS_PROFILE", h_profile, &
143  "Profile of initial layer thicknesses.", units="m", scale=us%m_to_Z, &
144  fail_if_missing=.true.)
145  call get_param(param_file, mdl, "NL_THICKNESS_PERT_AMP", pert_amp, &
146  "Amplitude of finite scale perturbations as fraction of depth.", &
147  units="nondim", default=0.)
148  call get_param(param_file, mdl, "NL_THICKNESS_NOISE_AMP", h_noise, &
149  "Amplitude of noise to scale layer by.", units="nondim", default=0.)
150 
151  ! e0 is the notional position of interfaces
152  e0(1) = 0. ! The surface
153  do k=1,nz
154  e0(k+1) = e0(k) - h_profile(k)
155  enddo
156 
157  do j=js,je ; do i=is,ie
158  e_interface = -g%bathyT(i,j)
159  do k=nz,2,-1
160  h(i,j,k) = gv%Z_to_H * (e0(k) - e_interface) ! Nominal thickness
161  x=(g%geoLonT(i,j)-g%west_lon)/g%len_lon
162  y=(g%geoLatT(i,j)-g%south_lat)/g%len_lat
163  r1=sqrt((x-0.7)**2+(y-0.2)**2)
164  r2=sqrt((x-0.3)**2+(y-0.25)**2)
165  h(i,j,k) = h(i,j,k) + pert_amp * (e0(k) - e0(nz+1)) * gv%Z_to_H * &
166  (spike(r1,0.15)-spike(r2,0.15)) ! Prescribed perturbation
167  if (h_noise /= 0.) then
168  rns = initializerandomnumberstream( int( 4096*(x + (y+1.)) ) )
169  call getrandomnumbers(rns, noise) ! x will be in (0,1)
170  noise = h_noise * 2. * ( noise - 0.5 ) ! range -h_noise to h_noise
171  h(i,j,k) = ( 1. + noise ) * h(i,j,k)
172  endif
173  h(i,j,k) = max( gv%Angstrom_H, h(i,j,k) ) ! Limit to non-negative
174  e_interface = e_interface + gv%H_to_Z * h(i,j,k) ! Actual position of upper interface
175  enddo
176  h(i,j,1) = gv%Z_to_H * (e0(1) - e_interface) ! Nominal thickness
177  h(i,j,1) = max( gv%Angstrom_H, h(i,j,1) ) ! Limit to non-negative
178  enddo ; enddo
179 
180 end subroutine neverland_initialize_thickness
181 
182 end module neverland_initialization
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_sponge::set_up_sponge_field
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
Definition: MOM_sponge.F90:214
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
neverland_initialization::neverland_initialize_thickness
subroutine, public neverland_initialize_thickness(h, G, GV, US, param_file, eqn_of_state, P_ref)
This subroutine initializes layer thicknesses for the Neverland test case, by finding the depths of i...
Definition: Neverland_initialization.F90:114
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_get_input::directories
Container for paths and parameter file names.
Definition: MOM_get_input.F90:20
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
neverland_initialization::cosbell
real function cosbell(x, L)
Returns the value of a cosine-bell function evaluated at x/L.
Definition: Neverland_initialization.F90:90
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
neverland_initialization
Initialization for the "Neverland" configuration.
Definition: Neverland_initialization.F90:2
mom_eos::eos_type
A control structure for the equation of state.
Definition: MOM_EOS.F90:91
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
neverland_initialization::neverland_initialize_topography
subroutine, public neverland_initialize_topography(D, G, param_file, max_depth)
This subroutine sets up the Neverland test case topography.
Definition: Neverland_initialization.F90:36
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.F90:2
mom_eos::calculate_density_derivs
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Definition: MOM_EOS.F90:70
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
neverland_initialization::spike
real function spike(x, L)
Returns the value of a sin-spike function evaluated at x/L.
Definition: Neverland_initialization.F90:100
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_sponge::initialize_sponge
subroutine, public initialize_sponge(Iresttime, int_height, G, param_file, CS, GV, Iresttime_i_mean, int_height_i_mean)
This subroutine determines the number of points which are within sponges in this computational domain...
Definition: MOM_sponge.F90:90
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60