MOM6
Functions/Subroutines | Variables
adjustment_initialization Module Reference

Detailed Description

Configures the model for the geostrophic adjustment test case.

Functions/Subroutines

subroutine, public adjustment_initialize_thickness (h, G, GV, US, param_file, just_read_params)
 Initializes the layer thicknesses in the adjustment test case. More...
 
subroutine, public adjustment_initialize_temperature_salinity (T, S, h, G, GV, param_file, eqn_of_state, just_read_params)
 Initialization of temperature and salinity in the adjustment test case. More...
 

Variables

character(len=40) mdl = "adjustment_initialization"
 This module's name.
 

Function/Subroutine Documentation

◆ adjustment_initialize_temperature_salinity()

subroutine, public adjustment_initialization::adjustment_initialize_temperature_salinity ( real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)  T,
real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)  S,
real, dimension(szi_(g),szj_(g), szk_(g)), intent(in)  h,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
type(eos_type), pointer  eqn_of_state,
logical, intent(in), optional  just_read_params 
)

Initialization of temperature and salinity in the adjustment test case.

Parameters
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[out]tThe temperature that is being initialized.
[out]sThe salinity that is being initialized.
[in]hThe model thicknesses [H ~> m or kg m-2].
[in]param_fileA structure indicating the open file to parse for model parameter values.
eqn_of_stateEquation of state.
[in]just_read_paramsIf present and true, this call will only read parameters without changing T & S.

Definition at line 198 of file adjustment_initialization.F90.

