ROMS
Loading...
Searching...
No Matches
sum_grad_mod Module Reference

Functions/Subroutines

subroutine, public sum_grad (ng, tile, linp, lout)
 
subroutine sum_grad_tile (ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, linp, lout, tl_t_obc, tl_u_obc, tl_v_obc, tl_ubar_obc, tl_vbar_obc, tl_zeta_obc, tl_ustr, tl_vstr, tl_tflux, tl_t, tl_u, tl_v, tl_zeta)
 

Function/Subroutine Documentation

◆ sum_grad()

subroutine, public sum_grad_mod::sum_grad ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout )

Definition at line 25 of file sum_grad.F.

26!***********************************************************************
27!
28 USE mod_param
29# ifdef ADJUST_BOUNDARY
30 USE mod_boundary
31# endif
32# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
33 USE mod_forces
34# endif
35 USE mod_ocean
36!
37! Imported variable declarations.
38!
39 integer, intent(in) :: ng, tile, Linp, Lout
40!
41! Local variable declarations.
42!
43# include "tile.h"
44!
45 CALL sum_grad_tile (ng, tile, &
46 & lbi, ubi, lbj, ubj, lbij, ubij, &
47 & imins, imaxs, jmins, jmaxs, &
48 & linp, lout, &
49# ifdef ADJUST_BOUNDARY
50# ifdef SOLVE3D
51 & boundary(ng) % tl_t_obc, &
52 & boundary(ng) % tl_u_obc, &
53 & boundary(ng) % tl_v_obc, &
54# endif
55 & boundary(ng) % tl_ubar_obc, &
56 & boundary(ng) % tl_vbar_obc, &
57 & boundary(ng) % tl_zeta_obc, &
58# endif
59# ifdef ADJUST_WSTRESS
60 & forces(ng) % tl_ustr, &
61 & forces(ng) % tl_vstr, &
62# endif
63# ifdef SOLVE3D
64# ifdef ADJUST_STFLUX
65 & forces(ng) % tl_tflux, &
66# endif
67 & ocean(ng) % tl_t, &
68 & ocean(ng) % tl_u, &
69 & ocean(ng) % tl_v, &
70# else
71 & ocean(ng) % tl_ubar, &
72 & ocean(ng) % tl_vbar, &
73# endif
74 & ocean(ng) % tl_zeta)
75 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351

References mod_boundary::boundary, mod_forces::forces, mod_ocean::ocean, and sum_grad_tile().

Referenced by i4dvar_mod::analysis(), convolve_mod::error_covariance(), and rbl4dvar_mod::increment().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ sum_grad_tile()

subroutine sum_grad_mod::sum_grad_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
integer, intent(in) lout,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) tl_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_zeta )
private

Definition at line 79 of file sum_grad.F.

