ROMS
Loading...
Searching...
No Matches
ad_t3dmix2_s.h
Go to the documentation of this file.
1 MODULE ad_t3dmix2_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This subroutine computes adjoint horizontal harmonic mixing of !
11! tracers along S-coordinate levels surfaces. !
12! !
13! BASIC STATE variables needed: t, Hz !
14! !
15!=======================================================================
16!
17 implicit none
18!
19 PRIVATE
20 PUBLIC ad_t3dmix2
21!
22 CONTAINS
23!
24!***********************************************************************
25 SUBROUTINE ad_t3dmix2 (ng, tile)
26!***********************************************************************
27!
28 USE mod_param
29#ifdef TS_MIX_CLIMA
30 USE mod_clima
31#endif
32#ifdef DIAGNOSTICS_TS
33!! USE mod_diags
34#endif
35 USE mod_grid
36 USE mod_mixing
37 USE mod_ocean
38 USE mod_stepping
39!
40! Imported variable declarations.
41!
42 integer, intent(in) :: ng, tile
43!
44! Local variable declarations.
45!
46 character (len=*), parameter :: MyFile = &
47 & __FILE__
48!
49#include "tile.h"
50!
51#ifdef PROFILE
52 CALL wclock_on (ng, iadm, 24, __line__, myfile)
53#endif
54 CALL ad_t3dmix2_s_tile (ng, tile, &
55 & lbi, ubi, lbj, ubj, &
56 & imins, imaxs, jmins, jmaxs, &
57 & nrhs(ng), nstp(ng), nnew(ng), &
58#ifdef MASKING
59 & grid(ng) % umask, &
60 & grid(ng) % vmask, &
61#endif
62#ifdef WET_DRY_NOT_YET
63 & grid(ng) % umask_wet, &
64 & grid(ng) % vmask_wet, &
65#endif
66 & grid(ng) % Hz, &
67 & grid(ng) % ad_Hz, &
68 & grid(ng) % pmon_u, &
69 & grid(ng) % pnom_v, &
70 & grid(ng) % pm, &
71 & grid(ng) % pn, &
72#ifdef DIFF_3DCOEF
73 & mixing(ng) % diff3d_r, &
74#else
75 & mixing(ng) % diff2, &
76#endif
77#ifdef TS_MIX_CLIMA
78 & clima(ng) % tclm, &
79#endif
80#ifdef DIAGNOSTICS_TS
81!! & DIAGS(ng) % DiaTwrk, &
82#endif
83 & ocean(ng) % t, &
84 & ocean(ng) % ad_t)
85#ifdef PROFILE
86 CALL wclock_off (ng, iadm, 24, __line__, myfile)
87#endif
88!
89 RETURN
90 END SUBROUTINE ad_t3dmix2
91!
92!***********************************************************************
93 SUBROUTINE ad_t3dmix2_s_tile (ng, tile, &
94 & LBi, UBi, LBj, UBj, &
95 & IminS, ImaxS, JminS, JmaxS, &
96 & nrhs, nstp, nnew, &
97#ifdef MASKING
98 & umask, vmask, &
99#endif
100#ifdef WET_DRY_NOT_YET
101 & umask_wet, vmask_wet, &
102#endif
103 & Hz, ad_Hz, &
104 & pmon_u, pnom_v, pm, pn, &
105#ifdef DIFF_3DCOEF
106 & diff3d_r, &
107#else
108 & diff2, &
109#endif
110#ifdef TS_MIX_CLIMA
111 & tclm, &
112#endif
113#ifdef DIAGNOSTICS_TS
114!! & DiaTwrk, &
115#endif
116 & t, ad_t)
117!***********************************************************************
118!
119 USE mod_param
120 USE mod_scalars
121!
122! Imported variable declarations.
123!
124 integer, intent(in) :: ng, tile
125 integer, intent(in) :: LBi, UBi, LBj, UBj
126 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
127 integer, intent(in) :: nrhs, nstp, nnew
128
129#ifdef ASSUMED_SHAPE
130# ifdef MASKING
131 real(r8), intent(in) :: umask(LBi:,LBj:)
132 real(r8), intent(in) :: vmask(LBi:,LBj:)
133# endif
134# ifdef WET_DRY_NOT_YET
135 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
136 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
137# endif
138# ifdef DIFF_3DCOEF
139 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
140# else
141 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
142# endif
143 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
144 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
145 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
146 real(r8), intent(in) :: pm(LBi:,LBj:)
147 real(r8), intent(in) :: pn(LBi:,LBj:)
148# ifdef TS_MIX_CLIMA
149 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
150# endif
151# ifdef DIAGNOSTICS_TS
152!! real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
153# endif
154 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
155
156 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
157 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
158#else
159# ifdef MASKING
160 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
161 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
162# endif
163# ifdef WET_DRY_NOT_YET
164 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
165 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
166# endif
167# ifdef DIFF_3DCOEF
168 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
169# else
170 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
171# endif
172 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
173 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
174 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
175 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
177# ifdef TS_MIX_CLIMA
178 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
179# endif
180# ifdef DIAGNOSTICS_TS
181!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
182!! & NDT)
183# endif
184 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
185
186 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
187 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
188#endif
189!
190! Local variable declarations.
191!
192 integer :: i, itrc, j, k
193
194 real(r8) :: cff, cff1
195 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4, ad_cff, ad_cff1
196
197 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FE
198 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FX
199
200#include "set_bounds.h"
201!
202!-----------------------------------------------------------------------
203! Initialize adjoint private variables.
204!-----------------------------------------------------------------------
205!
206 ad_cff=0.0_r8
207 ad_cff1=0.0_r8
208
209 ad_fe(imins:imaxs,jmins:jmaxs)=0.0_r8
210 ad_fx(imins:imaxs,jmins:jmaxs)=0.0_r8
211!
212!----------------------------------------------------------------------
213! Compute adjoint horizontal harmonic diffusion along constant
214! S-surfaces.
215!----------------------------------------------------------------------
216!
217 DO itrc=1,nt(ng)
218 DO k=1,n(ng)
219!
220! Time-step harmonic, S-surfaces diffusion term (m Tunits).
221!
222 DO j=jstr,jend
223 DO i=istr,iend
224#ifdef DIAGNOSTICS_TS
225!! DiaTwrk(i,j,k,itrc,iThdif)=cff
226#endif
227!^ tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+tl_cff
228!^
229 ad_cff=ad_cff+ad_t(i,j,k,nnew,itrc)
230!^ tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
231!^ & (tl_FX(i+1,j)-tl_FX(i,j)+ &
232!^ & tl_FE(i,j+1)-tl_FE(i,j))
233!^
234 adfac=dt(ng)*pm(i,j)*pn(i,j)*ad_cff
235 ad_fx(i ,j)=ad_fx(i ,j)-adfac
236 ad_fx(i+1,j)=ad_fx(i+1,j)+adfac
237 ad_fe(i,j )=ad_fe(i,j )-adfac
238 ad_fe(i,j+1)=ad_fe(i,j+1)+adfac
239 ad_cff=0.0_r8
240 END DO
241 END DO
242!
243! Compute XI- and ETA-components of diffusive tracer flux (T m3/s).
244!
245 DO j=jstr,jend+1
246 DO i=istr,iend
247#ifdef DIFF_3DCOEF
248 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
249 & pnom_v(i,j)
250#else
251 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
252 & pnom_v(i,j)
253#endif
254#ifdef WET_DRY_NOT_YET
255 fe(i,j)=fe(i,j)*vmask_wet(i,j)
256#endif
257#ifdef MASKING
258!^ tl_FE(i,j)=tl_FE(i,j)*vmask(i,j)
259!^
260 ad_fe(i,j)=ad_fe(i,j)*vmask(i,j)
261#endif
262#if defined TS_MIX_STABILITY
263!^ tl_FE(i,j)=cff* &
264!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
265!^ & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
266!^ & t(i,j-1,k,nrhs,itrc))+ &
267!^ & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
268!^ & t(i,j-1,k,nstp,itrc)))+ &
269!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
270!^ & (0.75_r8*(tl_t(i,j ,k,nrhs,itrc)- &
271!^ & tl_t(i,j-1,k,nrhs,itrc))+ &
272!^ & 0.25_r8*(tl_t(i,j ,k,nstp,itrc)- &
273!^ & tl_t(i,j-1,k,nstp,itrc))))
274!^
275 adfac=cff*ad_fe(i,j)
276 adfac1=adfac*(0.75_r8*(t(i,j ,k,nrhs,itrc)- &
277 & t(i,j-1,k,nrhs,itrc))+ &
278 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
279 & t(i,j-1,k,nstp,itrc)))
280 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
281 adfac3=adfac2*0.75_r8
282 adfac4=adfac2*0.25_r8
283 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
284 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
285 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac3
286 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac3
287 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)-adfac4
288 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+adfac4
289 ad_fe(i,j)=0.0_r8
290#elif defined TS_MIX_CLIMA
291 IF (ltracerclm(itrc,ng)) THEN
292!^ tl_FE(i,j)=cff* &
293!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
294!^ & ((t(i,j ,k,nrhs,itrc)- &
295!^ & tclm(i,j ,k,itrc))- &
296!^ & (t(i,j-1,k,nrhs,itrc)- &
297!^ & tclm(i,j-1,k,itrc)))+ &
298!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
299!^ & (tl_t(i,j ,k,nrhs,itrc)- &
300!^ & tl_t(i,j-1,k,nrhs,itrc)))
301!^
302 adfac=cff*ad_fe(i,j)
303 adfac1=adfac*((t(i,j ,k,nrhs,itrc)- &
304 & tclm(i,j ,k,itrc))- &
305 & (t(i,j-1,k,nrhs,itrc)- &
306 & tclm(i,j-1,k,itrc)))
307 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
308 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
309 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
310 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
311 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
312 ad_fe(i,j)=0.0_r8
313 ELSE
314!^ tl_FE(i,j)=cff* &
315!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
316!^ & (t(i,j ,k,nrhs,itrc)- &
317!^ & t(i,j-1,k,nrhs,itrc))+ &
318!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
319!^ & (tl_t(i,j ,k,nrhs,itrc)- &
320!^ & tl_t(i,j-1,k,nrhs,itrc)))
321!^
322 adfac=cff*ad_fe(i,j)
323 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
324 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
325 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
326 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
327 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
328 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
329 ad_fe(i,j)=0.0_r8
330 END IF
331#else
332!^ tl_FE(i,j)=cff* &
333!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
334!^ & (t(i,j ,k,nrhs,itrc)- &
335!^ & t(i,j-1,k,nrhs,itrc))+ &
336!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
337!^ & (tl_t(i,j ,k,nrhs,itrc)- &
338!^ & tl_t(i,j-1,k,nrhs,itrc)))
339!^
340 adfac=cff*ad_fe(i,j)
341 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
342 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
343 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
344 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
345 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
346 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
347 ad_fe(i,j)=0.0_r8
348#endif
349 END DO
350 END DO
351 DO j=jstr,jend
352 DO i=istr,iend+1
353#ifdef DIFF_3DCOEF
354 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
355 & pmon_u(i,j)
356#else
357 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
358 & pmon_u(i,j)
359#endif
360#ifdef WET_DRY_NOT_YET
361 fx(i,j)=fx(i,j)*umask_wet(i,j)
362#endif
363#ifdef MASKING
364!^ tl_FX(i,j)=tl_FX(i,j)*umask(i,j)
365!^
366 ad_fx(i,j)=ad_fx(i,j)*umask(i,j)
367#endif
368#if defined TS_MIX_STABILITY
369!^ tl_FX(i,j)=cff* &
370!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
371!^ & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
372!^ & t(i-1,j,k,nrhs,itrc))+ &
373!^ & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
374!^ & t(i-1,j,k,nstp,itrc)))+ &
375!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
376!^ & (0.75_r8*(tl_t(i ,j,k,nrhs,itrc)- &
377!^ & tl_t(i-1,j,k,nrhs,itrc))+ &
378!^ & 0.25_r8*(tl_t(i ,j,k,nstp,itrc)- &
379!^ & tl_t(i-1,j,k,nstp,itrc))))
380!^
381 adfac=cff*ad_fx(i,j)
382 adfac1=adfac*(0.75_r8*(t(i ,j,k,nrhs,itrc)- &
383 & t(i-1,j,k,nrhs,itrc))+ &
384 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
385 & t(i-1,j,k,nstp,itrc)))
386 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
387 adfac3=adfac2*0.75_r8
388 adfac4=adfac2*0.25_r8
389 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
390 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
391 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac3
392 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac3
393 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)-adfac4
394 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+adfac4
395 ad_fx(i,j)=0.0_r8
396#elif defined TS_MIX_CLIMA
397 IF (ltracerclm(itrc,ng)) THEN
398!^ tl_FX(i,j)=cff* &
399!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
400!^ & ((t(i ,j,k,nrhs,itrc)- &
401!^ & tclm(i ,j,k,itrc))- &
402!^ & (t(i-1,j,k,nrhs,itrc)- &
403!^ & tclm(i-1,j,k,itrc)))+ &
404!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
405!^ & (tl_t(i ,j,k,nrhs,itrc)- &
406!^ & tl_t(i-1,j,k,nrhs,itrc)))
407!^
408 adfac=cff*ad_fx(i,j)
409 adfac1=adfac*((t(i ,j,k,nrhs,itrc)- &
410 & tclm(i ,j,k,itrc))- &
411 & (t(i-1,j,k,nrhs,itrc)- &
412 & tclm(i-1,j,k,itrc)))
413 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
414 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
415 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
416 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
417 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
418 ad_fx(i,j)=0.0_r8
419 ELSE
420!^ tl_FX(i,j)=cff* &
421!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
422!^ & (t(i ,j,k,nrhs,itrc)- &
423!^ & t(i-1,j,k,nrhs,itrc))+ &
424!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
425!^ & (tl_t(i ,j,k,nrhs,itrc)- &
426!^ & tl_t(i-1,j,k,nrhs,itrc)))
427!^
428 adfac=cff*ad_fx(i,j)
429 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
430 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
431 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
432 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
433 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
434 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
435 ad_fx(i,j)=0.0_r8
436 END IF
437#else
438!^ tl_FX(i,j)=cff* &
439!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
440!^ & (t(i ,j,k,nrhs,itrc)- &
441!^ & t(i-1,j,k,nrhs,itrc))+ &
442!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
443!^ & (tl_t(i ,j,k,nrhs,itrc)- &
444!^ & tl_t(i-1,j,k,nrhs,itrc)))
445!^
446 adfac=cff*ad_fx(i,j)
447 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
448 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
449 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
450 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
451 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
452 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
453 ad_fx(i,j)=0.0_r8
454#endif
455 END DO
456 END DO
457 END DO
458 END DO
459!
460 RETURN
461 END SUBROUTINE ad_t3dmix2_s_tile
462
463 END MODULE ad_t3dmix2_mod
subroutine ad_t3dmix2_s_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, hz, ad_hz, pmon_u, pnom_v, pm, pn, diff3d_r, diff2, tclm, t, ad_t)
subroutine, public ad_t3dmix2(ng, tile)
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter iadm
Definition mod_param.F:665
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable ltracerclm
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3