198  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
199  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
200  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: t !< The temperature that is being initialized.
201  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: s !< The salinity that is being initialized.
202  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(in) :: h !< The model thicknesses [H ~> m or kg m-2].
203  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file to
204  !! parse for model parameter values.
205  type(eos_type), pointer :: eqn_of_state !< Equation of state.
206  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
207  !! only read parameters without changing T & S.
208 
209  integer :: i, j, k, is, ie, js, je, nz
210  real :: x, y, yy
211  integer :: index_bay_z
212  real :: s_ref, t_ref ! Reference salinity and temerature within
213  ! surface layer
214  real :: s_range, t_range ! Range of salinities and temperatures over the
215  ! vertical
216  real :: xi0, xi1, dsdz, delta_s, delta_s_strat
217  real :: adjustment_width, adjustment_deltas
218  real :: front_wave_amp, front_wave_length, front_wave_asym
219  real :: eta1d(szk_(g)+1)
220  logical :: just_read ! If true, just read parameters but set nothing.
221  character(len=20) :: verticalcoordinate
222 
223  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
224 
225  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
226 
227  ! Parameters used by main model initialization
228  call get_param(param_file, mdl, "S_REF", s_ref, 'Reference salinity', &
229  default=35.0, units='1e-3', do_not_log=just_read)
230  call get_param(param_file, mdl,"T_REF",t_ref,'Reference temperature', units='C', &
231  fail_if_missing=.not.just_read, do_not_log=just_read)
232  call get_param(param_file, mdl,"S_RANGE",s_range,'Initial salinity range', units='1e-3', &
233  default=2.0, do_not_log=just_read)
234  call get_param(param_file, mdl,"T_RANGE",t_range,'Initial temperature range', units='C', &
235  default=0.0, do_not_log=just_read)
236  ! Parameters specific to this experiment configuration BUT logged in previous s/r
237  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE",verticalcoordinate, &
238  default=default_coordinate_mode, do_not_log=just_read)
239  call get_param(param_file, mdl,"ADJUSTMENT_WIDTH", adjustment_width, &
240  fail_if_missing=.not.just_read, do_not_log=.true.)
241  call get_param(param_file, mdl,"ADJUSTMENT_DELTAS", adjustment_deltas, &
242  fail_if_missing=.not.just_read, do_not_log=.true.)
243  call get_param(param_file, mdl,"DELTA_S_STRAT", delta_s_strat, &
244  fail_if_missing=.not.just_read, do_not_log=.true.)
245  call get_param(param_file, mdl,"FRONT_WAVE_AMP", front_wave_amp, default=0., &
246  do_not_log=.true.)
247  call get_param(param_file, mdl,"FRONT_WAVE_LENGTH",front_wave_length, &
248  default=0., do_not_log=.true.)
249  call get_param(param_file, mdl,"FRONT_WAVE_ASYM", front_wave_asym, default=0., &
250  do_not_log=.true.)
251 
252  if (just_read) return ! All run-time parameters have been read, so return.
253 
254  t(:,:,:) = 0.0
255  s(:,:,:) = 0.0
256 
257  ! Linear salinity profile
258  select case ( coordinatemode(verticalcoordinate) )
259 
260  case ( regridding_zstar, regridding_sigma )
261  dsdz = -delta_s_strat / g%max_depth
262  do j=js,je ; do i=is,ie
263  eta1d(nz+1) = -g%bathyT(i,j)
264  do k=nz,1,-1
265  eta1d(k) = eta1d(k+1) + h(i,j,k)*gv%H_to_Z
266  enddo
267  if (front_wave_length /= 0.) then
268  y = ( 0.125 + g%geoLatT(i,j) / front_wave_length ) * ( 4. * acos(0.) )
269  yy = 2. * ( g%geoLatT(i,j) - 0.5 * g%len_lat ) / front_wave_length
270  yy = min(1.0, yy); yy = max(-1.0, yy)
271  yy = yy * 2. * acos( 0. )
272  y = front_wave_amp*sin(y) + front_wave_asym*sin(yy)
273  else
274  y = 0.
275  endif
276  x = ( ( g%geoLonT(i,j) - 0.5 * g%len_lon ) + y ) / adjustment_width
277  x = min(1.0, x); x = max(-1.0, x)
278  x = x * acos( 0. )
279  delta_s = adjustment_deltas * 0.5 * (1. - sin( x ) )
280  do k=1,nz
281  s(i,j,k) = s_ref + delta_s + 0.5 * ( eta1d(k)+eta1d(k+1) ) * dsdz
282  x = abs(s(i,j,k) - 0.5*real(nz-1)/real(nz)*s_range)/s_range*real(2*nz)
283  x = 1. - min(1., x)
284  t(i,j,k) = x
285  enddo
286  ! x = GV%H_to_Z*sum(T(i,j,:)*h(i,j,:))
287  ! T(i,j,:) = (T(i,j,:) / x) * (G%max_depth*1.5/real(nz))
288  enddo ; enddo
289 
290  case ( regridding_layer, regridding_rho )
291  do k = 1,nz
292  s(:,:,k) = s_ref + s_range * ( (real(k)-0.5) / real( nz ) )
293  ! x = abs(S(1,1,k) - 0.5*real(nz-1)/real(nz)*S_range)/S_range*real(2*nz)
294  ! x = 1.-min(1., x)
295  ! T(:,:,k) = x
296  enddo
297 
298  case default
299  call mom_error(fatal,"adjustment_initialize_temperature_salinity: "// &
300  "Unrecognized i.c. setup - set ADJUSTMENT_IC")
301 
302  end select
303 

◆ adjustment_initialize_thickness()

subroutine, public adjustment_initialization::adjustment_initialize_thickness ( real, dimension(szi_(g),szj_(g),szk_(gv)), intent(out)  h,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
logical, intent(in), optional  just_read_params 
)

Initializes the layer thicknesses in the adjustment test case.

Parameters
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[out]hThe thickness that is being initialized [H ~> m or kg m-2].
[in]param_fileA structure indicating the open file to parse for model parameter values.
[in]just_read_paramsIf present and true, this call will only read parameters without changing h.

Definition at line 36 of file adjustment_initialization.F90.

