21 use coupler_types_mod,
only : coupler_type_set_data, ind_csurf
24 implicit none ;
private
26 #include <MOM_memory.h>
28 public register_dyed_obc_tracer, initialize_dyed_obc_tracer
29 public dyed_obc_tracer_column_physics, dyed_obc_tracer_end
34 logical :: coupled_tracers = .false.
35 character(len=200) :: tracer_ic_file
36 type(time_type),
pointer :: time => null()
38 real,
pointer :: tr(:,:,:,:) => null()
40 integer,
allocatable,
dimension(:) :: ind_tr
53 function register_dyed_obc_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
63 character(len=80) :: name, longname
65 #include "version_variable.h"
66 character(len=40) :: mdl =
"dyed_obc_tracer"
67 character(len=200) :: inputdir
68 character(len=48) :: flux_units
70 real,
pointer :: tr_ptr(:,:,:) => null()
71 logical :: register_dyed_obc_tracer
72 integer :: isd, ied, jsd, jed, nz, m
73 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
75 if (
associated(cs))
then
76 call mom_error(warning,
"dyed_obc_register_tracer called with an "// &
77 "associated control structure.")
84 call get_param(param_file, mdl,
"NUM_DYE_TRACERS", cs%ntr, &
85 "The number of dye tracers in this run. Each tracer "//&
86 "should have a separate boundary segment.", default=0)
87 allocate(cs%ind_tr(cs%ntr))
88 allocate(cs%tr_desc(cs%ntr))
90 call get_param(param_file, mdl,
"dyed_obc_TRACER_IC_FILE", cs%tracer_IC_file, &
91 "The name of a file from which to read the initial "//&
92 "conditions for the dyed_obc tracers, or blank to initialize "//&
93 "them internally.", default=
" ")
94 if (len_trim(cs%tracer_IC_file) >= 1)
then
95 call get_param(param_file, mdl,
"INPUTDIR", inputdir, default=
".")
96 inputdir = slasher(inputdir)
97 cs%tracer_IC_file = trim(inputdir)//trim(cs%tracer_IC_file)
98 call log_param(param_file, mdl,
"INPUTDIR/dyed_obc_TRACER_IC_FILE", &
102 allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr)) ; cs%tr(:,:,:,:) = 0.0
105 write(name,
'("dye_",I2.2)') m
106 write(longname,
'("Concentration of dyed_obc Tracer ",I2.2)') m
107 cs%tr_desc(m) = var_desc(name, units=
"kg kg-1", longname=longname, caller=mdl)
108 if (gv%Boussinesq)
then ; flux_units =
"kg kg-1 m3 s-1"
109 else ; flux_units =
"kg s-1" ;
endif
113 tr_ptr => cs%tr(:,:,:,m)
115 call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
116 name=name, longname=longname, units=
"kg kg-1", &
117 registry_diags=.true., flux_units=flux_units, &
118 restart_cs=restart_cs)
123 if (cs%coupled_tracers) &
124 cs%ind_tr(m) = aof_set_coupler_flux(trim(name)//
'_flux', &
125 flux_type=
' ', implementation=
' ', caller=
"register_dyed_obc_tracer")
129 cs%restart_CSp => restart_cs
130 register_dyed_obc_tracer = .true.
131 end function register_dyed_obc_tracer
134 subroutine initialize_dyed_obc_tracer(restart, day, G, GV, h, diag, OBC, CS)
137 logical,
intent(in) :: restart
139 type(time_type),
target,
intent(in) :: day
140 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
141 type(
diag_ctrl),
target,
intent(in) :: diag
147 real,
allocatable :: temp(:,:,:)
148 real,
pointer,
dimension(:,:,:) :: &
149 obc_tr1_u => null(), &
153 character(len=24) :: name
154 character(len=72) :: longname
155 character(len=48) :: units
156 character(len=48) :: flux_units
158 real,
pointer :: tr_ptr(:,:,:) => null()
161 real :: e(szk_(g)+1), e_top, e_bot, d_tr
162 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
163 integer :: isdb, iedb, jsdb, jedb
165 if (.not.
associated(cs))
return
166 if (cs%ntr < 1)
return
167 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
168 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
169 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
170 h_neglect = gv%H_subroundoff
175 if (.not.restart)
then
176 if (len_trim(cs%tracer_IC_file) >= 1)
then
178 if (.not.
file_exists(cs%tracer_IC_file, g%Domain)) &
179 call mom_error(fatal,
"dyed_obc_initialize_tracer: Unable to open "// &
182 call query_vardesc(cs%tr_desc(m), name, caller=
"initialize_dyed_obc_tracer")
183 call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
187 do k=1,nz ;
do j=js,je ;
do i=is,ie
189 enddo ;
enddo ;
enddo
194 end subroutine initialize_dyed_obc_tracer
202 subroutine dyed_obc_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, CS, &
203 evap_CFL_limit, minimum_forcing_depth)
206 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
208 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
210 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
214 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
218 type(
forcing),
intent(in) :: fluxes
220 real,
intent(in) :: dt
223 real,
optional,
intent(in) :: evap_cfl_limit
225 real,
optional,
intent(in) :: minimum_forcing_depth
230 real :: c1(szi_(g),szk_(g))
231 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work
232 integer :: i, j, k, is, ie, js, je, nz, m
233 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
235 if (.not.
associated(cs))
return
236 if (cs%ntr < 1)
return
238 if (
present(evap_cfl_limit) .and.
present(minimum_forcing_depth))
then
240 do k=1,nz ;
do j=js,je ;
do i=is,ie
241 h_work(i,j,k) = h_old(i,j,k)
242 enddo ;
enddo ;
enddo
243 call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m) , dt, fluxes, h_work, &
244 evap_cfl_limit, minimum_forcing_depth)
245 if (nz > 1)
call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
249 if (nz > 1)
call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
253 end subroutine dyed_obc_tracer_column_physics
256 subroutine dyed_obc_tracer_end(CS)
261 if (
associated(cs))
then
262 if (
associated(cs%tr))
deallocate(cs%tr)
266 end subroutine dyed_obc_tracer_end