ROMS
Loading...
Searching...
No Matches
tl_set_data.F File Reference
#include "cppdefs.h"
#include "tile.h"
#include "set_bounds.h"
Include dependency graph for tl_set_data.F:

Go to the source code of this file.

Functions/Subroutines

subroutine tl_set_data (ng, tile)
 
subroutine tl_set_data_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs)
 

Function/Subroutine Documentation

◆ tl_set_data()

subroutine tl_set_data ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 3 of file tl_set_data.F.

4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This subroutine processes forcing, boundary, climatology, and !
13! other input data. It time-interpolates between snapshots. !
14! !
15!=======================================================================
16!
17 USE mod_param
18!
19! Imported variable declarations.
20!
21 integer, intent(in) :: ng, tile
22!
23! Local variable declarations.
24!
25 character (len=*), parameter :: MyFile = &
26 & __FILE__
27!
28# include "tile.h"
29!
30# ifdef PROFILE
31 CALL wclock_on (ng, itlm, 4, __line__, myfile)
32# endif
33 CALL tl_set_data_tile (ng, tile, &
34 & lbi, ubi, lbj, ubj, &
35 & imins, imaxs, jmins, jmaxs)
36# ifdef PROFILE
37 CALL wclock_off (ng, itlm, 4, __line__, myfile)
38# endif
39!
40 RETURN
integer, parameter itlm
Definition mod_param.F:663
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
subroutine tl_set_data_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs)
Definition tl_set_data.F:47

References mod_param::itlm, tl_set_data_tile(), wclock_off(), and wclock_on().

Referenced by tl_main3d(), and roms_kernel_mod::tlm_initial().

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

◆ tl_set_data_tile()

subroutine tl_set_data_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) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs )

Definition at line 44 of file tl_set_data.F.

