MOM6
RGC_initialization.F90
Go to the documentation of this file.
2 !***********************************************************************
3 !* GNU General Public License *
4 !* This file is a part of MOM. *
5 !* *
6 !* MOM is free software; you can redistribute it and/or modify it and *
7 !* are expected to follow the terms of the GNU General Public License *
8 !* as published by the Free Software Foundation; either version 2 of *
9 !* the License, or (at your option) any later version. *
10 !* *
11 !* MOM is distributed in the hope that it will be useful, but WITHOUT *
12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY *
13 !* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public *
14 !* License for more details. *
15 !* *
16 !* For the full text of the GNU General Public License, *
17 !* write to: Free Software Foundation, Inc., *
18 !* 675 Mass Ave, Cambridge, MA 02139, USA. *
19 !* or see: http://www.gnu.org/licenses/gpl.html *
20 !* By Elizabeth Yankovsky, May 2018 *
21 !***********************************************************************
22 
25 use mom_domains, only : pass_var
27 use mom_error_handler, only : mom_mesg, mom_error, fatal, is_root_pe, warning
29 use mom_get_input, only : directories
30 use mom_grid, only : ocean_grid_type
31 use mom_io, only : file_exists, read_data
32 use mom_io, only : slasher
39 implicit none ; private
40 
41 #include <MOM_memory.h>
42 
43 character(len=40) :: mod = "RGC_initialization" ! This module's name.
45 
46 contains
47 
48 !> Sets up the the inverse restoration time, and the values towards which the interface heights,
49 !! velocities and tracers should be restored within the sponges for the RGC test case.
50 subroutine rgc_initialize_sponges(G, GV, US, tv, u, v, PF, use_ALE, CSp, ACSp)
51  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
52  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
53  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
54  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers
55  !! to any available thermodynamic
56  !! fields, potential temperature and
57  !! salinity or mixed layer density.
58  !! Absent fields have NULL ptrs.
59  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
60  target, intent(in) :: u !< Array with the u velocity [L T-1 ~> m s-1]
61  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
62  target, intent(in) :: v !< Array with the v velocity [L T-1 ~> m s-1]
63  type(param_file_type), intent(in) :: pf !< A structure indicating the
64  !! open file to parse for model
65  !! parameter values.
66  logical, intent(in) :: use_ale !< If true, indicates model is in ALE mode
67  type(sponge_cs), pointer :: csp !< Layer-mode sponge structure
68  type(ale_sponge_cs), pointer :: acsp !< ALE-mode sponge structure
69 
70 ! Local variables
71  real :: t(szi_(g),szj_(g),szk_(g)) ! A temporary array for temp
72  real :: s(szi_(g),szj_(g),szk_(g)) ! A temporary array for salt
73  real :: u1(szib_(g),szj_(g),szk_(g)) ! A temporary array for u [L T-1 ~> m s-1]
74  real :: v1(szi_(g),szjb_(g),szk_(g)) ! A temporary array for v [L T-1 ~> m s-1]
75  real :: rho(szi_(g),szj_(g),szk_(g)) ! A temporary array for RHO
76  real :: tmp(szi_(g),szj_(g)) ! A temporary array for tracers.
77  real :: h(szi_(g),szj_(g),szk_(g)) ! A temporary array for thickness at h points
78  real :: idamp(szi_(g),szj_(g)) ! The inverse damping rate at h points [s-1].
79  real :: tnudg ! Nudging time scale, days
80  real :: pres(szi_(g)) ! An array of the reference pressure, in Pa
81  real :: e0(szk_(g)+1) ! The resting interface heights, in m, usually !
82  ! negative because it is positive upward. !
83  real :: eta(szi_(g),szj_(g),szk_(g)+1) ! A temporary array for eta.
84  ! positive upward, in m.
85  logical :: sponge_uv ! Nudge velocities (u and v) towards zero
86  real :: min_depth, dummy1, z, delta_h
87  real :: damp, rho_dummy, min_thickness, rho_tmp, xi0
88  real :: lenlat, lenlon, lensponge
89  character(len=40) :: filename, state_file
90  character(len=40) :: temp_var, salt_var, eta_var, inputdir, h_var
91 
92  character(len=40) :: mod = "RGC_initialize_sponges" ! This subroutine's name.
93  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, iscb, iecb, jscb, jecb
94 
95  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
96  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
97  iscb = g%iscB ; iecb = g%iecB; jscb = g%jscB ; jecb = g%jecB
98 
99  call get_param(pf,mod,"MIN_THICKNESS",min_thickness,'Minimum layer thickness',units='m',default=1.e-3)
100 
101  call get_param(pf, mod, "RGC_TNUDG", tnudg, 'Nudging time scale for sponge layers (days)', default=0.0)
102 
103  call get_param(pf, mod, "LENLAT", lenlat, &
104  "The latitudinal or y-direction length of the domain", &
105  fail_if_missing=.true., do_not_log=.true.)
106 
107  call get_param(pf, mod, "LENLON", lenlon, &
108  "The longitudinal or x-direction length of the domain", &
109  fail_if_missing=.true., do_not_log=.true.)
110 
111  call get_param(pf, mod, "LENSPONGE", lensponge, &
112  "The length of the sponge layer (km).", &
113  default=10.0)
114 
115  call get_param(pf, mod, "SPONGE_UV", sponge_uv, &
116  "Nudge velocities (u and v) towards zero in the sponge layer.", &
117  default=.false., do_not_log=.true.)
118 
119  t(:,:,:) = 0.0 ; s(:,:,:) = 0.0 ; idamp(:,:) = 0.0; rho(:,:,:) = 0.0
120 
121  call get_param(pf, mod, "MINIMUM_DEPTH", min_depth, &
122  "The minimum depth of the ocean.", units="m", default=0.0)
123 
124  if (associated(csp)) call mom_error(fatal, &
125  "RGC_initialize_sponges called with an associated control structure.")
126  if (associated(acsp)) call mom_error(fatal, &
127  "RGC_initialize_sponges called with an associated ALE-sponge control structure.")
128 
129  ! Here the inverse damping time, in s-1, is set. Set Idamp to 0 !
130  ! wherever there is no sponge, and the subroutines that are called !
131  ! will automatically set up the sponges only where Idamp is positive!
132  ! and mask2dT is 1.
133 
134  do i=is,ie ; do j=js,je
135  if (g%geoLonT(i,j) <= lensponge) then
136  dummy1 = -(g%geoLonT(i,j))/lensponge + 1.0
137  !damp = 1.0/TNUDG * max(0.0,dummy1)
138  damp = 0.0
139  !write(*,*)'1st, G%geoLonT(i,j), damp',G%geoLonT(i,j), damp
140 
141  elseif (g%geoLonT(i,j) >= (lenlon - lensponge) .AND. g%geoLonT(i,j) <= lenlon) then
142 
143 ! 1 / day
144  dummy1=(g%geoLonT(i,j)-(lenlon - lensponge))/(lensponge)
145  damp = (1.0/tnudg) * max(0.0,dummy1)
146 
147  else ; damp=0.0
148  endif
149 
150 ! convert to 1 / seconds
151  if (g%bathyT(i,j) > min_depth) then
152  idamp(i,j) = damp/86400.0
153  else ; idamp(i,j) = 0.0 ; endif
154  enddo ; enddo
155 
156 
157  ! 1) Read eta, salt and temp from IC file
158  call get_param(pf, mod, "INPUTDIR", inputdir, default=".")
159  inputdir = slasher(inputdir)
160  ! GM: get two different files, one with temp and one with salt values
161  ! this is work around to avoid having wrong values near the surface
162  ! because of the FIT_SALINITY option. To get salt values right in the
163  ! sponge, FIT_SALINITY=False. The oposite is true for temp. One can
164  ! combined the *correct* temp and salt values in one file instead.
165  call get_param(pf, mod, "RGC_SPONGE_FILE", state_file, &
166  "The name of the file with temps., salts. and interfaces to \n"// &
167  " damp toward.", fail_if_missing=.true.)
168  call get_param(pf, mod, "SPONGE_PTEMP_VAR", temp_var, &
169  "The name of the potential temperature variable in \n"//&
170  "SPONGE_STATE_FILE.", default="Temp")
171  call get_param(pf, mod, "SPONGE_SALT_VAR", salt_var, &
172  "The name of the salinity variable in \n"//&
173  "SPONGE_STATE_FILE.", default="Salt")
174  call get_param(pf, mod, "SPONGE_ETA_VAR", eta_var, &
175  "The name of the interface height variable in \n"//&
176  "SPONGE_STATE_FILE.", default="eta")
177  call get_param(pf, mod, "SPONGE_H_VAR", h_var, &
178  "The name of the layer thickness variable in \n"//&
179  "SPONGE_STATE_FILE.", default="h")
180 
181  !read temp and eta
182  filename = trim(inputdir)//trim(state_file)
183  if (.not.file_exists(filename, g%Domain)) &
184  call mom_error(fatal, " RGC_initialize_sponges: Unable to open "//trim(filename))
185  call read_data(filename,temp_var,t(:,:,:), domain=g%Domain%mpp_domain)
186  call read_data(filename,salt_var,s(:,:,:), domain=g%Domain%mpp_domain)
187 
188  if (use_ale) then
189 
190  call read_data(filename,h_var,h(:,:,:), domain=g%Domain%mpp_domain)
191  call pass_var(h, g%domain)
192 
193  !call initialize_ALE_sponge(Idamp, h, nz, G, PF, ACSp)
194  call initialize_ale_sponge(idamp, g, pf, acsp, h, nz)
195 
196  ! The remaining calls to set_up_sponge_field can be in any order. !
197  if ( associated(tv%T) ) then
198  call set_up_ale_sponge_field(t,g,tv%T,acsp)
199  endif
200  if ( associated(tv%S) ) then
201  call set_up_ale_sponge_field(s,g,tv%S,acsp)
202  endif
203 
204  if (sponge_uv) then
205  u1(:,:,:) = 0.0; v1(:,:,:) = 0.0
206  call set_up_ale_sponge_vel_field(u1,v1,g,u,v,acsp)
207  endif
208 
209 
210  else ! layer mode
211 
212  !read eta
213  call read_data(filename,eta_var,eta(:,:,:), domain=g%Domain%mpp_domain)
214 
215  ! Set the inverse damping rates so that the model will know where to
216  ! apply the sponges, along with the interface heights.
217  call initialize_sponge(idamp, eta, g, pf, csp, gv)
218 
219  if ( gv%nkml>0 ) then
220  ! This call to set_up_sponge_ML_density registers the target values of the
221  ! mixed layer density, which is used in determining which layers can be
222  ! inflated without causing static instabilities.
223  do i=is-1,ie ; pres(i) = tv%P_Ref ; enddo
224 
225  do j=js,je
226  call calculate_density(t(:,j,1), s(:,j,1), pres, tmp(:,j), &
227  is, ie-is+1, tv%eqn_of_state, scale=us%kg_m3_to_R)
228  enddo
229 
230  call set_up_sponge_ml_density(tmp, g, csp)
231  endif
232 
233  ! Apply sponge in tracer fields
234  call set_up_sponge_field(t, tv%T, g, nz, csp)
235  call set_up_sponge_field(s, tv%S, g, nz, csp)
236 
237  endif
238 
239 end subroutine rgc_initialize_sponges
240 
241 end module rgc_initialization
rgc_initialization::mod
character(len=40) mod
Definition: RGC_initialization.F90:43
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
rgc_initialization
Definition: RGC_initialization.F90:1
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
rgc_initialization::rgc_initialize_sponges
subroutine, public rgc_initialize_sponges(G, GV, US, tv, u, v, PF, use_ALE, CSp, ACSp)
Sets up the the inverse restoration time, and the values towards which the interface heights,...
Definition: RGC_initialization.F90:51
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_ale_sponge::initialize_ale_sponge
Ddetermine the number of points which are within sponges in this computational domain.
Definition: MOM_ALE_sponge.F90:48
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_ale_sponge::ale_sponge_cs
ALE sponge control structure.
Definition: MOM_ALE_sponge.F90:84
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_sponge::set_up_sponge_ml_density
subroutine, public set_up_sponge_ml_density(sp_val, G, CS, sp_val_i_mean)
This subroutine stores the reference value for mixed layer density. It is handled differently from ot...
Definition: MOM_sponge.F90:277
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
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_ale_sponge
This module contains the routines used to apply sponge layers when using the ALE mode.
Definition: MOM_ALE_sponge.F90:11
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
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_ale_sponge::set_up_ale_sponge_field
Store the reference profile at h points for a variable.
Definition: MOM_ALE_sponge.F90:33
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
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_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_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
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_ale_sponge::set_up_ale_sponge_vel_field
This subroutine stores the reference profile at u and v points for a vector.
Definition: MOM_ALE_sponge.F90:39
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60