/[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.19 - (show annotations) (download)
Wed Apr 6 18:37:49 2005 UTC (19 years, 2 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57f_post, checkpoint57k_post, checkpoint57g_post, checkpoint57i_post, checkpoint57m_post, checkpoint57g_pre, checkpoint57h_post, checkpoint57h_done, checkpoint57j_post, checkpoint57h_pre, checkpoint57l_post
Changes since 1.18: +2 -2 lines
use baseTime as time origin ; DIFF_BASE_MULTIPLE replaces DIFFERENT_MULTIPLE

1 C $Header: /u/gcmpack/MITgcm/pkg/ecco/the_main_loop.F,v 1.18 2005/03/29 20:10:53 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 = NINT( (startTime-baseTime)/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 # if (defined (ALLOW_AUTODIFF_MONITOR))
380 C Include call to a dummy routine. Its adjoint will be
381 C called at the proper place in the adjoint code.
382 C The adjoint routine will print out adjoint values
383 C if requested. The location of the call is important,
384 C it has to be after the adjoint of the exchanges
385 C (DO_GTERM_BLOCKING_EXCHANGES).
386 CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )
387 # endif
388 #endif
389
390
391 #ifdef ALLOW_COST
392
393 c-- Accumulate time averages of temperature, salinity, and SSH.
394 call timer_start('COST_AVERAGESFIELDS [ECCO MAIN]', mythid)
395 call cost_averagesFields( mytime, mythid )
396 call timer_stop ('COST_AVERAGESFIELDS [ECCO MAIN]', mythid)
397 #ifdef ALLOW_COST_ATLANTIC
398 c-- Compute meridional heat transport
399 call timer_start('cost_atlantic [ECCO MAIN]', mythid)
400 call cost_atlantic( mytime, myiter,mythid )
401 call timer_stop ('cost_atlantic [ECCO MAIN]', mythid)
402 #endif
403 #endif /* ALLOW_COST */
404
405 #ifdef ALLOW_AUTODIFF_TAMC
406 c**************************************
407 #include "checkpoint_lev1_directives.h"
408 c**************************************
409 #endif
410
411 #ifndef DISABLE_DEBUGMODE
412 IF ( debugLevel .GE. debLevB )
413 & CALL DEBUG_CALL('EXF_GETFORCING',myThid)
414 #endif
415 CALL TIMER_START('EXF_GETFORCING [FORWARD_STEP]',mythid)
416 CALL EXF_GETFORCING( mytime, myiter, mythid )
417 CALL TIMER_STOP ('EXF_GETFORCING [FORWARD_STEP]',mythid)
418
419 #ifdef ALLOW_SEAICE
420 cph this simple runtime flag causes a lot of recomp.
421 IF ( useSEAICE ) THEN
422 #ifndef DISABLE_DEBUGMODE
423 IF ( debugLevel .GE. debLevB )
424 & CALL DEBUG_CALL('SEAICE_MODEL',myThid)
425 #endif
426 #ifdef ALLOW_AUTODIFF_TAMC
427 CADJ STORE area = comlev1, key = ikey_dynamics
428 #endif
429 CALL TIMER_START('SEAICE_MODEL [FORWARD_STEP]',myThid)
430 CALL SEAICE_MODEL( myTime, myIter, myThid )
431 CALL TIMER_STOP ('SEAICE_MODEL [FORWARD_STEP]',myThid)
432 #ifdef ALLOW_COST_ICE
433 CALL COST_ICE ( myTime, myIter, myThid )
434 #endif
435 ENDIF
436 #endif /* ALLOW_SEAICE */
437
438 #ifdef ALLOW_AUTODIFF_TAMC
439 # ifdef ALLOW_PTRACERS
440 cph this replaces _bibj storing of ptracer within thermodynamics
441 CADJ STORE ptracer = comlev1, key = ikey_dynamics
442 # endif
443 #endif
444
445 #ifdef ALLOW_EBM
446 IF ( useEBM ) THEN
447 # ifdef ALLOW_DEBUG
448 IF ( debugLevel .GE. debLevB )
449 & CALL DEBUG_CALL('EBM',myThid)
450 # endif
451 CALL TIMER_START('EBM [FORWARD_STEP]',mythid)
452 CALL EBM_DRIVER ( myTime, myIter, myThid )
453 CALL TIMER_STOP ('EBM [FORWARD_STEP]',mythid)
454 ENDIF
455 #endif
456
457 C-- Step forward fields and calculate time tendency terms.
458
459 #ifdef ALLOW_DEBUG
460 IF ( debugLevel .GE. debLevB )
461 & CALL DEBUG_CALL('DO_ATMOSPHERIC_PHYS',myThid)
462 #endif
463 CALL TIMER_START('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
464 CALL DO_ATMOSPHERIC_PHYS( myTime, myIter, myThid )
465 CALL TIMER_STOP ('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
466
467 #ifndef ALLOW_OFFLINE
468 #ifdef ALLOW_DEBUG
469 IF ( debugLevel .GE. debLevB )
470 & CALL DEBUG_CALL('DO_OCEANIC_PHYS',myThid)
471 #endif
472 CALL TIMER_START('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
473 CALL DO_OCEANIC_PHYS( myTime, myIter, myThid )
474 CALL TIMER_STOP ('DO_OCEANIC_PHYS [FORWARD_STEP]',mythid)
475 #endif
476
477 #ifdef ALLOW_AUTODIFF_TAMC
478 cph needed to be moved here from do_oceanic_physics
479 cph to be visible down the road
480 c
481 CADJ STORE surfaceForcingS = comlev1, key = ikey_dynamics
482 CADJ STORE surfaceForcingT = comlev1, key = ikey_dynamics
483 CADJ STORE surfaceForcingTice = comlev1, key = ikey_dynamics
484 ctest(
485 CADJ STORE IVDConvCount = comlev1, key = ikey_dynamics
486 ctest)
487 # ifdef ALLOW_PTRACERS
488 CADJ STORE surfaceForcingPtr = comlev1, key = ikey_dynamics
489 # endif
490 c
491 # ifdef ALLOW_GMREDI
492 CADJ STORE Kwx = comlev1, key = ikey_dynamics
493 CADJ STORE Kwy = comlev1, key = ikey_dynamics
494 CADJ STORE Kwz = comlev1, key = ikey_dynamics
495 # ifdef GM_BOLUS_ADVEC
496 CADJ STORE GM_PsiX = comlev1, key = ikey_dynamics
497 CADJ STORE GM_PsiY = comlev1, key = ikey_dynamics
498 # endif
499 # endif
500 c
501 # ifdef ALLOW_KPP
502 CADJ STORE KPPghat = comlev1, key = ikey_dynamics
503 CADJ STORE KPPfrac = comlev1, key = ikey_dynamics
504 CADJ STORE KPPdiffKzS = comlev1, key = ikey_dynamics
505 CADJ STORE KPPdiffKzT = comlev1, key = ikey_dynamics
506 # endif
507 #endif /* ALLOW_AUTODIFF_TAMC */
508
509
510 IF ( .NOT.staggerTimeStep ) THEN
511 #ifdef ALLOW_DEBUG
512 IF ( debugLevel .GE. debLevB )
513 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
514 #endif
515 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
516 CALL THERMODYNAMICS( myTime, myIter, myThid )
517 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
518 C-- if not staggerTimeStep: end
519 ENDIF
520
521 C-- Step forward fields and calculate time tendency terms.
522 #ifndef ALLOW_OFFLINE
523 #ifndef ALLOW_AUTODIFF_TAMC
524 IF ( momStepping ) THEN
525 #endif
526 #ifdef ALLOW_DEBUG
527 IF ( debugLevel .GE. debLevB )
528 & CALL DEBUG_CALL('DYNAMICS',myThid)
529 #endif
530 CALL TIMER_START('DYNAMICS [FORWARD_STEP]',mythid)
531 CALL DYNAMICS( myTime, myIter, myThid )
532 CALL TIMER_STOP ('DYNAMICS [FORWARD_STEP]',mythid)
533 #ifndef ALLOW_AUTODIFF_TAMC
534 ENDIF
535 #endif
536 #endif
537
538 #ifdef ALLOW_NONHYDROSTATIC
539 C-- Step forward W field in N-H algorithm
540 IF ( momStepping .AND. nonHydrostatic ) THEN
541 #ifdef ALLOW_DEBUG
542 IF ( debugLevel .GE. debLevB )
543 & CALL DEBUG_CALL('CALC_GW',myThid)
544 #endif
545 CALL TIMER_START('CALC_GW [FORWARD_STEP]',myThid)
546 CALL CALC_GW(myThid)
547 CALL TIMER_STOP ('CALC_GW [FORWARD_STEP]',myThid)
548 ENDIF
549 #endif
550
551 C-- Update time-counter
552 myIter = nIter0 + iLoop
553 myTime = startTime + deltaTClock * float(iLoop)
554
555 C-- Update geometric factors:
556 #ifdef NONLIN_FRSURF
557 C- update hfacC,W,S and recip_hFac according to etaH(n+1) :
558 IF ( nonlinFreeSurf.GT.0) THEN
559 IF ( select_rStar.GT.0 ) THEN
560 CALL TIMER_START('UPDATE_R_STAR [FORWARD_STEP]',myThid)
561 CALL UPDATE_R_STAR( myTime, myIter, myThid )
562 CALL TIMER_STOP ('UPDATE_R_STAR [FORWARD_STEP]',myThid)
563 ELSE
564 CALL TIMER_START('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
565 CALL UPDATE_SURF_DR( myTime, myIter, myThid )
566 CALL TIMER_STOP ('UPDATE_SURF_DR [FORWARD_STEP]',myThid)
567 ENDIF
568 ENDIF
569 C- update also CG2D matrix (and preconditioner)
570 IF ( momStepping .AND. nonlinFreeSurf.GT.2 ) THEN
571 CALL TIMER_START('UPDATE_CG2D [FORWARD_STEP]',myThid)
572 CALL UPDATE_CG2D( myTime, myIter, myThid )
573 CALL TIMER_STOP ('UPDATE_CG2D [FORWARD_STEP]',myThid)
574 ENDIF
575 #endif
576
577 C-- Apply Filters to u*,v* before SOLVE_FOR_PRESSURE
578 #ifdef ALLOW_SHAP_FILT
579 IF (useSHAP_FILT .AND. shap_filt_uvStar) THEN
580 CALL TIMER_START('SHAP_FILT [FORWARD_STEP]',myThid)
581 IF (implicDiv2Dflow.LT.1.) THEN
582 C-- Explicit+Implicit part of the Barotropic Flow Divergence
583 C => Filtering of uVel,vVel is necessary
584 CALL SHAP_FILT_APPLY_UV( uVel,vVel,
585 & myTime, myIter, myThid )
586 ENDIF
587 CALL SHAP_FILT_APPLY_UV( gU,gV,myTime,myIter,myThid)
588 CALL TIMER_STOP ('SHAP_FILT [FORWARD_STEP]',myThid)
589 ENDIF
590 #endif
591 #ifdef ALLOW_ZONAL_FILT
592 IF (useZONAL_FILT .AND. zonal_filt_uvStar) THEN
593 CALL TIMER_START('ZONAL_FILT_APPLY [FORWARD_STEP]',myThid)
594 IF (implicDiv2Dflow.LT.1.) THEN
595 C-- Explicit+Implicit part of the Barotropic Flow Divergence
596 C => Filtering of uVel,vVel is necessary
597 CALL ZONAL_FILT_APPLY_UV( uVel, vVel, myThid )
598 ENDIF
599 CALL ZONAL_FILT_APPLY_UV( gU, gV, myThid )
600 CALL TIMER_STOP ('ZONAL_FILT_APPLY [FORWARD_STEP]',myThid)
601 ENDIF
602 #endif
603
604 C-- Solve elliptic equation(s).
605 C Two-dimensional only for conventional hydrostatic or
606 C three-dimensional for non-hydrostatic and/or IGW scheme.
607 #ifndef ALLOW_OFFLINE
608 IF ( momStepping ) THEN
609 CALL TIMER_START('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
610 CALL SOLVE_FOR_PRESSURE(myTime, myIter, myThid)
611 CALL TIMER_STOP ('SOLVE_FOR_PRESSURE [FORWARD_STEP]',myThid)
612 ENDIF
613 #endif
614
615 C-- Correct divergence in flow field and cycle time-stepping momentum
616 c IF ( momStepping ) THEN
617 #ifndef ALLOW_OFFLINE
618 CALL TIMER_START('UV_CORRECTION_STEP [FORWARD_STEP]',myThid)
619 CALL MOMENTUM_CORRECTION_STEP(myTime, myIter, myThid)
620 CALL TIMER_STOP ('UV_CORRECTION_STEP [FORWARD_STEP]',myThid)
621 #endif
622 c ENDIF
623
624 #ifdef EXACT_CONSERV
625 IF (exactConserv) THEN
626 C-- Update etaH(n+1) :
627 CALL TIMER_START('UPDATE_ETAH [FORWARD_STEP]',mythid)
628 CALL UPDATE_ETAH( myTime, myIter, myThid )
629 CALL TIMER_STOP ('UPDATE_ETAH [FORWARD_STEP]',mythid)
630 ENDIF
631 #endif /* EXACT_CONSERV */
632
633 #ifdef NONLIN_FRSURF
634 IF ( select_rStar.NE.0 ) THEN
635 C-- r* : compute the future level thickness according to etaH(n+1)
636 CALL TIMER_START('CALC_R_STAR [FORWARD_STEP]',mythid)
637 CALL CALC_R_STAR(etaH, myTime, myIter, myThid )
638 CALL TIMER_STOP ('CALC_R_STAR [FORWARD_STEP]',mythid)
639 ELSEIF ( nonlinFreeSurf.GT.0) THEN
640 C-- compute the future surface level thickness according to etaH(n+1)
641 CALL TIMER_START('CALC_SURF_DR [FORWARD_STEP]',mythid)
642 CALL CALC_SURF_DR(etaH, myTime, myIter, myThid )
643 CALL TIMER_STOP ('CALC_SURF_DR [FORWARD_STEP]',mythid)
644 ENDIF
645 #endif /* NONLIN_FRSURF */
646
647 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
648 IF ( staggerTimeStep ) THEN
649 C-- do exchanges of U,V (needed for multiDim) when using stagger time-step :
650 #ifdef ALLOW_DEBUG
651 IF ( debugLevel .GE. debLevB )
652 & CALL DEBUG_CALL('DO_STAGGER_FIELDS_EXCH.',myThid)
653 #endif
654 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
655 CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )
656 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
657
658 #ifdef ALLOW_DEBUG
659 IF ( debugLevel .GE. debLevB )
660 & CALL DEBUG_CALL('THERMODYNAMICS',myThid)
661 #endif
662 CALL TIMER_START('THERMODYNAMICS [FORWARD_STEP]',mythid)
663 CALL THERMODYNAMICS( myTime, myIter, myThid )
664 CALL TIMER_STOP ('THERMODYNAMICS [FORWARD_STEP]',mythid)
665
666 C-- if staggerTimeStep: end
667 ENDIF
668 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
669
670 #ifdef ALLOW_AUTODIFF_TAMC
671 cph This is needed because convective_adjustment calls
672 cph find_rho which may use pressure()
673 CADJ STORE totphihyd = comlev1, key = ikey_dynamics
674 #endif
675 C-- Cycle time-stepping Tracers arrays (T,S,+pTracers)
676 CALL TIMER_START('TS_CORRECTION_STEP [FORWARD_STEP]',myThid)
677 CALL TRACERS_CORRECTION_STEP(myTime, myIter, myThid)
678 CALL TIMER_STOP ('TS_CORRECTION_STEP [FORWARD_STEP]',myThid)
679
680 C-- Do "blocking" sends and receives for tendency "overlap" terms
681 c CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
682 c CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )
683 c CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
684
685 C-- Do "blocking" sends and receives for field "overlap" terms
686 CALL TIMER_START('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
687 CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
688 CALL TIMER_STOP ('BLOCKING_EXCHANGES [FORWARD_STEP]',myThid)
689
690 #ifdef ALLOW_FLT
691 C-- Calculate float trajectories
692 IF (useFLT) THEN
693 CALL TIMER_START('FLOATS [FORWARD_STEP]',myThid)
694 CALL FLT_MAIN(myIter,myTime, myThid)
695 CALL TIMER_STOP ('FLOATS [FORWARD_STEP]',myThid)
696 ENDIF
697 #endif
698
699 #ifdef ALLOW_AUTODIFF_TAMC
700 CALL AUTODIFF_INADMODE_SET( myThid )
701 #endif
702
703 C-- State-variables statistics (time-aver, diagnostics ...)
704 CALL TIMER_START('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
705 CALL DO_STATEVARS_DIAGS( myTime, myIter, myThid )
706 CALL TIMER_STOP ('DO_STATEVARS_DIAGS [FORWARD_STEP]',myThid)
707
708 #ifndef ALLOW_OFFLINE
709 #ifdef ALLOW_MONITOR
710 C-- Check status of solution (statistics, cfl, etc...)
711 CALL TIMER_START('MONITOR [FORWARD_STEP]',myThid)
712 CALL MONITOR( myIter, myTime, myThid )
713 CALL TIMER_STOP ('MONITOR [FORWARD_STEP]',myThid)
714 #endif /* ALLOW_MONITOR */
715 #endif
716
717 C-- Do IO if needed.
718 #ifdef ALLOW_OFFLINE
719 CALL TIMER_START('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
720 CALL OFFLINE_MODEL_IO( myTime, myIter, myThid )
721 CALL TIMER_STOP ('OFFLINE_MODEL_IO [FORWARD_STEP]',myThid)
722 #else
723 CALL TIMER_START('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
724 CALL DO_THE_MODEL_IO( myTime, myIter, myThid )
725 CALL TIMER_STOP ('DO_THE_MODEL_IO [FORWARD_STEP]',myThid)
726 #endif
727
728 C-- Save state for restarts
729 CALL TIMER_START('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
730 CALL PACKAGES_WRITE_PICKUP(
731 I .FALSE., myTime, myIter, myThid )
732 #ifndef ALLOW_OFFLINE
733 CALL WRITE_CHECKPOINT(
734 I .FALSE., myTime, myIter, myThid )
735 #endif
736 CALL TIMER_STOP ('WRITE_CHECKPOINT [FORWARD_STEP]',myThid)
737
738 #ifdef ALLOW_AUTODIFF_TAMC
739 #ifdef ALLOW_TAMC_CHECKPOINTING
740 endif
741 enddo
742 endif
743 enddo
744 #ifndef AUTODIFF_2_LEVEL_CHECKPOINT
745 endif
746 enddo
747 #endif
748 #else
749 enddo
750 #endif
751
752 #else
753 enddo
754 #endif /* ALLOW_AUTODIFF_TAMC */
755
756 _BARRIER
757 call timer_stop ('ECCO MAIN LOOP', mythid)
758
759 call timer_start('ECCO SPIN-DOWN', mythid)
760
761 #ifdef ALLOW_COST
762
763 #ifdef ALLOW_DIVIDED_ADJOINT
764 CADJ STORE mytime = onetape
765 #endif
766 c-- Accumulate time averages of temperature, salinity, and SSH.
767 #ifndef DISABLE_DEBUGMODE
768 IF ( debugLevel .GE. debLevB )
769 & CALL DEBUG_CALL('cost_averagesfields',myThid)
770 #endif
771 call timer_start('cost_averagesfields [ECCO SPIN-DOWN]', mythid)
772 call cost_averagesfields( mytime, mythid )
773 call timer_stop ('cost_averagesfields [ECCO SPIN-DOWN]', mythid)
774 #ifdef ALLOW_DIVIDED_ADJOINT
775 c**************************************
776 #include "cost_averages_bar_directives.h"
777 c**************************************
778 #endif
779
780 #ifdef ALLOW_COST_ATLANTIC
781 c-- Compute meridional heat transport
782 #ifndef DISABLE_DEBUGMODE
783 IF ( debugLevel .GE. debLevB )
784 & CALL DEBUG_CALL('cost_atlantic',myThid)
785 #endif
786 call timer_start('cost_atlantic [ECCO SPIN-DOWN]', mythid)
787 call cost_atlantic( mytime, myiter,mythid )
788 call timer_stop ('cost_atlantic [ECCO SPIN-DOWN]', mythid)
789 #endif
790
791 c-- Compute the cost function contribution of the boundary forcing,
792 c-- i.e. heat flux, salt flux, zonal and meridional wind stress.
793 #ifndef DISABLE_DEBUGMODE
794 IF ( debugLevel .GE. debLevB )
795 & CALL DEBUG_CALL('cost_forcing',myThid)
796 #endif
797 call timer_start('cost_forcing [ECCO SPIN-DOWN]', mythid)
798 call cost_forcing( myiter, mytime, mythid )
799 call timer_stop ('cost_forcing [ECCO SPIN-DOWN]', mythid)
800
801 c-- Compute cost function contribution of Temperature and Salinity.
802 #ifndef DISABLE_DEBUGMODE
803 IF ( debugLevel .GE. debLevB )
804 & CALL DEBUG_CALL('cost_hyd',myThid)
805 #endif
806 call timer_start('cost_hyd [ECCO SPIN-DOWN]', mythid)
807 call cost_hyd( myiter, mytime, mythid )
808 call timer_stop ('cost_hyd [ECCO SPIN-DOWN]', mythid)
809
810 #ifdef ALLOW_OBCS_COST_CONTRIBUTION
811 #ifndef DISABLE_DEBUGMODE
812 IF ( debugLevel .GE. debLevB )
813 & CALL DEBUG_CALL('cost_obcs',myThid)
814 #endif
815 call timer_start('cost_obcs [ECCO SPIN-DOWN]', mythid)
816 call cost_obcs( myiter, mytime, mythid )
817 call timer_stop ('cost_obcs [ECCO SPIN-DOWN]', mythid)
818 #endif
819
820 #ifdef ALLOW_CURMTR_COST_CONTRIBUTION
821 #ifndef DISABLE_DEBUGMODE
822 IF ( debugLevel .GE. debLevB )
823 & CALL DEBUG_CALL('cost_curmtr',myThid)
824 #endif
825 call timer_start('cost_curmtr [ECCO SPIN-DOWN]', mythid)
826 call cost_curmtr( myiter, mytime, mythid )
827 call timer_stop ('cost_curmtr [ECCO SPIN-DOWN]', mythid)
828 #endif
829
830 c-- Compute cost function contribution of SSH.
831 #ifdef ALLOW_SSH_COST_CONTRIBUTION
832 #ifndef DISABLE_DEBUGMODE
833 IF ( debugLevel .GE. debLevB )
834 & CALL DEBUG_CALL('cost_ssh',myThid)
835 #endif
836 call timer_start('cost_ssh [ECCO SPIN-DOWN]', mythid)
837 call cost_ssh( myiter, mytime, mythid )
838 call timer_stop ('cost_ssh [ECCO SPIN-DOWN]', mythid)
839 #endif
840
841 c-- Compute cost function contribution of drifter's velocities.
842 #ifdef ALLOW_DRIFTER_COST_CONTRIBUTION
843 #ifndef DISABLE_DEBUGMODE
844 IF ( debugLevel .GE. debLevB )
845 & CALL DEBUG_CALL('cost_drifter',myThid)
846 #endif
847 call timer_start('cost_drifter [ECCO SPIN-DOWN]', mythid)
848 call cost_drifter( myiter, mytime, mythid )
849 call timer_stop ('cost_drifter [ECCO SPIN-DOWN]', mythid)
850 #endif
851
852 c-- Compute cost function contribution of wind stress observations.
853 #ifdef ALLOW_SCAT_COST_CONTRIBUTION
854 #ifndef DISABLE_DEBUGMODE
855 IF ( debugLevel .GE. debLevB )
856 & CALL DEBUG_CALL('cost_scat',myThid)
857 #endif
858 call timer_start('cost_scat [ECCO SPIN-DOWN]', mythid)
859 call cost_scat( myiter, mytime, mythid )
860 call timer_stop ('cost_scat [ECCO SPIN-DOWN]', mythid)
861 #endif
862
863 c-- Compute cost function contribution of wind stress observations.
864 #ifdef ALLOW_MEAN_HFLUX_COST_CONTRIBUTION
865 call timer_start('cost_mean_heatflux [ECCO SPIN-DOWN]', mythid)
866 call cost_mean_heatflux( myiter, mytime, mythid )
867 call timer_stop ('cost_mean_heatflux [ECCO SPIN-DOWN]', mythid)
868 #endif
869
870 c-- Compute cost function contribution of wind stress observations.
871 #ifdef ALLOW_MEAN_SFLUX_COST_CONTRIBUTION
872 call timer_start('cost_mean_saltflux [ECCO SPIN-DOWN]', mythid)
873 call cost_mean_saltflux( myiter, mytime, mythid )
874 call timer_stop ('cost_mean_saltflux [ECCO SPIN-DOWN]', mythid)
875 #endif
876
877 c-- Compute cost function contribution of drift between the first
878 c and the last year.
879 #ifdef ALLOW_DRIFT_COST_CONTRIBUTION
880 #ifndef DISABLE_DEBUGMODE
881 IF ( debugLevel .GE. debLevB )
882 & CALL DEBUG_CALL('cost_drift',myThid)
883 #endif
884 call timer_start('cost_drift [ECCO SPIN-DOWN]', mythid)
885 call cost_drift( myiter, mytime, mythid )
886 call timer_stop ('cost_drift [ECCO SPIN-DOWN]', mythid)
887 #endif
888 #ifdef ALLOW_DRIFTW_COST_CONTRIBUTION
889 #ifndef DISABLE_DEBUGMODE
890 IF ( debugLevel .GE. debLevB )
891 & CALL DEBUG_CALL('cost_driftw',myThid)
892 #endif
893 call timer_start('cost_driftw [ECCO SPIN-DOWN]', mythid)
894 call cost_driftw( myiter, mytime, mythid )
895 call timer_stop ('cost_driftw [ECCO SPIN-DOWN]', mythid)
896 #endif
897 _BARRIER
898
899 c-- Compute initial vs. final T/S deviation
900 #ifdef ALLOW_COST_INI_FIN
901 call timer_start('cost_ini_fin [ECCO SPIN-DOWN]', mythid)
902 call cost_theta_ini_fin( myiter, mytime, mythid )
903 call cost_salt_ini_fin( myiter, mytime, mythid )
904 call timer_stop ('cost_ini_fin [ECCO SPIN-DOWN]', mythid)
905 #endif
906 _BARRIER
907
908 c-- Eddy stress penalty term
909 #ifdef ALLOW_COST_TAU_EDDY
910 call timer_start('cost_tau_eddy [ECCO SPIN-DOWN]', mythid)
911 call cost_tau_eddy( mythid )
912 call timer_stop ('cost_tau_eddy [ECCO SPIN-DOWN]', mythid)
913 #endif
914
915 c-- Sum all cost function contributions.
916 #ifndef DISABLE_DEBUGMODE
917 IF ( debugLevel .GE. debLevB )
918 & CALL DEBUG_CALL('cost_final',myThid)
919 #endif
920 call timer_start('COST_FINAL [ECCO SPIN-DOWN]', mythid)
921 call ecco_cost_final( mythid )
922 call timer_stop ('COST_FINAL [ECCO SPIN-DOWN]', mythid)
923
924 #endif /* ALLOW_COST */
925
926 call timer_stop ('ECCO SPIN-DOWN', mythid)
927
928 #ifndef DISABLE_DEBUGMODE
929 IF ( debugLevel .GE. debLevB )
930 & CALL DEBUG_LEAVE('THE_MAIN_LOOP',myThid)
931 #endif
932
933 return
934 end
935

  ViewVC Help
Powered by ViewVC 1.1.22