MOM6
BFB_initialization.F90
Go to the documentation of this file.
1 !> Initialization of the boundary-forced-basing configuration
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
8 use mom_get_input, only : directories
9 use mom_grid, only : ocean_grid_type
16 implicit none ; private
17 
18 #include <MOM_memory.h>
19 
20 public bfb_set_coord
22 
23 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
24 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
25 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
26 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
27 
28 !> Unsafe model variable
29 !! \todo Remove this module variable
30 logical :: first_call = .true.
31 
32 contains
33 
34 !> This subroutine specifies the vertical coordinate in terms of temperature at the surface and at the bottom.
35 !! This case is set up in such a way that the temperature of the topmost layer is equal to the SST at the
36 !! southern edge of the domain. The temperatures are then converted to densities of the top and bottom layers
37 !! and linearly interpolated for the intermediate layers.
38 subroutine bfb_set_coord(Rlay, g_prime, GV, US, param_file, eqn_of_state)
39  real, dimension(NKMEM_), intent(out) :: rlay !< Layer potential density [R ~> kg m-3].
40  real, dimension(NKMEM_), intent(out) :: g_prime !< The reduced gravity at
41  !! each interface [L2 Z-1 T-2 ~> m s-2].
42  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
43  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
44  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
45  type(eos_type), pointer :: eqn_of_state !< Integer that selects the
46  !! equation of state.
47  ! Local variables
48  real :: drho_dt, sst_s, t_bot, rho_top, rho_bot
49  integer :: k, nz
50  character(len=40) :: mdl = "BFB_set_coord" ! This subroutine's name.
51 
52  call get_param(param_file, mdl, "DRHO_DT", drho_dt, &
53  "Rate of change of density with temperature.", &
54  units="kg m-3 K-1", default=-0.2, scale=us%kg_m3_to_R)
55  call get_param(param_file, mdl, "SST_S", sst_s, &
56  "SST at the suothern edge of the domain.", units="C", default=20.0)
57  call get_param(param_file, mdl, "T_BOT", t_bot, &
58  "Bottom Temp", units="C", default=5.0)
59  rho_top = gv%Rho0 + drho_dt*sst_s
60  rho_bot = gv%Rho0 + drho_dt*t_bot
61  nz = gv%ke
62 
63  do k = 1,nz
64  rlay(k) = (rho_bot - rho_top)/(nz-1)*real(k-1) + rho_top
65  if (k >1) then
66  g_prime(k) = (rlay(k) - rlay(k-1)) * gv%g_Earth / (gv%Rho0)
67  else
68  g_prime(k) = gv%g_Earth
69  endif
70  !Rlay(:) = 0.0
71  !g_prime(:) = 0.0
72  enddo
73 
74  if (first_call) call write_bfb_log(param_file)
75 
76 end subroutine bfb_set_coord
77 
78 !> This subroutine sets up the sponges for the southern bouundary of the domain. Maximum damping occurs
79 !! within 2 degrees lat of the boundary. The damping linearly decreases northward over the next 2 degrees.
80 subroutine bfb_initialize_sponges_southonly(G, GV, US, use_temperature, tv, param_file, CSp, h)
81  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
82  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
83  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
84  logical, intent(in) :: use_temperature !< If true, temperature and salinity are used as
85  !! state variables.
86  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
87  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
88  type(sponge_cs), pointer :: csp !< A pointer to the sponge control structure
89  real, dimension(NIMEM_, NJMEM_, NKMEM_), &
90  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
91 
92  ! Local variables
93  real :: eta(szi_(g),szj_(g),szk_(g)+1) ! A temporary array for eta, in depth units [Z ~> m].
94  real :: idamp(szi_(g),szj_(g)) ! The inverse damping rate [s-1].
95  real :: h0(szk_(g)) ! Resting layer thicknesses in depth units [Z ~> m].
96  real :: min_depth ! The minimum ocean depth in depth units [Z ~> m].
97  real :: damp, e_dense, damp_new, slat, wlon, lenlat, lenlon, nlat
98  character(len=40) :: mdl = "BFB_initialize_sponges_southonly" ! This subroutine's name.
99  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
100 
101  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
102  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
103 
104  eta(:,:,:) = 0.0 ; idamp(:,:) = 0.0
105 
106 ! Here the inverse damping time [s-1], is set. Set Idamp to 0 !
107 ! wherever there is no sponge, and the subroutines that are called !
108 ! will automatically set up the sponges only where Idamp is positive!
109 ! and mask2dT is 1. !
110 
111 ! Set up sponges for DOME configuration
112  call get_param(param_file, mdl, "MINIMUM_DEPTH", min_depth, &
113  "The minimum depth of the ocean.", units="m", default=0.0, scale=us%m_to_Z)
114 
115  call get_param(param_file, mdl, "SOUTHLAT", slat, &
116  "The southern latitude of the domain.", units="degrees")
117  call get_param(param_file, mdl, "LENLAT", lenlat, &
118  "The latitudinal length of the domain.", units="degrees")
119  call get_param(param_file, mdl, "WESTLON", wlon, &
120  "The western longitude of the domain.", units="degrees", default=0.0)
121  call get_param(param_file, mdl, "LENLON", lenlon, &
122  "The longitudinal length of the domain.", units="degrees")
123  nlat = slat + lenlat
124  do k=1,nz ; h0(k) = -g%max_depth * real(k-1) / real(nz) ; enddo
125 
126  ! Use for meridional thickness profile initialization
127 ! do k=1,nz ; H0(k) = -G%max_depth * real(k-1) / real(nz-1) ; enddo
128 
129  do i=is,ie; do j=js,je
130  if (g%geoLatT(i,j) < slat+2.0) then ; damp = 1.0
131  elseif (g%geoLatT(i,j) < slat+4.0) then
132  damp_new = 1.0*(slat+4.0-g%geoLatT(i,j))/2.0
133  else ; damp = 0.0
134  endif
135 
136  ! These will be streched inside of apply_sponge, so they can be in
137  ! depth space for Boussinesq or non-Boussinesq models.
138 
139  ! This section is used for uniform thickness initialization
140  do k = 1,nz; eta(i,j,k) = h0(k); enddo
141 
142  ! The below section is used for meridional temperature profile thickness initiation
143  ! do k = 1,nz; eta(i,j,k) = H0(k); enddo
144  ! if (G%geoLatT(i,j) > 40.0) then
145  ! do k = 1,nz
146  ! eta(i,j,k) = -G%Angstrom_Z*(k-1)
147  ! enddo
148  ! elseif (G%geoLatT(i,j) > 20.0) then
149  ! do k = 1,nz
150  ! eta(i,j,k) = min(H0(k) + (G%geoLatT(i,j) - 20.0)*(G%max_depth - nz*G%Angstrom_Z)/20.0, &
151  ! -(k-1)*G%Angstrom_Z)
152  ! enddo
153  ! endif
154  eta(i,j,nz+1) = -g%max_depth
155 
156  if (g%bathyT(i,j) > min_depth) then
157  idamp(i,j) = damp/86400.0
158  else ; idamp(i,j) = 0.0 ; endif
159  enddo ; enddo
160 
161 ! This call sets up the damping rates and interface heights.
162 ! This sets the inverse damping timescale fields in the sponges. !
163  call initialize_sponge(idamp, eta, g, param_file, csp, gv)
164 
165 ! Now register all of the fields which are damped in the sponge. !
166 ! By default, momentum is advected vertically within the sponge, but !
167 ! momentum is typically not damped within the sponge. !
168 
169  if (first_call) call write_bfb_log(param_file)
170 
172 
173 !> Write output about the parameter values being used.
174 subroutine write_bfb_log(param_file)
175  type(param_file_type), intent(in) :: param_file !< A structure indicating the
176  !! open file to parse for model
177  !! parameter values.
178 
179 ! This include declares and sets the variable "version".
180 #include "version_variable.h"
181  character(len=40) :: mdl = "BFB_initialization" ! This module's name.
182 
183  call log_version(param_file, mdl, version)
184  first_call = .false.
185 
186 end subroutine write_bfb_log
187 
188 end module bfb_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
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
bfb_initialization::write_bfb_log
subroutine write_bfb_log(param_file)
Write output about the parameter values being used.
Definition: BFB_initialization.F90:175
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
bfb_initialization::bfb_initialize_sponges_southonly
subroutine, public bfb_initialize_sponges_southonly(G, GV, US, use_temperature, tv, param_file, CSp, h)
This subroutine sets up the sponges for the southern bouundary of the domain. Maximum damping occurs ...
Definition: BFB_initialization.F90:81
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
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
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
bfb_initialization::bfb_set_coord
subroutine, public bfb_set_coord(Rlay, g_prime, GV, US, param_file, eqn_of_state)
This subroutine specifies the vertical coordinate in terms of temperature at the surface and at the b...
Definition: BFB_initialization.F90:39
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
bfb_initialization
Initialization of the boundary-forced-basing configuration.
Definition: BFB_initialization.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
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
bfb_initialization::first_call
logical first_call
Unsafe model variable.
Definition: BFB_initialization.F90:30
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_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