36  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
37  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
38  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
39  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
40  intent(out) :: h !< The thickness that is being initialized [H ~> m or kg m-2].
41  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
42  !! to parse for model parameter values.
43  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
44  !! only read parameters without changing h.
45  ! Local variables
46  real :: e0(szk_(g)+1) ! The resting interface heights, in depth units [Z ~> m], usually
47  ! negative because it is positive upward.
48  real :: eta1d(szk_(g)+1)! Interface height relative to the sea surface
49  ! positive upward, in depth units [Z ~> m].
50  real :: drho_ds ! The partial derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].
51  ! In this subroutine it is hard coded at 1.0 kg m-3 ppt-1.
52  real :: x, y, yy
53  real :: delta_s_strat, dsdz, delta_s, s_ref
54  real :: min_thickness, adjustment_width, adjustment_delta
55  real :: adjustment_deltas
56  real :: front_wave_amp, front_wave_length, front_wave_asym
57  real :: target_values(szk_(g)+1) ! Target densities or density anomalies [R ~> kg m-3]
58  logical :: just_read ! If true, just read parameters but set nothing.
59  character(len=20) :: verticalcoordinate
60 ! This include declares and sets the variable "version".
61 #include "version_variable.h"
62  integer :: i, j, k, is, ie, js, je, nz
63 
64  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
65 
66  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
67 
68  if (.not.just_read) &
69  call mom_mesg("initialize_thickness_uniform: setting thickness")
70 
71  ! Parameters used by main model initialization
72  if (.not.just_read) call log_version(param_file, mdl, version, "")
73  call get_param(param_file, mdl, "S_REF", s_ref, 'Reference salinity', &
74  default=35.0, units='1e-3', do_not_log=just_read)
75  call get_param(param_file, mdl,"MIN_THICKNESS",min_thickness,'Minimum layer thickness', &
76  default=1.0e-3, units='m', scale=us%m_to_Z, do_not_log=just_read)
77 
78  ! Parameters specific to this experiment configuration
79  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE",verticalcoordinate, &
80  default=default_coordinate_mode, do_not_log=just_read)
81  call get_param(param_file, mdl,"ADJUSTMENT_WIDTH",adjustment_width, &
82  "Width of frontal zone", &
83  units="same as x,y", fail_if_missing=.not.just_read, do_not_log=just_read)
84  call get_param(param_file, mdl,"DELTA_S_STRAT",delta_s_strat, &
85  "Top-to-bottom salinity difference of stratification", &
86  units="1e-3", fail_if_missing=.not.just_read, do_not_log=just_read)
87  call get_param(param_file, mdl,"ADJUSTMENT_DELTAS",adjustment_deltas, &
88  "Salinity difference across front", &
89  units="1e-3", fail_if_missing=.not.just_read, do_not_log=just_read)
90  call get_param(param_file, mdl,"FRONT_WAVE_AMP",front_wave_amp, &
91  "Amplitude of trans-frontal wave perturbation", &
92  units="same as x,y", default=0., do_not_log=just_read)
93  call get_param(param_file, mdl,"FRONT_WAVE_LENGTH",front_wave_length, &
94  "Wave-length of trans-frontal wave perturbation", &
95  units="same as x,y", default=0., do_not_log=just_read)
96  call get_param(param_file, mdl,"FRONT_WAVE_ASYM",front_wave_asym, &
97  "Amplitude of frontal asymmetric perturbation", &
98  units="same as x,y", default=0., do_not_log=just_read)
99 
100  if (just_read) return ! All run-time parameters have been read, so return.
101 
102  ! WARNING: this routine specifies the interface heights so that the last layer
103  ! is vanished, even at maximum depth. In order to have a uniform
104  ! layer distribution, use this line of code within the loop:
105  ! e0(k) = -G%max_depth * real(k-1) / real(nz)
106  ! To obtain a thickness distribution where the last layer is
107  ! vanished and the other thicknesses uniformly distributed, use:
108  ! e0(k) = -G%max_depth * real(k-1) / real(nz-1)
109 
110  dsdz = -delta_s_strat / g%max_depth
111 
112  select case ( coordinatemode(verticalcoordinate) )
113 
114  case ( regridding_layer, regridding_rho )
115  drho_ds = 1.0 * us%kg_m3_to_R
116  if (delta_s_strat /= 0.) then
117  ! This was previously coded ambiguously.
118  adjustment_delta = (adjustment_deltas / delta_s_strat) * g%max_depth
119  do k=1,nz+1
120  e0(k) = adjustment_delta - (g%max_depth + 2*adjustment_delta) * (real(k-1) / real(nz))
121  enddo
122  else
123  adjustment_delta = 2.*g%max_depth
124  do k=1,nz+1
125  e0(k) = -g%max_depth * (real(k-1) / real(nz))
126  enddo
127  endif
128  if (nz > 1) then
129  target_values(1) = ( gv%Rlay(1) + 0.5*(gv%Rlay(1)-gv%Rlay(2)) )
130  target_values(nz+1) = ( gv%Rlay(nz) + 0.5*(gv%Rlay(nz)-gv%Rlay(nz-1)) )
131  else ! This might not be needed, but it avoids segmentation faults if nz=1.
132  target_values(1) = 0.0
133  target_values(nz+1) = 2.0 * gv%Rlay(1)
134  endif
135  do k = 2,nz
136  target_values(k) = target_values(k-1) + ( gv%Rlay(nz) - gv%Rlay(1) ) / (nz-1)
137  enddo
138  target_values(:) = target_values(:) - 1000.*us%kg_m3_to_R
139  do j=js,je ; do i=is,ie
140  if (front_wave_length /= 0.) then
141  y = ( 0.125 + g%geoLatT(i,j) / front_wave_length ) * ( 4. * acos(0.) )
142  yy = 2. * ( g%geoLatT(i,j) - 0.5 * g%len_lat ) / adjustment_width
143  yy = min(1.0, yy); yy = max(-1.0, yy)
144  yy = yy * 2. * acos( 0. )
145  y = front_wave_amp*sin(y) + front_wave_asym*sin(yy)
146  else
147  y = 0.
148  endif
149  x = ( ( g%geoLonT(i,j) - 0.5 * g%len_lon ) + y ) / adjustment_width
150  x = min(1.0, x); x = max(-1.0, x)
151  x = x * acos( 0. )
152  delta_s = adjustment_deltas * 0.5 * (1. - sin( x ) )
153  do k=2,nz
154  if (drho_ds*dsdz /= 0.) then
155  eta1d(k) = ( target_values(k) - drho_ds*( s_ref + delta_s ) ) / (drho_ds*dsdz)
156  else
157  eta1d(k) = e0(k) - (0.5*adjustment_delta) * sin( x )
158  endif
159  eta1d(k) = max( eta1d(k), -g%max_depth )
160  eta1d(k) = min( eta1d(k), 0. )
161  enddo
162  eta1d(1) = 0.; eta1d(nz+1) = -g%max_depth
163  do k=nz,1,-1
164  if (eta1d(k) > 0.) then
165  eta1d(k) = max( eta1d(k+1) + min_thickness, 0. )
166  h(i,j,k) = gv%Z_to_H * max( eta1d(k) - eta1d(k+1), min_thickness )
167  elseif (eta1d(k) <= (eta1d(k+1) + min_thickness)) then
168  eta1d(k) = eta1d(k+1) + min_thickness
169  h(i,j,k) = gv%Z_to_H * min_thickness
170  else
171  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
172  endif
173  enddo
174  enddo ; enddo
175 
176  case ( regridding_zstar, regridding_sigma )
177  do k=1,nz+1
178  eta1d(k) = -g%max_depth * (real(k-1) / real(nz))
179  eta1d(k) = max(min(eta1d(k), 0.), -g%max_depth)
180  enddo
181  do j=js,je ; do i=is,ie
182  do k=nz,1,-1
183  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
184  enddo
185  enddo ; enddo
186 
187  case default
188  call mom_error(fatal,"adjustment_initialize_thickness: "// &
189  "Unrecognized i.c. setup - set ADJUSTMENT_IC")
190 
191  end select
192