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
11 use mom_domains,
only : root_pe, to_all, scalar_pair, cgrid_ne, agrid
27 use mom_eos,
only : int_specific_vol_dp
28 use mom_ale,
only : ale_remap_scalar
30 implicit none ;
private
32 #include <MOM_memory.h>
34 public :: mom_initialize_tracer_from_z
41 character(len=40) :: mdl =
"MOM_tracer_initialization_from_Z"
46 subroutine mom_initialize_tracer_from_z(h, tr, G, GV, US, PF, src_file, src_var_nam, &
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
103 call calltree_enter(trim(mdl)//
"(), MOM_state_initialization.F90")
105 call get_param(pf, mdl,
"Z_INIT_HOMOGENIZE", homog, &
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
186 call calltree_leave(trim(mdl)//
'()')
187 call cpu_clock_end(id_clock_routine)
189 end subroutine mom_initialize_tracer_from_z