Go to the documentation of this file.
7 use mom_coms,
only : max_across_pes, min_across_pes
8 use mom_cpu_clock,
only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
30 implicit none ;
private
32 #include <MOM_memory.h>
41 character(len=40) ::
mdl =
"MOM_tracer_initialization_from_Z"
47 src_var_unit_conversion, src_var_record, homogenize, &
48 useALEremapping, remappingScheme, src_var_gridspec )
52 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
54 real,
dimension(:,:,:),
pointer :: tr
56 character(len=*),
intent(in) :: src_file
57 character(len=*),
intent(in) :: src_var_nam
58 real,
optional,
intent(in) :: src_var_unit_conversion
59 integer,
optional,
intent(in) :: src_var_record
60 logical,
optional,
intent(in) :: homogenize
61 logical,
optional,
intent(in) :: usealeremapping
62 character(len=*),
optional,
intent(in) :: remappingscheme
63 character(len=*),
optional,
intent(in) :: src_var_gridspec
66 real :: land_fill = 0.0
67 character(len=200) :: inputdir
68 character(len=200) :: mesg
71 character(len=10) :: remapscheme
72 logical :: homog,useale
75 # include "version_variable.h"
76 character(len=40) ::
mdl =
"MOM_initialize_tracers_from_Z"
78 integer :: is, ie, js, je, nz
79 integer :: isd, ied, jsd, jed
80 integer :: i, j, k, kd
81 real,
allocatable,
dimension(:,:,:),
target :: tr_z, mask_z
82 real,
allocatable,
dimension(:),
target :: z_edges_in, z_in
85 real,
dimension(:,:,:),
allocatable :: hsrc
86 real,
dimension(:),
allocatable :: h1
87 real :: ztopofcell, zbottomofcell, z_bathy
92 integer :: id_clock_routine, id_clock_ale
93 logical :: reentrant_x, tripolar_n
95 id_clock_routine = cpu_clock_id(
'(Initialize tracer from Z)', grain=clock_routine)
96 id_clock_ale = cpu_clock_id(
'(Initialize tracer from Z) ALE', grain=clock_loop)
98 call cpu_clock_begin(id_clock_routine)
100 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
101 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
106 "If True, then horizontally homogenize the interpolated "//&
107 "initial conditions.", default=.false.)
108 call get_param(pf,
mdl,
"Z_INIT_ALE_REMAPPING", useale, &
109 "If True, then remap straight to model coordinate from file.",&
111 call get_param(pf,
mdl,
"Z_INIT_REMAPPING_SCHEME", remapscheme, &
112 "The remapping scheme to use if using Z_INIT_ALE_REMAPPING is True.", &
117 reentrant_x = .false. ;
call get_param(pf,
mdl,
"REENTRANT_X", reentrant_x,default=.true.)
118 tripolar_n = .false. ;
call get_param(pf,
mdl,
"TRIPOLAR_N", tripolar_n, default=.false.)
120 if (
PRESENT(homogenize)) homog=homogenize
121 if (
PRESENT(usealeremapping)) useale=usealeremapping
122 if (
PRESENT(remappingscheme)) remapscheme=remappingscheme
124 if (
PRESENT(src_var_record)) recnum = src_var_record
126 if (
PRESENT(src_var_unit_conversion)) convert = src_var_unit_conversion
129 g, tr_z, mask_z, z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, &
130 homog, m_to_z=us%m_to_Z)
132 kd =
size(z_edges_in,1)-1
139 call cpu_clock_begin(id_clock_ale)
142 allocate( hsrc(isd:ied,jsd:jed,kd) )
143 call initialize_remapping( remapcs, remapscheme, boundary_extrapolation=.false. )
146 do j = js, je ;
do i = is, ie
147 if (g%mask2dT(i,j)>0.)
then
149 ztopofcell = 0. ; zbottomofcell = 0. ; npoints = 0
150 z_bathy = g%bathyT(i,j)
152 if (mask_z(i,j,k) > 0.)
then
153 zbottomofcell = -min( z_edges_in(k+1), z_bathy )
155 zbottomofcell = -z_bathy
157 h1(k) = ztopofcell - zbottomofcell
158 if (h1(k)>0.) npoints = npoints + 1
159 ztopofcell = zbottomofcell
161 h1(kd) = h1(kd) + ( ztopofcell + z_bathy )
165 hsrc(i,j,:) = gv%Z_to_H * h1(:)
168 call ale_remap_scalar(remapcs, g, gv, kd, hsrc, tr_z, h, tr, all_cells=.false. )
174 call mystats(tr(:,:,k), missing_value, is, ie, js, je, k,
'Tracer from ALE()')
176 call cpu_clock_end(id_clock_ale)
180 do k=1,nz ;
do j=js,je ;
do i=is,ie
181 if (tr(i,j,k) == missing_value)
then
184 enddo ;
enddo ;
enddo
187 call cpu_clock_end(id_clock_routine)
subroutine, public setverticalgridaxes(Rlay, GV, scale)
This sets the coordinate data for the "layer mode" of the isopycnal model.
subroutine, public int_specific_vol_dp(T, S, p_t, p_b, alpha_ref, HI, EOS, dza, intp_dza, intx_dza, inty_dza, halo_size, bathyP, dP_tiny, useMassWghtInterp)
Calls the appropriate subroutine to alculate analytical and nearly-analytical integrals in pressure a...
Regridding control structure.
subroutine, public remapping_core_h(CS, n0, h0, u0, n1, h1, u1, h_neglect, h_neglect_edge)
Remaps column of values u0 on grid h0 to grid h1 assuming the top edge is aligned.
subroutine, public ale_remap_scalar(CS, G, GV, nk_src, h_src, s_src, h_dst, s_dst, all_cells, old_remap)
Remaps a single scalar between grids described by thicknesses h_src and h_dst. h_dst must be dimensio...
Extrapolate and interpolate data.
Provides a transparent vertical ocean grid type and supporting routines.
An overloaded interface to log version information about modules.
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Handy functions for manipulating strings.
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Provides subroutines for quantities specific to the equation of state.
Do a halo update on an array.
Container for remapping parameters.
A structure that can be parsed to read and document run-time parameters.
An overloaded interface to read and log the values of various types of parameters.
Initializes hydrography from z-coordinate climatology files.
Horizontal interpolation.
character(len=len(input_string)) function, public uppercase(input_string)
Return a string in which all uppercase letters have been replaced by their lowercase counterparts.
Describes various unit conversion factors.
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Provides column-wise vertical remapping functions.
Do a halo update on a pair of arrays representing the two components of a vector.
A control structure for the equation of state.
integer, parameter, public to_all
A flag for passing in all directions.
Describes the vertical ocean grid, including unit conversion factors.
Describes the decomposed MOM domain and has routines for communications across PEs.
This module contains the main regridding routines.
Provides transparent structures with groups of MOM6 variables and supporting routines.
Wraps the MPP cpu clock functions.
The MOM6 facility to parse input files for runtime parameters.
subroutine, public calltree_leave(mesg)
Writes a message about leaving a subroutine if call tree reporting is active.
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Provides the ocean grid type.
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
character(len=40) mdl
This module's name.
subroutine, public mystats(array, missing, is, ie, js, je, k, mesg)
Write to the terminal some basic statistics about the k-th level of an array.
Provides checksumming functions for debugging.
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...
subroutine, public mom_initialize_tracer_from_z(h, tr, G, GV, US, PF, src_file, src_var_nam, src_var_unit_conversion, src_var_record, homogenize, useALEremapping, remappingScheme, src_var_gridspec)
Initializes a tracer from a z-space data file.
subroutine, public calltree_waypoint(mesg, n)
Writes a message about reaching a milestone if call tree reporting is active.
Generates vertical grids as part of the ALE algorithm.
An overloaded interface to log the values of various types of parameters.
subroutine, public initialize_remapping(CS, remapping_scheme, boundary_extrapolation, check_reconstruction, check_remapping, force_bounds_in_subcell, answers_2018)
Constructor for remapping control structure.
Routines for error handling and I/O management.
subroutine, public calltree_enter(mesg, n)
Writes a message about entering a subroutine if call tree reporting is active.
logical function, public ispointincell(G, i, j, x, y)
Returns true if the coordinates (x,y) are within the h-cell (i,j)
Ocean grid type. See mom_grid for details.
Calculates density of sea water from T, S and P.
An overloaded interface to read various types of parameters.