/[MITgcm]/MITgcm/pkg/ecco/the_main_loop.F
ViewVC logotype

Contents of /MITgcm/pkg/ecco/the_main_loop.F

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


Revision 1.16 - (show annotations) (download)
Thu Feb 10 05:33:22 2005 UTC (19 years, 4 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57d_post, eckpoint57e_pre
Changes since 1.15: +9 -1 lines
Add hooks for inAdExact

1 C $Header: /u/gcmpack/MITgcm/pkg/ecco/the_main_loop.F,v 1.15 2005/01/23 00:40:25 heimbach Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 #ifdef ALLOW_OBCS
8 # include "OBCS_OPTIONS.h"
9 #endif
10 #ifdef ALLOW_SEAICE
11 # include "SEAICE_OPTIONS.h"
12 #endif
13 #ifdef ALLOW_GMREDI
14 # include "GMREDI_OPTIONS.h"
15 #endif
16
17 subroutine the_main_loop( myTime, myIter, mythid )
18
19 c ==================================================================
20 c SUBROUTINE the_main_loop
21 c ==================================================================
22 c
23 c o Run the ocean model and evaluate the specified cost function.
24 c
25 c *the_main_loop* is the top-level 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*. In order to do a so called
31 c checkpointing during the adjoint calculation and to account for the
32 c typical data involved in oceanographic applications a call tree
33 c that is divided into yearly, monthly, daily, and step parts can
34 c be used.
35 c
36 c This routine is to be used in conjuction with the MITgcmuv release
37 c checkpoint 24.
38 c
39 c started: Christian Eckert eckert@mit.edu 30-Jun-1999
40 c
41 c changed: Christian Eckert eckert@mit.edu 14-Jul-1999
42 c
43 c - The call to mapping was moved to initialise_varia,
44 c since this routine has to be called before
45 c ini_predictor.
46 c
47 c Christian Eckert eckert@mit.edu 11-Feb-2000
48 c
49 c - Restructured the code in order to create a package
50 c for the MITgcmUV.
51 c
52 c Patrick Heimbach heimbach@mit.edu 3-Jun-2000
53 c - corrected computation of ikey_dynamics and
54 c added computation of ikey_dynamics for the case
55 c undef ALLOW_TAMC_CHECKPOINTING
56 c
57 c Patrick Heimbach heimbach@mit.edu 6-Jun-2000
58 c - corrected initialisation of comlev1 common blocks
59 c
60 c Dimitris Menemenlis menemenlis@jpl.nasa.gov 26-Feb-2003
61 c - modifications for pkg/seaice
62 c
63 c ==================================================================
64 c SUBROUTINE the_main_loop
65 c ==================================================================
66
67 implicit none
68
69 c == global variables ==
70
71 #include "SIZE.h"
72 #include "EEPARAMS.h"
73 #include "PARAMS.h"
74
75 c**************************************
76 #ifdef ALLOW_AUTODIFF_TAMC
77
78 c These includes are needed for
79 c AD-checkpointing.
80 c They provide the fields to be stored.
81
82 # include "GRID.h"
83 # include "DYNVARS.h"
84 # include "FFIELDS.h"
85 # include "EOS.h"
86 # include "GAD.h"
87
88 # ifdef ALLOW_CD_CODE
89 # include "CD_CODE_VARS.h"
90 # endif
91 # ifdef ALLOW_PTRACERS
92 # include "PTRACERS_SIZE.h"
93 # include "PTRACERS.h"
94 # endif
95 # ifdef ALLOW_NONHYDROSTATIC
96 # include "CG3D.h"
97 # endif
98 # ifdef EXACT_CONSERV
99 # include "SURFACE.h"
100 # endif
101 # ifdef ALLOW_OBCS
102 # include "OBCS.h"
103 # endif
104 # ifdef ALLOW_EXF
105 # include "exf_fields.h"
106 # include "exf_clim_fields.h"
107 # ifdef ALLOW_BULKFORMULAE
108 # include "exf_constants.h"
109 # endif
110 # endif /* ALLOW_EXF */
111 # ifdef ALLOW_SEAICE
112 # include "SEAICE.h"
113 # endif
114 # ifdef ALLOW_KPP
115 # include "KPP.h"
116 # endif
117 # ifdef ALLOW_GMREDI
118 # include "GMREDI.h"
119 # endif
120 # ifdef ALLOW_DIVIDED_ADJOINT_MPI
121 # include "mpif.h"
122 # endif
123
124 # include "tamc.h"
125 # include "ctrl.h"
126 # include "ctrl_dummy.h"
127 # include "cost.h"
128 # include "ecco_cost.h"
129
130 #endif /* ALLOW_AUTODIFF_TAMC */
131 c**************************************
132
133 c == routine arguments ==
134 c note: under the multi-threaded model myiter and
135 c mytime are local variables passed around as routine
136 c arguments. Although this is fiddly it saves the need to
137 c impose additional synchronisation points when they are
138 c updated.
139 c myiter - iteration counter for this thread
140 c mytime - time counter for this thread
141 c mythid - thread number for this instance of the routine.
142 integer mythid
143 integer myiter
144 _RL mytime
145
146 c == local variables ==
147
148 integer bi,bj
149 integer iloop
150 integer mydate(4)
151 #ifdef ALLOW_SNAPSHOTS
152 character yprefix*3
153 #endif
154
155 #ifdef ALLOW_TAMC_CHECKPOINTING
156 integer ilev_1
157 integer ilev_2
158 integer ilev_3
159 integer max_lev2
160 integer max_lev3
161 #endif
162
163 c-- == end of interface ==
164
165 #ifndef DISABLE_DEBUGMODE
166 IF ( debugLevel .GE. debLevB )
167 & CALL DEBUG_ENTER('THE_MAIN_LOOP',myThid)
168 #endif
169
170 #ifdef ALLOW_AUTODIFF_TAMC
171 c-- Initialize storage for the initialisations.
172 CADJ INIT tapelev_ini_bibj_k = USER
173 CADJ INIT tapelev_init = USER
174 #ifdef AUTODIFF_2_LEVEL_CHECKPOINT
175 CADJ INIT tapelev2 = USER
176 #else
177 CADJ INIT tapelev3 = USER
178 #endif
179 # ifdef ALLOW_DIVIDED_ADJOINT
180 CADJ INIT onetape = user
181 cphCADJ INIT onetape = common, 1
182 cph We want to avoid common blocks except in the inner loop.
183 cph Reason: the active write and consecutive read may occur
184 cph in separate model executions for which the info
185 cph in common blocks are lost.
186 cph Thus, we can only store real values (no integers)
187 cph because we only have active file handling to real available.
188 # endif
189 # ifdef ALLOW_TAMC_CHECKPOINTING
190 ikey_dynamics = 1
191 # endif
192 #endif /* ALLOW_AUTODIFF_TAMC */
193
194 CALL TIMER_START('ECCO SPIN-UP', mythid)
195
196 c-- Get the current date.
197 call CAL_TIMESTAMP( myiter, mytime, mydate, mythid )
198
199 C-- Set initial conditions (variable arrays)
200 #ifndef DISABLE_DEBUGMODE
201 IF ( debugLevel .GE. debLevB )
202 & CALL DEBUG_CALL('INITIALISE_VARIA',myThid)
203 #endif
204 CALL TIMER_START('INITIALISE_VARIA [THE_MAIN_LOOP]', mythid)
205 CALL INITIALISE_VARIA( mythid )
206 CALL TIMER_STOP ('INITIALISE_VARIA [THE_MAIN_LOOP]', mythid)
207
208 #ifdef ALLOW_MONITOR
209 #ifdef ALLOW_DEBUG
210 IF (debugMode) CALL DEBUG_CALL('MONITOR',myThid)
211 #endif
212 C-- Check status of solution (statistics, cfl, etc...)
213 CALL TIMER_START('MONITOR [THE_MAIN_LOOP]', mythid)
214 CALL MONITOR( myIter, myTime, myThid )
215 CALL TIMER_STOP ('MONITOR [THE_MAIN_LOOP]', mythid)
216 #endif /* ALLOW_MONITOR */
217
218 C-- Do IO if needed (Dump for start state).
219 #ifdef ALLOW_DEBUG
220 IF (debugMode) CALL DEBUG_CALL('DO_THE_MODEL_IO',myThid)
221 #endif
222
223 #ifdef ALLOW_OFFLINE
224 CALL TIMER_START('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
225 CALL OFFLINE_MODEL_IO( myTime, myIter, myThid )
226 CALL TIMER_STOP ('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
227 #else
228 CALL TIMER_START('DO_THE_MODEL_IO [THE_MAIN_LOOP]', mythid)
229 CALL DO_THE_MODEL_IO( myTime, myIter, mythid )
230 CALL TIMER_STOP ('DO_THE_MODEL_IO [THE_MAIN_LOOP]', mythid)
231 #endif
232
233 call timer_stop ('ECCO SPIN-UP', mythid)
234 _BARRIER
235
236 c-- Do the model integration.
237 call timer_start('ECCO MAIN LOOP',mythid)
238
239 c >>>>>>>>>>>>>>>>>>>>>>>>>>> LOOP <<<<<<<<<<<<<<<<<<<<<<<<<<<<
240 c >>>>>>>>>>>>>>>>>>>>>>>>>>> STARTS <<<<<<<<<<<<<<<<<<<<<<<<<<<<
241
242 #ifdef ALLOW_AUTODIFF_TAMC
243 #ifdef ALLOW_TAMC_CHECKPOINTING
244
245 max_lev3=nTimeSteps/(nchklev_1*nchklev_2)+1
246 max_lev2=nTimeSteps/nchklev_1+1
247
248 c**************************************
249 #ifdef ALLOW_DIVIDED_ADJOINT
250 CADJ loop = divided
251 #endif
252 c**************************************
253
254 #ifndef AUTODIFF_2_LEVEL_CHECKPOINT
255
256 do ilev_3 = 1,nchklev_3
257 if(ilev_3.le.max_lev3) then
258 c**************************************
259 #include "checkpoint_lev3_directives.h"
260 c**************************************
261
262 c-- Initialise storage for the middle loop.
263 CADJ INIT tapelev2 = USER
264
265 #endif /* AUTODIFF_2_LEVEL_CHECKPOINT */
266
267 do ilev_2 = 1,nchklev_2
268 if(ilev_2.le.max_lev2) then
269 c**************************************
270 #include "checkpoint_lev2_directives.h"
271 c**************************************
272
273 c**************************************
274 #ifdef ALLOW_AUTODIFF_TAMC
275 c-- Initialize storage for the innermost loop.
276 c-- Always check common block sizes for the checkpointing!
277 c--
278 CADJ INIT comlev1 = COMMON,nchklev_1
279 CADJ INIT comlev1_bibj = COMMON,nchklev_1*nsx*nsy*nthreads_chkpt
280 CADJ INIT comlev1_bibj_k = COMMON,nchklev_1*nsx*nsy*nr*nthreads_chkpt
281 c--
282 # ifdef ALLOW_KPP
283 CADJ INIT comlev1_kpp = COMMON,nchklev_1*nsx*nsy
284 CADJ INIT comlev1_kpp_k = COMMON,nchklev_1*nsx*nsy*nr
285 # endif /* ALLOW_KPP */
286 c--
287 # ifdef ALLOW_GMREDI
288 CADJ INIT comlev1_gmredi_k_gad
289 CADJ & = COMMON,nchklev_1*nsx*nsy*nr*nthreads_chkpt*maxpass
290 # endif /* ALLOW_GMREDI */
291 c--
292 # ifdef ALLOW_PTRACERS
293 CADJ INIT comlev1_bibj_ptracers = COMMON,
294 CADJ & nchklev_1*nsx*nsy*nthreads_chkpt*PTRACERS_num
295 # endif /* ALLOW_PTRACERS */
296 c--
297 # ifndef DISABLE_MULTIDIM_ADVECTION
298 CADJ INIT comlev1_bibj_k_gad
299 CADJ & = COMMON,nchklev_1*nsx*nsy*nr*nthreads_chkpt*maxpass
300 CADJ INIT comlev1_bibj_k_gad_pass
301 CADJ & = COMMON,nchklev_1*nsx*nsy*nr*nthreads_chkpt*maxpass*maxcube
302 # endif /* DISABLE_MULTIDIM_ADVECTION */
303 c--
304 # if (defined (ALLOW_EXF) && defined (ALLOW_BULKFORMULAE))
305 CADJ INIT comlev1_exf_1
306 CADJ & = COMMON,nchklev_1*snx*nsx*sny*nsy*nthreads_chkpt
307 CADJ INIT comlev1_exf_2
308 CADJ & = COMMON,niter_bulk*nchklev_1*snx*nsx*sny*nsy*nthreads_chkpt
309 # endif
310 c--
311 # ifdef ALLOW_SEAICE
312 # ifdef SEAICE_ALLOW_DYNAMICS
313 CADJ INIT comlev1_lsr = COMMON,nchklev_1*2
314 # endif
315 # ifdef SEAICE_MULTILEVEL
316 CADJ INIT comlev1_multdim
317 CADJ & = COMMON,nchklev_1*nsx*nsy*nthreads_chkpt*multdim
318 # endif
319 # endif /* ALLOW_SEAICE */
320 c--
321 #endif /* ALLOW_AUTODIFF_TAMC */
322 c**************************************
323
324 do ilev_1 = 1,nchklev_1
325
326 c-- The if-statement below introduces a some flexibility in the
327 c-- choice of the 3-tupel ( nchklev_1, nchklev_2, nchklev_3 ).
328 c--
329 c-- Requirement: nchklev_1*nchklev_2*nchklev_3 .ge. nTimeSteps .
330
331 iloop = (ilev_2 - 1)*nchklev_1 + ilev_1
332 #ifndef AUTODIFF_2_LEVEL_CHECKPOINT
333 & + (ilev_3 - 1)*nchklev_2*nchklev_1
334 #endif
335
336 if ( iloop .le. nTimeSteps ) then
337
338 #else /* ALLOW_TAMC_CHECKPOINTING undefined */
339 c-- Initialise storage for the reference trajectory without TAMC check-
340 c-- pointing.
341 CADJ INIT history = USER
342 CADJ INIT comlev1_bibj = COMMON,nchklev_0*nsx*nsy*nthreads_chkpt
343 CADJ INIT comlev1_bibj_k = COMMON,nchklev_0*nsx*nsy*nr*nthreads_chkpt
344 CADJ INIT comlev1_kpp = COMMON,nchklev_0*nsx*nsy
345
346 c-- Check the choice of the checkpointing parameters in relation
347 c-- to nTimeSteps: (nchklev_0 .ge. nTimeSteps)
348 if (nchklev_0 .lt. nTimeSteps) then
349 print*
350 print*, ' the_main_loop: ',
351 & 'TAMC checkpointing parameter nchklev_0 = ',
352 & nchklev_0
353 print*, ' is not consistent with nTimeSteps = ',
354 & nTimeSteps
355 stop ' ... stopped in the_main_loop.'
356 endif
357
358 do iloop = 1, nTimeSteps
359
360 #endif /* ALLOW_TAMC_CHECKPOINTING */
361
362 #else /* ALLOW_AUTODIFF_TAMC undefined */
363 c-- Start the main loop of ecco_Objfunc. Automatic differentiation is
364 c-- NOT enabled.
365 do iloop = 1, nTimeSteps
366 #endif /* ALLOW_AUTODIFF_TAMC */
367
368 #ifdef ALLOW_TAMC_CHECKPOINTING
369 nIter0 = INT( startTime/deltaTClock )
370 ikey_dynamics = ilev_1
371 #endif
372
373 c-- Set the model iteration counter and the model time.
374 myiter = nIter0 + (iloop-1)
375 mytime = startTime + float(iloop-1)*deltaTclock
376
377 #ifdef ALLOW_AUTODIFF_TAMC
378 CALL AUTODIFF_INADMODE_UNSET( myThid )
379 #endif
380
381 #ifdef ALLOW_COST
382
383 c-- Accumulate time averages of temperature, salinity, and SSH.
384 call timer_start('COST_AVERAGESFIELDS [ECCO MAIN]', mythid)
385 call cost_averagesFields( mytime, mythid )
386 call timer_stop ('COST_AVERAGESFIELDS [ECCO MAIN]', mythid)
387 #ifdef ALLOW_COST_ATLANTIC
388 c-- Compute meridional heat transport
389 call timer_start('cost_atlantic [ECCO MAIN]', mythid)
390 call cost_atlantic( mytime, myiter,mythid )
391 call timer_stop ('cost_atlantic [ECCO MAIN]', mythid)
392 #endif
393 #endif /* ALLOW_COST */
394
395 #ifdef ALLOW_AUTODIFF_TAMC
396 c**************************************
397 #include "checkpoint_lev1_directives.h"
398 c**************************************
399 #endif
400
401 #ifndef DISABLE_DEBUGMODE
402 IF ( debugLevel .GE. debLevB )
403 & CALL DEBUG_CALL('EXF_GETFORCING',myThid)
404 #endif
405 CALL TIMER_START('EXF_GETFORCING [FORWARD_STEP]',mythid)
406 CALL EXF_GETFORCING( mytime, myiter, mythid )
407 CALL TIMER_STOP ('EXF_GETFORCING [FORWARD_STEP]',mythid)
408
409 #ifdef ALLOW_SEAICE
410 cph this simple runtime flag causes a lot of recomp.
411 IF ( useSEAICE ) THEN
412 #ifndef DISABLE_DEBUGMODE
413 IF ( debugLevel .GE. debLevB )
414 & CALL DEBUG_CALL('SEAICE_MODEL',myThid)
415 #endif
416 #ifdef ALLOW_AUTODIFF_TAMC
417 CADJ STORE area = comlev1, key = ikey_dynamics
418 #endif
419 CALL TIMER_START('SEAICE_MODEL [FORWARD_STEP]',myThid)
420 CALL SEAICE_MODEL( myTime, myIter, myThid )
421 CALL TIMER_STOP ('SEAICE_MODEL [FORWARD_STEP]',myThid)
422 #ifdef ALLOW_COST_ICE
423 CALL COST_ICE ( myTime, myIter, myThid )
424 #endif
425 ENDIF
426 #endif /* ALLOW_SEAICE */
427
428 #ifdef ALLOW_AUTODIFF_TAMC
429 # ifdef ALLOW_PTRACERS
430 cph this replaces _bibj storing of ptracer within thermodynamics
431 CADJ STORE ptracer = comlev1, key = ikey_dynamics
432 # endif
433 #endif
434
435 #if (defined (ALLOW_AUTODIFF_TAMC) && \
436 defined (ALLOW_AUTODIFF_MONITOR))
437 C Include call to a dummy routine. Its adjoint will be
438 C called at the proper place in the adjoint code.
439 C The adjoint routine will print out adjoint values
440 C if requested. The location of the call is important,
441 C it has to be after the adjoint of the exchanges
442 C (DO_GTERM_BLOCKING_EXCHANGES).
443 CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )
444 #endif
445
446 #ifdef ALLOW_EBM
447 IF ( useEBM ) THEN
448 # ifdef ALLOW_DEBUG
449 IF ( debugLevel .GE. debLevB )
450 & CALL DEBUG_CALL('EBM',myThid)
451 # endif
452 CALL TIMER_START('EBM [FORWARD_STEP]',mythid)
453 CALL EBM_DRIVER ( myTime, myIter, myThid )
454 CALL TIMER_STOP ('EBM [FORWARD_STEP]',mythid)
455 ENDIF
456 #endif
457
458 C-- Step forward fields and calculate time tendency terms.
459
460 #ifdef ALLOW_DEBUG
461 IF ( debugLevel .GE. debLevB )
462 & CALL DEBUG_CALL('DO_ATMOSPHERIC_PHYS',myThid)
463 #endif
464 CALL TIMER_START('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
465 CALL DO_ATMOSPHERIC_PHYS( myTime, myIter, myThid )
466 CALL TIMER_STOP ('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
467
468 #ifndef ALLOW_OFFLINE
469 #ifdef ALLOW_DEBUG
470 IF ( debugLevel .GE. debLevB )
471 & CALL DEBUG_CALL('DO_OCEANIC_PHYS',myThid)
472 #endif
473 CALL TIMER_START('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
474 CALL DO_OCEANIC_PHYS( myTime, myIter, myThid )
475 CALL TIMER_STOP ('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
476 #endif
477
478 #ifdef ALLOW_AUTODIFF_TAMC
479 cph needed to be moved here from do_oceanic_physics
480 cph to be visible down the road
481 c
482 CADJ STORE surfaceForcingS = comlev1, key = ikey_dynamics
483 CADJ STORE surfaceForcingT = comlev1, key = ikey_dynamics
484 CADJ STORE surfaceForcingTice = comlev1, key = ikey_dynamics
485 ctest(
486 CADJ STORE IVDConvCount = comlev1, key = ikey_dynamics
487 ctest)
488 # ifdef ALLOW_PTRACERS
489 CADJ STORE surfaceForcingPtr = comlev1, key = ikey_dynamics
490 # endif
491 c
492 # ifdef ALLOW_GMREDI
493 CADJ STORE Kwx = comlev1, key = ikey_dynamics
494 CADJ STORE Kwy = comlev1, key = ikey_dynamics
495 CADJ STORE Kwz = comlev1, key = ikey_dynamics
496 # ifdef GM_BOLUS_ADVEC
497 CADJ STORE GM_PsiX = comlev1, key = ikey_dynamics
498 CADJ STORE GM_PsiY = comlev1, key = ikey_dynamics
499 # endif
500 # endif
501 c
502 # ifdef ALLOW_KPP
503 CADJ STORE KPPghat = comlev1, key = ikey_dynamics
504 CADJ STORE KPPfrac = comlev1, key = ikey_dynamics
505 CADJ STORE KPPdiffKzS = comlev1, key = ikey_dynamics
506 CADJ STORE KPPdiffKzT = comlev1, key = ikey_dynamics
507 # endif
508 #endif /* ALLOW_AUTODIFF_TAMC */
509
510
511 IF ( .NOT.staggerTimeStep ) THEN
512 #ifdef ALLOW_DEBUG
513 IF ( debugLevel .GE. debLevB )
514 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
515 #endif
516 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
517 CALL THERMODYNAMICS( myTime, myIter, myThid )
518 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
519 C-- if not staggerTimeStep: end
520 ENDIF
521
522 C-- Step forward fields and calculate time tendency terms.
523 #ifndef ALLOW_OFFLINE
524 #ifndef ALLOW_AUTODIFF_TAMC
525 IF ( momStepping ) THEN
526 #endif
527 #ifdef ALLOW_DEBUG
528 IF ( debugLevel .GE. debLevB )
529 & CALL DEBUG_CALL('DYNAMICS',myThid)
530 #endif
531 CALL TIMER_START('DYNAMICS [FORWARD_STEP]',mythid)
532 CALL DYNAMICS( myTime, myIter, myThid )
533 CALL TIMER_STOP ('DYNAMICS [FORWARD_STEP]',mythid)
534 #ifndef ALLOW_AUTODIFF_TAMC
535 ENDIF
536 #endif
537 #endif
538
539 #ifdef ALLOW_NONHYDROSTATIC
540 C-- Step forward W field in N-H algorithm
541 IF ( momStepping .AND. nonHydrostatic ) THEN
542 #ifdef ALLOW_DEBUG
543 IF ( debugLevel .GE. debLevB )
544 & CALL DEBUG_CALL('CALC_GW',myThid)
545 #endif
546 CALL TIMER_START('CALC_GW [FORWARD_STEP]',myThid)
547 CALL CALC_GW(myThid)
548 CALL TIMER_STOP ('CALC_GW [FORWARD_STEP]',myThid)
549 ENDIF
550 #endif
551
552 C-- Update time-counter
553 myIter = nIter0 + iLoop
554 myTime = startTime + deltaTClock * float(iLoop)
555
556 C-- Update geometric factors:
557 #ifdef NONLIN_FRSURF
558 C- update hfacC,W,S and recip_hFac according to etaH(n+1) :
559 IF ( nonlinFreeSurf.GT.0) THEN
560 IF ( select_rStar.GT.0 ) THEN
561 CALL TIMER_START('UPDATE_R_STAR [FORWARD_STEP]',myThid)
562 CALL UPDATE_R_STAR( myTime, myIter, myThid )
563 CALL TIMER_STOP ('UPDATE_R_STAR [FORWARD_STEP]',myThid)
564 ELSE
565 CALL TIMER_START('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
566 CALL UPDATE_SURF_DR( myTime, myIter, myThid )
567 CALL TIMER_STOP ('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
568 ENDIF
569 ENDIF
570 C- update also CG2D matrix (and preconditioner)
571 IF ( momStepping .AND. nonlinFreeSurf.GT.2 ) THEN
572 CALL TIMER_START('UPDATE_CG2D [FORWARD_STEP]',myThid)
573 CALL UPDATE_CG2D( myTime, myIter, myThid )
574 CALL TIMER_STOP ('UPDATE_CG2D [FORWARD_STEP]',myThid)
575 ENDIF
576 #endif
577
578 C-- Apply Filters to u*,v* before SOLVE_FOR_PRESSURE
579 #ifdef ALLOW_SHAP_FILT
580 IF (useSHAP_FILT .AND. shap_filt_uvStar) THEN
581 CALL TIMER_START('SHAP_FILT [FORWARD_STEP]',myThid)
582 IF (implicDiv2Dflow.LT.1.) THEN
583 C-- Explicit+Implicit part of the Barotropic Flow Divergence
584 C => Filtering of uVel,vVel is necessary
585 CALL SHAP_FILT_APPLY_UV( uVel,vVel,
586 & myTime, myIter, myThid )
587 ENDIF
588 CALL SHAP_FILT_APPLY_UV( gU,gV,myTime,myIter,myThid)
589 CALL TIMER_STOP ('SHAP_FILT [FORWARD_STEP]',myThid)
590 ENDIF
591 #endif
592 #ifdef ALLOW_ZONAL_FILT
593 IF (useZONAL_FILT .AND. zonal_filt_uvStar) THEN
594 CALL TIMER_START('ZONAL_FILT_APPLY [FORWARD_STEP]',myThid)
595 IF (implicDiv2Dflow.LT.1.) THEN
596 C-- Explicit+Implicit part of the Barotropic Flow Divergence
597 C => Filtering of uVel,vVel is necessary
598 CALL ZONAL_FILT_APPLY_UV( uVel, vVel, myThid )
599 ENDIF
600 CALL ZONAL_FILT_APPLY_UV( gU, gV, myThid )
601 CALL TIMER_STOP ('ZONAL_FILT_APPLY [FORWARD_STEP]',myThid)
602 ENDIF
603 #endif
604
605 C-- Solve elliptic equation(s).
606 C Two-dimensional only for conventional hydrostatic or
607 C three-dimensional for non-hydrostatic and/or IGW scheme.
608 #ifndef ALLOW_OFFLINE
609 IF ( momStepping ) THEN
610 CALL TIMER_START('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
611 CALL SOLVE_FOR_PRESSURE(myTime, myIter, myThid)
612 CALL TIMER_STOP ('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
613 ENDIF
614 #endif
615
616 C-- Correct divergence in flow field and cycle time-stepping momentum
617 c IF ( momStepping ) THEN
618 #ifndef ALLOW_OFFLINE
619 CALL TIMER_START('UV_CORRECTION_STEP [FORWARD_STEP]',myThid)
620 CALL MOMENTUM_CORRECTION_STEP(myTime, myIter, myThid)
621 CALL TIMER_STOP ('UV_CORRECTION_STEP [FORWARD_STEP]',myThid)
622 #endif
623 c ENDIF
624
625 #ifdef EXACT_CONSERV
626 IF (exactConserv) THEN
627 C-- Update etaH(n+1) :
628 CALL TIMER_START('UPDATE_ETAH [FORWARD_STEP]',mythid)
629 CALL UPDATE_ETAH( myTime, myIter, myThid )
630 CALL TIMER_STOP ('UPDATE_ETAH [FORWARD_STEP]',mythid)
631 ENDIF
632 #endif /* EXACT_CONSERV */
633
634 #ifdef NONLIN_FRSURF
635 IF ( select_rStar.NE.0 ) THEN
636 C-- r* : compute the future level thickness according to etaH(n+1)
637 CALL TIMER_START('CALC_R_STAR [FORWARD_STEP]',mythid)
638 CALL CALC_R_STAR(etaH, myTime, myIter, myThid )
639 CALL TIMER_STOP ('CALC_R_STAR [FORWARD_STEP]',mythid)
640 ELSEIF ( nonlinFreeSurf.GT.0) THEN
641 C-- compute the future surface level thickness according to etaH(n+1)
642 CALL TIMER_START('CALC_SURF_DR [FORWARD_STEP]',mythid)
643 CALL CALC_SURF_DR(etaH, myTime, myIter, myThid )
644 CALL TIMER_STOP ('CALC_SURF_DR [FORWARD_STEP]',mythid)
645 ENDIF
646 #endif /* NONLIN_FRSURF */
647
648 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
649 IF ( staggerTimeStep ) THEN
650 C-- do exchanges of U,V (needed for multiDim) when using stagger time-step :
651 #ifdef ALLOW_DEBUG
652 IF ( debugLevel .GE. debLevB )
653 & CALL DEBUG_CALL('DO_STAGGER_FIELDS_EXCH.',myThid)
654 #endif
655 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
656 CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )
657 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
658
659 #ifdef ALLOW_DEBUG
660 IF ( debugLevel .GE. debLevB )
661 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
662 #endif
663 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
664 CALL THERMODYNAMICS( myTime, myIter, myThid )
665 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
666
667 C-- if staggerTimeStep: end
668 ENDIF
669 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
670
671 #ifdef ALLOW_AUTODIFF_TAMC
672 cph This is needed because convective_adjustment calls
673 cph find_rho which may use pressure()
674 CADJ STORE totphihyd = comlev1, key = ikey_dynamics
675 #endif
676 C-- Cycle time-stepping Tracers arrays (T,S,+pTracers)
677 CALL TIMER_START('TS_CORRECTION_STEP [FORWARD_STEP]',myThid)
678 CALL TRACERS_CORRECTION_STEP(myTime, myIter, myThid)
679 CALL TIMER_STOP ('TS_CORRECTION_STEP [FORWARD_STEP]',myThid)
680
681 C-- Do "blocking" sends and receives for tendency "overlap" terms
682 c CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
683 c CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )
684 c CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
685
686 C-- Do "blocking" sends and receives for field "overlap" terms
687 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
688 CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
689 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
690
691 #ifdef ALLOW_FLT
692 C-- Calculate float trajectories
693 IF (useFLT) THEN
694 CALL TIMER_START('FLOATS [FORWARD_STEP]',myThid)
695 CALL FLT_MAIN(myIter,myTime, myThid)
696 CALL TIMER_STOP ('FLOATS [FORWARD_STEP]',myThid)
697 ENDIF
698 #endif
699
700 #ifdef ALLOW_AUTODIFF_TAMC
701 CALL AUTODIFF_INADMODE_SET( myThid )
702 #endif
703
704 C-- State-variables statistics (time-aver, diagnostics ...)
705 CALL TIMER_START('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
706 CALL DO_STATEVARS_DIAGS( myTime, myIter, myThid )
707 CALL TIMER_STOP ('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
708
709 #ifndef ALLOW_OFFLINE
710 #ifdef ALLOW_MONITOR
711 C-- Check status of solution (statistics, cfl, etc...)
712 CALL TIMER_START('MONITOR [FORWARD_STEP]',myThid)
713 CALL MONITOR( myIter, myTime, myThid )
714 CALL TIMER_STOP ('MONITOR [FORWARD_STEP]',myThid)
715 #endif /* ALLOW_MONITOR */
716 #endif
717
718 C-- Do IO if needed.
719 #ifdef ALLOW_OFFLINE
720 CALL TIMER_START('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
721 CALL OFFLINE_MODEL_IO( myTime, myIter, myThid )
722 CALL TIMER_STOP ('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
723 #else
724 CALL TIMER_START('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
725 CALL DO_THE_MODEL_IO( myTime, myIter, myThid )
726 CALL TIMER_STOP ('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
727 #endif
728
729 C-- Save state for restarts
730 CALL TIMER_START('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
731 CALL PACKAGES_WRITE_PICKUP(
732 I .FALSE., myTime, myIter, myThid )
733 #ifndef ALLOW_OFFLINE
734 CALL WRITE_CHECKPOINT(
735 I .FALSE., myTime, myIter, myThid )
736 #endif
737 CALL TIMER_STOP ('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
738
739 #ifdef ALLOW_AUTODIFF_TAMC
740 #ifdef ALLOW_TAMC_CHECKPOINTING
741 endif
742 enddo
743 endif
744 enddo
745 #ifndef AUTODIFF_2_LEVEL_CHECKPOINT
746 endif
747 enddo
748 #endif
749 #else
750 enddo
751 #endif
752
753 #else
754 enddo
755 #endif /* ALLOW_AUTODIFF_TAMC */
756
757 _BARRIER
758 call timer_stop ('ECCO MAIN LOOP', mythid)
759
760 call timer_start('ECCO SPIN-DOWN', mythid)
761
762 #ifdef ALLOW_COST
763
764 #ifdef ALLOW_DIVIDED_ADJOINT
765 CADJ STORE mytime = onetape
766 #endif
767 c-- Accumulate time averages of temperature, salinity, and SSH.
768 #ifndef DISABLE_DEBUGMODE
769 IF ( debugLevel .GE. debLevB )
770 & CALL DEBUG_CALL('cost_averagesfields',myThid)
771 #endif
772 call timer_start('cost_averagesfields [ECCO SPIN-DOWN]', mythid)
773 call cost_averagesfields( mytime, mythid )
774 call timer_stop ('cost_averagesfields [ECCO SPIN-DOWN]', mythid)
775 #ifdef ALLOW_DIVIDED_ADJOINT
776 c**************************************
777 #include "cost_averages_bar_directives.h"
778 c**************************************
779 #endif
780
781 #ifdef ALLOW_COST_ATLANTIC
782 c-- Compute meridional heat transport
783 #ifndef DISABLE_DEBUGMODE
784 IF ( debugLevel .GE. debLevB )
785 & CALL DEBUG_CALL('cost_atlantic',myThid)
786 #endif
787 call timer_start('cost_atlantic [ECCO SPIN-DOWN]', mythid)
788 call cost_atlantic( mytime, myiter,mythid )
789 call timer_stop ('cost_atlantic [ECCO SPIN-DOWN]', mythid)
790 #endif
791
792 c-- Compute the cost function contribution of the boundary forcing,
793 c-- i.e. heat flux, salt flux, zonal and meridional wind stress.
794 #ifndef DISABLE_DEBUGMODE
795 IF ( debugLevel .GE. debLevB )
796 & CALL DEBUG_CALL('cost_forcing',myThid)
797 #endif
798 call timer_start('cost_forcing [ECCO SPIN-DOWN]', mythid)
799 call cost_forcing( myiter, mytime, mythid )
800 call timer_stop ('cost_forcing [ECCO SPIN-DOWN]', mythid)
801
802 c-- Compute cost function contribution of Temperature and Salinity.
803 #ifndef DISABLE_DEBUGMODE
804 IF ( debugLevel .GE. debLevB )
805 & CALL DEBUG_CALL('cost_hyd',myThid)
806 #endif
807 call timer_start('cost_hyd [ECCO SPIN-DOWN]', mythid)
808 call cost_hyd( myiter, mytime, mythid )
809 call timer_stop ('cost_hyd [ECCO SPIN-DOWN]', mythid)
810
811 #ifdef ALLOW_OBCS_COST_CONTRIBUTION
812 #ifndef DISABLE_DEBUGMODE
813 IF ( debugLevel .GE. debLevB )
814 & CALL DEBUG_CALL('cost_obcs',myThid)
815 #endif
816 call timer_start('cost_obcs [ECCO SPIN-DOWN]', mythid)
817 call cost_obcs( myiter, mytime, mythid )
818 call timer_stop ('cost_obcs [ECCO SPIN-DOWN]', mythid)
819 #endif
820
821 #ifdef ALLOW_CURMTR_COST_CONTRIBUTION
822 #ifndef DISABLE_DEBUGMODE
823 IF ( debugLevel .GE. debLevB )
824 & CALL DEBUG_CALL('cost_curmtr',myThid)
825 #endif
826 call timer_start('cost_curmtr [ECCO SPIN-DOWN]', mythid)
827 call cost_curmtr( myiter, mytime, mythid )
828 call timer_stop ('cost_curmtr [ECCO SPIN-DOWN]', mythid)
829 #endif
830
831 c-- Compute cost function contribution of SSH.
832 #ifdef ALLOW_SSH_COST_CONTRIBUTION
833 #ifndef DISABLE_DEBUGMODE
834 IF ( debugLevel .GE. debLevB )
835 & CALL DEBUG_CALL('cost_ssh',myThid)
836 #endif
837 call timer_start('cost_ssh [ECCO SPIN-DOWN]', mythid)
838 call cost_ssh( myiter, mytime, mythid )
839 call timer_stop ('cost_ssh [ECCO SPIN-DOWN]', mythid)
840 #endif
841
842 c-- Compute cost function contribution of drifter's velocities.
843 #ifdef ALLOW_DRIFTER_COST_CONTRIBUTION
844 #ifndef DISABLE_DEBUGMODE
845 IF ( debugLevel .GE. debLevB )
846 & CALL DEBUG_CALL('cost_drifter',myThid)
847 #endif
848 call timer_start('cost_drifter [ECCO SPIN-DOWN]', mythid)
849 call cost_drifter( myiter, mytime, mythid )
850 call timer_stop ('cost_drifter [ECCO SPIN-DOWN]', mythid)
851 #endif
852
853 c-- Compute cost function contribution of wind stress observations.
854 #ifdef ALLOW_SCAT_COST_CONTRIBUTION
855 #ifndef DISABLE_DEBUGMODE
856 IF ( debugLevel .GE. debLevB )
857 & CALL DEBUG_CALL('cost_scat',myThid)
858 #endif
859 call timer_start('cost_scat [ECCO SPIN-DOWN]', mythid)
860 call cost_scat( myiter, mytime, mythid )
861 call timer_stop ('cost_scat [ECCO SPIN-DOWN]', mythid)
862 #endif
863
864 c-- Compute cost function contribution of wind stress observations.
865 #ifdef ALLOW_MEAN_HFLUX_COST_CONTRIBUTION
866 call timer_start('cost_mean_heatflux [ECCO SPIN-DOWN]', mythid)
867 call cost_mean_heatflux( myiter, mytime, mythid )
868 call timer_stop ('cost_mean_heatflux [ECCO SPIN-DOWN]', mythid)
869 #endif
870
871 c-- Compute cost function contribution of wind stress observations.
872 #ifdef ALLOW_MEAN_SFLUX_COST_CONTRIBUTION
873 call timer_start('cost_mean_saltflux [ECCO SPIN-DOWN]', mythid)
874 call cost_mean_saltflux( myiter, mytime, mythid )
875 call timer_stop ('cost_mean_saltflux [ECCO SPIN-DOWN]', mythid)
876 #endif
877
878 c-- Compute cost function contribution of drift between the first
879 c and the last year.
880 #ifdef ALLOW_DRIFT_COST_CONTRIBUTION
881 #ifndef DISABLE_DEBUGMODE
882 IF ( debugLevel .GE. debLevB )
883 & CALL DEBUG_CALL('cost_drift',myThid)
884 #endif
885 call timer_start('cost_drift [ECCO SPIN-DOWN]', mythid)
886 call cost_drift( myiter, mytime, mythid )
887 call timer_stop ('cost_drift [ECCO SPIN-DOWN]', mythid)
888 #endif
889 #ifdef ALLOW_DRIFTW_COST_CONTRIBUTION
890 #ifndef DISABLE_DEBUGMODE
891 IF ( debugLevel .GE. debLevB )
892 & CALL DEBUG_CALL('cost_driftw',myThid)
893 #endif
894 call timer_start('cost_driftw [ECCO SPIN-DOWN]', mythid)
895 call cost_driftw( myiter, mytime, mythid )
896 call timer_stop ('cost_driftw [ECCO SPIN-DOWN]', mythid)
897 #endif
898 _BARRIER
899
900 c-- Compute initial vs. final T/S deviation
901 #ifdef ALLOW_COST_INI_FIN
902 call timer_start('cost_ini_fin [ECCO SPIN-DOWN]', mythid)
903 call cost_theta_ini_fin( myiter, mytime, mythid )
904 call cost_salt_ini_fin( myiter, mytime, mythid )
905 call timer_stop ('cost_ini_fin [ECCO SPIN-DOWN]', mythid)
906 #endif
907 _BARRIER
908
909 c-- Sum all cost function contributions.
910 #ifndef DISABLE_DEBUGMODE
911 IF ( debugLevel .GE. debLevB )
912 & CALL DEBUG_CALL('cost_final',myThid)
913 #endif
914 call timer_start('COST_FINAL [ECCO SPIN-DOWN]', mythid)
915 call ecco_cost_final( mythid )
916 call timer_stop ('COST_FINAL [ECCO SPIN-DOWN]', mythid)
917
918 #endif /* ALLOW_COST */
919
920 call timer_stop ('ECCO SPIN-DOWN', mythid)
921
922 #ifndef DISABLE_DEBUGMODE
923 IF ( debugLevel .GE. debLevB )
924 & CALL DEBUG_LEAVE('THE_MAIN_LOOP',myThid)
925 #endif
926
927 return
928 end
929

  ViewVC Help
Powered by ViewVC 1.1.22