102!***********************************************************************
103!
104 USE mod_param
105 USE mod_ncparam
106# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS || \
107 defined adjust_boundary
108 USE mod_scalars
109# endif
110# ifdef ADJUST_BOUNDARY
111 USE mod_boundary
112# endif
113!
114! Imported variable declarations.
115!
116 integer, intent(in) :: ng, tile
117 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
118 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
119 integer, intent(in) :: Linp, Lout
120!
121# ifdef ASSUMED_SHAPE
122# ifdef ADJUST_BOUNDARY
123# ifdef SOLVE3D
124 real(r8), intent(inout) :: tl_t_obc(LBij:,:,:,:,:,:)
125 real(r8), intent(inout) :: tl_u_obc(LBij:,:,:,:,:)
126 real(r8), intent(inout) :: tl_v_obc(LBij:,:,:,:,:)
127# endif
128 real(r8), intent(inout) :: tl_ubar_obc(LBij:,:,:,:)
129 real(r8), intent(inout) :: tl_vbar_obc(LBij:,:,:,:)
130 real(r8), intent(inout) :: tl_zeta_obc(LBij:,:,:,:)
131# endif
132# ifdef ADJUST_WSTRESS
133 real(r8), intent(inout) :: tl_ustr(LBi:,LBj:,:,:)
134 real(r8), intent(inout) :: tl_vstr(LBi:,LBj:,:,:)
135# endif
136# ifdef SOLVE3D
137# ifdef ADJUST_STFLUX
138 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
139# endif
140 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
141 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
142 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
143# else
144 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
145 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
146# endif
147 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
148# else
149# ifdef ADJUST_BOUNDARY
150# ifdef SOLVE3D
151 real(r8), intent(inout) :: tl_t_obc(LBij:UBij,N(ng),4, &
152 & Nbrec(ng),2,NT(ng))
153 real(r8), intent(inout) :: tl_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
154 real(r8), intent(inout) :: tl_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
155# endif
156 real(r8), intent(inout) :: tl_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
157 real(r8), intent(inout) :: tl_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
158 real(r8), intent(inout) :: tl_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
159# endif
160# ifdef ADJUST_WSTRESS
161 real(r8), intent(inout) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
162 real(r8), intent(inout) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
163# endif
164# ifdef SOLVE3D
165# ifdef ADJUST_STFLUX
166 real(r8), intent(inout) :: tl_tflux(LBi:UBi,LBj:UBj, &
167 & Nfrec(ng),2,NT(ng))
168# endif
169 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
170 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
171 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
172# else
173 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
174 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
175# endif
176 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
177# endif
178!
179! Local variable declarations.
180!
181 integer :: i, ib, ir, j, k
182# ifdef SOLVE3D
183 integer :: itrc
184# endif
185
186# include "set_bounds.h"
187!
188!-----------------------------------------------------------------------
189! Sum of the background cost function gradients (v-space).
190!-----------------------------------------------------------------------
191!
192! Free-surface.
193!
194 DO j=jstrr,jendr
195 DO i=istrr,iendr
196 tl_zeta(i,j,lout)=tl_zeta(i,j,linp)+ &
197 & tl_zeta(i,j,lout)
198 END DO
199 END DO
200
201# ifdef ADJUST_BOUNDARY
202!
203! Free-surface open boundaries.
204!
205 IF (any(lobc(:,isfsur,ng))) THEN
206 DO ir=1,nbrec(ng)
207 IF ((lobc(iwest,isfsur,ng)).and. &
208 & domain(ng)%Western_Edge(tile)) THEN
209 ib=iwest
210 DO j=jstr,jend
211 tl_zeta_obc(j,ib,ir,lout)=tl_zeta_obc(j,ib,ir,linp)+ &
212 & tl_zeta_obc(j,ib,ir,lout)
213 END DO
214 END IF
215 IF ((lobc(ieast,isfsur,ng)).and. &
216 & domain(ng)%Eastern_Edge(tile)) THEN
217 ib=ieast
218 DO j=jstr,jend
219 tl_zeta_obc(j,ib,ir,lout)=tl_zeta_obc(j,ib,ir,linp)+ &
220 & tl_zeta_obc(j,ib,ir,lout)
221 END DO
222 END IF
223 IF ((lobc(isouth,isfsur,ng)).and. &
224 & domain(ng)%Southern_Edge(tile)) THEN
225 ib=isouth
226 DO i=istr,iend
227 tl_zeta_obc(i,ib,ir,lout)=tl_zeta_obc(i,ib,ir,linp)+ &
228 & tl_zeta_obc(i,ib,ir,lout)
229 END DO
230 END IF
231 IF ((lobc(inorth,isfsur,ng)).and. &
232 & domain(ng)%Northern_Edge(tile)) THEN
233 ib=inorth
234 DO i=istr,iend
235 tl_zeta_obc(i,ib,ir,lout)=tl_zeta_obc(i,ib,ir,linp)+ &
236 & tl_zeta_obc(i,ib,ir,lout)
237 END DO
238 END IF
239 END DO
240 END IF
241# endif
242
243# ifndef SOLVE3D
244!
245! 2D U-momentum component.
246!
247 DO j=jstrr,jendr
248 DO i=istr,iendr
249 tl_ubar(i,j,lout)=tl_ubar(i,j,linp)+ &
250 & tl_ubar(i,j,lout)
251 END DO
252 END DO
253# endif
254
255# ifdef ADJUST_BOUNDARY
256!
257! 2D U-momentum open boundaries.
258!
259 IF (any(lobc(:,isubar,ng))) THEN
260 DO ir=1,nbrec(ng)
261 IF ((lobc(iwest,isubar,ng)).and. &
262 & domain(ng)%Western_Edge(tile)) THEN
263 ib=iwest
264 DO j=jstr,jend
265 tl_ubar_obc(j,ib,ir,lout)=tl_ubar_obc(j,ib,ir,linp)+ &
266 & tl_ubar_obc(j,ib,ir,lout)
267 END DO
268 END IF
269 IF ((lobc(ieast,isubar,ng)).and. &
270 & domain(ng)%Eastern_Edge(tile)) THEN
271 ib=ieast
272 DO j=jstr,jend
273 tl_ubar_obc(j,ib,ir,lout)=tl_ubar_obc(j,ib,ir,linp)+ &
274 & tl_ubar_obc(j,ib,ir,lout)
275 END DO
276 END IF
277 IF ((lobc(isouth,isubar,ng)).and. &
278 & domain(ng)%Southern_Edge(tile)) THEN
279 ib=isouth
280 DO i=istru,iend
281 tl_ubar_obc(i,ib,ir,lout)=tl_ubar_obc(i,ib,ir,linp)+ &
282 & tl_ubar_obc(i,ib,ir,lout)
283 END DO
284 END IF
285 IF ((lobc(inorth,isubar,ng)).and. &
286 & domain(ng)%Northern_Edge(tile)) THEN
287 ib=inorth
288 DO i=istru,iend
289 tl_ubar_obc(i,ib,ir,lout)=tl_ubar_obc(i,ib,ir,linp)+ &
290 & tl_ubar_obc(i,ib,ir,lout)
291 END DO
292 END IF
293 END DO
294 END IF
295# endif
296
297# ifndef SOLVE3D
298!
299! 2D V-momentum.
300!
301 DO j=jstr,jendr
302 DO i=istrr,iendr
303 tl_vbar(i,j,lout)=tl_vbar(i,j,linp)+ &
304 & tl_vbar(i,j,lout)
305 END DO
306 END DO
307# endif
308
309# ifdef ADJUST_BOUNDARY
310!
311! 2D V-momentum open boundaries.
312!
313 IF (any(lobc(:,isvbar,ng))) THEN
314 DO ir=1,nbrec(ng)
315 IF ((lobc(iwest,isvbar,ng)).and. &
316 & domain(ng)%Western_Edge(tile)) THEN
317 ib=iwest
318 DO j=jstrv,jend
319 tl_vbar_obc(j,ib,ir,lout)=tl_vbar_obc(j,ib,ir,linp)+ &
320 & tl_vbar_obc(j,ib,ir,lout)
321 END DO
322 END IF
323 IF ((lobc(ieast,isvbar,ng)).and. &
324 & domain(ng)%Eastern_Edge(tile)) THEN
325 ib=ieast
326 DO j=jstrv,jend
327 tl_vbar_obc(j,ib,ir,lout)=tl_vbar_obc(j,ib,ir,linp)+ &
328 & tl_vbar_obc(j,ib,ir,lout)
329 END DO
330 END IF
331 IF ((lobc(isouth,isvbar,ng)).and. &
332 & domain(ng)%Southern_Edge(tile)) THEN
333 ib=isouth
334 DO i=istr,iend
335 tl_vbar_obc(i,ib,ir,lout)=tl_vbar_obc(i,ib,ir,linp)+ &
336 & tl_vbar_obc(i,ib,ir,lout)
337 END DO
338 END IF
339 IF ((lobc(inorth,isvbar,ng)).and. &
340 & domain(ng)%Northern_Edge(tile)) THEN
341 ib=inorth
342 DO i=istr,iend
343 tl_vbar_obc(i,ib,ir,lout)=tl_vbar_obc(i,ib,ir,linp)+ &
344 & tl_vbar_obc(i,ib,ir,lout)
345 END DO
346 END IF
347 END DO
348 END IF
349# endif
350
351# ifdef ADJUST_WSTRESS
352!
353! Surface momentum stress.
354!
355 DO k=1,nfrec(ng)
356 DO j=jstrr,jendr
357 DO i=istr,iendr
358 tl_ustr(i,j,k,lout)=tl_ustr(i,j,k,linp)+ &
359 & tl_ustr(i,j,k,lout)
360 END DO
361 END DO
362 DO j=jstr,jendr
363 DO i=istrr,iendr
364 tl_vstr(i,j,k,lout)=tl_vstr(i,j,k,linp)+ &
365 & tl_vstr(i,j,k,lout)
366 END DO
367 END DO
368 END DO
369# endif
370
371# ifdef SOLVE3D
372!
373! 3D U-momentum component.
374!
375 DO k=1,n(ng)
376 DO j=jstrr,jendr
377 DO i=istr,iendr
378 tl_u(i,j,k,lout)=tl_u(i,j,k,linp)+ &
379 & tl_u(i,j,k,lout)
380 END DO
381 END DO
382 END DO
383
384# ifdef ADJUST_BOUNDARY
385!
386! 3D U-momentum open boundaries.
387!
388 IF (any(lobc(:,isuvel,ng))) THEN
389 DO ir=1,nbrec(ng)
390 IF ((lobc(iwest,isuvel,ng)).and. &
391 & domain(ng)%Western_Edge(tile)) THEN
392 ib=iwest
393 DO k=1,n(ng)
394 DO j=jstr,jend
395 tl_u_obc(j,k,ib,ir,lout)=tl_u_obc(j,k,ib,ir,linp)+ &
396 & tl_u_obc(j,k,ib,ir,lout)
397 END DO
398 END DO
399 END IF
400 IF ((lobc(ieast,isuvel,ng)).and. &
401 & domain(ng)%Eastern_Edge(tile)) THEN
402 ib=ieast
403 DO k=1,n(ng)
404 DO j=jstr,jend
405 tl_u_obc(j,k,ib,ir,lout)=tl_u_obc(j,k,ib,ir,linp)+ &
406 & tl_u_obc(j,k,ib,ir,lout)
407 END DO
408 END DO
409 END IF
410 IF ((lobc(isouth,isuvel,ng)).and. &
411 & domain(ng)%Southern_Edge(tile)) THEN
412 ib=isouth
413 DO k=1,n(ng)
414 DO i=istru,iend
415 tl_u_obc(i,k,ib,ir,lout)=tl_u_obc(i,k,ib,ir,linp)+ &
416 & tl_u_obc(i,k,ib,ir,lout)
417 END DO
418 END DO
419 END IF
420 IF ((lobc(inorth,isuvel,ng)).and. &
421 & domain(ng)%Northern_Edge(tile)) THEN
422 ib=inorth
423 DO k=1,n(ng)
424 DO i=istru,iend
425 tl_u_obc(i,k,ib,ir,lout)=tl_u_obc(i,k,ib,ir,linp)+ &
426 & tl_u_obc(i,k,ib,ir,lout)
427 END DO
428 END DO
429 END IF
430 END DO
431 END IF
432# endif
433!
434! 3D V-momentum component.
435!
436 DO k=1,n(ng)
437 DO j=jstr,jendr
438 DO i=istrr,iendr
439 tl_v(i,j,k,lout)=tl_v(i,j,k,linp)+ &
440 & tl_v(i,j,k,lout)
441 END DO
442 END DO
443 END DO
444
445# ifdef ADJUST_BOUNDARY
446!
447! 3D V-momentum open boundaries.
448!
449 IF (any(lobc(:,isvvel,ng))) THEN
450 DO ir=1,nbrec(ng)
451 IF ((lobc(iwest,isvvel,ng)).and. &
452 & domain(ng)%Western_Edge(tile)) THEN
453 ib=iwest
454 DO k=1,n(ng)
455 DO j=jstrv,jend
456 tl_v_obc(j,k,ib,ir,lout)=tl_v_obc(j,k,ib,ir,linp)+ &
457 & tl_v_obc(j,k,ib,ir,lout)
458 END DO
459 END DO
460 END IF
461 IF ((lobc(ieast,isvvel,ng)).and. &
462 & domain(ng)%Eastern_Edge(tile)) THEN
463 ib=ieast
464 DO k=1,n(ng)
465 DO j=jstrv,jend
466 tl_v_obc(j,k,ib,ir,lout)=tl_v_obc(j,k,ib,ir,linp)+ &
467 & tl_v_obc(j,k,ib,ir,lout)
468 END DO
469 END DO
470 END IF
471 IF ((lobc(isouth,isvvel,ng)).and. &
472 & domain(ng)%Southern_Edge(tile)) THEN
473 ib=isouth
474 DO k=1,n(ng)
475 DO i=istr,iend
476 tl_v_obc(i,k,ib,ir,lout)=tl_v_obc(i,k,ib,ir,linp)+ &
477 & tl_v_obc(i,k,ib,ir,lout)
478 END DO
479 END DO
480 END IF
481 IF ((lobc(inorth,isvvel,ng)).and. &
482 & domain(ng)%Northern_Edge(tile)) THEN
483 ib=inorth
484 DO k=1,n(ng)
485 DO i=istr,iend
486 tl_v_obc(i,k,ib,ir,lout)=tl_v_obc(i,k,ib,ir,linp)+ &
487 & tl_v_obc(i,k,ib,ir,lout)
488 END DO
489 END DO
490 END IF
491 END DO
492 END IF
493# endif
494!
495! Tracers.
496!
497 DO itrc=1,nt(ng)
498 DO k=1,n(ng)
499 DO j=jstrr,jendr
500 DO i=istrr,iendr
501 tl_t(i,j,k,lout,itrc)=tl_t(i,j,k,linp,itrc)+ &
502 & tl_t(i,j,k,lout,itrc)
503 END DO
504 END DO
505 END DO
506 END DO
507
508# ifdef ADJUST_BOUNDARY
509!
510! Tracers open boundaries.
511!
512 DO itrc=1,nt(ng)
513 IF (any(lobc(:,istvar(itrc),ng))) THEN
514 DO ir=1,nbrec(ng)
515 IF ((lobc(iwest,istvar(itrc),ng)).and. &
516 & domain(ng)%Western_Edge(tile)) THEN
517 ib=iwest
518 DO k=1,n(ng)
519 DO j=jstr,jend
520 tl_t_obc(j,k,ib,ir,lout,itrc)= &
521 & tl_t_obc(j,k,ib,ir,linp,itrc)+ &
522 & tl_t_obc(j,k,ib,ir,lout,itrc)
523 END DO
524 END DO
525 END IF
526 IF ((lobc(ieast,istvar(itrc),ng)).and. &
527 & domain(ng)%Eastern_Edge(tile)) THEN
528 ib=ieast
529 DO k=1,n(ng)
530 DO j=jstr,jend
531 tl_t_obc(j,k,ib,ir,lout,itrc)= &
532 & tl_t_obc(j,k,ib,ir,linp,itrc)+ &
533 & tl_t_obc(j,k,ib,ir,lout,itrc)
534 END DO
535 END DO
536 END IF
537 IF ((lobc(isouth,istvar(itrc),ng)).and. &
538 & domain(ng)%Southern_Edge(tile)) THEN
539 ib=isouth
540 DO k=1,n(ng)
541 DO i=istr,iend
542 tl_t_obc(i,k,ib,ir,lout,itrc)= &
543 & tl_t_obc(i,k,ib,ir,linp,itrc)+ &
544 & tl_t_obc(i,k,ib,ir,lout,itrc)
545 END DO
546 END DO
547 END IF
548 IF ((lobc(inorth,istvar(itrc),ng)).and. &
549 & domain(ng)%Northern_Edge(tile)) THEN
550 ib=inorth
551 DO k=1,n(ng)
552 DO i=istr,iend
553 tl_t_obc(i,k,ib,ir,lout,itrc)= &
554 & tl_t_obc(i,k,ib,ir,linp,itrc)+ &
555 & tl_t_obc(i,k,ib,ir,lout,itrc)
556 END DO
557 END DO
558 END IF
559 END DO
560 END IF
561 END DO
562# endif
563# ifdef ADJUST_STFLUX
564!
565! Surface tracers flux.
566!
567 DO itrc=1,nt(ng)
568 IF (lstflux(itrc,ng)) THEN
569 DO k=1,nfrec(ng)
570 DO j=jstrr,jendr
571 DO i=istrr,iendr
572 tl_tflux(i,j,k,lout,itrc)=tl_tflux(i,j,k,linp,itrc)+ &
573 & tl_tflux(i,j,k,lout,itrc)
574 END DO
575 END DO
576 END DO
577 END IF
578 END DO
579# endif
580# endif
581
582 RETURN
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isubar
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable nt
Definition mod_param.F:489
logical, dimension(:,:,:), allocatable lobc
integer, parameter iwest
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nbrec

References mod_param::domain, mod_scalars::ieast, mod_scalars::inorth, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::iwest, mod_scalars::lobc, and mod_scalars::lstflux.

Referenced by sum_grad().

Here is the caller graph for this function: