ROMS
Loading...
Searching...
No Matches
prsgrd44.h
Go to the documentation of this file.
1# undef NEUMANN
2 MODULE prsgrd_mod
3!
4!git $Id$
5!=======================================================================
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md Hernan G. Arango !
9!========================================== Alexander F. Shchepetkin ===
10! !
11! This subroutine evaluates the baroclinic, hydrostatic pressure !
12! gradient term using a finite-volume pressure Jacobian scheme. !
13! The scheme is based on local, conservative, limited-oscillation !
14! vertical quartic polynomial reconstruction of density field and !
15! subsequent projection fits of the derivatives of density into !
16! isosurface of vertical coordinate. The monotonicity constraint !
17! uses a PPM-style limitting algorithm with a power-law slope !
18! reconciliation step. !
19! !
20! isosurface of vertical coordinate. !
21! !
22! The pressure gradient terms (m4/s2) are loaded into right-hand- !
23! side arrays "ru" and "rv". !
24! !
25! Reference: !
26! !
27! Shchepetkin A.F and J.C. McWilliams, 2003: A method for !
28! computing horizontal pressure gradient force in an ocean !
29! model with non-aligned vertical coordinate, JGR, 108, !
30! 1-34. !
31! !
32!=======================================================================
33!
34 implicit none
35!
36 PRIVATE
37 PUBLIC :: prsgrd
38!
39 CONTAINS
40!
41!***********************************************************************
42 SUBROUTINE prsgrd (ng, tile)
43!***********************************************************************
44!
45 USE mod_param
46#ifdef DIAGNOSTICS
47 USE mod_diags
48#endif
49#ifdef ATM_PRESS
50 USE mod_forces
51#endif
52 USE mod_grid
53 USE mod_ocean
54 USE mod_stepping
55!
56! Imported variable declarations.
57!
58 integer, intent(in) :: ng, tile
59!
60! Local variable declarations.
61!
62 character (len=*), parameter :: MyFile = &
63 & __FILE__
64!
65#include "tile.h"
66!
67#ifdef PROFILE
68 CALL wclock_on (ng, inlm, 23, __line__, myfile)
69#endif
70 CALL prsgrd44_tile (ng, tile, &
71 & lbi, ubi, lbj, ubj, &
72 & imins, imaxs, jmins, jmaxs, &
73 & nrhs(ng), &
74#ifdef WET_DRY
75 & grid(ng)%umask_wet, &
76 & grid(ng)%vmask_wet, &
77#endif
78 & grid(ng) % Hz, &
79 & grid(ng) % om_v, &
80 & grid(ng) % on_u, &
81 & grid(ng) % z_w, &
82 & ocean(ng) % rho, &
83#ifdef TIDE_GENERATING_FORCES
84 & ocean(ng) % eq_tide, &
85#endif
86#ifdef WEC_VF
87 & ocean(ng) % zetat, &
88#endif
89#ifdef ATM_PRESS
90 & forces(ng) % Pair, &
91#endif
92#ifdef DIAGNOSTICS_UV
93 & diags(ng) % DiaRU, &
94 & diags(ng) % DiaRV, &
95#endif
96 & ocean(ng) % ru, &
97 & ocean(ng) % rv)
98#ifdef PROFILE
99 CALL wclock_off (ng, inlm, 23, __line__, myfile)
100#endif
101!
102 RETURN
103 END SUBROUTINE prsgrd
104!
105!***********************************************************************
106 SUBROUTINE prsgrd44_tile (ng, tile, &
107 & LBi, UBi, LBj, UBj, &
108 & IminS, ImaxS, JminS, JmaxS, &
109 & nrhs, &
110#ifdef WET_DRY
111 & umask_wet, vmask_wet, &
112#endif
113 & Hz, om_v, on_u, z_w, &
114 & rho, &
115#ifdef TIDE_GENERATING_FORCES
116 & eq_tide, &
117#endif
118#ifdef WEC_VF
119 & zetat, &
120#endif
121#ifdef ATM_PRESS
122 & Pair, &
123#endif
124#ifdef DIAGNOSTICS_UV
125 & DiaRU, DiaRV, &
126#endif
127 & ru, rv)
128!***********************************************************************
129!
130 USE mod_param
131 USE mod_scalars
132!
133! Imported variable declarations.
134!
135 integer, intent(in) :: ng, tile
136 integer, intent(in) :: LBi, UBi, LBj, UBj
137 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
138 integer, intent(in) :: nrhs
139!
140#ifdef ASSUMED_SHAPE
141# ifdef WET_DRY
142 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
143 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
144# endif
145 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
146 real(r8), intent(in) :: om_v(LBi:,LBj:)
147 real(r8), intent(in) :: on_u(LBi:,LBj:)
148 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
149 real(r8), intent(in) :: rho(LBi:,LBj:,:)
150# ifdef TIDE_GENERATING_FORCES
151 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
152# endif
153# ifdef WEC_VF
154 real(r8), intent(in) :: zetat(LBi:,LBj:)
155# endif
156# ifdef ATM_PRESS
157 real(r8), intent(in) :: Pair(LBi:,LBj:)
158# endif
159# ifdef DIAGNOSTICS_UV
160 real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
161 real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
162# endif
163 real(r8), intent(inout) :: ru(LBi:,LBj:,0:,:)
164 real(r8), intent(inout) :: rv(LBi:,LBj:,0:,:)
165#else
166# ifdef WET_DRY
167 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
168 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
169# endif
170 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
171 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
172 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
173 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
174 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
175# ifdef TIDE_GENERATING_FORCES
176 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
177# endif
178# ifdef WEC_VF
179 real(r8), intent(in) :: zetat(LBi:UBi,LBj:UBj)
180# endif
181# ifdef ATM_PRESS
182 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
183# endif
184# ifdef DIAGNOSTICS_UV
185 real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
186 real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
187# endif
188 real(r8), intent(inout) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
189 real(r8), intent(inout) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
190#endif
191!
192! Local variable declarations.
193!
194 integer :: i, j, k
195
196 real(r8), parameter :: eps = 1.0e-8_r8
197
198 real(r8) :: Ampl, Hdd, cff, cff1, cff2, cff3, cffL, cffR
199 real(r8) :: deltaL, deltaR, dh, delP, limtr, rr
200#ifdef ATM_PRESS
201 real(r8) :: OneAtm, fac
202#endif
203 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: P
204 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: r
205 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: d
206
207 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: FX
208
209 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
210 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: aL
211 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: aR
212 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dL
213 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dR
214 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: d1
215 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: r1
216
217#include "set_bounds.h"
218!
219!---------------------------------------------------------------------
220! Finite-volume pressure gradient force algorithm.
221!---------------------------------------------------------------------
222!
223#ifdef ATM_PRESS
224 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
225 fac=100.0_r8/g
226#endif
227 DO j=jstrv-1,jend
228 DO k=n(ng)-1,1,-1
229 DO i=istru-1,iend
230 fc(i,k)=1.0_r8/(hz(i,j,k+1)+hz(i,j,k))
231 r(i,j,k)=fc(i,k)*(rho(i,j,k+1)*hz(i,j,k )+ &
232 & rho(i,j,k )*hz(i,j,k+1))
233 d(i,j,k)=fc(i,k)*(rho(i,j,k+1)-rho(i,j,k))
234 END DO
235 END DO
236!
237! Parabolic WENO reconstruction of density field. Compute left and
238! right side limits aL and aR for the density assuming monotonized
239! quartic polynomial distributions within each grid box. Also
240! compute dL and dR which are used as a measure of quadratic
241! variation during subsquent WENO reconciliation of side limits.
242!
243 DO k=2,n(ng)-1
244 DO i=istru-1,iend
245 deltar=hz(i,j,k)*d(i,j,k )
246 deltal=hz(i,j,k)*d(i,j,k-1)
247 IF ((deltar*deltal).lt.0.0_r8) THEN
248 deltar=0.0_r8
249 deltal=0.0_r8
250 END IF
251 cff=hz(i,j,k-1)+2.0_r8*hz(i,j,k)+hz(i,j,k+1)
252 cffr=cff*d(i,j,k )
253 cffl=cff*d(i,j,k-1)
254 IF (abs(deltar).gt.abs(cffl)) deltar=cffl
255 IF (abs(deltal).gt.abs(cffr)) deltal=cffr
256 cff=(deltar-deltal)/(hz(i,j,k-1)+hz(i,j,k)+hz(i,j,k+1))
257 deltar=deltar-cff*hz(i,j,k+1)
258 deltal=deltal+cff*hz(i,j,k-1)
259 ar(i,k)=rho(i,j,k)+deltar
260 al(i,k)=rho(i,j,k)-deltal
261 dr(i,k)=(2.0_r8*deltar-deltal)**2
262 dl(i,k)=(2.0_r8*deltal-deltar)**2
263 END DO
264 END DO
265!
266 DO i=istru-1,iend
267 al(i,n(ng))=ar(i,n(ng)-1)
268 ar(i,n(ng))=2.0_r8*rho(i,j,n(ng))-al(i,n(ng))
269 dr(i,n(ng))=(2.0_r8*ar(i,n(ng))+al(i,n(ng))- &
270 & 3.0_r8*rho(i,j,n(ng)))**2
271 dl(i,n(ng))=(3.0_r8*rho(i,j,n(ng))- &
272 & 2.0_r8*al(i,n(ng))-ar(i,n(ng)))**2
273 ar(i,1)=al(i,2)
274 al(i,1)=2.0_r8*rho(i,j,1)-ar(i,1)
275 dr(i,1)=(2.0_r8*ar(i,1)+al(i,1)-3.0_r8*rho(i,j,1))**2
276 dl(i,1)=(3.0_r8*rho(i,j,1)-2.0_r8*al(i,1)-ar(i,1))**2
277 END DO
278!
279 DO k=1,n(ng)-1
280 DO i=istru-1,iend
281 deltal=max(dl(i,k ),eps)
282 deltar=max(dr(i,k+1),eps)
283 r1(i,k)=(deltar*ar(i,k)+deltal*al(i,k+1))/(deltar+deltal)
284 END DO
285 END DO
286!
287 DO i=istru-1,iend
288#ifdef NEUMANN
289 r1(i,n(ng))=1.5_r8*rho(i,j,n(ng))-0.5_r8*r1(i,n(ng)-1)
290 r1(i,0)=1.5_r8*rho(i,j,1)-0.5_r8*r1(i,1)
291#else
292 r1(i,n(ng))=2.0_r8*rho(i,j,n(ng))-r1(i,n(ng)-1)
293 r1(i,0)=2.0_r8*rho(i,j,1)-r1(i,1)
294#endif
295 END DO
296!
297! Power-law reconciliation step. It starts with the computation of
298! side limits dR and dL of the first derivative assuming parabolic
299! distributions within each grid box. In this version of the code,
300! before doing so (see "else" branch of 3-way switch below), in the
301! situation when interfacial deviations deltaR and deltaL differ by
302! more than a factor of two (hence monotonic parabolic fit becomes
303! impossible), the parabolic assumption is switched to power-law
304! function, such that its derivative is zero at one end and,
305! consequently, larger than that of (would be) limited parabolic
306! on the other end. The basic parabolic version of the code is
307! commented out, but left here for reference.
308!
309 DO k=1,n(ng)
310 DO i=istru-1,iend
311!! cff=2.0_r8/Hz(i,j,k)
312!! dR(i,k)=cff*(2.0_r8*r1(i,k)+r1(i,k-1)-3.0_r8*rho(i,j,k))
313!! dL(i,k)=cff*(3.0_r8*rho(i,j,k)-2.0_r8*r1(i,k-1)-r1(i,k))
314!! cff=r(i,j,k)-r(i,j,k-1)
315!! if (cff*dR(i,k).lt.0.0_r8) dR(i,k)=0.0_r8
316!! if (cff*dL(i,k).lt.0.0_r8) dL(i,k)=0.0_r8
317 deltar=r1(i,k)-rho(i,j,k)
318 deltal=rho(i,j,k)-r1(i,k-1)
319 cff=deltar*deltal
320 IF (cff.gt.eps) THEN
321 cff=(deltar+deltal)/cff
322 ELSE
323 cff=0.0_r8
324 END IF
325 cffl=cff*deltal
326 cffr=cff*deltar
327 IF (cffl.gt.3.0_r8) THEN
328 cffl=cffl*deltal
329 cffr=0.0_r8
330 ELSE IF (cffr.gt.3.0_r8) THEN
331 cffl=0.0_r8
332 cffr=cffr*deltar
333 ELSE
334 cffl=4.0_r8*deltal-2.0_r8*deltar
335 cffr=4.0_r8*deltar-2.0_r8*deltal
336 END IF
337 cff=1.0_r8/hz(i,j,k)
338 dr(i,k)=cff*cffr
339 dl(i,k)=cff*cffl
340 END DO
341 END DO
342!
343! Compute final value of derivative at each interface by reconciling
344! two side limits dR(k) and dL(k+1) coming from adjacent grid boxes.
345! The difference between these two also causes change of interfacial
346! value r(k) by Ampl. The commented code (left here for reference)
347! computes the exact value of Ampl assuming power law reconciliation
348! and solving associated quadratic equation. The code segment below
349! corresponds to Pade fit to exact solution, which avoids computation
350! of SQRT for the sake of computational efficiency.
351!
352 DO k=n(ng)-1,1,-1
353 DO i=istru-1,iend
354 d(i,j,k)=fc(i,k)*(hz(i,j,k+1)*dl(i,k+1)+hz(i,j,k)*dr(i,k))
355 cffr=8.0_r8*(dr(i,k )+2.0_r8*dl(i,k ))
356 cffl=8.0_r8*(dl(i,k+1)+2.0_r8*dr(i,k+1))
357 IF (abs(d(i,j,k)).gt.abs(cffr)) d(i,j,k)=cffr
358 IF (abs(d(i,j,k)).gt.abs(cffl)) d(i,j,k)=cffl
359 IF ((dl(i,k+1)-dr(i,k))* &
360 & (rho(i,j,k+1)-rho(i,j,k)).gt.0.0_r8) THEN
361 hdd=hz(i,j,k)*(d(i,j,k)-dr(i,k))
362 rr=rho(i,j,k)-r1(i,k-1)
363 ELSE
364 hdd=hz(i,j,k+1)*(dl(i,k+1)-d(i,j,k))
365 rr=r1(i,k+1)-rho(i,j,k+1)
366 END IF
367 rr=abs(rr)
368!! Ampl=0.4_r8*Hdd*rr
369!! Hdd=ABS(Hdd)
370!! cff=rr*(rr+0.16_r8*Hdd)
371!! if (cff.gt.eps) Ampl=Ampl/(rr+sqrt(cff))
372 ampl=0.2_r8*hdd*rr
373 hdd=abs(hdd)
374 cff=rr*rr+0.0763636363636363636_r8*hdd* &
375 & (rr+0.004329004329004329_r8*hdd)
376 IF (cff.gt.eps) THEN
377 ampl=ampl*(rr+0.0363636363636363636_r8*hdd)/cff
378 ELSE
379 ampl=0.0_r8
380 END IF
381 r(i,j,k)=r1(i,k)+ampl
382 END DO
383 END DO
384 DO i=istru-1,iend
385#ifdef NEUMANN
386 r(i,j,0)=1.5_r8*rho(i,j,1)-0.5_r8*r(i,j,1)
387 r(i,j,n(ng))=1.5_r8*rho(i,j,n(ng))-0.5_r8*r(i,j,n(ng)-1)
388 d(i,j,0)=0.0_r8
389 d(i,j,n(ng))=0.0_r8
390#else
391 r(i,j,0)=2.0_r8*rho(i,j,1)-r(i,j,1)
392 r(i,j,n(ng))=2.0_r8*rho(i,j,n(ng))-r(i,j,n(ng)-1)
393 d(i,j,0)=d(i,j,1)
394 d(i,j,n(ng))=d(i,j,n(ng)-1)
395#endif
396 END DO
397!
398! Compute pressure (P) and lateral pressure force (FX). Initialize
399! pressure at the free-surface as zero
400!
401 DO i=istru-1,iend
402 p(i,j,n(ng))=0.0_r8
403#ifdef WEC_VF
404 p(i,j,n(ng))=p(i,j,n(ng))+zetat(i,j)
405#endif
406#ifdef ATM_PRESS
407 p(i,j,n(ng))=p(i,j,n(ng))+fac*(pair(i,j)-oneatm)
408#endif
409#ifdef TIDE_GENERATING_FORCES
410 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
411#endif
412 END DO
413 cff3=1.0_r8/12.0_r8
414 DO k=n(ng),1,-1
415 DO i=istru-1,iend
416 p(i,j,k-1)=p(i,j,k)+hz(i,j,k)*rho(i,j,k)
417 fx(i,j,k)=0.5_r8*hz(i,j,k)* &
418 & (p(i,j,k)+p(i,j,k-1)+ &
419 & 0.2_r8*hz(i,j,k)* &
420 & (r(i,j,k)-r(i,j,k-1)- &
421 & cff3*hz(i,j,k)*(d(i,j,k)+d(i,j,k-1))))
422 END DO
423 END DO
424!
425! Compute net pressure gradient forces in the XI- and ETA-directions.
426!
427 IF (j.ge.jstr) THEN
428 DO i=istru,iend
429 fc(i,n(ng))=0.0_r8
430 END DO
431 cff=0.5_r8*g
432 cff1=g/rho0
433 cff2=1.0_r8/6.0_r8
434 cff3=1.0_r8/12.0_r8
435 DO k=n(ng),1,-1
436 DO i=istru,iend
437 dh=z_w(i,j,k-1)-z_w(i-1,j,k-1)
438 delp=p(i-1,j,k-1)-p(i,j,k-1)
439 rr=0.5_r8*dh*(r(i,j,k-1)+r(i-1,j,k-1)- &
440 & cff2*dh*(d(i,j,k-1)-d(i-1,j,k-1)))
441 limtr=2.0_r8*delp*rr
442 rr=rr*rr+delp*delp
443 IF (limtr.gt.eps*rr) THEN
444 limtr=limtr/rr
445 ELSE
446 limtr=0.0_r8
447 END IF
448 fc(i,k-1)=0.5_r8*dh* &
449 & (p(i,j,k-1)+p(i-1,j,k-1)+ &
450 & limtr*0.2_r8*dh* &
451 & (r(i,j,k-1)-r(i-1,j,k-1)- &
452 & cff3*dh*(d(i,j,k-1)+d(i-1,j,k-1))))
453 ru(i,j,k,nrhs)=(cff*(hz(i-1,j,k)+hz(i,j,k))* &
454 & (z_w(i-1,j,n(ng))-z_w(i,j,n(ng)))+ &
455 & cff1*(fx(i-1,j,k)-fx(i,j,k)+ &
456 & fc(i,k)-fc(i,k-1)))*on_u(i,j)
457#ifdef WET_DRY
458 ru(i,j,k,nrhs)=ru(i,j,k,nrhs)*umask_wet(i,j)
459#endif
460#ifdef DIAGNOSTICS_UV
461 diaru(i,j,k,nrhs,m3pgrd)=ru(i,j,k,nrhs)
462#endif
463 END DO
464 END DO
465 END IF
466!
467 IF (j.ge.jstrv) THEN
468 DO i=istr,iend
469 fc(i,n(ng))=0.0_r8
470 END DO
471 cff=0.5_r8*g
472 cff1=g/rho0
473 cff2=1.0_r8/6.0_r8
474 cff3=1.0_r8/12.0_r8
475 DO k=n(ng),1,-1
476 DO i=istr,iend
477 dh=z_w(i,j,k-1)-z_w(i,j-1,k-1)
478 delp=p(i,j-1,k-1)-p(i,j,k-1)
479 rr=0.5_r8*dh*(r(i,j,k-1)+r(i,j-1,k-1)- &
480 & cff2*dh*(d(i,j,k-1)-d(i,j-1,k-1)))
481 limtr=2.0_r8*delp*rr
482 rr=rr*rr+delp*delp
483 IF (limtr.gt.eps*rr) THEN
484 limtr=limtr/rr
485 ELSE
486 limtr=0.0_r8
487 END IF
488 fc(i,k-1)=0.5_r8*dh* &
489 & (p(i,j,k-1)+p(i,j-1,k-1)+ &
490 & limtr*0.2_r8*dh* &
491 & (r(i,j,k-1)-r(i,j-1,k-1)- &
492 & cff3*dh*(d(i,j,k-1)+d(i,j-1,k-1))))
493 rv(i,j,k,nrhs)=(cff*(hz(i,j-1,k)+hz(i,j,k))* &
494 & (z_w(i,j-1,n(ng))-z_w(i,j,n(ng)))+ &
495 & cff1*(fx(i,j-1,k)-fx(i,j,k)+ &
496 & fc(i,k)-fc(i,k-1)))*om_v(i,j)
497#ifdef WET_DRY
498 rv(i,j,k,nrhs)=rv(i,j,k,nrhs)*vmask_wet(i,j)
499#endif
500#ifdef DIAGNOSTICS_UV
501 diarv(i,j,k,nrhs,m3pgrd)=rv(i,j,k,nrhs)
502#endif
503 END DO
504 END DO
505 END IF
506 END DO
507!
508 RETURN
509 END SUBROUTINE prsgrd44_tile
510
511 END MODULE prsgrd_mod
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
real(dp) g
real(dp) rho0
integer m3pgrd
integer, dimension(:), allocatable nrhs
subroutine prsgrd44_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, umask_wet, vmask_wet, hz, om_v, on_u, z_w, rho, eq_tide, zetat, pair, diaru, diarv, ru, rv)
Definition prsgrd44.h:128
subroutine, public prsgrd(ng, tile)
Definition prsgrd31.h:36
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