30 use coupler_types_mod,
only : coupler_type_set_data, ind_csurf
33 implicit none ;
private
35 #include <MOM_memory.h>
38 public register_isomip_tracer, initialize_isomip_tracer
39 public isomip_tracer_column_physics, isomip_tracer_surface_state, isomip_tracer_end
41 integer,
parameter :: ntr = 1
45 logical :: coupled_tracers = .false.
46 character(len = 200) :: tracer_ic_file
47 type(time_type),
pointer :: time
49 real,
pointer :: tr(:,:,:,:) => null()
50 real :: land_val(ntr) = -1.0
53 integer,
dimension(NTR) :: ind_tr
67 function register_isomip_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
77 character(len=80) :: name, longname
79 #include "version_variable.h"
80 character(len=40) :: mdl =
"ISOMIP_tracer"
81 character(len=200) :: inputdir
82 character(len=48) :: flux_units
84 real,
pointer :: tr_ptr(:,:,:) => null()
85 logical :: register_isomip_tracer
86 integer :: isd, ied, jsd, jed, nz, m
87 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
89 if (
associated(cs))
then
90 call mom_error(warning,
"ISOMIP_register_tracer called with an "// &
91 "associated control structure.")
98 call get_param(param_file, mdl,
"ISOMIP_TRACER_IC_FILE", cs%tracer_IC_file, &
99 "The name of a file from which to read the initial "//&
100 "conditions for the ISOMIP tracers, or blank to initialize "//&
101 "them internally.", default=
" ")
102 if (len_trim(cs%tracer_IC_file) >= 1)
then
103 call get_param(param_file, mdl,
"INPUTDIR", inputdir, default=
".")
104 inputdir = slasher(inputdir)
105 cs%tracer_IC_file = trim(inputdir)//trim(cs%tracer_IC_file)
106 call log_param(param_file, mdl,
"INPUTDIR/ISOMIP_TRACER_IC_FILE", &
109 call get_param(param_file, mdl,
"SPONGE", cs%use_sponge, &
110 "If true, sponges may be applied anywhere in the domain. "//&
111 "The exact location and properties of those sponges are "//&
112 "specified from MOM_initialization.F90.", default=.false.)
114 allocate(cs%tr(isd:ied,jsd:jed,nz,ntr)) ; cs%tr(:,:,:,:) = 0.0
117 if (m < 10)
then ;
write(name,
'("tr_D",I1.1)') m
118 else ;
write(name,
'("tr_D",I2.2)') m ;
endif
119 write(longname,
'("Concentration of ISOMIP Tracer ",I2.2)') m
120 cs%tr_desc(m) = var_desc(name, units=
"kg kg-1", longname=longname, caller=mdl)
121 if (gv%Boussinesq)
then ; flux_units =
"kg kg-1 m3 s-1"
122 else ; flux_units =
"kg s-1" ;
endif
126 tr_ptr => cs%tr(:,:,:,m)
128 call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
129 name=name, longname=longname, units=
"kg kg-1", &
130 registry_diags=.true., flux_units=flux_units, &
131 restart_cs=restart_cs)
136 if (cs%coupled_tracers) &
137 cs%ind_tr(m) = aof_set_coupler_flux(trim(name)//
'_flux', &
138 flux_type=
' ', implementation=
' ', caller=
"register_ISOMIP_tracer")
142 register_isomip_tracer = .true.
143 end function register_isomip_tracer
147 subroutine initialize_isomip_tracer(restart, day, G, GV, h, diag, OBC, CS, &
152 logical,
intent(in) :: restart
154 type(time_type),
target,
intent(in) :: day
155 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
156 type(
diag_ctrl),
target,
intent(in) :: diag
167 real,
allocatable :: temp(:,:,:)
168 real,
pointer,
dimension(:,:,:) :: &
169 obc_tr1_u => null(), &
173 character(len=16) :: name
174 character(len=72) :: longname
175 character(len=48) :: units
176 character(len=48) :: flux_units
178 real,
pointer :: tr_ptr(:,:,:) => null()
184 real :: e(szk_(g)+1), e_top, e_bot, d_tr
185 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
186 integer :: isdb, iedb, jsdb, jedb
188 if (.not.
associated(cs))
return
189 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
190 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
191 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
192 h_neglect = gv%H_subroundoff
197 if (.not.restart)
then
198 if (len_trim(cs%tracer_IC_file) >= 1)
then
200 if (.not.
file_exists(cs%tracer_IC_file, g%Domain)) &
201 call mom_error(fatal,
"ISOMIP_initialize_tracer: Unable to open "// &
204 call query_vardesc(cs%tr_desc(m), name, caller=
"initialize_ISOMIP_tracer")
205 call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
209 do k=1,nz ;
do j=js,je ;
do i=is,ie
211 enddo ;
enddo ;
enddo
246 end subroutine initialize_isomip_tracer
250 subroutine isomip_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, CS, &
251 evap_CFL_limit, minimum_forcing_depth)
254 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
256 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
258 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
262 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
266 type(
forcing),
intent(in) :: fluxes
268 real,
intent(in) :: dt
271 real,
optional,
intent(in) :: evap_cfl_limit
273 real,
optional,
intent(in) :: minimum_forcing_depth
282 real :: c1(szi_(g),szk_(g))
283 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work
284 real :: melt(szi_(g),szj_(g))
286 character(len=256) :: mesg
287 integer :: i, j, k, is, ie, js, je, nz, m
288 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
290 if (.not.
associated(cs))
return
292 melt(:,:) = fluxes%iceshelf_melt
295 mmax = maxval(melt(is:ie,js:je))
296 call max_across_pes(mmax)
301 do j=js,je ;
do i=is,ie
302 if (melt(i,j) > 0.0)
then
303 cs%tr(i,j,1:2,m) = melt(i,j)/mmax
305 cs%tr(i,j,1:2,m) = 0.0
310 if (
present(evap_cfl_limit) .and.
present(minimum_forcing_depth))
then
312 do k=1,nz ;
do j=js,je ;
do i=is,ie
313 h_work(i,j,k) = h_old(i,j,k)
314 enddo ;
enddo ;
enddo
315 call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m) , dt, fluxes, h_work, &
316 evap_cfl_limit, minimum_forcing_depth)
317 call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
321 call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
325 end subroutine isomip_tracer_column_physics
330 subroutine isomip_tracer_surface_state(state, h, G, CS)
332 type(
surface),
intent(inout) :: state
334 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
342 integer :: m, is, ie, js, je, isd, ied, jsd, jed
343 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
344 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
346 if (.not.
associated(cs))
return
348 if (cs%coupled_tracers)
then
352 call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
353 state%tr_fields, idim=(/isd, is, ie, ied/), &
354 jdim=(/jsd, js, je, jed/) )
358 end subroutine isomip_tracer_surface_state
361 subroutine isomip_tracer_end(CS)
366 if (
associated(cs))
then
367 if (
associated(cs%tr))
deallocate(cs%tr)
370 end subroutine isomip_tracer_end