1 | #include "cppdefs.h"
|
---|
2 | #ifdef NONLINEAR
|
---|
3 | SUBROUTINE set_data (ng, tile)
|
---|
4 | !
|
---|
5 | !svn $Id: set_data.F 799 2008-10-20 20:38:55Z jcwarner $
|
---|
6 | !================================================== Hernan G. Arango ===
|
---|
7 | ! Copyright (c) 2002-2008 The ROMS/TOMS Group !
|
---|
8 | ! Licensed under a MIT/X style license !
|
---|
9 | ! See License_ROMS.txt !
|
---|
10 | !=======================================================================
|
---|
11 | ! !
|
---|
12 | ! This subroutine processes forcing, boundary, climatology, and !
|
---|
13 | ! assimilation 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 | # include "tile.h"
|
---|
26 | !
|
---|
27 | # ifdef PROFILE
|
---|
28 | CALL wclock_on (ng, iNLM, 4)
|
---|
29 | # endif
|
---|
30 | CALL set_data_tile (ng, tile, &
|
---|
31 | & LBi, UBi, LBj, UBj)
|
---|
32 | # ifdef PROFILE
|
---|
33 | CALL wclock_off (ng, iNLM, 4)
|
---|
34 | # endif
|
---|
35 | RETURN
|
---|
36 | END SUBROUTINE set_data
|
---|
37 | !
|
---|
38 | !***********************************************************************
|
---|
39 | SUBROUTINE set_data_tile (ng, tile, &
|
---|
40 | & LBi, UBi, LBj, UBj)
|
---|
41 | !***********************************************************************
|
---|
42 | !
|
---|
43 | USE mod_param
|
---|
44 | USE mod_boundary
|
---|
45 | # ifdef CLIMATOLOGY
|
---|
46 | USE mod_clima
|
---|
47 | # endif
|
---|
48 | USE mod_forces
|
---|
49 | USE mod_grid
|
---|
50 | USE mod_mixing
|
---|
51 | USE mod_ncparam
|
---|
52 | USE mod_ocean
|
---|
53 | USE mod_stepping
|
---|
54 | # if defined ASSIMILATION || defined NUDGING
|
---|
55 | USE mod_obs
|
---|
56 | # endif
|
---|
57 | USE mod_scalars
|
---|
58 | # if defined TS_PSOURCE || defined UV_PSOURCE || defined Q_PSOURCE
|
---|
59 | USE mod_sources
|
---|
60 | # endif
|
---|
61 | !
|
---|
62 | # ifdef ANALYTICAL
|
---|
63 | USE analytical_mod
|
---|
64 | # endif
|
---|
65 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
66 | USE exchange_2d_mod
|
---|
67 | # endif
|
---|
68 | # ifdef SOLVE3D
|
---|
69 | # if defined NUDGING_UVsur || defined ASSIMILATION_UVsur
|
---|
70 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
71 | USE exchange_3d_mod, ONLY : exchange_r3d_tile
|
---|
72 | # endif
|
---|
73 | # endif
|
---|
74 | # endif
|
---|
75 | # if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
|
---|
76 | USE frc_adjust_mod, ONLY : frc_NLadjust
|
---|
77 | # endif
|
---|
78 | # ifdef DISTRIBUTE
|
---|
79 | USE mp_exchange_mod, ONLY : mp_exchange2d
|
---|
80 | # ifdef SOLVE3D
|
---|
81 | USE mp_exchange_mod, ONLY : mp_exchange3d
|
---|
82 | # endif
|
---|
83 | # endif
|
---|
84 | USE set_2dfld_mod
|
---|
85 | # ifdef SOLVE3D
|
---|
86 | USE set_3dfld_mod
|
---|
87 | # endif
|
---|
88 | !
|
---|
89 | implicit none
|
---|
90 | !
|
---|
91 | ! Imported variable declarations.
|
---|
92 | !
|
---|
93 | integer, intent(in) :: ng, tile
|
---|
94 | integer, intent(in) :: LBi, UBi, LBj, UBj
|
---|
95 | !
|
---|
96 | ! Local variable declarations.
|
---|
97 | !
|
---|
98 | # ifdef DISTRIBUTE
|
---|
99 | # ifdef EW_PERIODIC
|
---|
100 | logical :: EWperiodic=.TRUE.
|
---|
101 | # else
|
---|
102 | logical :: EWperiodic=.FALSE.
|
---|
103 | # endif
|
---|
104 | # ifdef NS_PERIODIC
|
---|
105 | logical :: NSperiodic=.TRUE.
|
---|
106 | # else
|
---|
107 | logical :: NSperiodic=.FALSE.
|
---|
108 | # endif
|
---|
109 | # endif
|
---|
110 | logical :: update = .FALSE.
|
---|
111 |
|
---|
112 | # ifdef OBC
|
---|
113 | integer :: ILB, IUB, JLB, JUB
|
---|
114 | # endif
|
---|
115 | integer :: i, itrc, j, k, order
|
---|
116 |
|
---|
117 | real(r8) :: Zr, cff, cff1, cff2
|
---|
118 |
|
---|
119 | real(r8), dimension(PRIVATE_2D_SCRATCH_ARRAY) :: work1
|
---|
120 | real(r8), dimension(PRIVATE_2D_SCRATCH_ARRAY) :: work2
|
---|
121 |
|
---|
122 | # include "set_bounds.h"
|
---|
123 |
|
---|
124 | # ifdef SOLVE3D
|
---|
125 |
|
---|
126 | # if defined CLOUDS && !defined AIR_OCEAN
|
---|
127 | !
|
---|
128 | !-----------------------------------------------------------------------
|
---|
129 | ! Set cloud fraction (nondimensional). Notice that clouds are
|
---|
130 | ! processed first in case that they are used to adjust shortwave
|
---|
131 | ! radiation.
|
---|
132 | !-----------------------------------------------------------------------
|
---|
133 | !
|
---|
134 | # ifdef ANA_CLOUD
|
---|
135 | CALL ana_cloud (ng, tile, iNLM)
|
---|
136 | # else
|
---|
137 | CALL set_2dfld_tile (ng, tile, iNLM, idCfra, &
|
---|
138 | & LBi, UBi, LBj, UBj, &
|
---|
139 | & FORCES(ng)%cloudG, &
|
---|
140 | & FORCES(ng)%cloud, &
|
---|
141 | & update)
|
---|
142 | # endif
|
---|
143 | # endif
|
---|
144 |
|
---|
145 | # if (defined BULK_FLUXES || defined ECOSIM || \
|
---|
146 | (defined SHORTWAVE && defined ANA_SRFLUX && defined ALBEDO)) && \
|
---|
147 | !defined AIR_OCEAN
|
---|
148 | !
|
---|
149 | !-----------------------------------------------------------------------
|
---|
150 | ! Set surface air temperature (degC).
|
---|
151 | !-----------------------------------------------------------------------
|
---|
152 | !
|
---|
153 | # ifdef ANA_TAIR
|
---|
154 | CALL ana_tair (ng, tile, iNLM)
|
---|
155 | # else
|
---|
156 | CALL set_2dfld_tile (ng, tile, iNLM, idTair, &
|
---|
157 | & LBi, UBi, LBj, UBj, &
|
---|
158 | & FORCES(ng)%TairG, &
|
---|
159 | & FORCES(ng)%Tair, &
|
---|
160 | & update)
|
---|
161 | # endif
|
---|
162 | # endif
|
---|
163 |
|
---|
164 | # if (defined BULK_FLUXES || defined ECOSIM || \
|
---|
165 | (defined SHORTWAVE && defined ANA_SRFLUX && defined ALBEDO)) && \
|
---|
166 | !defined AIR_OCEAN
|
---|
167 | !
|
---|
168 | !-----------------------------------------------------------------------
|
---|
169 | ! Set surface air relative or specific humidity.
|
---|
170 | !-----------------------------------------------------------------------
|
---|
171 | !
|
---|
172 | # ifdef ANA_HUMIDITY
|
---|
173 | CALL ana_humid (ng, tile, iNLM)
|
---|
174 | # else
|
---|
175 | CALL set_2dfld_tile (ng, tile, iNLM, idQair, &
|
---|
176 | & LBi, UBi, LBj, UBj, &
|
---|
177 | & FORCES(ng)%HairG, &
|
---|
178 | & FORCES(ng)%Hair, &
|
---|
179 | & update)
|
---|
180 | # endif
|
---|
181 | # endif
|
---|
182 |
|
---|
183 | # if defined SHORTWAVE && !defined AIR_OCEAN
|
---|
184 | !
|
---|
185 | !-----------------------------------------------------------------------
|
---|
186 | ! Set kinematic surface solar shortwave radiation flux (degC m/s).
|
---|
187 | !-----------------------------------------------------------------------
|
---|
188 | !
|
---|
189 | # ifdef ANA_SRFLUX
|
---|
190 | CALL ana_srflux (ng, tile, iNLM)
|
---|
191 | # else
|
---|
192 | CALL set_2dfld_tile (ng, tile, iNLM, idSrad, &
|
---|
193 | & LBi, UBi, LBj, UBj, &
|
---|
194 | & FORCES(ng)%srflxG, &
|
---|
195 | & FORCES(ng)%srflx, &
|
---|
196 | & update)
|
---|
197 | # endif
|
---|
198 | # ifdef DIURNAL_SRFLUX
|
---|
199 | !
|
---|
200 | ! Modulate the averaged shortwave radiation flux by the local diurnal
|
---|
201 | ! cycle.
|
---|
202 | !
|
---|
203 | CALL ana_srflux (ng, tile, iNLM)
|
---|
204 | # endif
|
---|
205 | # endif
|
---|
206 |
|
---|
207 | # if defined BULK_FLUXES && !defined LONGWAVE && !defined LONGWAVE_OUT \
|
---|
208 | && !defined AIR_OCEAN
|
---|
209 | !
|
---|
210 | !-----------------------------------------------------------------------
|
---|
211 | ! Surface net longwave radiation (degC m/s).
|
---|
212 | !-----------------------------------------------------------------------
|
---|
213 | !
|
---|
214 | CALL set_2dfld_tile (ng, tile, iNLM, idLrad, &
|
---|
215 | & LBi, UBi, LBj, UBj, &
|
---|
216 | & FORCES(ng)%lrflxG, &
|
---|
217 | & FORCES(ng)%lrflx, &
|
---|
218 | & update)
|
---|
219 | # endif
|
---|
220 |
|
---|
221 | # if defined LONGWAVE_OUT && defined BULK_FLUXES && !defined AIR_OCEAN
|
---|
222 | !
|
---|
223 | !-----------------------------------------------------------------------
|
---|
224 | ! Surface downwelling longwave radiation (degC m/s).
|
---|
225 | !-----------------------------------------------------------------------
|
---|
226 | !
|
---|
227 | CALL set_2dfld_tile (ng, tile, iNLM, idLdwn, &
|
---|
228 | & LBi, UBi, LBj, UBj, &
|
---|
229 | & FORCES(ng)%lrflxG, &
|
---|
230 | & FORCES(ng)%lrflx, &
|
---|
231 | & update)
|
---|
232 | # endif
|
---|
233 |
|
---|
234 | # if (defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS) && \
|
---|
235 | !defined AIR_OCEAN
|
---|
236 | !
|
---|
237 | !-----------------------------------------------------------------------
|
---|
238 | ! Set surface air pressure (mb).
|
---|
239 | !-----------------------------------------------------------------------
|
---|
240 | !
|
---|
241 | # ifdef ANA_PAIR
|
---|
242 | CALL ana_pair (ng, tile, iNLM)
|
---|
243 | # else
|
---|
244 | CALL set_2dfld_tile (ng, tile, iNLM, idPair, &
|
---|
245 | & LBi, UBi, LBj, UBj, &
|
---|
246 | & FORCES(ng)%PairG, &
|
---|
247 | & FORCES(ng)%Pair, &
|
---|
248 | & update)
|
---|
249 | # endif
|
---|
250 | # endif
|
---|
251 | # if (defined BULK_FLUXES || defined ECOSIM) && !defined AIR_OCEAN
|
---|
252 | !
|
---|
253 | !-----------------------------------------------------------------------
|
---|
254 | ! Set surface winds (m/s).
|
---|
255 | !-----------------------------------------------------------------------
|
---|
256 | !
|
---|
257 | # ifdef ANA_WINDS
|
---|
258 | CALL ana_winds (ng, tile, iNLM)
|
---|
259 | # else
|
---|
260 | CALL set_2dfld_tile (ng, tile, iNLM, idUair, &
|
---|
261 | & LBi, UBi, LBj, UBj, &
|
---|
262 | & FORCES(ng)%UwindG, &
|
---|
263 | & FORCES(ng)%Uwind, &
|
---|
264 | & update)
|
---|
265 | CALL set_2dfld_tile (ng, tile, iNLM, idVair, &
|
---|
266 | & LBi, UBi, LBj, UBj, &
|
---|
267 | & FORCES(ng)%VwindG, &
|
---|
268 | & FORCES(ng)%Vwind, &
|
---|
269 | & update)
|
---|
270 | # ifdef CURVGRID
|
---|
271 | !
|
---|
272 | ! If input point surface winds or interpolated from coarse data, rotate
|
---|
273 | ! to curvilinear grid.
|
---|
274 | !
|
---|
275 | IF (.not.Linfo(1,idUair,ng).or. &
|
---|
276 | & (Iinfo(5,idUair,ng).ne.Lm(ng)+2).or. &
|
---|
277 | & (Iinfo(6,idUair,ng).ne.Mm(ng)+2)) THEN
|
---|
278 | DO j=JstrR,JendR
|
---|
279 | DO i=IstrR,IendR
|
---|
280 | cff1=FORCES(ng)%Uwind(i,j)*GRID(ng)%CosAngler(i,j)+ &
|
---|
281 | & FORCES(ng)%Vwind(i,j)*GRID(ng)%SinAngler(i,j)
|
---|
282 | cff2=FORCES(ng)%Vwind(i,j)*GRID(ng)%CosAngler(i,j)- &
|
---|
283 | & FORCES(ng)%Uwind(i,j)*GRID(ng)%SinAngler(i,j)
|
---|
284 | FORCES(ng)%Uwind(i,j)=cff1
|
---|
285 | FORCES(ng)%Vwind(i,j)=cff2
|
---|
286 | END DO
|
---|
287 | END DO
|
---|
288 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
289 | CALL exchange_r2d_tile (ng, tile, &
|
---|
290 | & LBi, UBi, LBj, UBj, &
|
---|
291 | & FORCES(ng)%UWind)
|
---|
292 | CALL exchange_r2d_tile (ng, tile, &
|
---|
293 | & LBi, UBi, LBj, UBj, &
|
---|
294 | & FORCES(ng)%VWind)
|
---|
295 | # endif
|
---|
296 | # ifdef DISTRIBUTE
|
---|
297 | CALL mp_exchange2d (ng, tile, iNLM, 2, &
|
---|
298 | & LBi, UBi, LBj, UBj, &
|
---|
299 | & NghostPoints, EWperiodic, NSperiodic, &
|
---|
300 | & FORCES(ng)%UWind, &
|
---|
301 | & FORCES(ng)%VWind)
|
---|
302 | # endif
|
---|
303 | END IF
|
---|
304 | # endif
|
---|
305 | # endif
|
---|
306 | # endif
|
---|
307 |
|
---|
308 | # if defined BULK_FLUXES && !defined AIR_OCEAN
|
---|
309 | !
|
---|
310 | !-----------------------------------------------------------------------
|
---|
311 | ! Set rain fall rate (kg/m2/s).
|
---|
312 | !-----------------------------------------------------------------------
|
---|
313 | !
|
---|
314 | # ifdef ANA_RAIN
|
---|
315 | CALL ana_rain (ng, tile, iNLM)
|
---|
316 | # else
|
---|
317 | CALL set_2dfld_tile (ng, tile, iNLM, idrain, &
|
---|
318 | & LBi, UBi, LBj, UBj, &
|
---|
319 | & FORCES(ng)%rainG, &
|
---|
320 | & FORCES(ng)%rain, &
|
---|
321 | & update)
|
---|
322 | # endif
|
---|
323 | # endif
|
---|
324 |
|
---|
325 | # ifndef BULK_FLUXES
|
---|
326 | !
|
---|
327 | !-----------------------------------------------------------------------
|
---|
328 | ! Set kinematic surface net heat flux (degC m/s).
|
---|
329 | !-----------------------------------------------------------------------
|
---|
330 | !
|
---|
331 | # ifdef ANA_STFLUX
|
---|
332 | CALL ana_stflux (ng, tile, iNLM, itemp)
|
---|
333 | # else
|
---|
334 | CALL set_2dfld_tile (ng, tile, iNLM, idTsur(itemp), &
|
---|
335 | & LBi, UBi, LBj, UBj, &
|
---|
336 | & FORCES(ng)%stflxG(:,:,:,itemp), &
|
---|
337 | & FORCES(ng)%stflx (:,:,itemp), &
|
---|
338 | & update)
|
---|
339 | # endif
|
---|
340 | # endif
|
---|
341 |
|
---|
342 | # ifdef QCORRECTION
|
---|
343 | !
|
---|
344 | !-----------------------------------------------------------------------
|
---|
345 | ! Set sea surface temperature (SST) and heat flux sensitivity to
|
---|
346 | ! SST (dQdSST) which are used for surface heat flux correction.
|
---|
347 | !-----------------------------------------------------------------------
|
---|
348 | !
|
---|
349 | # ifdef ANA_SST
|
---|
350 | CALL ana_sst (ng, tile, iNLM)
|
---|
351 | # else
|
---|
352 | CALL set_2dfld_tile (ng, tile, iNLM, idSSTc, &
|
---|
353 | & LBi, UBi, LBj, UBj, &
|
---|
354 | & FORCES(ng)%sstG, &
|
---|
355 | & FORCES(ng)%sst, &
|
---|
356 | & update)
|
---|
357 | CALL set_2dfld_tile (ng, tile, iNLM, iddQdT, &
|
---|
358 | & LBi, UBi, LBj, UBj, &
|
---|
359 | & FORCES(ng)%dqdtG, &
|
---|
360 | & FORCES(ng)%dqdt, &
|
---|
361 | & update)
|
---|
362 | # endif
|
---|
363 | # endif
|
---|
364 | !
|
---|
365 | !-----------------------------------------------------------------------
|
---|
366 | ! Set kinematic bottom net heat flux (degC m/s).
|
---|
367 | !-----------------------------------------------------------------------
|
---|
368 | !
|
---|
369 | # ifdef ANA_BTFLUX
|
---|
370 | CALL ana_btflux (ng, tile, iNLM, itemp)
|
---|
371 | # else
|
---|
372 | CALL set_2dfld_tile (ng, tile, iNLM, idTbot(itemp), &
|
---|
373 | & LBi, UBi, LBj, UBj, &
|
---|
374 | & FORCES(ng)%btflxG(:,:,:,itemp), &
|
---|
375 | & FORCES(ng)%btflx (:,:,itemp), &
|
---|
376 | & update)
|
---|
377 | # endif
|
---|
378 |
|
---|
379 | # ifdef SALINITY
|
---|
380 | !
|
---|
381 | !-----------------------------------------------------------------------
|
---|
382 | ! Set kinematic surface freshwater (E-P) flux (m/s).
|
---|
383 | !-----------------------------------------------------------------------
|
---|
384 | !
|
---|
385 | # ifdef ANA_SSFLUX
|
---|
386 | CALL ana_stflux (ng, tile, iNLM, isalt)
|
---|
387 | # else
|
---|
388 | # if !(defined EMINUSP || defined SRELAXATION)
|
---|
389 | CALL set_2dfld_tile (ng, tile, iNLM, idsfwf, &
|
---|
390 | & LBi, UBi, LBj, UBj, &
|
---|
391 | & FORCES(ng)%stflxG(:,:,:,isalt), &
|
---|
392 | & FORCES(ng)%stflx (:,:,isalt), &
|
---|
393 | & update)
|
---|
394 | # endif
|
---|
395 | # endif
|
---|
396 |
|
---|
397 | # if defined SCORRECTION || defined SRELAXATION
|
---|
398 | !
|
---|
399 | !-----------------------------------------------------------------------
|
---|
400 | ! Set surface salinity for freshwater flux correction.
|
---|
401 | !-----------------------------------------------------------------------
|
---|
402 | !
|
---|
403 | # ifdef ANA_SSS
|
---|
404 | CALL ana_sss (ng, tile, iNLM)
|
---|
405 | # else
|
---|
406 | CALL set_2dfld_tile (ng, tile, iNLM, idSSSc, &
|
---|
407 | & LBi, UBi, LBj, UBj, &
|
---|
408 | & FORCES(ng)%sssG, &
|
---|
409 | & FORCES(ng)%sss, &
|
---|
410 | & update)
|
---|
411 | # endif
|
---|
412 | # endif
|
---|
413 | !
|
---|
414 | !-----------------------------------------------------------------------
|
---|
415 | ! Set kinematic bottom salt flux (m/s).
|
---|
416 | !-----------------------------------------------------------------------
|
---|
417 | !
|
---|
418 | # ifdef ANA_BSFLUX
|
---|
419 | CALL ana_btflux (ng, tile, iNLM, isalt)
|
---|
420 | # else
|
---|
421 | CALL set_2dfld_tile (ng, tile, iNLM, idTbot(isalt), &
|
---|
422 | & LBi, UBi, LBj, UBj, &
|
---|
423 | & FORCES(ng)%btflxG(:,:,:,isalt), &
|
---|
424 | & FORCES(ng)%btflx (:,:,isalt), &
|
---|
425 | & update)
|
---|
426 | # endif
|
---|
427 | # endif
|
---|
428 |
|
---|
429 | # if defined SEDIMENT || defined BIOLOGY
|
---|
430 | !
|
---|
431 | !-----------------------------------------------------------------------
|
---|
432 | ! Set kinematic surface and bottom pasive tracer fluxes (T m/s).
|
---|
433 | !-----------------------------------------------------------------------
|
---|
434 | !
|
---|
435 | DO itrc=NAT+1,NT(ng)
|
---|
436 | # ifdef ANA_SPFLUX
|
---|
437 | CALL ana_stflux (ng, tile, iNLM, itrc)
|
---|
438 | # else
|
---|
439 | CALL set_2dfld_tile (ng, tile, iNLM, idTsur(itrc), &
|
---|
440 | & LBi, UBi, LBj, UBj, &
|
---|
441 | & FORCES(ng)%stflxG(:,:,:,itrc), &
|
---|
442 | & FORCES(ng)%stflx (:,:,itrc), &
|
---|
443 | & update)
|
---|
444 | # endif
|
---|
445 | # ifdef ANA_BPFLUX
|
---|
446 | CALL ana_btflux (ng, tile, iNLM, itrc)
|
---|
447 | # else
|
---|
448 | CALL set_2dfld_tile (ng, tile, iNLM, idTbot(itrc), &
|
---|
449 | & LBi, UBi, LBj, UBj, &
|
---|
450 | & FORCES(ng)%btflxG(:,:,:,itrc), &
|
---|
451 | & FORCES(ng)%btflx (:,:,itrc), &
|
---|
452 | & update)
|
---|
453 | # endif
|
---|
454 | END DO
|
---|
455 | # endif
|
---|
456 | # endif
|
---|
457 |
|
---|
458 | # ifndef AIR_OCEAN
|
---|
459 | # ifndef BULK_FLUXES
|
---|
460 | !
|
---|
461 | !-----------------------------------------------------------------------
|
---|
462 | ! Set kinematic surface momentum flux (m2/s2).
|
---|
463 | !-----------------------------------------------------------------------
|
---|
464 | !
|
---|
465 | # ifdef ANA_SMFLUX
|
---|
466 | CALL ana_smflux (ng, tile, iNLM)
|
---|
467 | # else
|
---|
468 | CALL set_2dfld_tile (ng, tile, iNLM, idUsms, &
|
---|
469 | & LBi, UBi, LBj, UBj, &
|
---|
470 | & FORCES(ng)%sustrG, &
|
---|
471 | & FORCES(ng)%sustr, &
|
---|
472 | & update)
|
---|
473 | CALL set_2dfld_tile (ng, tile, iNLM, idVsms, &
|
---|
474 | & LBi, UBi, LBj, UBj, &
|
---|
475 | & FORCES(ng)%svstrG, &
|
---|
476 | & FORCES(ng)%svstr, &
|
---|
477 | & update)
|
---|
478 | # ifdef CURVGRID
|
---|
479 | !
|
---|
480 | ! If input point wind stress, rotate to curvilinear grid. Notice
|
---|
481 | ! that rotation is done at RHO-points. It does not matter.
|
---|
482 | !
|
---|
483 | IF (.not.Linfo(1,idUsms,ng).or. &
|
---|
484 | & (Iinfo(5,idUsms,ng).ne.Lm(ng)+1).or. &
|
---|
485 | & (Iinfo(6,idUsms,ng).ne.Mm(ng)+2)) THEN
|
---|
486 | DO j=JstrR,JendR
|
---|
487 | DO i=IstrR,IendR
|
---|
488 | cff1=FORCES(ng)%sustr(i,j)*GRID(ng)%CosAngler(i,j)+ &
|
---|
489 | & FORCES(ng)%svstr(i,j)*GRID(ng)%SinAngler(i,j)
|
---|
490 | cff2=FORCES(ng)%svstr(i,j)*GRID(ng)%CosAngler(i,j)- &
|
---|
491 | & FORCES(ng)%sustr(i,j)*GRID(ng)%SinAngler(i,j)
|
---|
492 | FORCES(ng)%sustr(i,j)=cff1
|
---|
493 | FORCES(ng)%svstr(i,j)=cff2
|
---|
494 | END DO
|
---|
495 | END DO
|
---|
496 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
497 | CALL exchange_u2d_tile (ng, tile, &
|
---|
498 | & LBi, UBi, LBj, UBj, &
|
---|
499 | & FORCES(ng)%sustr)
|
---|
500 | CALL exchange_v2d_tile (ng, tile, &
|
---|
501 | & LBi, UBi, LBj, UBj, &
|
---|
502 | & FORCES(ng)%svstr)
|
---|
503 | # endif
|
---|
504 | # ifdef DISTRIBUTE
|
---|
505 | CALL mp_exchange2d (ng, tile, iNLM, 2, &
|
---|
506 | & LBi, UBi, LBj, UBj, &
|
---|
507 | & NghostPoints, EWperiodic, NSperiodic, &
|
---|
508 | & FORCES(ng)%sustr, &
|
---|
509 | & FORCES(ng)%svstr)
|
---|
510 | # endif
|
---|
511 | END IF
|
---|
512 | # endif
|
---|
513 | # endif
|
---|
514 | # endif
|
---|
515 | # endif
|
---|
516 |
|
---|
517 | # ifdef WAVE_DATA
|
---|
518 | !
|
---|
519 | !-----------------------------------------------------------------------
|
---|
520 | ! Set surface wind-induced wave amplitude, direction and period.
|
---|
521 | !-----------------------------------------------------------------------
|
---|
522 | !
|
---|
523 | # ifdef ANA_WWAVE
|
---|
524 | CALL ana_wwave (ng, tile, iNLM)
|
---|
525 | # else
|
---|
526 | # ifdef WAVES_DIR
|
---|
527 | CALL set_2dfld_tile (ng, tile, iNLM, idWdir, &
|
---|
528 | & LBi, UBi, LBj, UBj, &
|
---|
529 | & FORCES(ng)%DwaveG, &
|
---|
530 | & FORCES(ng)%Dwave, &
|
---|
531 | & update)
|
---|
532 | # ifdef CURVGRID
|
---|
533 | !
|
---|
534 | ! If input point-data, rotate direction to curvilinear coordinates.
|
---|
535 | !
|
---|
536 | IF (.not.Linfo(1,idWdir,ng).or. &
|
---|
537 | & (Iinfo(5,idWdir,ng).ne.Lm(ng)+2).or. &
|
---|
538 | & (Iinfo(6,idWdir,ng).ne.Mm(ng)+2)) THEN
|
---|
539 | DO j=JstrR,JendR
|
---|
540 | DO i=IstrR,IendR
|
---|
541 | FORCES(ng)%Dwave(i,j)=FORCES(ng)%Dwave(i,j)- &
|
---|
542 | & GRID(ng)%angler(i,j)
|
---|
543 | END DO
|
---|
544 | END DO
|
---|
545 | END IF
|
---|
546 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
547 | CALL exchange_r2d_tile (ng, tile, &
|
---|
548 | & LBi, UBi, LBj, UBj, &
|
---|
549 | & FORCES(ng)%Dwave)
|
---|
550 | # endif
|
---|
551 | # ifdef DISTRIBUTE
|
---|
552 | CALL mp_exchange2d (ng, tile, iNLM, 1, &
|
---|
553 | & LBi, UBi, LBj, UBj, &
|
---|
554 | & NghostPoints, EWperiodic, NSperiodic, &
|
---|
555 | & FORCES(ng)%Dwave)
|
---|
556 | # endif
|
---|
557 | # endif
|
---|
558 | # endif
|
---|
559 | # ifdef WAVES_HEIGHT
|
---|
560 | CALL set_2dfld_tile (ng, tile, iNLM, idWamp, &
|
---|
561 | & LBi, UBi, LBj, UBj, &
|
---|
562 | & FORCES(ng)%HwaveG, &
|
---|
563 | & FORCES(ng)%Hwave, &
|
---|
564 | & update)
|
---|
565 | # endif
|
---|
566 | # ifdef WAVES_LENGTH
|
---|
567 | CALL set_2dfld_tile (ng, tile, iNLM, idWlen, &
|
---|
568 | & LBi, UBi, LBj, UBj, &
|
---|
569 | & FORCES(ng)%LwaveG, &
|
---|
570 | & FORCES(ng)%Lwave, &
|
---|
571 | & update)
|
---|
572 | # endif
|
---|
573 | # ifdef WAVES_TOP_PERIOD
|
---|
574 | CALL set_2dfld_tile (ng, tile, iNLM, idWptp, &
|
---|
575 | & LBi, UBi, LBj, UBj, &
|
---|
576 | & FORCES(ng)%Pwave_topG, &
|
---|
577 | & FORCES(ng)%Pwave_top, &
|
---|
578 | & update)
|
---|
579 | # endif
|
---|
580 | # ifdef WAVES_BOT_PERIOD
|
---|
581 | CALL set_2dfld_tile (ng, tile, iNLM, idWpbt, &
|
---|
582 | & LBi, UBi, LBj, UBj, &
|
---|
583 | & FORCES(ng)%Pwave_botG, &
|
---|
584 | & FORCES(ng)%Pwave_bot, &
|
---|
585 | & update)
|
---|
586 | # endif
|
---|
587 | # if defined WAVES_UB
|
---|
588 | CALL set_2dfld_tile (ng, tile, iNLM, idWorb, &
|
---|
589 | & LBi, UBi, LBj, UBj, &
|
---|
590 | & FORCES(ng)%Ub_swanG, &
|
---|
591 | & FORCES(ng)%Ub_swan, &
|
---|
592 | & update)
|
---|
593 | # endif
|
---|
594 | # if defined TKE_WAVEDISS
|
---|
595 | CALL set_2dfld_tile (ng, tile, iNLM, idWdis, &
|
---|
596 | & LBi, UBi, LBj, UBj, &
|
---|
597 | & FORCES(ng)%Wave_dissipG, &
|
---|
598 | & FORCES(ng)%Wave_dissip, &
|
---|
599 | & update)
|
---|
600 | # endif
|
---|
601 | # if defined SVENDSEN_ROLLER
|
---|
602 | CALL set_2dfld_tile (ng, tile, iNLM, idWbrk, &
|
---|
603 | & LBi, UBi, LBj, UBj, &
|
---|
604 | & FORCES(ng)%Wave_breakG, &
|
---|
605 | & FORCES(ng)%Wave_break, &
|
---|
606 | & update)
|
---|
607 | # endif
|
---|
608 | # endif
|
---|
609 | # endif
|
---|
610 |
|
---|
611 | # if defined ECOSIM && defined SOLVE3D
|
---|
612 | !
|
---|
613 | !-----------------------------------------------------------------------
|
---|
614 | ! Compute spectral irradiance and cosine of average zenith angle of
|
---|
615 | ! downwelling spectral photons.
|
---|
616 | !-----------------------------------------------------------------------
|
---|
617 | !
|
---|
618 | CALL ana_specir (ng, tile, iNLM)
|
---|
619 | # endif
|
---|
620 |
|
---|
621 | # ifdef ANA_SPINNING
|
---|
622 | !
|
---|
623 | !-----------------------------------------------------------------------
|
---|
624 | ! Set time-varying rotation force (centripetal accelerations) for
|
---|
625 | ! polar coordinate grids.
|
---|
626 | !-----------------------------------------------------------------------
|
---|
627 | !
|
---|
628 | CALL ana_spinning (ng, tile, iNLM)
|
---|
629 | # endif
|
---|
630 |
|
---|
631 | # if defined UV_PSOURCE || defined TS_PSOURCE || defined Q_PSOURCE
|
---|
632 | !
|
---|
633 | !-----------------------------------------------------------------------
|
---|
634 | ! Set point Sources/Sinks (river runoff).
|
---|
635 | !-----------------------------------------------------------------------
|
---|
636 | !
|
---|
637 | IF (SOUTH_WEST_TEST) THEN
|
---|
638 | # ifdef ANA_PSOURCE
|
---|
639 | CALL ana_psource (ng, tile, iNLM)
|
---|
640 | # else
|
---|
641 | # if defined UV_PSOURCE || defined Q_PSOURCE
|
---|
642 | CALL set_ngfld (ng, iNLM, idRtra, 1, Nsrc(ng), 1, &
|
---|
643 | & 1, Nsrc(ng), 1, &
|
---|
644 | & SOURCES(ng) % QbarG, &
|
---|
645 | & SOURCES(ng) % Qbar, &
|
---|
646 | & update)
|
---|
647 | # ifdef SOLVE3D
|
---|
648 | DO k=1,N(ng)
|
---|
649 | DO i=1,Nsrc(ng)
|
---|
650 | SOURCES(ng)%Qsrc(i,k)=SOURCES(ng)%Qbar(i)* &
|
---|
651 | & SOURCES(ng)%Qshape(i,k)
|
---|
652 | END DO
|
---|
653 | END DO
|
---|
654 | # endif
|
---|
655 | # endif
|
---|
656 | # if defined TS_PSOURCE && defined SOLVE3D
|
---|
657 | DO itrc=1,NT(ng)
|
---|
658 | IF (SOURCES(ng)%Ltracer(itrc)) THEN
|
---|
659 | CALL set_ngfld (ng, iNLM, idRtrc(itrc), 1, Nsrc(ng), N(ng), &
|
---|
660 | & 1, Nsrc(ng), N(ng), &
|
---|
661 | & SOURCES(ng) % TsrcG(1,1,1,itrc), &
|
---|
662 | & SOURCES(ng) % Tsrc(1,1,itrc), &
|
---|
663 | & update)
|
---|
664 | IF (.not.update) THEN
|
---|
665 | DO i=1,Nsrc(ng)
|
---|
666 | SOURCES(ng) % Lsrc(i,itrc)=.FALSE.
|
---|
667 | END DO
|
---|
668 | exit_flag=0
|
---|
669 | END IF
|
---|
670 | END IF
|
---|
671 | END DO
|
---|
672 | # endif
|
---|
673 | # endif
|
---|
674 | END IF
|
---|
675 | # endif
|
---|
676 |
|
---|
677 | # ifdef OBC
|
---|
678 | # ifdef REFINED_GRID
|
---|
679 | IF (ng.eq.1) THEN
|
---|
680 | # endif
|
---|
681 | !
|
---|
682 | !-----------------------------------------------------------------------
|
---|
683 | ! Set open boundary conditions fields.
|
---|
684 | !-----------------------------------------------------------------------
|
---|
685 | !
|
---|
686 | ! Lower and upper bounds for nontiled arrays.
|
---|
687 | !
|
---|
688 | ILB=LOWER_BOUND_I
|
---|
689 | IUB=UPPER_BOUND_I
|
---|
690 | JLB=LOWER_BOUND_J
|
---|
691 | JUB=UPPER_BOUND_J
|
---|
692 |
|
---|
693 | # ifdef ANA_FSOBC
|
---|
694 | CALL ana_fsobc (ng, tile, iNLM)
|
---|
695 | # else
|
---|
696 | IF (SOUTH_WEST_TEST) THEN
|
---|
697 | # ifdef WEST_FSOBC
|
---|
698 | CALL set_ngfld (ng, iNLM, idZbry(iwest), JLB, JUB, 1, &
|
---|
699 | & 0, Mm(ng)+1, 1, &
|
---|
700 | & BOUNDARY(ng) % zetaG_west(JLB,1), &
|
---|
701 | & BOUNDARY(ng) % zeta_west(JLB), &
|
---|
702 | & update)
|
---|
703 | # endif
|
---|
704 | # ifdef EAST_FSOBC
|
---|
705 | CALL set_ngfld (ng, iNLM, idZbry(ieast), JLB, JUB, 1, &
|
---|
706 | & 0, Mm(ng)+1, 1, &
|
---|
707 | & BOUNDARY(ng) % zetaG_east(JLB,1), &
|
---|
708 | & BOUNDARY(ng) % zeta_east(JLB), &
|
---|
709 | & update)
|
---|
710 | # endif
|
---|
711 | # ifdef SOUTH_FSOBC
|
---|
712 | CALL set_ngfld (ng, iNLM, idZbry(isouth), ILB, IUB, 1, &
|
---|
713 | & 0, Lm(ng)+1 ,1, &
|
---|
714 | & BOUNDARY(ng) % zetaG_south(ILB,1), &
|
---|
715 | & BOUNDARY(ng) % zeta_south(ILB), &
|
---|
716 | & update)
|
---|
717 | # endif
|
---|
718 | # ifdef NORTH_FSOBC
|
---|
719 | CALL set_ngfld (ng, iNLM, idZbry(inorth), ILB, IUB, 1, &
|
---|
720 | & 0, Lm(ng)+1, 1, &
|
---|
721 | & BOUNDARY(ng) % zetaG_north(ILB,1), &
|
---|
722 | & BOUNDARY(ng) % zeta_north(ILB), &
|
---|
723 | & update)
|
---|
724 | # endif
|
---|
725 | END IF
|
---|
726 | # endif
|
---|
727 | # if defined WET_DRY
|
---|
728 | !
|
---|
729 | !-----------------------------------------------------------------------
|
---|
730 | ! Ensure that water level on boundary cells is above bed elevation.
|
---|
731 | !-----------------------------------------------------------------------
|
---|
732 | IF (SOUTH_WEST_TEST) THEN
|
---|
733 | cff=Dcrit(ng)
|
---|
734 | # ifdef WEST_FSOBC
|
---|
735 | DO j=JLB,JUB
|
---|
736 | IF (BOUNDARY(ng)%zeta_west(j).le. &
|
---|
737 | & (Dcrit(ng)-GRID(ng)%h(0,j))) THEN
|
---|
738 | BOUNDARY(ng)%zeta_west(j)=cff-GRID(ng)%h(0,j)
|
---|
739 | END IF
|
---|
740 | END DO
|
---|
741 | # endif
|
---|
742 | # ifdef EAST_FSOBC
|
---|
743 | DO j=JLB,JUB
|
---|
744 | IF (BOUNDARY(ng)%zeta_east(j).le. &
|
---|
745 | & (Dcrit(ng)-GRID(ng)%h(Lm(ng)+1))) THEN
|
---|
746 | BOUNDARY(ng)%zeta_est(j)=cff-GRID(ng)%h(Lm(ng)+1,j)
|
---|
747 | END IF
|
---|
748 | END DO
|
---|
749 | # endif
|
---|
750 | # ifdef SOUTH_FSOBC
|
---|
751 | DO i=ILB,IUB
|
---|
752 | IF (BOUNDARY(ng)%zeta_south(i).le. &
|
---|
753 | & (Dcrit(ng)-GRID(ng)%h(i,0))) THEN
|
---|
754 | BOUNDARY(ng)%zeta_south(i)=cff-GRID(ng)%h(i,0)
|
---|
755 | END IF
|
---|
756 | END DO
|
---|
757 | # endif
|
---|
758 | # ifdef NORTH_FSOBC
|
---|
759 | DO i=ILB,IUB
|
---|
760 | IF (BOUNDARY(ng)%zeta_north(i).le. &
|
---|
761 | & (Dcrit(ng)-GRID(ng)%h(i,Mm(ng)+1))) THEN
|
---|
762 | BOUNDARY(ng)%zeta_north(i)=cff-GRID(ng)%h(i,Mm(ng)+1)
|
---|
763 | END IF
|
---|
764 | END DO
|
---|
765 | # endif
|
---|
766 | END IF
|
---|
767 | # endif
|
---|
768 | # ifdef ANA_M2OBC
|
---|
769 | CALL ana_m2obc (ng, tile, iNLM)
|
---|
770 | # else
|
---|
771 | IF (SOUTH_WEST_TEST) THEN
|
---|
772 | # ifdef WEST_M2OBC
|
---|
773 | CALL set_ngfld (ng, iNLM, idU2bc(iwest), JLB, JUB, 1, &
|
---|
774 | & 0, Mm(ng)+1, 1, &
|
---|
775 | & BOUNDARY(ng) % ubarG_west(JLB,1), &
|
---|
776 | & BOUNDARY(ng) % ubar_west(JLB), &
|
---|
777 | & update)
|
---|
778 | CALL set_ngfld (ng, iNLM, idV2bc(iwest), JLB, JUB, 1, &
|
---|
779 | & 1, Mm(ng)+1, 1, &
|
---|
780 | & BOUNDARY(ng) % vbarG_west(JLB,1), &
|
---|
781 | & BOUNDARY(ng) % vbar_west(JLB), &
|
---|
782 | & update)
|
---|
783 | # endif
|
---|
784 | # ifdef EAST_M2OBC
|
---|
785 | CALL set_ngfld (ng, iNLM, idU2bc(ieast), JLB, JUB, 1, &
|
---|
786 | & 0, Mm(ng)+1, 1, &
|
---|
787 | & BOUNDARY(ng) % ubarG_east(JLB,1), &
|
---|
788 | & BOUNDARY(ng) % ubar_east(JLB), &
|
---|
789 | & update)
|
---|
790 | CALL set_ngfld (ng, iNLM, idV2bc(ieast), JLB, JUB, 1, &
|
---|
791 | & 1, Mm(ng)+1, 1, &
|
---|
792 | & BOUNDARY(ng) % vbarG_east(JLB,1), &
|
---|
793 | & BOUNDARY(ng) % vbar_east(JLB), &
|
---|
794 | & update)
|
---|
795 | # endif
|
---|
796 | # ifdef SOUTH_M2OBC
|
---|
797 | CALL set_ngfld (ng, iNLM, idU2bc(isouth), ILB, IUB, 1, &
|
---|
798 | & 1, Lm(ng)+1, 1, &
|
---|
799 | & BOUNDARY(ng) % ubarG_south(ILB,1), &
|
---|
800 | & BOUNDARY(ng) % ubar_south(ILB), &
|
---|
801 | & update)
|
---|
802 | CALL set_ngfld (ng, iNLM, idV2bc(isouth), ILB, IUB, 1, &
|
---|
803 | & 0, Lm(ng)+1, 1, &
|
---|
804 | & BOUNDARY(ng) % vbarG_south(ILB,1), &
|
---|
805 | & BOUNDARY(ng) % vbar_south(ILB), &
|
---|
806 | & update)
|
---|
807 | # endif
|
---|
808 | # ifdef NORTH_M2OBC
|
---|
809 | CALL set_ngfld (ng, iNLM, idU2bc(inorth), ILB, IUB, 1, &
|
---|
810 | & 1, Lm(ng)+1, 1, &
|
---|
811 | & BOUNDARY(ng) % ubarG_north(ILB,1), &
|
---|
812 | & BOUNDARY(ng) % ubar_north(ILB), &
|
---|
813 | & update)
|
---|
814 | CALL set_ngfld (ng, iNLM, idV2bc(inorth), ILB, IUB, 1, &
|
---|
815 | & 0, Lm(ng)+1, 1, &
|
---|
816 | & BOUNDARY(ng) % vbarG_north(ILB,1), &
|
---|
817 | & BOUNDARY(ng) % vbar_north(ILB), &
|
---|
818 | & update)
|
---|
819 | # endif
|
---|
820 | END IF
|
---|
821 | # endif
|
---|
822 | # ifdef SOLVE3D
|
---|
823 | # ifdef ANA_M3OBC
|
---|
824 | CALL ana_m3obc (ng, tile, iNLM)
|
---|
825 | # else
|
---|
826 | IF (SOUTH_WEST_TEST) THEN
|
---|
827 | # ifdef WEST_M3OBC
|
---|
828 | CALL set_ngfld (ng, iNLM, idU3bc(iwest), JLB, JUB, N(ng), &
|
---|
829 | & 0, Mm(ng)+1, N(ng), &
|
---|
830 | & BOUNDARY(ng) % uG_west(JLB,1,1), &
|
---|
831 | & BOUNDARY(ng) % u_west(JLB,1), &
|
---|
832 | & update)
|
---|
833 | CALL set_ngfld (ng, iNLM, idV3bc(iwest), JLB, JUB, N(ng), &
|
---|
834 | & 1, Mm(ng)+1, N(ng), &
|
---|
835 | & BOUNDARY(ng) % vG_west(JLB,1,1), &
|
---|
836 | & BOUNDARY(ng) % v_west(JLB,1), &
|
---|
837 | & update)
|
---|
838 | # endif
|
---|
839 | # ifdef EAST_M3OBC
|
---|
840 | CALL set_ngfld (ng, iNLM, idU3bc(ieast), JLB, JUB, N(ng), &
|
---|
841 | & 0, Mm(ng)+1, N(ng), &
|
---|
842 | & BOUNDARY(ng) % uG_east(JLB,1,1), &
|
---|
843 | & BOUNDARY(ng) % u_east(JLB,1), &
|
---|
844 | & update)
|
---|
845 | CALL set_ngfld (ng, iNLM, idV3bc(ieast), JLB, JUB, N(ng), &
|
---|
846 | & 1, Mm(ng)+1, N(ng), &
|
---|
847 | & BOUNDARY(ng) % vG_east(JLB,1,1), &
|
---|
848 | & BOUNDARY(ng) % v_east(JLB,1), &
|
---|
849 | & update)
|
---|
850 | # endif
|
---|
851 | # ifdef SOUTH_M3OBC
|
---|
852 | CALL set_ngfld (ng, iNLM, idU3bc(isouth), ILB, IUB, N(ng), &
|
---|
853 | & 1, Lm(ng)+1, N(ng), &
|
---|
854 | & BOUNDARY(ng) % uG_south(ILB,1,1), &
|
---|
855 | & BOUNDARY(ng) % u_south(ILB,1), &
|
---|
856 | & update)
|
---|
857 | CALL set_ngfld (ng, iNLM, idV3bc(isouth), ILB, IUB, N(ng), &
|
---|
858 | & 0, Lm(ng)+1, N(ng), &
|
---|
859 | & BOUNDARY(ng) % vG_south(ILB,1,1), &
|
---|
860 | & BOUNDARY(ng) % v_south(ILB,1), &
|
---|
861 | & update)
|
---|
862 | # endif
|
---|
863 | # ifdef NORTH_M3OBC
|
---|
864 | CALL set_ngfld (ng, iNLM, idU3bc(inorth), ILB, IUB, N(ng), &
|
---|
865 | & 1, Lm(ng)+1, N(ng), &
|
---|
866 | & BOUNDARY(ng) % uG_north(ILB,1,1), &
|
---|
867 | & BOUNDARY(ng) % u_north(ILB,1), &
|
---|
868 | & update)
|
---|
869 | CALL set_ngfld (ng, iNLM, idV3bc(inorth), ILB, IUB, N(ng), &
|
---|
870 | & 0, Lm(ng)+1, N(ng), &
|
---|
871 | & BOUNDARY(ng) % vG_north(ILB,1,1), &
|
---|
872 | & BOUNDARY(ng) % v_north(ILB,1), &
|
---|
873 | & update)
|
---|
874 | # endif
|
---|
875 | END IF
|
---|
876 | # endif
|
---|
877 | # ifdef ANA_TOBC
|
---|
878 | CALL ana_tobc (ng, tile, iNLM)
|
---|
879 | # else
|
---|
880 | IF (SOUTH_WEST_TEST) THEN
|
---|
881 | DO itrc=1,NT(ng)
|
---|
882 | # ifdef WEST_TOBC
|
---|
883 | CALL set_ngfld (ng, iNLM, idTbry(iwest,itrc), &
|
---|
884 | & JLB, JUB, N(ng), 0, Mm(ng)+1, N(ng), &
|
---|
885 | & BOUNDARY(ng) % tG_west(JLB,1,1,itrc), &
|
---|
886 | & BOUNDARY(ng) % t_west(JLB,1,itrc), &
|
---|
887 | & update)
|
---|
888 | # endif
|
---|
889 | # ifdef EAST_TOBC
|
---|
890 | CALL set_ngfld (ng, iNLM, idTbry(ieast,itrc), &
|
---|
891 | & JLB, JUB, N(ng), 0, Mm(ng)+1, N(ng), &
|
---|
892 | & BOUNDARY(ng) % tG_east(JLB,1,1,itrc), &
|
---|
893 | & BOUNDARY(ng) % t_east(JLB,1,itrc), &
|
---|
894 | & update)
|
---|
895 | # endif
|
---|
896 | # ifdef SOUTH_TOBC
|
---|
897 | CALL set_ngfld (ng, iNLM, idTbry(isouth,itrc), &
|
---|
898 | & ILB, IUB, N(ng), 0, Lm(ng)+1, N(ng), &
|
---|
899 | & BOUNDARY(ng) % tG_south(ILB,1,1,itrc), &
|
---|
900 | & BOUNDARY(ng) % t_south(ILB,1,itrc), &
|
---|
901 | & update)
|
---|
902 | # endif
|
---|
903 | # ifdef NORTH_TOBC
|
---|
904 | CALL set_ngfld (ng, iNLM, idTbry(inorth,itrc), &
|
---|
905 | & ILB, IUB, N(ng), 0, Lm(ng)+1, N(ng), &
|
---|
906 | & BOUNDARY(ng) % tG_north(ILB,1,1,itrc), &
|
---|
907 | & BOUNDARY(ng) % t_north(ILB,1,itrc), &
|
---|
908 | & update)
|
---|
909 | # endif
|
---|
910 | END DO
|
---|
911 | END IF
|
---|
912 | # endif
|
---|
913 | # endif
|
---|
914 | # ifdef REFINED_GRID
|
---|
915 | END IF
|
---|
916 | # endif
|
---|
917 | # endif
|
---|
918 |
|
---|
919 | # ifdef ZCLIMATOLOGY
|
---|
920 | !
|
---|
921 | !-----------------------------------------------------------------------
|
---|
922 | ! Set sea surface height climatology (m).
|
---|
923 | !-----------------------------------------------------------------------
|
---|
924 | !
|
---|
925 | # ifdef ANA_SSH
|
---|
926 | CALL ana_ssh (ng, tile, iNLM)
|
---|
927 | # else
|
---|
928 | CALL set_2dfld_tile (ng, tile, iNLM, idSSHc, &
|
---|
929 | & LBi, UBi, LBj, UBj, &
|
---|
930 | & CLIMA(ng)%sshG, &
|
---|
931 | & CLIMA(ng)%ssh, &
|
---|
932 | & update)
|
---|
933 | # endif
|
---|
934 | # endif
|
---|
935 |
|
---|
936 | # ifdef M2CLIMATOLOGY
|
---|
937 | !
|
---|
938 | !-----------------------------------------------------------------------
|
---|
939 | ! Set 2D momentum climatology (m/s).
|
---|
940 | !-----------------------------------------------------------------------
|
---|
941 | !
|
---|
942 | # ifdef ANA_M2CLIMA
|
---|
943 | CALL ana_m2clima (ng, tile, iNLM)
|
---|
944 | # else
|
---|
945 | CALL set_2dfld_tile (ng, tile, iNLM, idUbcl, &
|
---|
946 | & LBi, UBi, LBj, UBj, &
|
---|
947 | & CLIMA(ng)%ubarclmG, &
|
---|
948 | & CLIMA(ng)%ubarclm, &
|
---|
949 | & update)
|
---|
950 | CALL set_2dfld_tile (ng, tile, iNLM, idVbcl, &
|
---|
951 | & LBi, UBi, LBj, UBj, &
|
---|
952 | & CLIMA(ng)%vbarclmG, &
|
---|
953 | & CLIMA(ng)%vbarclm, &
|
---|
954 | & update)
|
---|
955 | # endif
|
---|
956 | # endif
|
---|
957 |
|
---|
958 | # if defined SOLVE3D && defined TCLIMATOLOGY
|
---|
959 | # if defined REFINED_GRID
|
---|
960 | IF (ng.eq.1) THEN
|
---|
961 | # endif
|
---|
962 | !
|
---|
963 | !-----------------------------------------------------------------------
|
---|
964 | ! Set tracer climatology.
|
---|
965 | !-----------------------------------------------------------------------
|
---|
966 | !
|
---|
967 | # ifdef ANA_TCLIMA
|
---|
968 | CALL ana_tclima (ng, tile, iNLM)
|
---|
969 | # else
|
---|
970 | DO itrc=1,NAT
|
---|
971 | CALL set_3dfld_tile (ng, tile, iNLM, idTclm(itrc), &
|
---|
972 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
973 | & CLIMA(ng)%tclmG(:,:,:,:,itrc), &
|
---|
974 | & CLIMA(ng)%tclm (:,:,:,itrc), &
|
---|
975 | & update)
|
---|
976 | END DO
|
---|
977 | # endif
|
---|
978 | # if defined REFINED_GRID
|
---|
979 | END IF
|
---|
980 | # endif
|
---|
981 | # endif
|
---|
982 |
|
---|
983 | # if defined SOLVE3D && defined M3CLIMATOLOGY
|
---|
984 | !
|
---|
985 | !-----------------------------------------------------------------------
|
---|
986 | ! Set 3D momentum climatology (m/s).
|
---|
987 | !-----------------------------------------------------------------------
|
---|
988 | !
|
---|
989 | # ifdef ANA_M3CLIMA
|
---|
990 | CALL ana_m3clima (ng, tile, iNLM)
|
---|
991 | # else
|
---|
992 | CALL set_3dfld_tile (ng, tile, iNLM, idUclm, &
|
---|
993 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
994 | & CLIMA(ng)%uclmG, &
|
---|
995 | & CLIMA(ng)%uclm, &
|
---|
996 | & update)
|
---|
997 | CALL set_3dfld_tile (ng, tile, iNLM, idVclm, &
|
---|
998 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
999 | & CLIMA(ng)%vclmG, &
|
---|
1000 | & CLIMA(ng)%vclm, &
|
---|
1001 | & update)
|
---|
1002 | # endif
|
---|
1003 | # endif
|
---|
1004 |
|
---|
1005 | # if defined NUDGING_SSH
|
---|
1006 | !
|
---|
1007 | !-----------------------------------------------------------------------
|
---|
1008 | ! Set sea surface height observations and error variance.
|
---|
1009 | !-----------------------------------------------------------------------
|
---|
1010 | !
|
---|
1011 | IF (assi_SSH(ng)) THEN
|
---|
1012 | CALL set_2dfld_tile (ng, tile, iNLM, idSSHo, &
|
---|
1013 | & LBi, UBi, LBj, UBj, &
|
---|
1014 | & OBS(ng)%SSHdat, &
|
---|
1015 | & OBS(ng)%SSHobs, &
|
---|
1016 | & update)
|
---|
1017 | CALL set_2dfld_tile (ng, tile, iNLM, idSSHe, &
|
---|
1018 | & LBi, UBi, LBj, UBj, &
|
---|
1019 | & OBS(ng)%EdatSSH, &
|
---|
1020 | & OBS(ng)%EobsSSH, &
|
---|
1021 | & update)
|
---|
1022 | IF (.not.update.and.SOUTH_WEST_TEST) THEN
|
---|
1023 | update_SSH(ng)=.FALSE.
|
---|
1024 | exit_flag=0
|
---|
1025 | END IF
|
---|
1026 | END IF
|
---|
1027 | # endif
|
---|
1028 |
|
---|
1029 | # ifdef SOLVE3D
|
---|
1030 | # if defined NUDGING_SST || defined ASSIMILATION_SST
|
---|
1031 | !
|
---|
1032 | !-----------------------------------------------------------------------
|
---|
1033 | ! Set sea surface temperature observations and error variance.
|
---|
1034 | !-----------------------------------------------------------------------
|
---|
1035 | !
|
---|
1036 | IF (assi_SST(ng)) THEN
|
---|
1037 | # ifdef NUDGING_SST
|
---|
1038 | CALL set_2dfld_tile (ng, tile, iNLM, idSSTo, &
|
---|
1039 | & LBi, UBi, LBj, UBj, &
|
---|
1040 | & OBS(ng)%SSTdat, &
|
---|
1041 | & OBS(ng)%SSTobs, &
|
---|
1042 | & update)
|
---|
1043 | CALL set_2dfld_tile (ng, tile, iNLM, idSSTe, &
|
---|
1044 | & LBi, UBi, LBj, UBj, &
|
---|
1045 | & OBS(ng)%EdatSST, &
|
---|
1046 | & OBS(ng)%EobsSST, &
|
---|
1047 | & update)
|
---|
1048 | IF (.not.update.and.SOUTH_WEST_TEST) THEN
|
---|
1049 | update_SST(ng)=.FALSE.
|
---|
1050 | update_T(itemp,ng)=.FALSE.
|
---|
1051 | exit_flag=0
|
---|
1052 | END IF
|
---|
1053 | # endif
|
---|
1054 | !
|
---|
1055 | ! Extend sea surface temperature and associated error variance using
|
---|
1056 | ! provided basis function polynomials.
|
---|
1057 | !
|
---|
1058 | IF (extend_SST(ng).and.update_SST(ng)) THEN
|
---|
1059 | IF (SOUTH_WEST_TEST) THEN
|
---|
1060 | update_SST(ng)=.FALSE.
|
---|
1061 | update_T(itemp,ng)=.TRUE.
|
---|
1062 | # ifdef ASSIMILATION_SST
|
---|
1063 | tTobs(1,itemp,ng)=Vtime(1,idSSTo,ng)
|
---|
1064 | tsTobs(itemp,ng)=Vtime(1,idSSTo,ng)*day2sec
|
---|
1065 | Finfo(7,idSSTo,ng)=tsTobs(itemp,ng)
|
---|
1066 | Finfo(7,idSSTe,ng)=tsTobs(itemp,ng)
|
---|
1067 | EobsTmin(itemp,ng)=Finfo(8,idSSTe,ng)
|
---|
1068 | EobsTmax(itemp,ng)=Finfo(9,idSSTe,ng)
|
---|
1069 | # endif
|
---|
1070 | END IF
|
---|
1071 | DO k=1,N(ng)
|
---|
1072 | DO j=JstrR,JendR
|
---|
1073 | DO i=IstrR,IendR
|
---|
1074 | Zr=GRID(ng)%z_r(i,j,k)/GRID(ng)%h(i,j)
|
---|
1075 | cff=perr_SST(npSST(ng),ng)
|
---|
1076 | cff1=pcoef_SST(npSST(ng),ng)
|
---|
1077 | DO order=npSST(ng)-1,0,-1
|
---|
1078 | cff=Zr*cff+perr_SST(order,ng)
|
---|
1079 | cff1=Zr*cff+pcoef_SST(order,ng)
|
---|
1080 | END DO
|
---|
1081 | OBS(ng)%EobsT(i,j,k,itemp)=MIN(1.0_r8,cff* &
|
---|
1082 | & OBS(ng)%EobsSST(i,j))
|
---|
1083 | OBS(ng)%Tobs(i,j,k,itemp)=cff1*OBS(ng)%SSTobs(i,j)
|
---|
1084 | END DO
|
---|
1085 | END DO
|
---|
1086 | END DO
|
---|
1087 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
1088 | CALL exchange_r3d_tile (ng, tile, &
|
---|
1089 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1090 | & OBS(ng)%EobsT(:,:,:,itemp))
|
---|
1091 | CALL exchange_r3d_tile (ng, tile, &
|
---|
1092 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1093 | & OBS(ng)%Tobs(:,:,:,itemp))
|
---|
1094 | # endif
|
---|
1095 | # ifdef DISTRIBUTE
|
---|
1096 | CALL mp_exchange3d (ng, tile, iNLM, 2, &
|
---|
1097 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1098 | & NghostPoints, EWperiodic, NSperiodic, &
|
---|
1099 | & OBS(ng)%EobsT(:,:,:,itemp), &
|
---|
1100 | & OBS(ng)%Tobs(:,:,:,itemp))
|
---|
1101 | # endif
|
---|
1102 | END IF
|
---|
1103 | END IF
|
---|
1104 | # endif
|
---|
1105 |
|
---|
1106 | # if defined NUDGING_T
|
---|
1107 | !
|
---|
1108 | !-----------------------------------------------------------------------
|
---|
1109 | ! Set tracers observations and error variance.
|
---|
1110 | !-----------------------------------------------------------------------
|
---|
1111 | !
|
---|
1112 | DO itrc=1,NAT
|
---|
1113 | IF (assi_T(itrc,ng)) THEN
|
---|
1114 | CALL set_3dfld_tile (ng, tile, iNLM, idTobs(itrc), &
|
---|
1115 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1116 | & OBS(ng)%Tdat(:,:,:,:,itrc), &
|
---|
1117 | & OBS(ng)%Tobs(:,:,:,itrc), &
|
---|
1118 | & update)
|
---|
1119 | CALL set_3dfld_tile (ng, tile, iNLM, idTerr(itrc), &
|
---|
1120 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1121 | & OBS(ng)%EdatT(:,:,:,:,itrc), &
|
---|
1122 | & OBS(ng)%EobsT(:,:,:,itrc), &
|
---|
1123 | & update)
|
---|
1124 | IF (.not.update.and.SOUTH_WEST_TEST) THEN
|
---|
1125 | update_T(itrc,ng)=.FALSE.
|
---|
1126 | exit_flag=0
|
---|
1127 | END IF
|
---|
1128 | END IF
|
---|
1129 | END DO
|
---|
1130 | # endif
|
---|
1131 |
|
---|
1132 | # if defined NUDGING_UVsur || defined ASSIMILATION_UVsur
|
---|
1133 | !
|
---|
1134 | !-----------------------------------------------------------------------
|
---|
1135 | ! Set surface current observations and error variance.
|
---|
1136 | !-----------------------------------------------------------------------
|
---|
1137 | !
|
---|
1138 | IF (assi_UVsur(ng)) THEN
|
---|
1139 | # ifdef NUDGING_UVsur
|
---|
1140 | CALL set_2dfld_tile (ng, tile, iNLM, idUsur, &
|
---|
1141 | & LBi, UBi, LBj, UBj, &
|
---|
1142 | & OBS(ng)%Usurdat, &
|
---|
1143 | & OBS(ng)%Usur, &
|
---|
1144 | & update)
|
---|
1145 | CALL set_2dfld_tile (ng, tile, iNLM, idVsur, &
|
---|
1146 | & LBi, UBi, LBj, UBj, &
|
---|
1147 | & OBS(ng)%Vsurdat, &
|
---|
1148 | & OBS(ng)%Vsur, &
|
---|
1149 | & update)
|
---|
1150 | CALL set_2dfld_tile (ng, tile, iNLM, idUVse, &
|
---|
1151 | & LBi, UBi, LBj, UBj, &
|
---|
1152 | & OBS(ng)%EdatVsur, &
|
---|
1153 | & OBS(ng)%EobsVsur, &
|
---|
1154 | & update)
|
---|
1155 | IF (.not.update.and.SOUTH_WEST_TEST) THEN
|
---|
1156 | update_UVsur(ng)=.FALSE.
|
---|
1157 | update_UV(ng)=.FALSE.
|
---|
1158 | exit_flag=0
|
---|
1159 | END IF
|
---|
1160 | # endif
|
---|
1161 | !
|
---|
1162 | ! Extend surface currents observations and associated error variance
|
---|
1163 | ! using provided basis function polynomials.
|
---|
1164 | !
|
---|
1165 | IF (extend_UV(ng).and.update_UVsur(ng)) THEN
|
---|
1166 | IF (SOUTH_WEST_TEST) THEN
|
---|
1167 | update_UVsur(ng)=.FALSE.
|
---|
1168 | update_UV(ng)=.TRUE.
|
---|
1169 | # ifdef ASSIMILATION_UVsur
|
---|
1170 | tVobs(1,ng)=Vtime(1,idVsur,ng)
|
---|
1171 | tsVobs(ng)=Vtime(1,idVsur,ng)*day2sec
|
---|
1172 | Finfo(7,idUsur,ng)=tsVobs(ng)
|
---|
1173 | Finfo(7,idVsur,ng)=tsVobs(ng)
|
---|
1174 | Finfo(7,idUVse,ng)=tsVobs(ng)
|
---|
1175 | EobsUVmin(ng)=Finfo(8,idUVse,ng)
|
---|
1176 | EobsUVmax(ng)=Finfo(9,idUVse,ng)
|
---|
1177 | # endif
|
---|
1178 | END IF
|
---|
1179 | DO k=1,N(ng)
|
---|
1180 | DO j=JstrR,JendR
|
---|
1181 | DO i=IstrR,IendR
|
---|
1182 | Zr=GRID(ng)%z_r(i,j,k)/GRID(ng)%h(i,j)
|
---|
1183 | cff=perr_V(npUV(ng),ng)
|
---|
1184 | DO order=npUV(ng)-1,0,-1
|
---|
1185 | cff=Zr*cff+perr_V(order,ng)
|
---|
1186 | END DO
|
---|
1187 | OBS(ng)%EobsUV(i,j,k)=MIN(1.0_r8,cff+ &
|
---|
1188 | & OBS(ng)%EobsVsur(i,j))
|
---|
1189 | END DO
|
---|
1190 | END DO
|
---|
1191 | DO j=JstrV-1,Jend
|
---|
1192 | DO i=IstrU-1,Iend
|
---|
1193 | Zr=GRID(ng)%z_r(i,j,k)/GRID(ng)%h(i,j)
|
---|
1194 | cff1=pcoef_U(npUV(ng),ng)
|
---|
1195 | cff2=pcoef_V(npUV(ng),ng)
|
---|
1196 | DO order=npUV(ng)-1,0,-1
|
---|
1197 | cff1=Zr*cff1+pcoef_U(order,ng)
|
---|
1198 | cff2=Zr*cff2+pcoef_V(order,ng)
|
---|
1199 | END DO
|
---|
1200 | work1(i,j)=cff1*OBS(ng)%Usur(i,j)- &
|
---|
1201 | & cff2*OBS(ng)%Vsur(i,j)
|
---|
1202 | work2(i,j)=cff2*OBS(ng)%Usur(i,j)+ &
|
---|
1203 | & cff1*OBS(ng)%Vsur(i,j)
|
---|
1204 | END DO
|
---|
1205 | END DO
|
---|
1206 | DO j=Jstr,Jend
|
---|
1207 | DO i=IstrU,Iend
|
---|
1208 | OBS(ng)%Uobs(i,j,k)=0.5_r8*(work1(i-1,j)+work1(i,j))
|
---|
1209 | END DO
|
---|
1210 | IF (j.ge.JstrV) THEN
|
---|
1211 | DO i=Istr,Iend
|
---|
1212 | OBS(ng)%Vobs(i,j,k)=0.5_r8*(work2(i,j-1)+work2(i,j))
|
---|
1213 | END DO
|
---|
1214 | END IF
|
---|
1215 | END DO
|
---|
1216 | END DO
|
---|
1217 | # if defined EW_PERIODIC || defined NS_PERIODIC
|
---|
1218 | CALL exchange_r3d_tile (ng, tile, &
|
---|
1219 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1220 | & OBS(ng)%EobsUV)
|
---|
1221 | # endif
|
---|
1222 | # ifdef DISTRIBUTE
|
---|
1223 | CALL mp_exchange3d (ng, tile, iNLM, 1, &
|
---|
1224 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1225 | & NghostPoints, EWperiodic, NSperiodic, &
|
---|
1226 | & OBS(ng)%EobsUV)
|
---|
1227 | # endif
|
---|
1228 | END IF
|
---|
1229 | END IF
|
---|
1230 | # endif
|
---|
1231 |
|
---|
1232 | # ifdef NUDGING_UV
|
---|
1233 | !
|
---|
1234 | !-----------------------------------------------------------------------
|
---|
1235 | ! Set horizontal current observations and error variance.
|
---|
1236 | !-----------------------------------------------------------------------
|
---|
1237 | !
|
---|
1238 | IF (assi_UV(ng)) THEN
|
---|
1239 | CALL set_3dfld_tile (ng, tile, iNLM, idUobs, &
|
---|
1240 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1241 | & OBS(ng)%Udat, &
|
---|
1242 | & OBS(ng)%Uobs, &
|
---|
1243 | & update)
|
---|
1244 | CALL set_3dfld_tile (ng, tile, iNLM, idVobs, &
|
---|
1245 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1246 | & OBS(ng)%Vdat, &
|
---|
1247 | & OBS(ng)%Vobs, &
|
---|
1248 | & update)
|
---|
1249 | CALL set_3dfld_tile (ng, tile, iNLM, idUVer, &
|
---|
1250 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1251 | & OBS(ng)%EdatUV, &
|
---|
1252 | & OBS(ng)%EobsUV, &
|
---|
1253 | & update)
|
---|
1254 | IF (.not.update.and.SOUTH_WEST_TEST) THEN
|
---|
1255 | update_UV(ng)=.FALSE.
|
---|
1256 | exit_flag=0
|
---|
1257 | END IF
|
---|
1258 | END IF
|
---|
1259 | # endif
|
---|
1260 | # endif
|
---|
1261 |
|
---|
1262 | # if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
|
---|
1263 | !
|
---|
1264 | !-----------------------------------------------------------------------
|
---|
1265 | ! Adjust surface forcing, add 4DVAR increments
|
---|
1266 | !-----------------------------------------------------------------------
|
---|
1267 | !
|
---|
1268 | CALL frc_NLadjust (ng, tile, Lfinp(ng))
|
---|
1269 | # endif
|
---|
1270 |
|
---|
1271 | # if defined W4DPSAS || defined NLM_OUTER
|
---|
1272 | !
|
---|
1273 | !-----------------------------------------------------------------------
|
---|
1274 | ! Set weak contraint forcing.
|
---|
1275 | !-----------------------------------------------------------------------
|
---|
1276 | !
|
---|
1277 | IF (FrequentImpulse) THEN
|
---|
1278 | !
|
---|
1279 | ! Set free-surface forcing.
|
---|
1280 | !
|
---|
1281 | CALL set_2dfld_tile (ng, tile, iNLM, idFsur, &
|
---|
1282 | & LBi, UBi, LBj, UBj, &
|
---|
1283 | & OCEAN(ng)%zetaG, &
|
---|
1284 | & OCEAN(ng)%f_zeta, &
|
---|
1285 | & update)
|
---|
1286 | !
|
---|
1287 | ! Set 2D momentum forcing.
|
---|
1288 | !
|
---|
1289 | CALL set_2dfld_tile (ng, tile, iNLM, idUbar, &
|
---|
1290 | & LBi, UBi, LBj, UBj, &
|
---|
1291 | & OCEAN(ng)%ubarG, &
|
---|
1292 | & OCEAN(ng)%f_ubar, &
|
---|
1293 | & update)
|
---|
1294 | CALL set_2dfld_tile (ng, tile, iNLM, idVbar, &
|
---|
1295 | & LBi, UBi, LBj, UBj, &
|
---|
1296 | & OCEAN(ng)%vbarG, &
|
---|
1297 | & OCEAN(ng)%f_vbar, &
|
---|
1298 | & update)
|
---|
1299 |
|
---|
1300 | # ifdef SOLVE3D
|
---|
1301 | !
|
---|
1302 | ! Set 3D momentum.
|
---|
1303 | !
|
---|
1304 | CALL set_3dfld_tile (ng, tile, iNLM, idUvel, &
|
---|
1305 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1306 | & OCEAN(ng)%uG, &
|
---|
1307 | & OCEAN(ng)%f_u, &
|
---|
1308 | & update)
|
---|
1309 | CALL set_3dfld_tile (ng, tile, iNLM, idVvel, &
|
---|
1310 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1311 | & OCEAN(ng)%vG, &
|
---|
1312 | & OCEAN(ng)%f_v, &
|
---|
1313 | & update)
|
---|
1314 | !
|
---|
1315 | ! Set 3D tracers.
|
---|
1316 | !
|
---|
1317 | DO itrc=1,NT(ng)
|
---|
1318 | CALL set_3dfld_tile (ng, tile, iNLM, idTvar(itrc), &
|
---|
1319 | & LBi, UBi, LBj, UBj, 1, N(ng), &
|
---|
1320 | & OCEAN(ng)%tG(:,:,:,:,itrc), &
|
---|
1321 | & OCEAN(ng)%f_t(:,:,:,itrc), &
|
---|
1322 | & update)
|
---|
1323 | END DO
|
---|
1324 | # endif
|
---|
1325 | END IF
|
---|
1326 | # endif
|
---|
1327 |
|
---|
1328 | RETURN
|
---|
1329 | END SUBROUTINE set_data_tile
|
---|
1330 | #else
|
---|
1331 | SUBROUTINE set_data
|
---|
1332 | RETURN
|
---|
1333 | END SUBROUTINE set_data
|
---|
1334 | #endif
|
---|