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

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

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

revision 1.103 by heimbach, Mon Sep 27 18:00:19 2004 UTC revision 1.120 by jmc, Sun Sep 11 18:17:20 2005 UTC
# Line 4  C $Name$ Line 4  C $Name$
4  #include "PACKAGES_CONFIG.h"  #include "PACKAGES_CONFIG.h"
5  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
6    
 #ifdef ALLOW_GCHEM  
 # include "GCHEM_OPTIONS.h"  
 #endif  
7  #ifdef ALLOW_OFFLINE  #ifdef ALLOW_OFFLINE
8  # include "OFFLINE_OPTIONS.h"  # include "OFFLINE_OPTIONS.h"
9  #endif  #endif
10    #ifdef ALLOW_GMREDI
11    # include "GMREDI_OPTIONS.h"
12    #endif
13    
14  CBOP  CBOP
15  C     !ROUTINE: FORWARD_STEP  C     !ROUTINE: FORWARD_STEP
# Line 86  C     == Global variables == Line 85  C     == Global variables ==
85  # ifdef EXACT_CONSERV  # ifdef EXACT_CONSERV
86  #  include "SURFACE.h"  #  include "SURFACE.h"
87  # endif  # endif
88    # ifdef ALLOW_KPP
89    #  include "KPP.h"
90    # endif
91    # ifdef ALLOW_GMREDI
92    #  include "GMREDI.h"
93    # endif
94  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
95    
96  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
# Line 114  CEOP Line 119  CEOP
119    
120  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
121  C--   Reset the model iteration counter and the model time.  C--   Reset the model iteration counter and the model time.
122        myiter = nIter0 + (iloop-1)        myIter = nIter0 + (iloop-1)
123        mytime = startTime + float(iloop-1)*deltaTclock        myTime = startTime + float(iloop-1)*deltaTclock
 #endif  
   
 #if (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_AUTODIFF_MONITOR))  
 C     Include call to a dummy routine. Its adjoint will be  
 C     called at the proper place in the adjoint code.  
 C     The adjoint routine will print out adjoint values  
 C     if requested. The location of the call is important,  
 C     it has to be after the adjoint of the exchanges  
 C     (DO_GTERM_BLOCKING_EXCHANGES).  
       CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )  
 cph   I've commented this line since it may conflict with MITgcm's adjoint  
 cph   However, need to check whether that's still consistent  
 cph   with the ecco-branch (it should).  
 cph      CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )  
