21 use coupler_types_mod,
only : coupler_type_set_data, ind_csurf
24 implicit none ;
private
26 #include <MOM_memory.h>
28 public user_register_tracer_example, user_initialize_tracer, user_tracer_stock
29 public tracer_column_physics, user_tracer_surface_state, user_tracer_example_end
31 integer,
parameter :: ntr = 1
35 logical :: coupled_tracers = .false.
36 character(len=200) :: tracer_ic_file
38 type(time_type),
pointer :: time => null()
40 real,
pointer :: tr(:,:,:,:) => null()
41 real :: land_val(ntr) = -1.0
44 integer,
dimension(NTR) :: ind_tr
55 function user_register_tracer_example(HI, GV, param_file, CS, tr_Reg, restart_CS)
67 character(len=80) :: name, longname
69 #include "version_variable.h"
70 character(len=40) :: mdl =
"tracer_example"
71 character(len=200) :: inputdir
72 character(len=48) :: flux_units
74 real,
pointer :: tr_ptr(:,:,:) => null()
75 logical :: user_register_tracer_example
76 integer :: isd, ied, jsd, jed, nz, m
77 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
79 if (
associated(cs))
then
80 call mom_error(warning,
"USER_register_tracer_example called with an "// &
81 "associated control structure.")
88 call get_param(param_file, mdl,
"TRACER_EXAMPLE_IC_FILE", cs%tracer_IC_file, &
89 "The name of a file from which to read the initial "//&
90 "conditions for the DOME tracers, or blank to initialize "//&
91 "them internally.", default=
" ")
92 if (len_trim(cs%tracer_IC_file) >= 1)
then
93 call get_param(param_file, mdl,
"INPUTDIR", inputdir, default=
".")
94 cs%tracer_IC_file = trim(slasher(inputdir))//trim(cs%tracer_IC_file)
95 call log_param(param_file, mdl,
"INPUTDIR/TRACER_EXAMPLE_IC_FILE", &
98 call get_param(param_file, mdl,
"SPONGE", cs%use_sponge, &
99 "If true, sponges may be applied anywhere in the domain. "//&
100 "The exact location and properties of those sponges are "//&
101 "specified from MOM_initialization.F90.", default=.false.)
103 allocate(cs%tr(isd:ied,jsd:jed,nz,ntr)) ; cs%tr(:,:,:,:) = 0.0
106 if (m < 10)
then ;
write(name,
'("tr",I1.1)') m
107 else ;
write(name,
'("tr",I2.2)') m ;
endif
108 write(longname,
'("Concentration of Tracer ",I2.2)') m
109 cs%tr_desc(m) = var_desc(name, units=
"kg kg-1", longname=longname, caller=mdl)
112 if (gv%Boussinesq)
then ; flux_units =
"kg kg-1 m3 s-1"
113 else ; flux_units =
"kg s-1" ;
endif
117 tr_ptr => cs%tr(:,:,:,m)
119 call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
120 name=name, longname=longname, units=
"kg kg-1", &
121 registry_diags=.true., flux_units=flux_units, &
122 restart_cs=restart_cs)
127 if (cs%coupled_tracers) &
128 cs%ind_tr(m) = aof_set_coupler_flux(trim(name)//
'_flux', &
129 flux_type=
' ', implementation=
' ', caller=
"USER_register_tracer_example")
133 user_register_tracer_example = .true.
134 end function user_register_tracer_example
138 subroutine user_initialize_tracer(restart, day, G, GV, h, diag, OBC, CS, &
140 logical,
intent(in) :: restart
142 type(time_type),
target,
intent(in) :: day
145 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
147 type(
diag_ctrl),
target,
intent(in) :: diag
158 real,
allocatable :: temp(:,:,:)
159 character(len=32) :: name
160 character(len=72) :: longname
161 character(len=48) :: units
162 character(len=48) :: flux_units
164 real,
pointer :: tr_ptr(:,:,:) => null()
168 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
169 integer :: isdb, iedb, jsdb, jedb, lntr
171 if (.not.
associated(cs))
return
172 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
173 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
174 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
180 if (.not.restart)
then
181 if (len_trim(cs%tracer_IC_file) >= 1)
then
183 if (.not.
file_exists(cs%tracer_IC_file, g%Domain)) &
184 call mom_error(fatal,
"USER_initialize_tracer: Unable to open "// &
187 call query_vardesc(cs%tr_desc(m), name, caller=
"USER_initialize_tracer")
188 call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
192 do k=1,nz ;
do j=js,je ;
do i=is,ie
193 cs%tr(i,j,k,m) = 1.0e-20
194 enddo ;
enddo ;
enddo
200 dist2 = (g%Rad_Earth * pi / 180.0)**2 * &
201 (g%geoLatT(i,j) - 40.0) * (g%geoLatT(i,j) - 40.0)
202 tr_y = 0.5*exp(-dist2/(1.0e5*1.0e5))
204 do k=1,nz ;
do i=is,ie
206 cs%tr(i,j,k,1) = cs%tr(i,j,k,1) + tr_y
212 if ( cs%use_sponge )
then
217 if (.not.
associated(sponge_csp)) &
218 call mom_error(fatal,
"USER_initialize_tracer: "// &
219 "The pointer to sponge_CSp must be associated if SPONGE is defined.")
221 allocate(temp(g%isd:g%ied,g%jsd:g%jed,nz))
222 do k=1,nz ;
do j=js,je ;
do i=is,ie
223 if (g%geoLatT(i,j) > 700.0 .and. (k > nz/2))
then
228 enddo ;
enddo ;
enddo
234 tr_ptr => cs%tr(:,:,:,m)
235 call set_up_sponge_field(temp, tr_ptr, g, nz, sponge_csp)
240 if (
associated(obc))
then
241 call query_vardesc(cs%tr_desc(1), name, caller=
"USER_initialize_tracer")
242 if (obc%specified_v_BCs_exist_globally)
then
250 call query_vardesc(cs%tr_desc(m), name, caller=
"USER_initialize_tracer")
255 end subroutine user_initialize_tracer
262 subroutine tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, CS)
265 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
267 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
269 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
273 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
277 type(
forcing),
intent(in) :: fluxes
279 real,
intent(in) :: dt
284 real :: hold0(szi_(g))
287 real :: c1(szi_(g),szk_(g))
292 integer :: i, j, k, is, ie, js, je, nz, m
312 data trdc / 1.0,0.0,0.0 /
313 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
315 if (.not.
associated(cs))
return
316 h_neglect = gv%H_subroundoff
323 hold0(i) = h_old(i,j,1)
330 b_denom_1 = h_old(i,j,1) + ea(i,j,1) + h_neglect
331 b1(i) = 1.0 / (b_denom_1 + eb(i,j,1))
333 d1(i) = trdc(1) * (b_denom_1 * b1(i)) + (1.0 - trdc(1))
335 cs%tr(i,j,1,m) = b1(i)*(hold0(i)*cs%tr(i,j,1,m) + trdc(3)*eb(i,j,1))
339 do k=2,nz ;
do i=is,ie
340 c1(i,k) = trdc(1) * eb(i,j,k-1) * b1(i)
341 b_denom_1 = h_old(i,j,k) + d1(i)*ea(i,j,k) + h_neglect
342 b1(i) = 1.0 / (b_denom_1 + eb(i,j,k))
343 d1(i) = trdc(1) * (b_denom_1 * b1(i)) + (1.0 - trdc(1))
345 cs%tr(i,j,k,m) = b1(i) * (h_old(i,j,k)*cs%tr(i,j,k,m) + &
346 ea(i,j,k)*(trdc(1)*cs%tr(i,j,k-1,m)+trdc(2)) + &
350 do m=1,ntr ;
do k=nz-1,1,-1 ;
do i=is,ie
351 cs%tr(i,j,k,m) = cs%tr(i,j,k,m) + c1(i,k+1)*cs%tr(i,j,k+1,m)
352 enddo ;
enddo ;
enddo
355 end subroutine tracer_column_physics
360 function user_tracer_stock(h, stocks, G, GV, CS, names, units, stock_index)
363 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
365 real,
dimension(:),
intent(out) :: stocks
369 character(len=*),
dimension(:),
intent(out) :: names
370 character(len=*),
dimension(:),
intent(out) :: units
371 integer,
optional,
intent(in) :: stock_index
373 integer :: user_tracer_stock
377 integer :: i, j, k, is, ie, js, je, nz, m
378 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
380 user_tracer_stock = 0
381 if (.not.
associated(cs))
return
383 if (
present(stock_index))
then ;
if (stock_index > 0)
then
391 call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller=
"USER_tracer_stock")
392 units(m) = trim(units(m))//
" kg"
394 do k=1,nz ;
do j=js,je ;
do i=is,ie
395 stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
396 (g%mask2dT(i,j) * g%areaT(i,j) * h(i,j,k))
397 enddo ;
enddo ;
enddo
398 stocks(m) = gv%H_to_kg_m2 * stocks(m)
400 user_tracer_stock = ntr
402 end function user_tracer_stock
406 subroutine user_tracer_surface_state(state, h, G, CS)
408 type(
surface),
intent(inout) :: state
410 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
418 integer :: m, is, ie, js, je, isd, ied, jsd, jed
419 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
420 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
422 if (.not.
associated(cs))
return
424 if (cs%coupled_tracers)
then
428 call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
429 state%tr_fields, idim=(/isd, is, ie, ied/), &
430 jdim=(/jsd, js, je, jed/) )
434 end subroutine user_tracer_surface_state
437 subroutine user_tracer_example_end(CS)
442 if (
associated(cs))
then
443 if (
associated(cs%tr))
deallocate(cs%tr)
446 end subroutine user_tracer_example_end