47!***********************************************************************
48!
49 USE mod_param
50 USE mod_boundary
51 USE mod_clima
52# if defined FORWARD_READ && defined SOLVE3D
53 USE mod_coupling
54# endif
55 USE mod_forces
56 USE mod_grid
57 USE mod_mixing
58 USE mod_ncparam
59 USE mod_ocean
60 USE mod_stepping
61 USE mod_scalars
62 USE mod_sources
63!
64# ifdef ANALYTICAL
66# endif
69# ifdef SOLVE3D
71# endif
72# ifdef DISTRIBUTE
73# if defined WET_DRY
74 USE distribute_mod, ONLY : mp_boundary
75# endif
77# ifdef SOLVE3D
79# endif
80# endif
81 USE strings_mod, ONLY : founderror
82!
83 implicit none
84!
85! Imported variable declarations.
86!
87 integer, intent(in) :: ng, tile
88 integer, intent(in) :: LBi, UBi, LBj, UBj
89 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
90!
91! Local variable declarations.
92!
93 logical :: SetBC
94 logical :: update = .false.
95# if defined WET_DRY
96 logical :: bry_update
97# endif
98!
99 integer :: ILB, IUB, JLB, JUB
100 integer :: i, ic, itrc, j, k, my_tile
101!
102 real(r8) :: cff, cff1, cff2
103!
104 character (len=*), parameter :: MyFile = &
105 & __FILE__//", tl_set_data_tile"
106
107# include "set_bounds.h"
108!
109! Lower and upper bounds for nontiled (global values) boundary arrays.
110!
111 my_tile=-1 ! for global values
112 ilb=bounds(ng)%LBi(my_tile)
113 iub=bounds(ng)%UBi(my_tile)
114 jlb=bounds(ng)%LBj(my_tile)
115 jub=bounds(ng)%UBj(my_tile)
116
117# ifdef SOLVE3D
118
119# ifdef CLOUDS
120!
121!-----------------------------------------------------------------------
122! Set cloud fraction (nondimensional). Notice that clouds are
123! processed first in case that they are used to adjust shortwave
124! radiation.
125!-----------------------------------------------------------------------
126!
127# ifdef ANA_CLOUD
128 CALL ana_cloud (ng, tile, itlm)
129# else
130 CALL set_2dfld_tile (ng, tile, itlm, idcfra, &
131 & lbi, ubi, lbj, ubj, &
132 & forces(ng)%cloudG, &
133 & forces(ng)%cloud, &
134 & update)
135 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
136# endif
137# endif
138
139# if (defined BULK_FLUXES && !defined FORWARD_FLUXES) || \
140 defined ecosim || \
141 (defined shortwave && defined ana_srflux && defined albedo)
142!
143!-----------------------------------------------------------------------
144! Set surface air temperature (degC).
145!-----------------------------------------------------------------------
146!
147# ifdef ANA_TAIR
148 CALL ana_tair (ng, tile, itlm)
149# else
150 CALL set_2dfld_tile (ng, tile, itlm, idtair, &
151 & lbi, ubi, lbj, ubj, &
152 & forces(ng)%TairG, &
153 & forces(ng)%Tair, &
154 & update)
155 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
156# endif
157# endif
158
159# if (defined BULK_FLUXES && !defined FORWARD_FLUXES) || \
160 defined ecosim || \
161 (defined shortwave && defined ana_srflux && defined albedo)
162!
163!-----------------------------------------------------------------------
164! Set surface air relative or specific humidity.
165!-----------------------------------------------------------------------
166!
167# ifdef ANA_HUMIDITY
168 CALL ana_humid (ng, tile, itlm)
169# else
170 CALL set_2dfld_tile (ng, tile, itlm, idqair, &
171 & lbi, ubi, lbj, ubj, &
172 & forces(ng)%HairG, &
173 & forces(ng)%Hair, &
174 & update)
175 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
176# endif
177# endif
178
179# ifdef SHORTWAVE
180!
181!-----------------------------------------------------------------------
182! Set kinematic surface solar shortwave radiation flux (degC m/s).
183!-----------------------------------------------------------------------
184!
185# ifdef ANA_SRFLUX
186 CALL ana_srflux (ng, tile, itlm)
187# else
188 CALL set_2dfld_tile (ng, tile, itlm, idsrad, &
189 & lbi, ubi, lbj, ubj, &
190 & forces(ng)%srflxG, &
191 & forces(ng)%srflx, &
192 & update)
193 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
194# endif
195
196# if defined DIURNAL_SRFLUX && !defined FORWARD_FLUXES
197!
198! Modulate the averaged shortwave radiation flux by the local diurnal
199! cycle.
200!
201 CALL ana_srflux (ng, tile, itlm)
202# endif
203# endif
204
205# if (defined BULK_FLUXES && !defined FORWARD_FLUXES) && \
206 !defined LONGWAVE && !defined LONGWAVE_OUT
207!
208!-----------------------------------------------------------------------
209! Surface net longwave radiation (degC m/s).
210!-----------------------------------------------------------------------
211!
212 CALL set_2dfld_tile (ng, tile, itlm, idlrad, &
213 & lbi, ubi, lbj, ubj, &
214 & forces(ng)%lrflxG, &
215 & forces(ng)%lrflx, &
216 & update)
217 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
218# endif
219
220# if defined LONGWAVE_OUT && \
221 (defined bulk_fluxes && !defined FORWARD_FLUXES)
222!
223!-----------------------------------------------------------------------
224! Surface downwelling longwave radiation (degC m/s).
225!-----------------------------------------------------------------------
226!
227 CALL set_2dfld_tile (ng, tile, itlm, idldwn, &
228 & lbi, ubi, lbj, ubj, &
229 & forces(ng)%lrflxG, &
230 & forces(ng)%lrflx, &
231 & update)
232 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
233# endif
234
235# if (defined BULK_FLUXES && !defined FORWARD_FLUXES) || \
236 defined ecosim
237!
238!-----------------------------------------------------------------------
239! Set surface winds (m/s).
240!-----------------------------------------------------------------------
241!
242# ifdef ANA_WINDS
243 CALL ana_winds (ng, tile, itlm)
244# else
245 CALL set_2dfld_tile (ng, tile, itlm, iduair, &
246 & lbi, ubi, lbj, ubj, &
247 & forces(ng)%UwindG, &
248 & forces(ng)%Uwind, &
249 & update)
250 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
251
252 CALL set_2dfld_tile (ng, tile, itlm, idvair, &
253 & lbi, ubi, lbj, ubj, &
254 & forces(ng)%VwindG, &
255 & forces(ng)%Vwind, &
256 & update)
257 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
258
259# ifdef CURVGRID
260!
261! If input point surface winds or interpolated from coarse data, rotate
262! to curvilinear grid.
263!
264 IF (.not.linfo(1,iduair,ng).or. &
265 & (iinfo(5,iduair,ng).ne.lm(ng)+2).or. &
266 & (iinfo(6,iduair,ng).ne.mm(ng)+2)) THEN
267 DO j=jstrr,jendr
268 DO i=istrr,iendr
269 cff1=forces(ng)%Uwind(i,j)*grid(ng)%CosAngler(i,j)+ &
270 & forces(ng)%Vwind(i,j)*grid(ng)%SinAngler(i,j)
271 cff2=forces(ng)%Vwind(i,j)*grid(ng)%CosAngler(i,j)- &
272 & forces(ng)%Uwind(i,j)*grid(ng)%SinAngler(i,j)
273 forces(ng)%Uwind(i,j)=cff1
274 forces(ng)%Vwind(i,j)=cff2
275 END DO
276 END DO
277
278 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
279 CALL exchange_r2d_tile (ng, tile, &
280 & lbi, ubi, lbj, ubj, &
281 & forces(ng)%UWind)
282 CALL exchange_r2d_tile (ng, tile, &
283 & lbi, ubi, lbj, ubj, &
284 & forces(ng)%VWind)
285 END IF
286
287# ifdef DISTRIBUTE
288 CALL mp_exchange2d (ng, tile, itlm, 2, &
289 & lbi, ubi, lbj, ubj, &
290 & nghostpoints, &
291 & ewperiodic(ng), nsperiodic(ng), &
292 & forces(ng)%UWind, &
293 & forces(ng)%VWind)
294# endif
295 END IF
296# endif
297# endif
298# endif
299
300# if defined BULK_FLUXES && !defined FORWARD_FLUXES
301!
302!-----------------------------------------------------------------------
303! Set rain fall rate (kg/m2/s).
304!-----------------------------------------------------------------------
305!
306# ifdef ANA_RAIN
307 CALL ana_rain (ng, tile, itlm)
308# else
309 CALL set_2dfld_tile (ng, tile, itlm, idrain, &
310 & lbi, ubi, lbj, ubj, &
311 & forces(ng)%rainG, &
312 & forces(ng)%rain, &
313 & update)
314 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
315# endif
316# endif
317
318# if !defined BULK_FLUXES || defined FORWARD_FLUXES
319!
320!-----------------------------------------------------------------------
321! Set kinematic surface net heat flux (degC m/s).
322!-----------------------------------------------------------------------
323!
324# ifdef ANA_STFLUX
325 CALL ana_stflux (ng, tile, itlm, itemp)
326# else
327 CALL set_2dfld_tile (ng, tile, itlm, idtsur(itemp), &
328 & lbi, ubi, lbj, ubj, &
329 & forces(ng)%stfluxG(:,:,:,itemp), &
330 & forces(ng)%stflux (:,:,itemp), &
331 & update)
332 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
333# endif
334# endif
335
336# ifdef QCORRECTION
337!
338!-----------------------------------------------------------------------
339! Set sea surface temperature (SST) and heat flux sensitivity to
340! SST (dQdSST) which are used for surface heat flux correction.
341!-----------------------------------------------------------------------
342!
343# ifdef ANA_SST
344 CALL ana_sst (ng, tile, itlm)
345# else
346 CALL set_2dfld_tile (ng, tile, itlm, idsstc, &
347 & lbi, ubi, lbj, ubj, &
348 & forces(ng)%sstG, &
349 & forces(ng)%sst, &
350 & update)
351 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
352# endif
353!
354# ifdef ANA_DQDSST
355 CALL ana_dqdsst (ng, tile, itlm)
356# else
357 CALL set_2dfld_tile (ng, tile, itlm, iddqdt, &
358 & lbi, ubi, lbj, ubj, &
359 & forces(ng)%dqdtG, &
360 & forces(ng)%dqdt, &
361 & update)
362 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
363# endif
364# endif
365!
366!-----------------------------------------------------------------------
367! Set kinematic bottom net heat flux (degC m/s).
368!-----------------------------------------------------------------------
369!
370# ifdef ANA_BTFLUX
371 CALL ana_btflux (ng, tile, itlm, itemp)
372# else
373 CALL set_2dfld_tile (ng, tile, itlm, idtbot(itemp), &
374 & lbi, ubi, lbj, ubj, &
375 & forces(ng)%btfluxG(:,:,:,itemp), &
376 & forces(ng)%btflux (:,:,itemp), &
377 & update)
378 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
379# endif
380
381# ifdef SALINITY
382!
383!-----------------------------------------------------------------------
384! Set kinematic surface freshwater (E-P) flux (m/s).
385!-----------------------------------------------------------------------
386!
387# ifdef ANA_SSFLUX
388 CALL ana_stflux (ng, tile, itlm, isalt)
389# else
390# if !(defined EMINUSP || defined FORWARD_FLUXES || \
391 defined frc_coupling || defined srelaxation)
392 CALL set_2dfld_tile (ng, tile, itlm, idsfwf, &
393 & lbi, ubi, lbj, ubj, &
394 & forces(ng)%stfluxG(:,:,:,isalt), &
395 & forces(ng)%stflux (:,:,isalt), &
396 & update)
397 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
398
399# elif (defined EMINUSP || defined FORWARD_FLUXES || \
400 defined frc_coupling)
401 CALL set_2dfld_tile (ng, tile, itlm, idempf, &
402 & lbi, ubi, lbj, ubj, &
403 & forces(ng)%stfluxG(:,:,:,isalt), &
404 & forces(ng)%stflux (:,:,isalt), &
405 & update)
406 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
407# endif
408# endif
409
410# if defined SCORRECTION || defined SRELAXATION
411!
412!-----------------------------------------------------------------------
413! Set surface salinity for freshwater flux correction.
414!-----------------------------------------------------------------------
415!
416# ifdef ANA_SSS
417 CALL ana_sss (ng, tile, itlm)
418# else
419 CALL set_2dfld_tile (ng, tile, itlm, idsssc, &
420 & lbi, ubi, lbj, ubj, &
421 & forces(ng)%sssG, &
422 & forces(ng)%sss, &
423 & update)
424 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
425# endif
426# endif
427!
428!-----------------------------------------------------------------------
429! Set kinematic bottom salt flux (m/s).
430!-----------------------------------------------------------------------
431!
432# ifdef ANA_BSFLUX
433 CALL ana_btflux (ng, tile, itlm, isalt)
434# else
435 CALL set_2dfld_tile (ng, tile, itlm, idtbot(isalt), &
436 & lbi, ubi, lbj, ubj, &
437 & forces(ng)%btfluxG(:,:,:,isalt), &
438 & forces(ng)%btflux (:,:,isalt), &
439 & update)
440 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
441# endif
442# endif
443
444# if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE
445!
446!-----------------------------------------------------------------------
447! Set kinematic surface and bottom passive tracer fluxes (T m/s).
448!-----------------------------------------------------------------------
449!
450 DO itrc=nat+1,nt(ng)
451# ifdef ANA_SPFLUX
452 CALL ana_stflux (ng, tile, itlm, itrc)
453# else
454 CALL set_2dfld_tile (ng, tile, itlm, idtsur(itrc), &
455 & lbi, ubi, lbj, ubj, &
456 & forces(ng)%stfluxG(:,:,:,itrc), &
457 & forces(ng)%stflux (:,:,itrc), &
458 & update)
459 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
460# endif
461
462# ifdef ANA_BPFLUX
463 CALL ana_btflux (ng, tile, itlm, itrc)
464# else
465 CALL set_2dfld_tile (ng, tile, itlm, idtbot(itrc), &
466 & lbi, ubi, lbj, ubj, &
467 & forces(ng)%btfluxG(:,:,:,itrc), &
468 & forces(ng)%btflux (:,:,itrc), &
469 & update)
470 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
471# endif
472 END DO
473# endif
474# endif
475
476# if !defined BULK_FLUXES || defined FORWARD_FLUXES || \
477 defined frc_coupling
478!
479!-----------------------------------------------------------------------
480! Set kinematic surface momentum flux (m2/s2).
481!-----------------------------------------------------------------------
482!
483# ifdef ANA_SMFLUX
484 CALL ana_smflux (ng, tile, itlm)
485# else
486 CALL set_2dfld_tile (ng, tile, itlm, idusms, &
487 & lbi, ubi, lbj, ubj, &
488 & forces(ng)%sustrG, &
489 & forces(ng)%sustr, &
490 & update)
491 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
492!
493 CALL set_2dfld_tile (ng, tile, itlm, idvsms, &
494 & lbi, ubi, lbj, ubj, &
495 & forces(ng)%svstrG, &
496 & forces(ng)%svstr, &
497 & update)
498 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
499
500# ifdef CURVGRID
501!
502! If input point wind stress, rotate to curvilinear grid. Notice
503! that rotation is done at RHO-points. It does not matter.
504!
505 IF (.not.linfo(1,idusms,ng).or. &
506 & (iinfo(5,idusms,ng).ne.lm(ng)+1).or. &
507 & (iinfo(6,idusms,ng).ne.mm(ng)+2)) THEN
508 DO j=jstrr,jendr
509 DO i=istrr,iendr
510 cff1=forces(ng)%sustr(i,j)*grid(ng)%CosAngler(i,j)+ &
511 & forces(ng)%svstr(i,j)*grid(ng)%SinAngler(i,j)
512 cff2=forces(ng)%svstr(i,j)*grid(ng)%CosAngler(i,j)- &
513 & forces(ng)%sustr(i,j)*grid(ng)%SinAngler(i,j)
514 forces(ng)%sustr(i,j)=cff1
515 forces(ng)%svstr(i,j)=cff2
516 END DO
517 END DO
518
519 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
520 CALL exchange_u2d_tile (ng, tile, &
521 & lbi, ubi, lbj, ubj, &
522 & forces(ng)%sustr)
523 CALL exchange_v2d_tile (ng, tile, &
524 & lbi, ubi, lbj, ubj, &
525 & forces(ng)%svstr)
526 END IF
527
528# ifdef DISTRIBUTE
529 CALL mp_exchange2d (ng, tile, itlm, 2, &
530 & lbi, ubi, lbj, ubj, &
531 & nghostpoints, &
532 & ewperiodic(ng), nsperiodic(ng), &
533 & forces(ng)%sustr, &
534 & forces(ng)%svstr)
535# endif
536 END IF
537# endif
538# endif
539# endif
540
541# if (defined BULK_FLUXES && !defined FORWARD_FLUXES) || \
542 defined ecosim || defined atm_press
543!
544!-----------------------------------------------------------------------
545! Set surface air pressure (mb).
546!-----------------------------------------------------------------------
547!
548# ifdef ANA_PAIR
549 CALL ana_pair (ng, tile, itlm)
550# else
551 setbc=.true.
552! SetBC=.FALSE.
553 CALL set_2dfld_tile (ng, tile, itlm, idpair, &
554 & lbi, ubi, lbj, ubj, &
555 & forces(ng)%PairG, &
556 & forces(ng)%Pair, &
557 & update, setbc)
558 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
559# endif
560# endif
561
562# ifdef WAVE_DATA
563!
564!-----------------------------------------------------------------------
565! Set surface wind-induced wave amplitude, direction and period.
566!-----------------------------------------------------------------------
567!
568# ifdef ANA_WWAVE
569 CALL ana_wwave (ng, tile, itlm)
570# else
571# ifdef WAVES_DIR
572 CALL set_2dfld_tile (ng, tile, itlm, idwdir, &
573 & lbi, ubi, lbj, ubj, &
574 & forces(ng)%DwaveG, &
575 & forces(ng)%Dwave, &
576 & update)
577 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
578
579# ifdef CURVGRID
580!
581! If input point-data, rotate direction to curvilinear coordinates.
582!
583 IF (.not.linfo(1,idwdir,ng).or. &
584 & (iinfo(5,idwdir,ng).ne.lm(ng)+2).or. &
585 & (iinfo(6,idwdir,ng).ne.mm(ng)+2)) THEN
586 DO j=jstrr,jendr
587 DO i=istrr,iendr
588 forces(ng)%Dwave(i,j)=forces(ng)%Dwave(i,j)- &
589 & grid(ng)%angler(i,j)
590 END DO
591 END DO
592 END IF
593
594 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
595 CALL exchange_r2d_tile (ng, tile, &
596 & lbi, ubi, lbj, ubj, &
597 & forces(ng)%Dwave)
598 END IF
599
600# ifdef DISTRIBUTE
601 CALL mp_exchange2d (ng, tile, itlm, 1, &
602 & lbi, ubi, lbj, ubj, &
603 & nghostpoints, &
604 & ewperiodic(ng), nsperiodic(ng), &
605 & forces(ng)%Dwave)
606# endif
607# endif
608# endif
609
610# ifdef WAVES_HEIGHT
611!
612 CALL set_2dfld_tile (ng, tile, itlm, idwamp, &
613 & lbi, ubi, lbj, ubj, &
614 & forces(ng)%HwaveG, &
615 & forces(ng)%Hwave, &
616 & update)
617 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
618# endif
619
620# ifdef WAVES_LENGTH
621!
622 CALL set_2dfld_tile (ng, tile, itlm, idwlen, &
623 & lbi, ubi, lbj, ubj, &
624 & forces(ng)%LwaveG, &
625 & forces(ng)%Lwave, &
626 & update)
627 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
628# endif
629
630# ifdef WAVES_TOP_PERIOD
631!
632 CALL set_2dfld_tile (ng, tile, itlm, idwptp, &
633 & lbi, ubi, lbj, ubj, &
634 & forces(ng)%Pwave_topG, &
635 & forces(ng)%Pwave_top, &
636 & update)
637 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
638# endif
639
640# ifdef WAVES_BOT_PERIOD
641!
642 CALL set_2dfld_tile (ng, tile, itlm, idwpbt, &
643 & lbi, ubi, lbj, ubj, &
644 & forces(ng)%Pwave_botG, &
645 & forces(ng)%Pwave_bot, &
646 & update)
647 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
648# endif
649
650# if defined WAVES_UB
651!
652 CALL set_2dfld_tile (ng, tile, itlm, idworb, &
653 & lbi, ubi, lbj, ubj, &
654 & forces(ng)%Ub_swanG, &
655 & forces(ng)%Ub_swan, &
656 & update)
657 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
658# endif
659
660# if defined TKE_WAVEDISS
661!
662 CALL set_2dfld_tile (ng, tile, itlm, idwdis, &
663 & lbi, ubi, lbj, ubj, &
664 & forces(ng)%Wave_dissipG, &
665 & forces(ng)%Wave_dissip, &
666 & update)
667 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
668# endif
669
670# if defined SVENDSEN_ROLLER
671!
672 CALL set_2dfld_tile (ng, tile, itlm, idwbrk, &
673 & lbi, ubi, lbj, ubj, &
674 & forces(ng)%Wave_breakG, &
675 & forces(ng)%Wave_break, &
676 & update)
677 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
678# endif
679# endif
680# endif
681
682# if defined ECOSIM && defined SOLVE3D
683!
684!-----------------------------------------------------------------------
685! Compute spectral irradiance and cosine of average zenith angle of
686! downwelling spectral photons.
687!-----------------------------------------------------------------------
688!
689 CALL ana_specir (ng, tile, itlm)
690# endif
691
692# ifdef ANA_SPINNING
693!
694!-----------------------------------------------------------------------
695! Set time-varying rotation force (centripetal accelerations) for
696! polar coordinate grids.
697!-----------------------------------------------------------------------
698!
699 CALL ana_spinning (ng, tile, itlm)
700# endif
701!
702!-----------------------------------------------------------------------
703! Set point Sources/Sinks (river runoff).
704!-----------------------------------------------------------------------
705!
706# ifdef ANA_PSOURCE
707 IF (luvsrc(ng).or.lwsrc(ng).or.any(ltracersrc(:,ng))) THEN
708 CALL ana_psource (ng, tile itlm)
709 END IF
710# else
711 IF (domain(ng)%SouthWest_Test(tile)) THEN
712 IF (luvsrc(ng).or.lwsrc(ng)) THEN
713 CALL set_ngfld (ng, itlm, idrtra, 1, nsrc(ng), 1, &
714 & 1, nsrc(ng), 1, &
715 & sources(ng) % QbarG, &
716 & sources(ng) % Qbar, &
717 & update)
718 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
719
720# ifdef SOLVE3D
721 DO k=1,n(ng)
722 DO i=1,nsrc(ng)
723 sources(ng)%Qsrc(i,k)=sources(ng)%Qbar(i)* &
724 & sources(ng)%Qshape(i,k)
725 END DO
726 END DO
727# endif
728 END IF
729
730# ifdef SOLVE3D
731 DO itrc=1,nt(ng)
732 IF (ltracersrc(itrc,ng)) THEN
733 CALL set_ngfld (ng, itlm, idrtrc(itrc), 1, nsrc(ng), n(ng), &
734 & 1, nsrc(ng), n(ng), &
735 & sources(ng) % TsrcG(:,:,:,itrc), &
736 & sources(ng) % Tsrc(:,:,itrc), &
737 & update)
738 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
739 END IF
740 END DO
741# endif
742 END IF
743# endif
744!
745!-----------------------------------------------------------------------
746! Set open boundary conditions fields.
747!-----------------------------------------------------------------------
748!
749! Free-surface.
750!
751 IF (lprocessobc(ng)) THEN
752# ifdef ANA_FSOBC
753 CALL ana_fsobc (ng, tile, itlm)
754# else
755 IF (domain(ng)%SouthWest_Test(tile)) THEN
756 IF (tl_lbc(iwest,isfsur,ng)%acquire) THEN
757 CALL set_ngfld (ng, itlm, idzbry(iwest), jlb, jub, 1, &
758 & 0, mm(ng)+1, 1, &
759 & boundary(ng) % zetaG_west, &
760 & boundary(ng) % zeta_west, &
761 & update)
762 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
763 END IF
764
765 IF (tl_lbc(ieast,isfsur,ng)%acquire) THEN
766 CALL set_ngfld (ng, itlm, idzbry(ieast), jlb, jub, 1, &
767 & 0, mm(ng)+1, 1, &
768 & boundary(ng) % zetaG_east, &
769 & boundary(ng) % zeta_east, &
770 & update)
771 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
772 END IF
773
774 IF (tl_lbc(isouth,isfsur,ng)%acquire) THEN
775 CALL set_ngfld (ng, itlm, idzbry(isouth), ilb, iub, 1, &
776 & 0, lm(ng)+1 ,1, &
777 & boundary(ng) % zetaG_south, &
778 & boundary(ng) % zeta_south, &
779 & update)
780 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
781 END IF
782
783 IF (tl_lbc(inorth,isfsur,ng)%acquire) THEN
784 CALL set_ngfld (ng, itlm, idzbry(inorth), ilb, iub, 1, &
785 & 0, lm(ng)+1, 1, &
786 & boundary(ng) % zetaG_north, &
787 & boundary(ng) % zeta_north, &
788 & update)
789 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
790 END IF
791 END IF
792# endif
793
794# if defined WET_DRY
795!
796! Ensure that water level on boundary cells is above bed elevation.
797!
798 IF (tl_lbc(iwest,isfsur,ng)%acquire) THEN
799 bry_update=.false.
800 IF (domain(ng)%Western_Edge(tile)) THEN
801 DO j=jstrr,jendr
802 cff=dcrit(ng)-grid(ng)%h(0,j)
803 IF (boundary(ng)%zeta_west(j).le.cff) THEN
804 boundary(ng)%zeta_west(j)=cff
805 END IF
806 END DO
807 bry_update=.true.
808 END IF
809# ifdef DISTRIBUTE
810 CALL mp_boundary (ng, itlm, jstrr, jendr, jlb, jub, 1, 1, &
811 & bry_update, &
812 & boundary(ng)%zeta_west)
813# endif
814 END IF
815
816 IF (tl_lbc(ieast,isfsur,ng)%acquire) THEN
817 bry_update=.false.
818 IF (domain(ng)%Eastern_Edge(tile)) THEN
819 DO j=jstrr,jendr
820 cff=dcrit(ng)-grid(ng)%h(lm(ng)+1,j)
821 IF (boundary(ng)%zeta_east(j).le.cff) THEN
822 boundary(ng)%zeta_east(j)=cff
823 END IF
824 END DO
825 bry_update=.true.
826 END IF
827# ifdef DISTRIBUTE
828 CALL mp_boundary (ng, itlm, jstrr, jendr, jlb, jub, 1, 1, &
829 & bry_update, &
830 & boundary(ng)%zeta_east)
831# endif
832 END IF
833
834 IF (tl_lbc(isouth,isfsur,ng)%acquire) THEN
835 bry_update=.false.
836 IF (domain(ng)%Southern_Edge(tile)) THEN
837 DO i=istrr,iendr
838 cff=dcrit(ng)-grid(ng)%h(i,0)
839 IF (boundary(ng)%zeta_south(i).le.cff) THEN
840 boundary(ng)%zeta_south(i)=cff
841 END IF
842 END DO
843 bry_update=.true.
844 END IF
845# ifdef DISTRIBUTE
846 CALL mp_boundary (ng, itlm, istrr, iendr, ilb, iub, 1, 1, &
847 & bry_update, &
848 & boundary(ng)%zeta_south)
849# endif
850 END IF
851
852 IF (tl_lbc(inorth,isfsur,ng)%acquire) THEN
853 bry_update=.false.
854 IF (domain(ng)%Northern_Edge(tile)) THEN
855 DO i=istrr,iendr
856 cff=dcrit(ng)-grid(ng)%h(i,mm(ng)+1)
857 IF (boundary(ng)%zeta_north(i).le.cff) THEN
858 boundary(ng)%zeta_north(i)=cff
859 END IF
860 END DO
861 bry_update=.true.
862 END IF
863# ifdef DISTRIBUTE
864 CALL mp_boundary (ng, itlm, istrr, iendr, ilb, iub, 1, 1, &
865 & bry_update, &
866 & boundary(ng)%zeta_north)
867# endif
868 END IF
869# endif
870 END IF
871!
872! 2D momentum.
873!
874 IF (lprocessobc(ng)) THEN
875# ifdef ANA_M2OBC
876 CALL ana_m2obc (ng, tile, itlm)
877# else
878 IF (domain(ng)%SouthWest_Test(tile)) THEN
879 IF (tl_lbc(iwest,isubar,ng)%acquire) THEN
880 CALL set_ngfld (ng, itlm, idu2bc(iwest), jlb, jub, 1, &
881 & 0, mm(ng)+1, 1, &
882 & boundary(ng) % ubarG_west, &
883 & boundary(ng) % ubar_west, &
884 & update)
885 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
886 END IF
887
888 IF (tl_lbc(iwest,isvbar,ng)%acquire) THEN
889 CALL set_ngfld (ng, itlm, idv2bc(iwest), jlb, jub, 1, &
890 & 1, mm(ng)+1, 1, &
891 & boundary(ng) % vbarG_west, &
892 & boundary(ng) % vbar_west, &
893 & update)
894 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
895 END IF
896
897 IF (tl_lbc(ieast,isubar,ng)%acquire) THEN
898 CALL set_ngfld (ng, itlm, idu2bc(ieast), jlb, jub, 1, &
899 & 0, mm(ng)+1, 1, &
900 & boundary(ng) % ubarG_east, &
901 & boundary(ng) % ubar_east, &
902 & update)
903 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
904 END IF
905
906 IF (tl_lbc(ieast,isvbar,ng)%acquire) THEN
907 CALL set_ngfld (ng, itlm, idv2bc(ieast), jlb, jub, 1, &
908 & 1, mm(ng)+1, 1, &
909 & boundary(ng) % vbarG_east, &
910 & boundary(ng) % vbar_east, &
911 & update)
912 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
913 END IF
914
915 IF (tl_lbc(isouth,isubar,ng)%acquire) THEN
916 CALL set_ngfld (ng, itlm, idu2bc(isouth), ilb, iub, 1, &
917 & 1, lm(ng)+1, 1, &
918 & boundary(ng) % ubarG_south, &
919 & boundary(ng) % ubar_south, &
920 & update)
921 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
922 END IF
923
924 IF (tl_lbc(isouth,isvbar,ng)%acquire) THEN
925 CALL set_ngfld (ng, itlm, idv2bc(isouth), ilb, iub, 1, &
926 & 0, lm(ng)+1, 1, &
927 & boundary(ng) % vbarG_south, &
928 & boundary(ng) % vbar_south, &
929 & update)
930 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
931 END IF
932
933 IF (tl_lbc(inorth,isubar,ng)%acquire) THEN
934 CALL set_ngfld (ng, itlm, idu2bc(inorth), ilb, iub, 1, &
935 & 1, lm(ng)+1, 1, &
936 & boundary(ng) % ubarG_north, &
937 & boundary(ng) % ubar_north, &
938 & update)
939 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
940 END IF
941
942 IF (tl_lbc(inorth,isvbar,ng)%acquire) THEN
943 CALL set_ngfld (ng, itlm, idv2bc(inorth), ilb, iub, 1, &
944 & 0, lm(ng)+1, 1, &
945 & boundary(ng) % vbarG_north, &
946 & boundary(ng) % vbar_north, &
947 & update)
948 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
949 END IF
950 END IF
951# endif
952 END IF
953
954# ifdef SOLVE3D
955!
956! 3D momentum.
957!
958 IF (lprocessobc(ng)) THEN
959# ifdef ANA_M3OBC
960 CALL ana_m3obc (ng, tile, itlm)
961# else
962 IF (domain(ng)%SouthWest_Test(tile)) THEN
963 IF (tl_lbc(iwest,isuvel,ng)%acquire) THEN
964 CALL set_ngfld (ng, itlm, idu3bc(iwest), jlb, jub, n(ng), &
965 & 0, mm(ng)+1, n(ng), &
966 & boundary(ng) % uG_west, &
967 & boundary(ng) % u_west, &
968 & update)
969 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
970 END IF
971
972 IF (tl_lbc(iwest,isvvel,ng)%acquire) THEN
973 CALL set_ngfld (ng, itlm, idv3bc(iwest), jlb, jub, n(ng), &
974 & 1, mm(ng)+1, n(ng), &
975 & boundary(ng) % vG_west, &
976 & boundary(ng) % v_west, &
977 & update)
978 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
979 END IF
980
981 IF (tl_lbc(ieast,isuvel,ng)%acquire) THEN
982 CALL set_ngfld (ng, itlm, idu3bc(ieast), jlb, jub, n(ng), &
983 & 0, mm(ng)+1, n(ng), &
984 & boundary(ng) % uG_east, &
985 & boundary(ng) % u_east, &
986 & update)
987 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
988 END IF
989
990 IF (tl_lbc(ieast,isvvel,ng)%acquire) THEN
991 CALL set_ngfld (ng, itlm, idv3bc(ieast), jlb, jub, n(ng), &
992 & 1, mm(ng)+1, n(ng), &
993 & boundary(ng) % vG_east, &
994 & boundary(ng) % v_east, &
995 & update)
996 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
997 END IF
998
999 IF (tl_lbc(isouth,isuvel,ng)%acquire) THEN
1000 CALL set_ngfld (ng, itlm, idu3bc(isouth), ilb, iub, n(ng), &
1001 & 1, lm(ng)+1, n(ng), &
1002 & boundary(ng) % uG_south, &
1003 & boundary(ng) % u_south, &
1004 & update)
1005 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1006 END IF
1007
1008 IF (tl_lbc(isouth,isvvel,ng)%acquire) THEN
1009 CALL set_ngfld (ng, itlm, idv3bc(isouth), ilb, iub, n(ng), &
1010 & 0, lm(ng)+1, n(ng), &
1011 & boundary(ng) % vG_south, &
1012 & boundary(ng) % v_south, &
1013 & update)
1014 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1015 END IF
1016
1017 IF (tl_lbc(inorth,isuvel,ng)%acquire) THEN
1018 CALL set_ngfld (ng, itlm, idu3bc(inorth), ilb, iub, n(ng), &
1019 & 1, lm(ng)+1, n(ng), &
1020 & boundary(ng) % uG_north, &
1021 & boundary(ng) % u_north, &
1022 & update)
1023 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1024 END IF
1025
1026 IF (tl_lbc(inorth,isvvel,ng)%acquire) THEN
1027 CALL set_ngfld (ng, itlm, idv3bc(inorth), ilb, iub, n(ng), &
1028 & 0, lm(ng)+1, n(ng), &
1029 & boundary(ng) % vG_north, &
1030 & boundary(ng) % v_north, &
1031 & update)
1032 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1033 END IF
1034 END IF
1035# endif
1036 END IF
1037!
1038! Tracers.
1039!
1040 IF (lprocessobc(ng)) THEN
1041# ifdef ANA_TOBC
1042 CALL ana_tobc (ng, tile, itlm)
1043# else
1044 IF (domain(ng)%SouthWest_Test(tile)) THEN
1045 DO itrc=1,nt(ng)
1046 IF (tl_lbc(iwest,istvar(itrc),ng)%acquire) THEN
1047 CALL set_ngfld (ng, itlm, idtbry(iwest,itrc), &
1048 & jlb, jub, n(ng), 0, mm(ng)+1, n(ng), &
1049 & boundary(ng) % tG_west(:,:,:,itrc), &
1050 & boundary(ng) % t_west(:,:,itrc), &
1051 & update)
1053 & __line__, myfile)) RETURN
1054 END IF
1055
1056 IF (tl_lbc(ieast,istvar(itrc),ng)%acquire) THEN
1057 CALL set_ngfld (ng, itlm, idtbry(ieast,itrc), &
1058 & jlb, jub, n(ng), 0, mm(ng)+1, n(ng), &
1059 & boundary(ng) % tG_east(:,:,:,itrc), &
1060 & boundary(ng) % t_east(:,:,itrc), &
1061 & update)
1063 & __line__, myfile)) RETURN
1064 END IF
1065
1066 IF (tl_lbc(isouth,istvar(itrc),ng)%acquire) THEN
1067 CALL set_ngfld (ng, itlm, idtbry(isouth,itrc), &
1068 & ilb, iub, n(ng), 0, lm(ng)+1, n(ng), &
1069 & boundary(ng) % tG_south(:,:,:,itrc), &
1070 & boundary(ng) % t_south(:,:,itrc), &
1071 & update)
1073 & __line__, myfile)) RETURN
1074 END IF
1075
1076 IF (tl_lbc(inorth,istvar(itrc),ng)%acquire) THEN
1077 CALL set_ngfld (ng, itlm, idtbry(inorth,itrc), &
1078 & ilb, iub, n(ng), 0, lm(ng)+1, n(ng), &
1079 & boundary(ng) % tG_north(:,:,:,itrc), &
1080 & boundary(ng) % t_north(:,:,itrc), &
1081 & update)
1083 & __line__, myfile)) RETURN
1084 END IF
1085 END DO
1086 END IF
1087# endif
1088 END IF
1089# endif
1090!
1091!-----------------------------------------------------------------------
1092! Set sea surface height climatology (m).
1093!-----------------------------------------------------------------------
1094!
1095 IF (lsshclm(ng)) THEN
1096# ifdef ANA_SSH
1097 CALL ana_ssh (ng, tile, itlm)
1098# else
1099 CALL set_2dfld_tile (ng, tile, itlm, idsshc, &
1100 & lbi, ubi, lbj, ubj, &
1101 & clima(ng)%sshG, &
1102 & clima(ng)%ssh, &
1103 & update)
1104 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1105# endif
1106 END IF
1107!
1108!-----------------------------------------------------------------------
1109! Set 2D momentum climatology (m/s).
1110!-----------------------------------------------------------------------
1111!
1112 IF (lm2clm(ng)) THEN
1113# ifdef ANA_M2CLIMA
1114 CALL ana_m2clima (ng, tile, itlm)
1115# else
1116 CALL set_2dfld_tile (ng, tile, itlm, idubcl, &
1117 & lbi, ubi, lbj, ubj, &
1118 & clima(ng)%ubarclmG, &
1119 & clima(ng)%ubarclm, &
1120 & update)
1121 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1122!
1123 CALL set_2dfld_tile (ng, tile, itlm, idvbcl, &
1124 & lbi, ubi, lbj, ubj, &
1125 & clima(ng)%vbarclmG, &
1126 & clima(ng)%vbarclm, &
1127 & update)
1128 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1129# endif
1130 END IF
1131
1132# ifdef SOLVE3D
1133!
1134!-----------------------------------------------------------------------
1135! Set 3D momentum climatology (m/s).
1136!-----------------------------------------------------------------------
1137!
1138 IF (lm3clm(ng)) THEN
1139# ifdef ANA_M3CLIMA
1140 CALL ana_m3clima (ng, tile, itlm)
1141# else
1142 CALL set_3dfld_tile (ng, tile, itlm, iduclm, &
1143 & lbi, ubi, lbj, ubj, 1, n(ng), &
1144 & clima(ng)%uclmG, &
1145 & clima(ng)%uclm, &
1146 & update)
1147 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1148!
1149 CALL set_3dfld_tile (ng, tile, itlm, idvclm, &
1150 & lbi, ubi, lbj, ubj, 1, n(ng), &
1151 & clima(ng)%vclmG, &
1152 & clima(ng)%vclm, &
1153 & update)
1154 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1155# endif
1156 END IF
1157!
1158!-----------------------------------------------------------------------
1159! Set tracer climatology.
1160!-----------------------------------------------------------------------
1161!
1162# ifdef ANA_TCLIMA
1163 IF (any(ltracerclm(:,ng))) THEN
1164 CALL ana_tclima (ng, tile, itlm)
1165 END IF
1166# else
1167 ic=0
1168 DO itrc=1,nt(ng)
1169 IF (ltracerclm(itrc,ng)) THEN
1170 ic=ic+1
1171 CALL set_3dfld_tile (ng, tile, itlm, idtclm(itrc), &
1172 & lbi, ubi, lbj, ubj, 1, n(ng), &
1173 & clima(ng)%tclmG(:,:,:,:,ic), &
1174 & clima(ng)%tclm (:,:,:,ic), &
1175 & update)
1176 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1177 END IF
1178 END DO
1179# endif
1180# endif
1181
1182# ifdef FORWARD_READ
1183!
1184!-----------------------------------------------------------------------
1185! Set forward solution needed by Tangent/Adjoint models.
1186!-----------------------------------------------------------------------
1187!
1188! Set forward free-surface.
1189!
1190 DO k=1,3
1191 CALL set_2dfld_tile (ng, tile, itlm, idfsur, &
1192 & lbi, ubi, lbj, ubj, &
1193 & ocean(ng)%zetaG, &
1194 & ocean(ng)%zeta(:,:,k), &
1195 & update)
1196 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1197 END DO
1198
1199# ifdef SOLVE3D
1200 DO j=jstrr,jendr
1201 DO i=istrr,iendr
1202 coupling(ng)%Zt_avg1(i,j)=ocean(ng)%zeta(i,j,1)
1203 END DO
1204 END DO
1205# endif
1206!
1207! Set forward 2D momentum.
1208!
1209 DO k=1,3
1210 CALL set_2dfld_tile (ng, tile, itlm, idubar, &
1211 & lbi, ubi, lbj, ubj, &
1212 & ocean(ng)%ubarG, &
1213 & ocean(ng)%ubar(:,:,k), &
1214 & update)
1215 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1216
1217 CALL set_2dfld_tile (ng, tile, itlm, idvbar, &
1218 & lbi, ubi, lbj, ubj, &
1219 & ocean(ng)%vbarG, &
1220 & ocean(ng)%vbar(:,:,k), &
1221 & update)
1222 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1223 END DO
1224
1225# ifdef FORWARD_RHS
1226!
1227! Set forward variables associated with 2D right-hand-side terms.
1228!
1229 DO k=1,2
1230 CALL set_2dfld_tile (ng, tile, itlm, idrzet, &
1231 & lbi, ubi, lbj, ubj, &
1232 & ocean(ng)%rzetaG, &
1233 & ocean(ng)%rzeta(:,:,k), &
1234 & update)
1235 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1236!
1237 CALL set_2dfld_tile (ng, tile, itlm, idru2d, &
1238 & lbi, ubi, lbj, ubj, &
1239 & ocean(ng)%rubarG, &
1240 & ocean(ng)%rubar(:,:,k), &
1241 & update)
1242 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1243!
1244 CALL set_2dfld_tile (ng, tile, itlm, idrv2d, &
1245 & lbi, ubi, lbj, ubj, &
1246 & ocean(ng)%rvbarG, &
1247 & ocean(ng)%rvbar(:,:,k), &
1248 & update)
1249 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1250 END DO
1251# endif
1252
1253# ifdef SOLVE3D
1254!
1255! Set forward time-averaged barotropic fluxes.
1256!
1257 CALL set_2dfld_tile (ng, tile, itlm, idufx1, &
1258 & lbi, ubi, lbj, ubj, &
1259 & coupling(ng)%DU_avg1G, &
1260 & coupling(ng)%DU_avg1, &
1261 & update)
1262 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1263!
1264 CALL set_2dfld_tile (ng, tile, itlm, idufx2, &
1265 & lbi, ubi, lbj, ubj, &
1266 & coupling(ng)%DU_avg2G, &
1267 & coupling(ng)%DU_avg2, &
1268 & update)
1269 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1270!
1271 CALL set_2dfld_tile (ng, tile, itlm, idvfx1, &
1272 & lbi, ubi, lbj, ubj, &
1273 & coupling(ng)%DV_avg1G, &
1274 & coupling(ng)%DV_avg1, &
1275 & update)
1276 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1277!
1278 CALL set_2dfld_tile (ng, tile, itlm, idvfx2, &
1279 & lbi, ubi, lbj, ubj, &
1280 & coupling(ng)%DV_avg2G, &
1281 & coupling(ng)%DV_avg2, &
1282 & update)
1283 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1284!
1285! Set 3D momentum.
1286!
1287 DO k=1,2
1288 CALL set_3dfld_tile (ng, tile, itlm, iduvel, &
1289 & lbi, ubi, lbj, ubj, 1, n(ng), &
1290 & ocean(ng)%uG, &
1291 & ocean(ng)%u(:,:,:,k), &
1292 & update)
1293 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1294!
1295 CALL set_3dfld_tile (ng, tile, itlm, idvvel, &
1296 & lbi, ubi, lbj, ubj, 1, n(ng), &
1297 & ocean(ng)%vG, &
1298 & ocean(ng)%v(:,:,:,k), &
1299 & update)
1300 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1301 END DO
1302
1303# ifdef FORWARD_RHS
1304!
1305! Set variables associated with 3D momentum right-hand-side terms.
1306!
1307 DO k=1,2
1308 CALL set_3dfld_tile (ng, tile, itlm, idru3d, &
1309 & lbi, ubi, lbj, ubj, 1, n(ng), &
1310 & ocean(ng)%ruG, &
1311 & ocean(ng)%ru(:,:,:,k), &
1312 & update)
1313 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1314!
1315 CALL set_3dfld_tile (ng, tile, itlm, idrv3d, &
1316 & lbi, ubi, lbj, ubj, 1, n(ng), &
1317 & ocean(ng)%rvG, &
1318 & ocean(ng)%rv(:,:,:,k), &
1319 & update)
1320 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1321 END DO
1322# endif
1323!
1324! Set 3D tracers.
1325!
1326 DO itrc=1,nt(ng)
1327 DO k=1,3
1328 CALL set_3dfld_tile (ng, tile, itlm, idtvar(itrc), &
1329 & lbi, ubi, lbj, ubj, 1, n(ng), &
1330 & ocean(ng)%tG(:,:,:,:,itrc), &
1331 & ocean(ng)%t(:,:,:,k,itrc), &
1332 & update)
1333 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1334 END DO
1335 END DO
1336
1337# ifdef FORWARD_MIXING
1338!
1339! Set forward vertical mixing variables.
1340!
1341 DO itrc=1,nat
1342 CALL set_3dfld_tile (ng, tile, itlm, iddiff(itrc), &
1343 & lbi, ubi, lbj, ubj, 0, n(ng), &
1344 & mixing(ng)%AktG(:,:,:,:,itrc), &
1345 & mixing(ng)%Akt(:,:,:,itrc), &
1346 & update)
1347 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1348 END DO
1349!
1350 CALL set_3dfld_tile (ng, tile, itlm, idvvis, &
1351 & lbi, ubi, lbj, ubj, 0, n(ng), &
1352 & mixing(ng)%AkvG, &
1353 & mixing(ng)%Akv, &
1354 & update)
1355 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1356# endif
1357
1358# if defined MY25_MIXING_NOT_YET || defined GLS_MIXING_NOT_YET
1359!
1360! Set forward turbulent kinetic energy.
1361!
1362 DO k=1,3
1363 CALL set_3dfld_tile (ng, tile, itlm, idmtke, &
1364 & lbi, ubi, lbj, ubj, 0, n(ng), &
1365 & mixing(ng)%tkeG, &
1366 & mixing(ng)%tke(:,:,:,k), &
1367 & update)
1368 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1369 END DO
1370!
1371! Set forward turbulent kinetic energy times length scale.
1372!
1373 DO k=1,3
1374 CALL set_3dfld_tile (ng, tile, itlm, idmtls, &
1375 & lbi, ubi, lbj, ubj, 0, n(ng), &
1376 & mixing(ng)%glsG, &
1377 & mixing(ng)%gls(:,:,:,k), &
1378 & update)
1379 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1380 END DO
1381!
1382! Set forward vertical mixing length scale.
1383!
1384 CALL set_3dfld_tile (ng, tile, itlm, idvmls, &
1385 & lbi, ubi, lbj, ubj, 0, n(ng), &
1386 & mixing(ng)%LscaleG, &
1387 & mixing(ng)%Lscale, &
1388 & update)
1389 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1390!
1391! Set forward vertical mixing coefficient for turbulent kinetic energy.
1392!
1393 CALL set_3dfld_tile (ng, tile, itlm, idvmkk, &
1394 & lbi, ubi, lbj, ubj, 0, n(ng), &
1395 & mixing(ng)%AkkG, &
1396 & mixing(ng)%Akk, &
1397 & update)
1398 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1399
1400# ifdef GLS_MIXING_NOT_YET
1401!
1402! Set forward vertical mixing coefficient for turbulent length scale.
1403!
1404 CALL set_3dfld_tile (ng, tile, itlm, idvmkp, &
1405 & lbi, ubi, lbj, ubj, 0, n(ng), &
1406 & mixing(ng)%AkpG, &
1407 & mixing(ng)%Akp, &
1408 & update)
1409 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1410# endif
1411# endif
1412
1413# ifdef LMD_MIXING_NOT_YET
1414!
1415! Set forward depth of surface oceanic boundary layer.
1416!
1417 CALL set_2dfld_tile (ng, tile, itlm, idhsbl, &
1418 & lbi, ubi, lbj, ubj, &
1419 & mixing(ng)%hsblG, &
1420 & mixing(ng)%hsbl, &
1421 & update)
1422 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1423# endif
1424
1425# ifdef LMD_BKPP_NOT_YET
1426!
1427! Set forward depth of bottom oceanic boundary layer.
1428!
1429 CALL set_2dfld_tile (ng, tile, itlm, idhbbl, &
1430 & lbi, ubi, lbj, ubj, &
1431 & mixing(ng)%hbblG, &
1432 & mixing(ng)%hbbl, &
1433 & update)
1434 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1435# endif
1436
1437# ifdef LMD_NONLOCAL_NOT_YET
1438!
1439! Set forward boundary layer nonlocal transport.
1440!
1441 DO itrc=1,nat
1442 CALL set_3dfld_tile (ng, tile, itlm, idghat(itrc), &
1443 & lbi, ubi, lbj, ubj, 0, n(ng), &
1444 & mixing(ng)%ghatsG(:,:,:,:,itrc), &
1445 & mixing(ng)%ghat(:,:,:,itrc), &
1446 & update)
1447 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1448 END DO
1449# endif
1450# endif
1451
1452# ifdef WEAK_CONSTRAINT
1453!
1454!-----------------------------------------------------------------------
1455! Set weak contraint frequent impulse forcing.
1456!-----------------------------------------------------------------------
1457!
1458 IF (frequentimpulse(ng)) THEN
1459!
1460! Set free-surface forcing.
1461!
1462 CALL set_2dfld_tile (ng, tile, itlm, idztlf, &
1463 & lbi, ubi, lbj, ubj, &
1464 & ocean(ng)%f_zetaG, &
1465 & ocean(ng)%f_zeta, &
1466 & update)
1467 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1468
1469# ifndef SOLVE3D
1470!
1471! Set 2D momentum forcing.
1472!
1473 CALL set_2dfld_tile (ng, tile, itlm, idubtf, &
1474 & lbi, ubi, lbj, ubj, &
1475 & ocean(ng)%f_ubarG, &
1476 & ocean(ng)%f_ubar, &
1477 & update)
1478 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1479
1480 CALL set_2dfld_tile (ng, tile, itlm, idvbtf, &
1481 & lbi, ubi, lbj, ubj, &
1482 & ocean(ng)%f_vbarG, &
1483 & ocean(ng)%f_vbar, &
1484 & update)
1485 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1486# else
1487!
1488! Set 3D momentum forcing.
1489!
1490 CALL set_3dfld_tile (ng, tile, itlm, idutlf, &
1491 & lbi, ubi, lbj, ubj, 1, n(ng), &
1492 & ocean(ng)%f_uG, &
1493 & ocean(ng)%f_u, &
1494 & update)
1495 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1496!
1497 CALL set_3dfld_tile (ng, tile, itlm, idvtlf, &
1498 & lbi, ubi, lbj, ubj, 1, n(ng), &
1499 & ocean(ng)%f_vG, &
1500 & ocean(ng)%f_v, &
1501 & update)
1502 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1503!
1504! Set 3D tracers forcing.
1505!
1506 DO itrc=1,nt(ng)
1507 CALL set_3dfld_tile (ng, tile, itlm, idttlf(itrc), &
1508 & lbi, ubi, lbj, ubj, 1, n(ng), &
1509 & ocean(ng)%f_tG(:,:,:,:,itrc), &
1510 & ocean(ng)%f_t(:,:,:,itrc), &
1511 & update)
1512 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1513 END DO
1514# endif
1515 END IF
1516# endif
1517# endif
1518!
1519 RETURN
subroutine ana_spinning(ng, tile, model)
Definition ana_spinning.h:3
subroutine ana_fsobc(ng, tile, model)
Definition ana_fsobc.h:3
subroutine ana_wwave(ng, tile, model)
Definition ana_wwave.h:3
subroutine ana_sst(ng, tile, model)
Definition ana_sst.h:3
subroutine ana_m3obc(ng, tile, model)
Definition ana_m3obc.h:3
subroutine ana_dqdsst(ng, tile, model)
Definition ana_dqdsst.h:3
subroutine ana_tclima(ng, tile, model)
Definition ana_tclima.h:3
subroutine ana_btflux(ng, tile, model, itrc)
Definition ana_btflux.h:3
subroutine ana_ssh(ng, tile, model)
Definition ana_ssh.h:3
subroutine ana_m2obc(ng, tile, model)
Definition ana_m2obc.h:3
subroutine ana_tobc(ng, tile, model)
Definition ana_tobc.h:3
subroutine ana_sss(ng, tile, model)
Definition ana_sss.h:3
subroutine ana_psource(ng, tile, model)
Definition ana_psource.h:3
subroutine ana_winds(ng, tile, model)
Definition ana_winds.h:3
subroutine ana_srflux(ng, tile, model)
Definition ana_srflux.h:3
subroutine ana_smflux(ng, tile, model)
Definition ana_smflux.h:3
subroutine ana_pair(ng, tile, model)
Definition ana_pair.h:3
subroutine ana_m2clima(ng, tile, model)
Definition ana_m2clima.h:3
subroutine ana_tair(ng, tile, model)
Definition ana_tair.h:3
subroutine ana_m3clima(ng, tile, model)
Definition ana_m3clima.h:3
subroutine ana_stflux(ng, tile, model, itrc)
Definition ana_stflux.h:3
subroutine ana_rain(ng, tile, model)
Definition ana_rain.h:3
subroutine ana_specir(ng, tile, model)
Definition ana_specir.h:3
subroutine ana_humid(ng, tile, model)
Definition ana_humid.h:3
subroutine ana_cloud(ng, tile, model)
Definition ana_cloud.h:3
subroutine mp_boundary(ng, model, imin, imax, lbi, ubi, lbk, ubk, update, a)
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_boundary), dimension(:), allocatable boundary
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_coupling), dimension(:), allocatable coupling
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer iddqdt
integer idvair
integer idvmls
integer idubtf
integer, dimension(2) iddiff
integer, dimension(4) idu3bc
integer idrtra
integer idwdis
integer idrv3d
integer idubar
integer, dimension(4) idzbry
integer idvvel
integer idhsbl
integer, dimension(4) idu2bc
integer idvsms
logical, dimension(:,:,:), allocatable linfo
integer idcfra
integer idwlen
integer isvvel
integer isvbar
integer idvclm
integer, dimension(:), allocatable idrtrc
integer idpair
integer, dimension(:), allocatable idtbot
integer idvfx2
integer, dimension(:), allocatable idtsur
integer idru2d
integer idvmkp
integer idempf
integer, dimension(:), allocatable idtclm
integer idfsur
integer idwbrk
integer, dimension(:), allocatable idtvar
integer, dimension(:), allocatable istvar
integer idhbbl
integer idvfx1
integer idldwn
integer isuvel
integer idufx2
integer idsssc
integer isfsur
integer iduclm
integer idztlf
integer idsfwf
integer idvbtf
integer, dimension(4) idv3bc
integer iduair
integer, dimension(:), allocatable idttlf
integer idmtke
integer iduvel
integer, dimension(2) idghat
integer idqair
integer idutlf
integer isubar
integer, dimension(:,:,:), allocatable iinfo
integer idubcl
integer, dimension(4) idv2bc
integer idlrad
integer idvtlf
integer idru3d
integer idusms
integer idvmkk
integer idwamp
integer idwdir
integer idvvis
integer idwptp
integer idrzet
integer idufx1
integer idrain
integer idvbcl
integer idsrad
integer idmtls
integer idsshc
integer idwpbt
integer idworb
integer idsstc
integer, dimension(:,:), allocatable idtbry
integer idtair
integer idrv2d
integer idvbar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nat
Definition mod_param.F:499
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, dimension(:), allocatable mm
Definition mod_param.F:456
logical, dimension(:), allocatable luvsrc
logical, dimension(:,:), allocatable ltracersrc
real(r8), dimension(:), allocatable dcrit
logical, dimension(:), allocatable lprocessobc
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable lm3clm
logical, dimension(:), allocatable lsshclm
logical, dimension(:), allocatable frequentimpulse
logical, dimension(:), allocatable lwsrc
integer exit_flag
integer isalt
integer itemp
integer, parameter isouth
logical, dimension(:), allocatable lm2clm
integer, parameter ieast
logical, dimension(:,:), allocatable ltracerclm
integer, parameter inorth
integer noerror
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine set_2dfld_tile(ng, tile, model, ifield, lbi, ubi, lbj, ubj, finp, fout, update, setbc)
Definition set_2dfld.F:25
subroutine set_3dfld_tile(ng, tile, model, ifield, lbi, ubi, lbj, ubj, lbk, ubk, finp, fout, update, setbc)
Definition set_3dfld.F:26
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine set_ngfld(ng, model, ifield, lbi, ubi, ubj, istr, iend, jrec, finp, fout, update)
Definition set_ngfld.F:5