124  #endif  #endif
125    
126  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 138  c************************************** Line 129  c**************************************
129  c**************************************  c**************************************
130  #endif  #endif
131    
132    C--   Switch on/off diagnostics for snap-shot output:
133    #ifdef ALLOW_DIAGNOSTICS
134          IF ( useDiagnostics ) THEN
135            CALL DIAGNOSTICS_SWITCH_ONOFF( myTime, myIter, myThid )
136          ENDIF
137    #endif
138    
139    C--   State-variables diagnostics
140          IF ( usediagnostics ) THEN
141            CALL TIMER_START('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)
142            CALL DO_STATEVARS_DIAGS( myTime, 0, myIter, myThid )
143            CALL TIMER_STOP ('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)
144          ENDIF
145    
146  C--   Call external forcing package  C--   Call external forcing package
147  #ifdef ALLOW_BULK_FORCE  #ifdef ALLOW_BULK_FORCE
148        IF ( useBulkForce ) THEN        IF ( useBulkForce ) THEN
# Line 192  c--   Add control vector for forcing and Line 197  c--   Add control vector for forcing and
197       &     CALL CTRL_MAP_FORCING (mythid)       &     CALL CTRL_MAP_FORCING (mythid)
198  #endif  #endif
199    
200    #if (defined (ALLOW_AUTODIFF_TAMC) && defined (ALLOW_AUTODIFF_MONITOR))
201    C     Include call to a dummy routine. Its adjoint will be
202    C     called at the proper place in the adjoint code.
203    C     The adjoint routine will print out adjoint values
204    C     if requested. The location of the call is important,
205    C     it has to be after the adjoint of the exchanges
206    C     (DO_GTERM_BLOCKING_EXCHANGES).
207          CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )
208    cph   I've commented this line since it may conflict with MITgcm's adjoint
209    cph   However, need to check whether that's still consistent
210    cph   with the ecco-branch (it should).
211    cph      CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
212    #endif
213    
214  # ifdef ALLOW_SEAICE  # ifdef ALLOW_SEAICE
215  C--   Call sea ice model to compute forcing/external data fields.  In  C--   Call sea ice model to compute forcing/external data fields.  In
216  C     addition to computing prognostic sea-ice variables and diagnosing the  C     addition to computing prognostic sea-ice variables and diagnosing the
# Line 227  CADJ STORE ptracer  = comlev1, key = ike Line 246  CADJ STORE ptracer  = comlev1, key = ike
246    
247  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
248  # ifdef ALLOW_GCHEM  # ifdef ALLOW_GCHEM
249            IF ( useGCHEM ) THEN
250    #ifdef ALLOW_DEBUG
251             IF ( debugLevel .GE. debLevB )
252         &        CALL DEBUG_CALL('GCHEM_FIELDS_LOAD',myThid)
253    #endif /* ALLOW_DEBUG */
254           CALL GCHEM_FIELDS_LOAD( mytime, myiter, mythid )           CALL GCHEM_FIELDS_LOAD( mytime, myiter, mythid )
255            ENDIF
256  # endif  # endif
257  #endif  #endif
258    
# Line 264  C--     Step forward fields and calculat Line 289  C--     Step forward fields and calculat
289         CALL DO_ATMOSPHERIC_PHYS( myTime, myIter, myThid )         CALL DO_ATMOSPHERIC_PHYS( myTime, myIter, myThid )
290         CALL TIMER_STOP ('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)         CALL TIMER_STOP ('DO_ATMOSPHERIC_PHYS [FORWARD_STEP]',mythid)
291    
292    #ifdef ALLOW_AUTODIFF_TAMC
293    CADJ STORE theta              = comlev1, key = ikey_dynamics
294    CADJ STORE salt               = comlev1, key = ikey_dynamics
295    CADJ STORE totphihyd          = comlev1, key = ikey_dynamics
296    CADJ STORE surfaceforcingtice = comlev1, key = ikey_dynamics
297    # ifdef ALLOW_KPP
298    CADJ STORE uvel               = comlev1, key = ikey_dynamics
299    CADJ STORE vvel               = comlev1, key = ikey_dynamics
300    # endif
301    # ifdef EXACT_CONSERV
302    CADJ STORE empmr              = comlev1, key = ikey_dynamics
303    CADJ STORE pmepr              = comlev1, key = ikey_dynamics
304    # endif
305    #endif /* ALLOW_AUTODIFF_TAMC */
306    
307  #ifndef ALLOW_OFFLINE  #ifndef ALLOW_OFFLINE
308  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
309         IF ( debugLevel .GE. debLevB )         IF ( debugLevel .GE. debLevB )
# Line 274  C--     Step forward fields and calculat Line 314  C--     Step forward fields and calculat
314         CALL TIMER_STOP ('DO_OCEANIC_PHYS     [FORWARD_STEP]',mythid)         CALL TIMER_STOP ('DO_OCEANIC_PHYS     [FORWARD_STEP]',mythid)
315  #endif  #endif
316    
317    #ifdef ALLOW_GCHEM
318    C     GCHEM package is an interface for any bio-geochemical or
319    C     ecosystem model you would like to include.
320    C     If GCHEM_SEPARATE_FORCING is not defined, you are
321    C     responsible for computing tendency terms for passive
322    C     tracers and storing them on a 3DxNumPtracers-array called
323    C     gchemTendency in GCHEM_CALC_TENDENCY. This tendency is then added
324    C     to gPtr in ptracers_forcing later-on.
325    C     If GCHEM_SEPARATE_FORCING is defined, you are reponsible for
326    C     UPDATING ptracers directly in GCHEM_FORCING_SEP. This amounts
327    C     to a completely separate time step that you have to implement
328    C     yourself (Eulerian seems to be fine in most cases).
329    CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
330    C     CAVEAT: Up to now, when GCHEM is turned on the field ptracerForcingSurf,
331    C     which is needed for KPP is not set properly. ptracerForcingSurf must
332    C     be treated differently depending on whether GCHEM_SEPARATE_FORCING
333    C     is define or not. TBD.
334    CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
335           IF ( useGCHEM ) THEN
336    #ifdef ALLOW_DEBUG
337            IF ( debugLevel .GE. debLevB )
338         &       CALL DEBUG_CALL('GCHEM_CALC_TENDENCY',myThid)
339    #endif
340            CALL TIMER_START('GCHEM_CALC_TENDENCY [FORWARD_STEP]',myThid)
341            CALL GCHEM_CALC_TENDENCY( myTime, myIter, myThid )
342            CALL TIMER_STOP ('GCHEM_CALC_TENDENCY [FORWARD_STEP]',myThid)
343           ENDIF
344    #endif /* ALLOW_GCHEM */
345    
346    #ifdef ALLOW_AUTODIFF_TAMC
347    cph needed to be moved here from do_oceanic_physics
348    cph to be visible down the road
349    c
350    CADJ STORE surfaceForcingS    = comlev1, key = ikey_dynamics
351    CADJ STORE surfaceForcingT    = comlev1, key = ikey_dynamics
352    CADJ STORE surfaceForcingTice = comlev1, key = ikey_dynamics
353    ctest(
354    CADJ STORE IVDConvCount       = comlev1, key = ikey_dynamics
355    ctest)
356    # ifdef ALLOW_PTRACERS
357    CADJ STORE surfaceForcingPtr  = comlev1, key = ikey_dynamics
358    # endif
359    c
360    # ifdef ALLOW_GMREDI
361    CADJ STORE Kwx                = comlev1, key = ikey_dynamics
362    CADJ STORE Kwy                = comlev1, key = ikey_dynamics
363    CADJ STORE Kwz                = comlev1, key = ikey_dynamics
364    #  ifdef GM_BOLUS_ADVEC
365    CADJ STORE GM_PsiX            = comlev1, key = ikey_dynamics
366    CADJ STORE GM_PsiY            = comlev1, key = ikey_dynamics
367    #  endif
368    # endif
369    c
370    # ifdef ALLOW_KPP
371    CADJ STORE KPPghat            = comlev1, key = ikey_dynamics
372    CADJ STORE KPPfrac            = comlev1, key = ikey_dynamics
373    CADJ STORE KPPdiffKzS         = comlev1, key = ikey_dynamics
374    CADJ STORE KPPdiffKzT         = comlev1, key = ikey_dynamics
375    # endif
376    #endif /* ALLOW_AUTODIFF_TAMC */
377    
378        IF ( .NOT.staggerTimeStep ) THEN        IF ( .NOT.staggerTimeStep ) THEN
379  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
380          IF ( debugLevel .GE. debLevB )          IF ( debugLevel .GE. debLevB )
# Line 319  C--   Step forward fields and calculate Line 420  C--   Step forward fields and calculate
420  #endif  #endif
421  #endif  #endif
422    
 #ifdef ALLOW_NONHYDROSTATIC  
 C--   Step forward W field in N-H algorithm  
       IF ( momStepping .AND. nonHydrostatic ) THEN  
 #ifdef ALLOW_DEBUG  
           IF ( debugLevel .GE. debLevB )  
      &     CALL DEBUG_CALL('CALC_GW',myThid)  
 #endif  
          CALL TIMER_START('CALC_GW          [FORWARD_STEP]',myThid)  
          CALL CALC_GW(myThid)  
          CALL TIMER_STOP ('CALC_GW          [FORWARD_STEP]',myThid)  
       ENDIF  
 #endif  
   
423  C--   Update time-counter  C--   Update time-counter
424        myIter = nIter0 + iLoop        myIter = nIter0 + iLoop
425        myTime = startTime + deltaTClock * float(iLoop)        myTime = startTime + deltaTClock * float(iLoop)
426    
427    #ifdef ALLOW_MNC
428    C     Update the default next iter for MNC
429          IF ( useMNC ) THEN
430             CALL MNC_CW_CITER_SETG( 1, 1, -1, myIter , myThid )
431    
432    C        TODO: Logic should be added here so that users can specify, on
433    C        a per-citer-group basis, when it is time to update the
434    C        "current" (and not just the "next") iteration
435           ENDIF
436    #endif
437    
438  C--   Update geometric factors:  C--   Update geometric factors:
439  #ifdef NONLIN_FRSURF  #ifdef NONLIN_FRSURF
440  C-    update hfacC,W,S and recip_hFac according to etaH(n+1) :  C-    update hfacC,W,S and recip_hFac according to etaH(n+1) :
# Line 439  C--   do exchanges of U,V (needed for mu Line 538  C--   do exchanges of U,V (needed for mu
538          CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )          CALL DO_STAGGER_FIELDS_EXCHANGES( myTime, myIter, myThid )
539          CALL TIMER_STOP ('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)          CALL TIMER_STOP ('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)
540    
541    C--   State-variables diagnostics
542           IF ( usediagnostics ) THEN
543            CALL TIMER_START('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)
544            CALL DO_STATEVARS_DIAGS( myTime, 1, myIter, myThid )
545            CALL TIMER_STOP ('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)
546           ENDIF
547    
548  #ifdef ALLOW_DEBUG  #ifdef ALLOW_DEBUG
549          IF ( debugLevel .GE. debLevB )          IF ( debugLevel .GE. debLevB )
550       &    CALL DEBUG_CALL('THERMODYNAMICS',myThid)       &    CALL DEBUG_CALL('THERMODYNAMICS',myThid)
# Line 461  C--   Cycle time-stepping Tracers arrays Line 567  C--   Cycle time-stepping Tracers arrays
567          CALL TRACERS_CORRECTION_STEP(myTime, myIter, myThid)          CALL TRACERS_CORRECTION_STEP(myTime, myIter, myThid)
568          CALL TIMER_STOP ('TS_CORRECTION_STEP  [FORWARD_STEP]',myThid)          CALL TIMER_STOP ('TS_CORRECTION_STEP  [FORWARD_STEP]',myThid)
569    
570    #ifdef ALLOW_GCHEM
571    C     Add separate timestepping of chemical/biological/forcing
572    C     of ptracers here in GCHEM_FORCING_SEP
573            IF ( useGCHEM ) THEN
574    #ifdef ALLOW_DEBUG
575             IF ( debugLevel .GE. debLevB )
576         &        CALL DEBUG_CALL('GCHEM_FORCING_SEP',myThid)
577    #endif /* ALLOW_DEBUG */
578             CALL TIMER_START('GCHEM_FORCING_SEP  [FORWARD_STEP]',myThid)
579             CALL GCHEM_FORCING_SEP( myTime,myIter,myThid )
580             CALL TIMER_STOP ('GCHEM_FORCING_SEP  [FORWARD_STEP]',myThid)
581            ENDIF  
582    #endif /* ALLOW_GCHEM */
583    
584  C--   Do "blocking" sends and receives for tendency "overlap" terms  C--   Do "blocking" sends and receives for tendency "overlap" terms
585  c     CALL TIMER_START('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)  c     CALL TIMER_START('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)
586  c     CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )  c     CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )
# Line 471  C--   Do "blocking" sends and receives f Line 591  C--   Do "blocking" sends and receives f
591        CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )        CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
592        CALL TIMER_STOP ('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)        CALL TIMER_STOP ('BLOCKING_EXCHANGES  [FORWARD_STEP]',myThid)
593    
 cswdptr -- add for seperate timestepping of chemical/biological/forcing  
 cswdptr    of ptracers ---  
 #ifdef ALLOW_GCHEM  
 ceh3 This is broken -- this ifdef should not be visible!  
 #ifdef PTRACERS_SEPARATE_FORCING  
 ceh3 needs an IF ( use GCHEM ) THEN  
         call GCHEM_FORCING_SEP( myTime,myIter,myThid )  
 #endif /* PTRACERS_SEPARATE_FORCING */  
 #endif /* ALLOW_GCHEM */  
 cswdptr -- end add ---  
594    
595  C AMM  C AMM
596  #ifdef ALLOW_GRIDALT  #ifdef ALLOW_GRIDALT
# Line 509  C--   Calculate float trajectories Line 619  C--   Calculate float trajectories
619        ENDIF        ENDIF
620  #endif  #endif
621    
622  C--   State-variables statistics (time-aver, diagnostics ...)  C--   State-variables time-averaging
623        CALL TIMER_START('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)        CALL TIMER_START('DO_STATEVARS_TAVE   [FORWARD_STEP]',myThid)
624        CALL DO_STATEVARS_DIAGS( myTime, myIter, myThid )        CALL DO_STATEVARS_TAVE( myTime, myIter, myThid )
625        CALL TIMER_STOP ('DO_STATEVARS_DIAGS  [FORWARD_STEP]',myThid)        CALL TIMER_STOP ('DO_STATEVARS_TAVE   [FORWARD_STEP]',myThid)
626    
627  #ifndef ALLOW_OFFLINE  #ifndef ALLOW_OFFLINE
628  #ifdef ALLOW_MONITOR  #ifdef ALLOW_MONITOR
# Line 526  C--   Check status of solution (statisti Line 636  C--   Check status of solution (statisti
636  #ifdef ALLOW_COST  #ifdef ALLOW_COST
637  C--     compare model with data and compute cost function  C--     compare model with data and compute cost function
638  C--     this is done after exchanges to allow interpolation  C--     this is done after exchanges to allow interpolation
 ceh3 perhaps needs an IF ( useCOST ) THEN  
 CADJ STORE theta, uvel, vvel = comlev1, key = ikey_dynamics  
639        CALL TIMER_START('COST_TILE           [FORWARD_STEP]',myThid)        CALL TIMER_START('COST_TILE           [FORWARD_STEP]',myThid)
640        CALL COST_TILE  ( mytime, myiter, myThid )        CALL COST_TILE  ( mytime, myiter, myThid )
641        CALL TIMER_STOP ('COST_TILE           [FORWARD_STEP]',myThid)        CALL TIMER_STOP ('COST_TILE           [FORWARD_STEP]',myThid)
Legend:
Removed from v.1.103  
changed lines
  Added in v.1.120
  ViewVC Help
Powered by ViewVC 1.1.22