MOM6
baroclinic_zone_initialization.F90
Go to the documentation of this file.
1 !> Initial conditions for an idealized baroclinic zone
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
8 use mom_grid, only : ocean_grid_type
11 
12 implicit none ; private
13 
14 #include <MOM_memory.h>
15 #include "version_variable.h"
16 
17 ! Private (module-wise) parameters
18 character(len=40) :: mdl = "baroclinic_zone_initialization" !< This module's name.
19 
21 
22 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
23 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
24 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
25 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
26 
27 contains
28 
29 !> Reads the parameters unique to this module
30 subroutine bcz_params(G, GV, US, param_file, S_ref, dSdz, delta_S, dSdx, T_ref, dTdz, &
31  delta_T, dTdx, L_zone, just_read_params)
32  type(ocean_grid_type), intent(in) :: G !< Grid structure
33  type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
34  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
35  type(param_file_type), intent(in) :: param_file !< Parameter file handle
36  real, intent(out) :: S_ref !< Reference salinity [ppt]
37  real, intent(out) :: dSdz !< Salinity stratification [ppt Z-1 ~> ppt m-1]
38  real, intent(out) :: delta_S !< Salinity difference across baroclinic zone [ppt]
39  real, intent(out) :: dSdx !< Linear salinity gradient [ppt m-1]
40  real, intent(out) :: T_ref !< Reference temperature [degC]
41  real, intent(out) :: dTdz !< Temperature stratification [degC Z-1 ~> degC m-1]
42  real, intent(out) :: delta_T !< Temperature difference across baroclinic zone [degC]
43  real, intent(out) :: dTdx !< Linear temperature gradient [degC m-1]
44  real, intent(out) :: L_zone !< Width of baroclinic zone [m]
45  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
46  !! only read parameters without changing h.
47 
48  logical :: just_read ! If true, just read parameters but set nothing.
49 
50  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
51 
52  if (.not.just_read) &
53  call log_version(param_file, mdl, version, 'Initialization of an analytic baroclinic zone')
54  call openparameterblock(param_file,'BCZIC')
55  call get_param(param_file, mdl, "S_REF", s_ref, 'Reference salinity', units='ppt', &
56  default=35., do_not_log=just_read)
57  call get_param(param_file, mdl, "DSDZ", dsdz, 'Salinity stratification', &
58  units='ppt/m', default=0.0, scale=us%Z_to_m, do_not_log=just_read)
59  call get_param(param_file, mdl,"DELTA_S",delta_s,'Salinity difference across baroclinic zone', &
60  units='ppt', default=0.0, do_not_log=just_read)
61  call get_param(param_file, mdl,"DSDX",dsdx,'Meridional salinity difference', &
62  units='ppt/'//trim(g%x_axis_units), default=0.0, do_not_log=just_read)
63  call get_param(param_file, mdl,"T_REF",t_ref,'Reference temperature',units='C', &
64  default=10., do_not_log=just_read)
65  call get_param(param_file, mdl, "DTDZ", dtdz, 'Temperature stratification', &
66  units='C/m', default=0.0, scale=us%Z_to_m, do_not_log=just_read)
67  call get_param(param_file, mdl,"DELTA_T",delta_t,'Temperature difference across baroclinic zone', &
68  units='C', default=0.0, do_not_log=just_read)
69  call get_param(param_file, mdl,"DTDX",dtdx,'Meridional temperature difference', &
70  units='C/'//trim(g%x_axis_units), default=0.0, do_not_log=just_read)
71  call get_param(param_file, mdl,"L_ZONE",l_zone,'Width of baroclinic zone', &
72  units=g%x_axis_units, default=0.5*g%len_lat, do_not_log=just_read)
73  call closeparameterblock(param_file)
74 
75 end subroutine bcz_params
76 
77 !> Initialization of temperature and salinity with the baroclinic zone initial conditions
78 subroutine baroclinic_zone_init_temperature_salinity(T, S, h, G, GV, US, param_file, &
79  just_read_params)
80  type(ocean_grid_type), intent(in) :: g !< Grid structure
81  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
82  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
83  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: t !< Potential temperature [degC]
84  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: s !< Salinity [ppt]
85  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(in) :: h !< The model thicknesses [H ~> m or kg m-2]
86  type(param_file_type), intent(in) :: param_file !< Parameter file handle
87  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
88  !! only read parameters without changing T & S.
89 
90  integer :: i, j, k, is, ie, js, je, nz
91  real :: t_ref, dtdz, dtdx, delta_t ! Parameters describing temperature distribution
92  real :: s_ref, dsdz, dsdx, delta_s ! Parameters describing salinity distribution
93  real :: l_zone ! Width of baroclinic zone
94  real :: zc, zi ! Depths in depth units [Z ~> m]
95  real :: x, xd, xs, y, yd, fn
96  real :: pi ! 3.1415926... calculated as 4*atan(1)
97  logical :: just_read ! If true, just read parameters but set nothing.
98 
99  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
100  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
101 
102  call bcz_params(g, gv, us, param_file, s_ref, dsdz, delta_s, dsdx, t_ref, dtdz, &
103  delta_t, dtdx, l_zone, just_read_params)
104 
105  if (just_read) return ! All run-time parameters have been read, so return.
106 
107  t(:,:,:) = 0.
108  s(:,:,:) = 0.
109  pi = 4.*atan(1.)
110 
111  do j = g%jsc,g%jec ; do i = g%isc,g%iec
112  zi = -g%bathyT(i,j)
113  x = g%geoLonT(i,j) - (g%west_lon + 0.5*g%len_lon) ! Relative to center of domain
114  xd = x / g%len_lon ! -1/2 < xd 1/2
115  y = g%geoLatT(i,j) - (g%south_lat + 0.5*g%len_lat) ! Relative to center of domain
116  yd = y / g%len_lat ! -1/2 < yd 1/2
117  if (l_zone/=0.) then
118  xs = min(1., max(-1., x/l_zone)) ! -1 < ys < 1
119  fn = sin((0.5*pi)*xs)
120  else
121  xs = sign(1., x) ! +/- 1
122  fn = xs
123  endif
124  do k = nz, 1, -1
125  zc = zi + 0.5*h(i,j,k)*gv%H_to_Z ! Position of middle of cell
126  zi = zi + h(i,j,k)*gv%H_to_Z ! Top interface position
127  t(i,j,k) = t_ref + dtdz * zc & ! Linear temperature stratification
128  + dtdx * x & ! Linear gradient
129  + delta_t * fn ! Smooth fn of width L_zone
130  s(i,j,k) = s_ref + dsdz * zc & ! Linear temperature stratification
131  + dsdx * x & ! Linear gradient
132  + delta_s * fn ! Smooth fn of width L_zone
133  enddo
134  enddo ; enddo
135 
137 
138 !> \namespace baroclinic_zone_initialization
139 !!
140 !! \section section_baroclinic_zone Description of the baroclinic zone initial conditions
141 !!
142 !! yada yada yada
143 
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
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
baroclinic_zone_initialization::mdl
character(len=40) mdl
This module's name.
Definition: baroclinic_zone_initialization.F90:18
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
baroclinic_zone_initialization::baroclinic_zone_init_temperature_salinity
subroutine, public baroclinic_zone_init_temperature_salinity(T, S, h, G, GV, US, param_file, just_read_params)
Initialization of temperature and salinity with the baroclinic zone initial conditions.
Definition: baroclinic_zone_initialization.F90:80
baroclinic_zone_initialization::bcz_params
subroutine bcz_params(G, GV, US, param_file, S_ref, dSdz, delta_S, dSdx, T_ref, dTdz, delta_T, dTdx, L_zone, just_read_params)
Reads the parameters unique to this module.
Definition: baroclinic_zone_initialization.F90:32
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_file_parser::openparameterblock
subroutine, public openparameterblock(CS, blockName, desc)
Tags blockName onto the end of the active parameter block name.
Definition: MOM_file_parser.F90:2015
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_file_parser::closeparameterblock
subroutine, public closeparameterblock(CS)
Remove the lowest level of recursion from the active block name.
Definition: MOM_file_parser.F90:2033
baroclinic_zone_initialization
Initial conditions for an idealized baroclinic zone.
Definition: baroclinic_zone_initialization.F90:2