References analytical_mod::ana_btflux(), analytical_mod::ana_cloud(), analytical_mod::ana_dqdsst(), analytical_mod::ana_fsobc(), analytical_mod::ana_humid(), analytical_mod::ana_m2clima(), analytical_mod::ana_m2obc(), analytical_mod::ana_m3clima(), analytical_mod::ana_m3obc(), analytical_mod::ana_pair(), analytical_mod::ana_psource(), analytical_mod::ana_rain(), analytical_mod::ana_smflux(), analytical_mod::ana_specir(), analytical_mod::ana_spinning(), analytical_mod::ana_srflux(), analytical_mod::ana_ssh(), analytical_mod::ana_sss(), analytical_mod::ana_sst(), analytical_mod::ana_stflux(), analytical_mod::ana_tair(), analytical_mod::ana_tclima(), analytical_mod::ana_tobc(), analytical_mod::ana_winds(), analytical_mod::ana_wwave(), mod_boundary::boundary, mod_param::bounds, mod_clima::clima, mod_coupling::coupling, mod_scalars::dcrit, mod_param::domain, mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), exchange_2d_mod::exchange_u2d_tile(), exchange_2d_mod::exchange_v2d_tile(), mod_scalars::exit_flag, mod_forces::forces, strings_mod::founderror(), mod_scalars::frequentimpulse, mod_grid::grid, mod_ncparam::idcfra, mod_ncparam::iddiff, mod_ncparam::iddqdt, mod_ncparam::idempf, mod_ncparam::idfsur, mod_ncparam::idghat, mod_ncparam::idhbbl, mod_ncparam::idhsbl, mod_ncparam::idldwn, mod_ncparam::idlrad, mod_ncparam::idmtke, mod_ncparam::idmtls, mod_ncparam::idpair, mod_ncparam::idqair, mod_ncparam::idrain, mod_ncparam::idrtra, mod_ncparam::idrtrc, mod_ncparam::idru2d, mod_ncparam::idru3d, mod_ncparam::idrv2d, mod_ncparam::idrv3d, mod_ncparam::idrzet, mod_ncparam::idsfwf, mod_ncparam::idsrad, mod_ncparam::idsshc, mod_ncparam::idsssc, mod_ncparam::idsstc, mod_ncparam::idtair, mod_ncparam::idtbot, mod_ncparam::idtbry, mod_ncparam::idtclm, mod_ncparam::idtsur, mod_ncparam::idttlf, mod_ncparam::idtvar, mod_ncparam::idu2bc, mod_ncparam::idu3bc, mod_ncparam::iduair, mod_ncparam::idubar, mod_ncparam::idubcl, mod_ncparam::idubtf, mod_ncparam::iduclm, mod_ncparam::idufx1, mod_ncparam::idufx2, mod_ncparam::idusms, mod_ncparam::idutlf, mod_ncparam::iduvel, mod_ncparam::idv2bc, mod_ncparam::idv3bc, mod_ncparam::idvair, mod_ncparam::idvbar, mod_ncparam::idvbcl, mod_ncparam::idvbtf, mod_ncparam::idvclm, mod_ncparam::idvfx1, mod_ncparam::idvfx2, mod_ncparam::idvmkk, mod_ncparam::idvmkp, mod_ncparam::idvmls, mod_ncparam::idvsms, mod_ncparam::idvtlf, mod_ncparam::idvvel, mod_ncparam::idvvis, mod_ncparam::idwamp, mod_ncparam::idwbrk, mod_ncparam::idwdir, mod_ncparam::idwdis, mod_ncparam::idwlen, mod_ncparam::idworb, mod_ncparam::idwpbt, mod_ncparam::idwptp, mod_ncparam::idzbry, mod_ncparam::idztlf, mod_scalars::ieast, mod_ncparam::iinfo, mod_scalars::inorth, mod_scalars::isalt, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_param::itlm, mod_scalars::iwest, mod_ncparam::linfo, mod_param::lm, mod_scalars::lm2clm, mod_scalars::lm3clm, mod_scalars::lprocessobc, mod_scalars::lsshclm, mod_scalars::ltracerclm, mod_scalars::ltracersrc, mod_scalars::luvsrc, mod_scalars::lwsrc, mod_mixing::mixing, mod_param::mm, distribute_mod::mp_boundary(), mp_exchange_mod::mp_exchange2d(), mp_exchange_mod::mp_exchange3d(), mod_param::n, mod_param::nat, mod_param::nghostpoints, mod_scalars::noerror, mod_scalars::nsperiodic, mod_sources::nsrc, mod_param::nt, mod_ocean::ocean, set_2dfld_mod::set_2dfld_tile(), set_3dfld_mod::set_3dfld_tile(), set_ngfld(), mod_sources::sources, and mod_param::tl_lbc.

Referenced by tl_set_data().

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