/[MITgcm]/MITgcm/model/src/forward_step.F
ViewVC logotype

Contents of /MITgcm/model/src/forward_step.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph

Revision 1.134 - (show annotations) (download)
Wed Mar 8 06:36:39 2006 UTC (18 years, 2 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint58b_post
Changes since 1.133: +93 -26 lines 
o Another overhaul of store dirs. for NLFS to eliminate
  "hidden" recomputations.
o TBD: "hidden" mom_vecinv recomp. in dynamics
1 C $Header: /u/gcmpack/MITgcm/model/src/forward_step.F,v 1.133 2006/03/06 18:25:49 heimbach Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 #ifdef ALLOW_OFFLINE
8 # include "OFFLINE_OPTIONS.h"
9 #endif
10 #ifdef ALLOW_GMREDI
11 # include "GMREDI_OPTIONS.h"
12 #endif
13
14 CBOP
15 C !ROUTINE: FORWARD_STEP
16 C !INTERFACE:
17 SUBROUTINE FORWARD_STEP( iloop, myTime, myIter, myThid )
18
19 C !DESCRIPTION: \bv
20 C *==================================================================
21 C | SUBROUTINE forward_step
22 C | o Run the ocean model and, optionally, evaluate a cost function.
23 C *==================================================================
24 C |
25 C | THE_MAIN_LOOP is the toplevel routine for the Tangent Linear and
26 C | Adjoint Model Compiler (TAMC). For this purpose the initialization
27 C | of the model was split into two parts. Those parameters that do
28 C | not depend on a specific model run are set in INITIALISE_FIXED,
29 C | whereas those that do depend on the specific realization are
30 C | initialized in INITIALISE_VARIA.
31 C |
32 C *==================================================================
33 C \ev
34
35 C !USES:
36 IMPLICIT NONE
37 C == Global variables ==
38 #include "SIZE.h"
39 #include "EEPARAMS.h"
40 #include "PARAMS.h"
41 #include "DYNVARS.h"
42
43 #ifdef ALLOW_MNC
44 #include "MNC_PARAMS.h"
45 EXTERNAL DIFFERENT_MULTIPLE
46 LOGICAL DIFFERENT_MULTIPLE
47 #endif
48
49 #ifdef HAVE_SIGREG
50 #include "SIGREG.h"
51 #endif
52
53 #ifdef ALLOW_SHAP_FILT
54 # include "SHAP_FILT.h"
55 #endif
56 #ifdef ALLOW_ZONAL_FILT
57 # include "ZONAL_FILT.h"
58 #endif
59 #ifdef COMPONENT_MODULE
60 # include "CPL_PARAMS.h"
61 #endif
62
63 #ifdef ALLOW_AUTODIFF_TAMC
64 # include "FFIELDS.h"
65
66 # include "tamc.h"
67 # include "ctrl.h"
68 # include "ctrl_dummy.h"
69 # include "cost.h"
70 # include "EOS.h"
71 # ifdef NONLIN_FRSURF
72 # include "GRID.h"
73 # endif
74 # ifdef ALLOW_EXF
75 # include "exf_fields.h"
76 # include "exf_clim_fields.h"
77 # ifdef ALLOW_BULKFORMULAE
78 # include "exf_constants.h"
79 # endif
80 # endif
81 # ifdef ALLOW_OBCS
82 # include "OBCS.h"
83 # endif
84 # ifdef ALLOW_PTRACERS
85 # include "PTRACERS_SIZE.h"
86 # include "PTRACERS.h"
87 # endif
88 # ifdef ALLOW_CD_CODE
89 # include "CD_CODE_VARS.h"
90 # endif
91 # ifdef ALLOW_EBM
92 # include "EBM.h"
93 # endif
94 # ifdef EXACT_CONSERV
95 # include "SURFACE.h"
96 # endif
97 # ifdef ALLOW_KPP
98 # include "KPP.h"
99 # endif
100 # ifdef ALLOW_GMREDI
101 # include "GMREDI.h"
102 # endif
103 #endif /* ALLOW_AUTODIFF_TAMC */
104
105 C !LOCAL VARIABLES:
106 C == Routine arguments ==
107 C note: under the multi-threaded model myiter and
108 C mytime are local variables passed around as routine
109 C arguments. Although this is fiddly it saves the need to
110 C impose additional synchronisation points when they are
111 C updated.
112 C myIter - iteration counter for this thread
113 C myTime - time counter for this thread
114 C myThid - thread number for this instance of the routine.
115 INTEGER iloop
116 INTEGER myThid
117 INTEGER myIter
118 _RL myTime
119
120 C == Local variables ==
121 #ifdef COMPONENT_MODULE
122 INTEGER myItP1
123 #endif
124 CEOP
125
126 #ifdef ALLOW_DEBUG
127 IF ( debugLevel .GE. debLevB )
128 & CALL DEBUG_ENTER('FORWARD_STEP',myThid)
129 #endif
130
131 #ifdef ALLOW_AUTODIFF_TAMC
132 C-- Reset the model iteration counter and the model time.
133 myIter = nIter0 + (iloop-1)
134 myTime = startTime + float(iloop-1)*deltaTclock
135 #endif
136
137 #ifdef ALLOW_AUTODIFF_TAMC
138 c**************************************
139 #include "checkpoint_lev1_directives.h"
140 c**************************************
141 #endif
142
143 C-- Switch on/off diagnostics for snap-shot output:
144 #ifdef ALLOW_DIAGNOSTICS
145 IF ( useDiagnostics ) THEN
146 CALL DIAGNOSTICS_SWITCH_ONOFF( myTime, myIter, myThid )
147 C-- State-variables diagnostics
148 CALL TIMER_START('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
149 CALL DO_STATEVARS_DIAGS( myTime, 0, myIter, myThid )
150 CALL TIMER_STOP ('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
151 ENDIF
152 #endif
153
154 C-- Call Bulk-Formulae forcing package
155 #ifdef ALLOW_BULK_FORCE
156 IF ( useBulkForce ) THEN
157 #ifdef ALLOW_DEBUG
158 IF ( debugLevel .GE. debLevB )
159 & CALL DEBUG_CALL('BULKF_FIELDS_LOAD',myThid)
160 #endif
161 CALL TIMER_START('BULKF_FORCING [FORWARD_STEP]',mythid)
162 C- load all forcing fields at current time
163 CALL BULKF_FIELDS_LOAD( myTime, myIter, myThid )
164 C- calculate qnet and empmr (and wind stress)
165 CALL BULKF_FORCING( myTime, myIter, myThid )
166 CALL TIMER_STOP ('BULKF_FORCING [FORWARD_STEP]',mythid)
167 ENDIF
168 #endif /* ALLOW_BULK_FORCE */
169
170 C-- Call external forcing package
171 # ifdef ALLOW_EXF
172 # ifdef ALLOW_DEBUG
173 IF ( debugLevel .GE. debLevB )
174 & CALL DEBUG_CALL('EXF_GETFORCING',myThid)
175 # endif
176 CALL TIMER_START('EXF_GETFORCING [FORWARD_STEP]',mythid)
177 CALL EXF_GETFORCING( mytime, myiter, mythid )
178 CALL TIMER_STOP ('EXF_GETFORCING [FORWARD_STEP]',mythid)
179 # else /* ALLOW_EXF undef */
180 cph The following IF-statement creates an additional dependency
181 cph for the forcing fields requiring additional storing.
182 cph Therefore, the IF-statement will be put between CPP-OPTIONS,
183 cph assuming that ALLOW_SEAICE has not yet been differentiated.
184 # if (defined (ALLOW_SEAICE) || defined (ALLOW_EBM))
185 IF ( .NOT. useSEAICE .AND. .NOT. useEBM ) THEN
186 # endif
187 #ifdef ALLOW_DEBUG
188 IF ( debugLevel .GE. debLevB )
189 & CALL DEBUG_CALL('EXTERNAL_FIELDS_LOAD',myThid)
190 #endif
191 CALL TIMER_START('EXTERNAL_FIELDS_LOAD[FORWARD_STEP]',mythid)
192 CALL EXTERNAL_FIELDS_LOAD( mytime, myiter, mythid )
193 CALL TIMER_STOP ('EXTERNAL_FIELDS_LOAD[FORWARD_STEP]',mythid)
194 # ifdef NONLIN_FRSURF
195 CADJ STORE SST = comlev1, key = ikey_dynamics
196 CADJ STORE SSS = comlev1, key = ikey_dynamics
197 # ifdef SHORTWAVE_HEATING
198 CADJ STORE Qsw = comlev1, key = ikey_dynamics
199 # endif
200 # endif
201 # if (defined (ALLOW_SEAICE) || defined (ALLOW_EBM))
202 ENDIF
203 # endif
204 # endif /* ALLOW_EXF */
205
206 #ifdef ALLOW_AUTODIFF
207 c-- Add control vector for forcing and parameter fields
208 if ( myiter .EQ. nIter0 )
209 & CALL CTRL_MAP_FORCING (mythid)
210 #endif
211
212 #if (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_AUTODIFF_MONITOR))
213 C Include call to a dummy routine. Its adjoint will be
214 C called at the proper place in the adjoint code.
215 C The adjoint routine will print out adjoint values
216 C if requested. The location of the call is important,
217 C it has to be after the adjoint of the exchanges
218 C (DO_GTERM_BLOCKING_EXCHANGES).
219 CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )
220 cph I've commented this line since it may conflict with MITgcm's adjoint
221 cph However, need to check whether that's still consistent
222 cph with the ecco-branch (it should).
223 cph CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
224 #endif
225
226 # ifdef ALLOW_SEAICE
227 C-- Call sea ice model to compute forcing/external data fields. In
228 C addition to computing prognostic sea-ice variables and diagnosing the
229 C forcing/external data fields that drive the ocean model, SEAICE_MODEL
230 C also sets theta to the freezing point under sea-ice. The implied
231 C surface heat flux is then stored in variable surfaceTendencyTice,
232 C which is needed by KPP package (kpp_calc.F and kpp_transport_t.F)
233 C to diagnose surface buoyancy fluxes and for the non-local transport
234 C term. Because this call precedes model thermodynamics, temperature
235 C under sea-ice may not be "exactly" at the freezing point by the time
236 C theta is dumped or time-averaged.
237 IF ( useSEAICE ) THEN
238 #ifdef ALLOW_DEBUG
239 IF ( debugLevel .GE. debLevB )
240 & CALL DEBUG_CALL('SEAICE_MODEL',myThid)
241 #endif
242 CALL TIMER_START('SEAICE_MODEL [FORWARD_STEP]',myThid)
243 CALL SEAICE_MODEL( myTime, myIter, myThid )
244 CALL TIMER_STOP ('SEAICE_MODEL [FORWARD_STEP]',myThid)
245 ENDIF
246 # endif /* ALLOW_SEAICE */
247
248 #ifdef ALLOW_AUTODIFF_TAMC
249 # ifdef ALLOW_PTRACERS
250 cph this replaces _bibj storing of ptracer within thermodynamics
251 CADJ STORE ptracer = comlev1, key = ikey_dynamics
252 # endif
253 #endif
254
255 #ifdef ALLOW_OFFLINE
256 call OFFLINE_FIELDS_LOAD( myTime, myIter, myThid )
257 #endif
258
259 #ifdef ALLOW_PTRACERS
260 # ifdef ALLOW_GCHEM
261 IF ( useGCHEM ) THEN
262 #ifdef ALLOW_DEBUG
263 IF ( debugLevel .GE. debLevB )
264 & CALL DEBUG_CALL('GCHEM_FIELDS_LOAD',myThid)
265 #endif /* ALLOW_DEBUG */
266 CALL GCHEM_FIELDS_LOAD( mytime, myiter, mythid )
267 ENDIF
268 # endif
269 #endif
270
271 #ifdef ALLOW_RBCS
272 IF ( useRBCS ) THEN
273 CALL RBCS_FIELDS_LOAD( mytime, myiter, mythid )
274 ENDIF
275 #endif
276
277 #ifdef COMPONENT_MODULE
278 IF ( useCoupler .AND. cpl_earlyExpImpCall ) THEN
279 C Post coupling data that I export.
280 C Read in coupling data that I import.
281 CALL TIMER_START('CPL_EXPORT-IMPORT [FORWARD_STEP]',myThid)
282 CALL CPL_EXPORT_MY_DATA( myIter, myTime, myThid )
283 CALL CPL_IMPORT_EXTERNAL_DATA( myIter, myTime, myThid )
284 CALL TIMER_STOP ('CPL_EXPORT-IMPORT [FORWARD_STEP]',myThid)
285 ENDIF
286 #endif /* COMPONENT_MODULE */
287
288 #ifdef ALLOW_EBM
289 IF ( useEBM ) THEN
290 # ifdef ALLOW_DEBUG
291 IF ( debugLevel .GE. debLevB )
292 & CALL DEBUG_CALL('EBM',myThid)
293 # endif
294 CALL TIMER_START('EBM [FORWARD_STEP]',mythid)
295 CALL EBM_DRIVER ( myTime, myIter, myThid )
296 CALL TIMER_STOP ('EBM [FORWARD_STEP]',mythid)
297 ENDIF
298 #endif
299
300 C-- Step forward fields and calculate time tendency terms.
301
302 #ifdef ALLOW_DEBUG
303 IF ( debugLevel .GE. debLevB )
304 & CALL DEBUG_CALL('DO_ATMOSPHERIC_PHYS',myThid)
305 #endif
306 CALL TIMER_START('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
307 CALL DO_ATMOSPHERIC_PHYS( myTime, myIter, myThid )
308 CALL TIMER_STOP ('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
309
310 #ifdef ALLOW_AUTODIFF_TAMC
311 # ifdef ALLOW_KPP
312 CADJ STORE uvel = comlev1, key = ikey_dynamics
313 CADJ STORE vvel = comlev1, key = ikey_dynamics
314 # endif
315 # ifdef EXACT_CONSERV
316 cphCADJ STORE empmr = comlev1, key = ikey_dynamics
317 cphCADJ STORE pmepr = comlev1, key = ikey_dynamics
318 # endif
319 # ifdef NONLIN_FRSURF
320 cph-test
321 CADJ STORE hFacC = comlev1, key = ikey_dynamics
322 # endif
323 #endif /* ALLOW_AUTODIFF_TAMC */
324
325 #ifndef ALLOW_OFFLINE
326 #ifdef ALLOW_DEBUG
327 IF ( debugLevel .GE. debLevB )
328 & CALL DEBUG_CALL('DO_OCEANIC_PHYS',myThid)
329 #endif
330 CALL TIMER_START('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
331 CALL DO_OCEANIC_PHYS( myTime, myIter, myThid )
332 CALL TIMER_STOP ('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
333 #ifdef ALLOW_AUTODIFF_TAMC
334 CADJ STORE EmPmR = comlev1, key = ikey_dynamics
335 # ifdef EXACT_CONSERV
336 CADJ STORE pmepr = comlev1, key = ikey_dynamics
337 # endif
338 #endif
339 #endif
340 C
341 #ifdef ALLOW_AUTODIFF_TAMC
342 # ifdef NONLIN_FRSURF
343 cph-test
344 CADJ STORE hFac_surfC = comlev1, key = ikey_dynamics
345 CADJ STORE hfac_surfs = comlev1, key = ikey_dynamics
346 CADJ STORE hfac_surfw = comlev1, key = ikey_dynamics
347 CADJ STORE hFacC, hFacS, hFacW
348 CADJ & = comlev1, key = ikey_dynamics
349 CADJ STORE recip_hFacC, recip_hFacS, recip_hFacW
350 CADJ & = comlev1, key = ikey_dynamics
351 c
352 CADJ STORE surfaceforcingu = comlev1, key = ikey_dynamics
353 CADJ STORE surfaceforcingv = comlev1, key = ikey_dynamics
354 # endif
355 #endif /* ALLOW_AUTODIFF_TAMC */
356
357 #ifdef ALLOW_GCHEM
358 C GCHEM package is an interface for any bio-geochemical or
359 C ecosystem model you would like to include.
360 C If GCHEM_SEPARATE_FORCING is not defined, you are
361 C responsible for computing tendency terms for passive
362 C tracers and storing them on a 3DxNumPtracers-array called
363 C gchemTendency in GCHEM_CALC_TENDENCY. This tendency is then added
364 C to gPtr in ptracers_forcing later-on.
365 C If GCHEM_SEPARATE_FORCING is defined, you are reponsible for
366 C UPDATING ptracers directly in GCHEM_FORCING_SEP. This amounts
367 C to a completely separate time step that you have to implement
368 C yourself (Eulerian seems to be fine in most cases).
369 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
370 C CAVEAT: Up to now, when GCHEM is turned on the field ptracerForcingSurf,
371 C which is needed for KPP is not set properly. ptracerForcingSurf must
372 C be treated differently depending on whether GCHEM_SEPARATE_FORCING
373 C is define or not. TBD.
374 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
375 IF ( useGCHEM ) THEN
376 #ifdef ALLOW_DEBUG
377 IF ( debugLevel .GE. debLevB )
378 & CALL DEBUG_CALL('GCHEM_CALC_TENDENCY',myThid)
379 #endif
380 CALL TIMER_START('GCHEM_CALC_TENDENCY [FORWARD_STEP]',myThid)
381 CALL GCHEM_CALC_TENDENCY( myTime, myIter, myThid )
382 CALL TIMER_STOP ('GCHEM_CALC_TENDENCY [FORWARD_STEP]',myThid)
383 ENDIF
384 #endif /* ALLOW_GCHEM */
385
386 #ifdef ALLOW_AUTODIFF_TAMC
387 cph needed to be moved here from do_oceanic_physics
388 cph to be visible down the road
389 c
390 CADJ STORE surfaceForcingS = comlev1, key = ikey_dynamics
391 CADJ STORE surfaceForcingT = comlev1, key = ikey_dynamics
392 CADJ STORE surfaceForcingTice = comlev1, key = ikey_dynamics
393 ctest(
394 CADJ STORE IVDConvCount = comlev1, key = ikey_dynamics
395 ctest)
396 # ifdef ALLOW_PTRACERS
397 CADJ STORE surfaceForcingPtr = comlev1, key = ikey_dynamics
398 # endif
399 c
400 # ifdef ALLOW_GMREDI
401 CADJ STORE Kwx = comlev1, key = ikey_dynamics
402 CADJ STORE Kwy = comlev1, key = ikey_dynamics
403 CADJ STORE Kwz = comlev1, key = ikey_dynamics
404 # ifdef GM_BOLUS_ADVEC
405 CADJ STORE GM_PsiX = comlev1, key = ikey_dynamics
406 CADJ STORE GM_PsiY = comlev1, key = ikey_dynamics
407 # endif
408 # endif
409 c
410 # ifdef ALLOW_KPP
411 CADJ STORE KPPghat = comlev1, key = ikey_dynamics
412 CADJ STORE KPPfrac = comlev1, key = ikey_dynamics
413 CADJ STORE KPPdiffKzS = comlev1, key = ikey_dynamics
414 CADJ STORE KPPdiffKzT = comlev1, key = ikey_dynamics
415 # endif
416 c
417 # ifdef NONLIN_FRSURF
418 CADJ STORE etaH = comlev1, key = ikey_dynamics
419 # ifdef ALLOW_CD_CODE
420 CADJ STORE etanm1 = comlev1, key = ikey_dynamics
421 # endif
422 # endif
423 #endif /* ALLOW_AUTODIFF_TAMC */
424
425 IF ( .NOT.staggerTimeStep ) THEN
426 #ifdef ALLOW_DEBUG
427 IF ( debugLevel .GE. debLevB )
428 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
429 #endif
430 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
431 CALL THERMODYNAMICS( myTime, myIter, myThid )
432 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
433 C-- if not staggerTimeStep: end
434 ENDIF
435 c #ifdef ALLOW_NONHYDROSTATIC
436 IF ( implicitIntGravWave ) THEN
437 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
438 CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )
439 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
440 ENDIF
441 c #endif
442
443 #ifdef COMPONENT_MODULE
444 IF ( useCoupler .AND. .NOT.cpl_earlyExpImpCall ) THEN
445 C Post coupling data that I export.
446 C Read in coupling data that I import.
447 myItP1 = myIter + 1
448 CALL TIMER_START('CPL_EXPORT-IMPORT [FORWARD_STEP]',myThid)
449 CALL CPL_EXPORT_MY_DATA( myItP1, myTime, myThid )
450 CALL CPL_IMPORT_EXTERNAL_DATA( myItP1, myTime, myThid )
451 CALL TIMER_STOP ('CPL_EXPORT-IMPORT [FORWARD_STEP]',myThid)
452 # ifndef ALLOW_AIM
453 IF ( useRealFreshWaterFlux ) THEN
454 CALL OCN_APPLY_IMPORT( .FALSE., myTime, myIter, myThid )
455 ENDIF
456 # endif
457 ENDIF
458 #endif /* COMPONENT_MODULE */
459
460 #ifdef ALLOW_AUTODIFF_TAMC
461 # ifdef NONLIN_FRSURF
462 CADJ STORE hFacC = comlev1, key = ikey_dynamics
463 CADJ STORE hFacS = comlev1, key = ikey_dynamics
464 CADJ STORE hFacW = comlev1, key = ikey_dynamics
465 CADJ STORE recip_hFacC = comlev1, key = ikey_dynamics
466 CADJ STORE recip_hFacS = comlev1, key = ikey_dynamics
467 CADJ STORE recip_hFacW = comlev1, key = ikey_dynamics
468 CADJ STORE etaN = comlev1, key = ikey_dynamics
469 # endif
470 #endif
471 C-- Step forward fields and calculate time tendency terms.
472 #ifndef ALLOW_OFFLINE
473 #ifndef ALLOW_AUTODIFF_TAMC
474 IF ( momStepping ) THEN
475 #endif
476 #ifdef ALLOW_DEBUG
477 IF ( debugLevel .GE. debLevB )
478 & CALL DEBUG_CALL('DYNAMICS',myThid)
479 #endif
480 CALL TIMER_START('DYNAMICS [FORWARD_STEP]',mythid)
481 CALL DYNAMICS( myTime, myIter, myThid )
482 CALL TIMER_STOP ('DYNAMICS [FORWARD_STEP]',mythid)
483 # ifndef ALLOW_AUTODIFF_TAMC
484 ENDIF
485 # endif
486 #endif
487 C
488 #ifdef ALLOW_AUTODIFF_TAMC
489 # ifdef NONLIN_FRSURF
490 cph-test
491 CADJ STORE gU, gV = comlev1, key = ikey_dynamics
492 # endif
493 #endif
494
495 C-- Update time-counter
496 myIter = nIter0 + iLoop
497 myTime = startTime + deltaTClock * float(iLoop)
498
499 #ifdef ALLOW_MNC
500 C Update the default next iter for MNC
501 IF ( useMNC ) THEN
502 CALL MNC_CW_CITER_SETG( 1, 1, -1, myIter , myThid )
503
504 C TODO: Logic should be added here so that users can specify, on
505 C a per-citer-group basis, when it is time to update the
506 C "current" (and not just the "next") iteration
507
508 C TODO: the following is just a temporary band-aid (mostly, for
509 C Baylor) until someone writes a routine that better handles time
510 C boundaries such as weeks, months, years, etc.
511 IF ( mnc_filefreq .GT. 0 ) THEN
512 IF (DIFFERENT_MULTIPLE(mnc_filefreq,myTime,deltaTClock))
513 & THEN
514 CALL MNC_CW_CITER_SETG( 1, 1, myIter, -1 , myThid )
515 ENDIF
516 ENDIF
517 ENDIF
518 #endif
519
520 C-- Update geometric factors:
521 #ifdef NONLIN_FRSURF
522 C- update hfacC,W,S and recip_hFac according to etaH(n+1) :
523 IF ( nonlinFreeSurf.GT.0) THEN
524 IF ( select_rStar.GT.0 ) THEN
525 # ifndef DISABLE_RSTAR_CODE
526 # ifdef ALLOW_AUTODIFF_TAMC
527 cph-test
528 CADJ STORE hFacC = comlev1, key = ikey_dynamics
529 CADJ STORE hFacS = comlev1, key = ikey_dynamics
530 CADJ STORE hFacW = comlev1, key = ikey_dynamics
531 CADJ STORE recip_hFacC = comlev1, key = ikey_dynamics
532 CADJ STORE recip_hFacS = comlev1, key = ikey_dynamics
533 CADJ STORE recip_hFacW = comlev1, key = ikey_dynamics
534 # endif
535 CALL TIMER_START('UPDATE_R_STAR [FORWARD_STEP]',myThid)
536 CALL UPDATE_R_STAR( myTime, myIter, myThid )
537 CALL TIMER_STOP ('UPDATE_R_STAR [FORWARD_STEP]',myThid)
538 # ifdef ALLOW_AUTODIFF_TAMC
539 cph-test
540 CADJ STORE hFacC = comlev1, key = ikey_dynamics
541 CADJ STORE hFacS = comlev1, key = ikey_dynamics
542 CADJ STORE hFacW = comlev1, key = ikey_dynamics
543 CADJ STORE recip_hFacC = comlev1, key = ikey_dynamics
544 CADJ STORE recip_hFacS = comlev1, key = ikey_dynamics
545 CADJ STORE recip_hFacW = comlev1, key = ikey_dynamics
546 # endif
547 # endif /* DISABLE_RSTAR_CODE */
548 ELSE
549 #ifdef ALLOW_AUTODIFF_TAMC
550 CADJ STORE hFac_surfC, hFac_surfS, hFac_surfW
551 CADJ & = comlev1, key = ikey_dynamics
552 #endif
553 CALL TIMER_START('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
554 CALL UPDATE_SURF_DR( myTime, myIter, myThid )
555 CALL TIMER_STOP ('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
556 ENDIF
557 ENDIF
558 # ifdef ALLOW_AUTODIFF_TAMC
559 cph-test
560 CADJ STORE hFacC = comlev1, key = ikey_dynamics
561 CADJ STORE hFacS = comlev1, key = ikey_dynamics
562 CADJ STORE hFacW = comlev1, key = ikey_dynamics
563 CADJ STORE recip_hFacC = comlev1, key = ikey_dynamics
564 CADJ STORE recip_hFacS = comlev1, key = ikey_dynamics
565 CADJ STORE recip_hFacW = comlev1, key = ikey_dynamics
566 # endif
567 C- update also CG2D matrix (and preconditioner)
568 IF ( momStepping .AND. nonlinFreeSurf.GT.2 ) THEN
569 CALL TIMER_START('UPDATE_CG2D [FORWARD_STEP]',myThid)
570 CALL UPDATE_CG2D( myTime, myIter, myThid )
571 CALL TIMER_STOP ('UPDATE_CG2D [FORWARD_STEP]',myThid)
572 ENDIF
573 #endif /* NONLIN_FRSURF */
574
575 C-- Apply Filters to u*,v* before SOLVE_FOR_PRESSURE
576 #ifdef ALLOW_SHAP_FILT
577 IF (useSHAP_FILT .AND. shap_filt_uvStar) THEN
578 CALL TIMER_START('SHAP_FILT_UV [FORWARD_STEP]',myThid)
579 IF (implicDiv2Dflow.LT.1.) THEN
580 C-- Explicit+Implicit part of the Barotropic Flow Divergence
581 C => Filtering of uVel,vVel is necessary
582 CALL SHAP_FILT_APPLY_UV( uVel,vVel,
583 & myTime, myIter, myThid )
584 ENDIF
585 CALL SHAP_FILT_APPLY_UV( gU,gV,myTime,myIter,myThid)
586 CALL TIMER_STOP ('SHAP_FILT_UV [FORWARD_STEP]',myThid)
587 ENDIF
588 #endif
589 #ifdef ALLOW_ZONAL_FILT
590 IF (useZONAL_FILT .AND. zonal_filt_uvStar) THEN
591 CALL TIMER_START('ZONAL_FILT_UV [FORWARD_STEP]',myThid)
592 IF (implicDiv2Dflow.LT.1.) THEN
593 C-- Explicit+Implicit part of the Barotropic Flow Divergence
594 C => Filtering of uVel,vVel is necessary
595 CALL ZONAL_FILT_APPLY_UV( uVel, vVel, myThid )
596 ENDIF
597 CALL ZONAL_FILT_APPLY_UV( gU, gV, myThid )
598 CALL TIMER_STOP ('ZONAL_FILT_UV [FORWARD_STEP]',myThid)
599 ENDIF
600 #endif
601
602 C-- Solve elliptic equation(s).
603 C Two-dimensional only for conventional hydrostatic or
604 C three-dimensional for non-hydrostatic and/or IGW scheme.
605 #ifndef ALLOW_OFFLINE
606 IF ( momStepping ) THEN
607 #ifdef ALLOW_AUTODIFF_TAMC
608 # ifdef NONLIN_FRSURF
609 CADJ STORE uvel, vvel
610 CADJ & = comlev1, key = ikey_dynamics
611 CADJ STORE empmr,hfacs,hfacw
612 CADJ & = comlev1, key = ikey_dynamics
613 # endif
614 #endif
615 CALL TIMER_START('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
616 CALL SOLVE_FOR_PRESSURE(myTime, myIter, myThid)
617 CALL TIMER_STOP ('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
618 ENDIF
619 #endif
620
621 C-- Correct divergence in flow field and cycle time-stepping momentum
622 c IF ( momStepping ) THEN
623 #ifndef ALLOW_OFFLINE
624 CALL TIMER_START('MOM_CORRECTION_STEP [FORWARD_STEP]',myThid)
625 CALL MOMENTUM_CORRECTION_STEP(myTime, myIter, myThid)
626 CALL TIMER_STOP ('MOM_CORRECTION_STEP [FORWARD_STEP]',myThid)
627 #endif
628 c ENDIF
629
630 #ifdef EXACT_CONSERV
631 IF (exactConserv) THEN
632 #ifdef ALLOW_AUTODIFF_TAMC
633 cph-test
634 cphCADJ STORE etaH = comlev1, key = ikey_dynamics
635 #endif
636 C-- Update etaH(n+1) :
637 CALL TIMER_START('UPDATE_ETAH [FORWARD_STEP]',mythid)
638 CALL UPDATE_ETAH( myTime, myIter, myThid )
639 CALL TIMER_STOP ('UPDATE_ETAH [FORWARD_STEP]',mythid)
640 ENDIF
641 #endif /* EXACT_CONSERV */
642
643 #ifdef NONLIN_FRSURF
644 IF ( select_rStar.NE.0 ) THEN
645 # ifndef DISABLE_RSTAR_CODE
646 C-- r* : compute the future level thickness according to etaH(n+1)
647 CALL TIMER_START('CALC_R_STAR [FORWARD_STEP]',mythid)
648 CALL CALC_R_STAR(etaH, myTime, myIter, myThid )
649 CALL TIMER_STOP ('CALC_R_STAR [FORWARD_STEP]',mythid)
650 # endif /* DISABLE_RSTAR_CODE */
651 ELSEIF ( nonlinFreeSurf.GT.0) THEN
652 C-- compute the future surface level thickness according to etaH(n+1)
653 # ifdef ALLOW_AUTODIFF_TAMC
654 CADJ STORE etaH = comlev1, key = ikey_dynamics
655 # endif
656 CALL TIMER_START('CALC_SURF_DR [FORWARD_STEP]',mythid)
657 CALL CALC_SURF_DR(etaH, myTime, myIter, myThid )
658 CALL TIMER_STOP ('CALC_SURF_DR [FORWARD_STEP]',mythid)
659 ENDIF
660 # ifdef ALLOW_AUTODIFF_TAMC
661 cph-test
662 CADJ STORE hFac_surfC = comlev1, key = ikey_dynamics
663 # endif
664 #endif /* NONLIN_FRSURF */
665
666 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
667 IF ( staggerTimeStep ) THEN
668 C-- do exchanges of U,V (needed for multiDim) when using stagger time-step :
669 #ifdef ALLOW_DEBUG
670 IF ( debugLevel .GE. debLevB )
671 & CALL DEBUG_CALL('DO_STAGGER_FIELDS_EXCH.',myThid)
672 #endif
673 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
674 CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )
675 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
676
677 #ifdef ALLOW_DIAGNOSTICS
678 C-- State-variables diagnostics
679 IF ( usediagnostics ) THEN
680 CALL TIMER_START('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
681 CALL DO_STATEVARS_DIAGS( myTime, 1, myIter, myThid )
682 CALL TIMER_STOP ('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
683 ENDIF
684 #endif
685
686 #ifdef ALLOW_DEBUG
687 IF ( debugLevel .GE. debLevB )
688 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
689 #endif
690 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
691 CALL THERMODYNAMICS( myTime, myIter, myThid )
692 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
693
694 C-- if staggerTimeStep: end
695 ENDIF
696 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
697
698 #ifdef ALLOW_AUTODIFF_TAMC
699 cph This is needed because convective_adjustment calls
700 cph find_rho which may use pressure()
701 CADJ STORE totphihyd = comlev1, key = ikey_dynamics
702 #endif
703 C-- Cycle time-stepping Tracers arrays (T,S,+pTracers)
704 CALL TIMER_START('TRC_CORRECTION_STEP [FORWARD_STEP]',myThid)
705 CALL TRACERS_CORRECTION_STEP(myTime, myIter, myThid)
706 CALL TIMER_STOP ('TRC_CORRECTION_STEP [FORWARD_STEP]',myThid)
707
708 #ifdef ALLOW_GCHEM
709 C Add separate timestepping of chemical/biological/forcing
710 C of ptracers here in GCHEM_FORCING_SEP
711 IF ( useGCHEM ) THEN
712 #ifdef ALLOW_DEBUG
713 IF ( debugLevel .GE. debLevB )
714 & CALL DEBUG_CALL('GCHEM_FORCING_SEP',myThid)
715 #endif /* ALLOW_DEBUG */
716 CALL TIMER_START('GCHEM_FORCING_SEP [FORWARD_STEP]',myThid)
717 CALL GCHEM_FORCING_SEP( myTime,myIter,myThid )
718 CALL TIMER_STOP ('GCHEM_FORCING_SEP [FORWARD_STEP]',myThid)
719 ENDIF
720 #endif /* ALLOW_GCHEM */
721
722 C-- Do "blocking" sends and receives for tendency "overlap" terms
723 c CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
724 c CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )
725 c CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
726
727 C-- Do "blocking" sends and receives for field "overlap" terms
728 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
729 CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
730 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
731
732 #ifdef ALLOW_DIAGNOSTICS
733 IF ( usediagnostics ) THEN
734 CALL TIMER_START('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
735 CALL DO_STATEVARS_DIAGS( myTime, 2, myIter, myThid )
736 CALL TIMER_STOP ('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
737 ENDIF
738 #endif
739
740 C AMM
741 #ifdef ALLOW_GRIDALT
742 if (useGRIDALT) then
743 CALL GRIDALT_UPDATE(myThid)
744 endif
745 #endif
746 C AMM
747
748 C AMM
749 #ifdef ALLOW_FIZHI
750 if( useFIZHI) then
751 CALL TIMER_START('FIZHI [FORWARD_STEP]',mythid)
752 CALL STEP_FIZHI_CORR ( myTime, myIter, myThid, dTtracerLev(1) )
753 CALL TIMER_STOP('FIZHI [FORWARD_STEP]',mythid)
754 endif
755 #endif
756 C AMM
757
758 #ifdef ALLOW_FLT
759 C-- Calculate float trajectories
760 IF (useFLT) THEN
761 CALL TIMER_START('FLOATS [FORWARD_STEP]',myThid)
762 CALL FLT_MAIN(myIter,myTime, myThid)
763 CALL TIMER_STOP ('FLOATS [FORWARD_STEP]',myThid)
764 ENDIF
765 #endif
766
767 #ifdef ALLOW_TIMEAVE
768 C-- State-variables time-averaging
769 CALL TIMER_START('DO_STATEVARS_TAVE [FORWARD_STEP]',myThid)
770 CALL DO_STATEVARS_TAVE( myTime, myIter, myThid )
771 CALL TIMER_STOP ('DO_STATEVARS_TAVE [FORWARD_STEP]',myThid)
772 #endif
773
774 #ifndef ALLOW_OFFLINE
775 #ifdef ALLOW_MONITOR
776 C-- Check status of solution (statistics, cfl, etc...)
777 CALL TIMER_START('MONITOR [FORWARD_STEP]',myThid)
778 CALL MONITOR( myIter, myTime, myThid )
779 CALL TIMER_STOP ('MONITOR [FORWARD_STEP]',myThid)
780 #endif /* ALLOW_MONITOR */
781 #endif
782
783 #ifdef ALLOW_COST
784 C-- compare model with data and compute cost function
785 C-- this is done after exchanges to allow interpolation
786 CALL TIMER_START('COST_TILE [FORWARD_STEP]',myThid)
787 CALL COST_TILE ( mytime, myiter, myThid )
788 CALL TIMER_STOP ('COST_TILE [FORWARD_STEP]',myThid)
789 #endif
790
791 C-- Do IO if needed.
792 #ifdef ALLOW_OFFLINE
793 CALL TIMER_START('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
794 CALL OFFLINE_MODEL_IO( myTime, myIter, myThid )
795 CALL TIMER_STOP ('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
796 #else
797 CALL TIMER_START('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
798 CALL DO_THE_MODEL_IO( myTime, myIter, myThid )
799 CALL TIMER_STOP ('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
800 #endif
801
802 #ifdef HAVE_SIGREG
803 IF ( useSIGREG ) THEN
804 IF ( i_got_signal .GT. 0 ) THEN
805 CALL PACKAGES_WRITE_PICKUP(
806 I .TRUE., myTime, myIter, myThid )
807 #ifndef ALLOW_OFFLINE
808 CALL WRITE_CHECKPOINT(
809 I .TRUE., myTime, myIter, myThid )
810 #endif
811 STOP 'Checkpoint completed -- killed by signal handler'
812 ENDIF
813 ENDIF
814 #endif
815
816 C-- Save state for restarts
817 CALL TIMER_START('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
818 CALL PACKAGES_WRITE_PICKUP(
819 I .FALSE., myTime, myIter, myThid )
820 #ifndef ALLOW_OFFLINE
821 CALL WRITE_CHECKPOINT(
822 I .FALSE., myTime, myIter, myThid )
823 #endif
824 CALL TIMER_STOP ('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
825
826 #ifdef ALLOW_DEBUG
827 IF ( debugLevel .GE. debLevB )
828 & CALL DEBUG_LEAVE('FORWARD_STEP',myThid)
829 #endif
830
831 RETURN
832 END
  ViewVC Help
Powered by ViewVC 1.1.22