/[MITgcm]/MITgcm_contrib/AITCZ/code/forward_step.F
ViewVC logotype

Annotation of /MITgcm_contrib/AITCZ/code/forward_step.F

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


Revision 1.1 - (hide annotations) (download)
Wed Aug 20 15:24:59 2003 UTC (21 years, 11 months ago) by czaja
Branch: MAIN
CVS Tags: HEAD
Initial creation of Arnaud's simple coupled simulation.

1 czaja 1.1 C $Header: /u/u0/gcmpack/models/MITgcmUV/model/src/forward_step.F,v 1.23 2001/09/27 20:12:10 heimbach Exp $
2     C $Name: release1_beta1 $
3    
4     #include "CPP_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: FORWARD_STEP
8     C !INTERFACE:
9     SUBROUTINE FORWARD_STEP( iloop, myTime, myIter, myThid )
10    
11     C !DESCRIPTION: \bv
12     C *==================================================================
13     C | SUBROUTINE forward_step
14     C | o Run the ocean model and, optionally, evaluate a cost function.
15     C *==================================================================
16     C |
17     C | THE_MAIN_LOOP is the toplevel routine for the Tangent Linear and
18     C | Adjoint Model Compiler (TAMC). For this purpose the initialization
19     C | of the model was split into two parts. Those parameters that do
20     C | not depend on a specific model run are set in INITIALISE_FIXED,
21     C | whereas those that do depend on the specific realization are
22     C | initialized in INITIALISE_VARIA.
23     C |
24     C *==================================================================
25     C \ev
26    
27     C !USES:
28     IMPLICIT NONE
29     C == Global variables ==
30     #include "SIZE.h"
31     #include "EEPARAMS.h"
32     #include "PARAMS.h"
33     #include "DYNVARS.h"
34     #include "FFIELDS.h"
35    
36     #ifdef ALLOW_NONHYDROSTATIC
37     #include "CG3D.h"
38     #endif
39    
40     #ifdef ALLOW_SHAP_FILT
41     #include "SHAP_FILT.h"
42     #endif
43    
44    
45     #ifdef ALLOW_AUTODIFF_TAMC
46     #include "tamc.h"
47     #include "ctrl.h"
48     #include "ctrl_dummy.h"
49     #include "cost.h"
50     #endif
51    
52    
53     C !LOCAL VARIABLES:
54     C == Routine arguments ==
55     C note: under the multi-threaded model myiter and
56     C mytime are local variables passed around as routine
57     C arguments. Although this is fiddly it saves the need to
58     C impose additional synchronisation points when they are
59     C updated.
60     C myiter - iteration counter for this thread
61     C mytime - time counter for this thread
62     C mythid - thread number for this instance of the routine.
63     integer iloop
64     integer mythid
65     integer myiter
66     _RL mytime
67     INTEGER bi,bj
68    
69     CEOP
70    
71    
72     #ifdef ALLOW_AUTODIFF_TAMC
73     C-- Reset the model iteration counter and the model time.
74     myiter = nIter0 + (iloop-1)
75     mytime = startTime + float(iloop-1)*deltaTclock
76     #endif
77    
78     #ifdef ALLOW_AUTODIFF_TAMC
79     C Include call to a dummy routine. Its adjoint will be
80     C called at the proper place in the adjoint code.
81     C The adjoint routine will print out adjoint values
82     C if requested. The location of the call is important,
83     C it has to be after the adjoint of the exchanges
84     C (DO_GTERM_BLOCKING_EXCHANGES).
85     CALL DUMMY_IN_STEPPING( myTime, myIter, myThid )
86     #endif
87    
88     #ifdef EXACT_CONSERV
89     IF (exactConserv) THEN
90     C-- Update etaH(n+1) :
91     DO bj=myByLo(myThid),myByHi(myThid)
92     DO bi=myBxLo(myThid),myBxHi(myThid)
93     CALL CALC_EXACT_ETA( .FALSE., bi,bj, uVel,vVel,
94     I startTime, nIter0, myThid )
95     ENDDO
96     ENDDO
97     IF (implicDiv2Dflow .NE. 1. _d 0 )
98     & _EXCH_XY_R8(etaH, myThid )
99     ENDIF
100     #endif /* EXACT_CONSERV */
101    
102     #ifdef NONLIN_FRSURF
103     C-- compute the future surface level thickness
104     C according to etaH(n+1)
105     IF ( nonlinFreeSurf.GT.0) THEN
106     CALL CALC_SURF_DR(etaH, myTime, myIter, myThid )
107     ENDIF
108     #endif /* NONLIN_FRSURF */
109    
110     C-- Load forcing/external data fields.
111     #ifdef INCLUDE_EXTERNAL_FORCING_PACKAGE
112     C NOTE, that although the exf package is part of the
113     C distribution, it is not currently maintained, i.e.
114     C exf is disabled by default in genmake.
115     #ifdef ALLOW_BULKFORMULAE
116     CADJ STORE theta = comlev1, key = ikey_dynamics
117     #endif
118     CALL EXF_GETFORCING( mytime, myiter, mythid )
119     #else
120     CALL TIMER_START('EXTERNAL_FIELDS_LOAD[THE_MAIN_LOOP]',mythid)
121     CALL EXTERNAL_FIELDS_LOAD( mytime, myiter, mythid )
122     CALL TIMER_STOP ('EXTERNAL_FIELDS_LOAD[THE_MAIN_LOOP]',mythid)
123     #endif
124    
125     C-- Step forward fields and calculate time tendency terms.
126     CALL TIMER_START('THERMODYNAMICS [THE_MAIN_LOOP]',mythid)
127     CALL THERMODYNAMICS( myTime, myIter, myThid )
128     CALL TIMER_STOP ('THERMODYNAMICS [THE_MAIN_LOOP]',mythid)
129    
130     #ifdef ALLOW_SHAP_FILT
131     IF (useSHAP_FILT .AND.
132     & staggerTimeStep .AND. shap_filt_TrStagg ) THEN
133     CALL TIMER_START('SHAP_FILT [THE_MAIN_LOOP]',myThid)
134     CALL SHAP_FILT_APPLY_TS( gT, gS, myTime, myIter, myThid )
135     CALL TIMER_STOP ('SHAP_FILT [THE_MAIN_LOOP]',myThid)
136     ENDIF
137     #endif
138    
139    
140     C-- Step forward fields and calculate time tendency terms.
141     IF ( momStepping ) THEN
142     CALL TIMER_START('DYNAMICS [THE_MAIN_LOOP]',mythid)
143     CALL DYNAMICS( myTime, myIter, myThid )
144     CALL TIMER_STOP ('DYNAMICS [THE_MAIN_LOOP]',mythid)
145     ENDIF
146    
147     #ifndef ALLOW_AUTODIFF_TAMC
148     #ifdef ALLOW_NONHYDROSTATIC
149     C-- Step forward W field in N-H algorithm
150     IF ( momStepping .AND. nonHydrostatic ) THEN
151     CALL TIMER_START('CALC_GW [THE_MAIN_LOOP]',myThid)
152     CALL CALC_GW(myThid)
153     CALL TIMER_STOP ('CALC_GW [THE_MAIN_LOOP]',myThid)
154     ENDIF
155     #endif
156     #endif /* ALLOW_AUTODIFF_TAMC */
157    
158     #ifdef NONLIN_FRSURF
159     C-- update hfacC,W,S and recip_hFac according to etaH(n+1) :
160     IF ( momStepping ) THEN
161     IF ( nonlinFreeSurf.GT.0) THEN
162     CALL UPDATE_SURF_DR( myTime, myIter, myThid )
163     ENDIF
164     C- update also CG2D matrix (and preconditioner)
165     IF ( nonlinFreeSurf.GT.2) THEN
166     CALL UPDATE_CG2D( myTime, myIter, myThid )
167     ENDIF
168     ENDIF
169     #endif
170    
171     C-- This block has been moved to the_correction_step(). We have
172     C left this code here to indicate where the filters used to be
173     C in the algorithm before JMC moved them to after the pressure
174     C solver.
175     C
176     C modif omp: zonal filter before elliptic solver
177     #ifdef ALLOW_SHAP_FILT
178     IF (useSHAP_FILT .AND. shap_filt_uvStar) THEN
179     CALL TIMER_START('SHAP_FILT [THE_MAIN_LOOP]',myThid)
180     CALL SHAP_FILT_APPLY_UV( gUnm1,gVnm1, myTime,myIter,myThid )
181     IF (implicDiv2Dflow.LT.1.) THEN
182     C-- Explicit+Implicit part of the Barotropic Flow Divergence
183     C => Filtering of uVel,vVel is necessary
184     CALL SHAP_FILT_APPLY_UV( uVel,vVel, myTime,myIter,myThid )
185     ENDIF
186     CALL TIMER_STOP ('SHAP_FILT [THE_MAIN_LOOP]',myThid)
187     ENDIF
188     #endif
189    
190     C
191     C modif omp: zonal filter before elliptic solver
192     #ifdef ALLOW_ZONAL_FILT
193     IF (zonal_filt_lat.LT.90.) THEN
194     CALL TIMER_START('ZONAL_FILT_APPLY [THE_MAIN_LOOP]',myThid)
195     CALL ZONAL_FILT_APPLY(
196     U gUnm1, gVnm1, gT, gS,
197     I myThid )
198     CALL TIMER_STOP ('ZONAL_FILT_APPLY [THE_MAIN_LOOP]',myThid)
199     ENDIF
200     #endif
201    
202     C-- Solve elliptic equation(s).
203     C Two-dimensional only for conventional hydrostatic or
204     C three-dimensional for non-hydrostatic and/or IGW scheme.
205     IF ( momStepping ) THEN
206     CALL TIMER_START('SOLVE_FOR_PRESSURE [THE_MAIN_LOOP]',myThid)
207     CALL SOLVE_FOR_PRESSURE( myThid )
208     CALL TIMER_STOP ('SOLVE_FOR_PRESSURE [THE_MAIN_LOOP]',myThid)
209     ENDIF
210    
211     C-- Correct divergence in flow field and cycle time-stepping
212     C arrays (for all fields) ; update time-counter
213     myIter = nIter0 + iLoop
214     myTime = startTime + deltaTClock * float(iLoop)
215     CALL TIMER_START('THE_CORRECTION_STEP [THE_MAIN_LOOP]',myThid)
216     CALL THE_CORRECTION_STEP(myTime, myIter, myThid)
217     CALL TIMER_STOP ('THE_CORRECTION_STEP [THE_MAIN_LOOP]',myThid)
218    
219     C-- Do "blocking" sends and receives for tendency "overlap" terms
220     c CALL TIMER_START('BLOCKING_EXCHANGES [THE_MAIN_LOOP]',myThid)
221     c CALL DO_GTERM_BLOCKING_EXCHANGES( myThid )
222     c CALL TIMER_STOP ('BLOCKING_EXCHANGES [THE_MAIN_LOOP]',myThid)
223    
224     C-- Do "blocking" sends and receives for field "overlap" terms
225     CALL TIMER_START('BLOCKING_EXCHANGES [THE_MAIN_LOOP]',myThid)
226     CALL DO_FIELDS_BLOCKING_EXCHANGES( myThid )
227     CALL TIMER_STOP ('BLOCKING_EXCHANGES [THE_MAIN_LOOP]',myThid)
228    
229     #ifdef ALLOW_FLT
230     C-- Calculate float trajectories
231     IF (useFLT) THEN
232     CALL TIMER_START('FLOATS [THE_MAIN_LOOP]',myThid)
233     CALL FLT_MAIN(myIter,myTime, myThid)
234     CALL TIMER_STOP ('FLOATS [THE_MAIN_LOOP]',myThid)
235     ENDIF
236     #endif
237    
238     #ifndef EXCLUDE_MONITOR
239     C-- Check status of solution (statistics, cfl, etc...)
240     CALL MONITOR( myIter, myTime, myThid )
241     #endif /* EXCLUDE_MONITOR */
242    
243     #ifndef ALLOW_AUTODIFF_TAMC
244     C-- Do IO if needed.
245     CALL TIMER_START('DO_THE_MODEL_IO [THE_MAIN_LOOP]',myThid)
246     CALL DO_THE_MODEL_IO( myTime, myIter, myThid )
247     CALL TIMER_STOP ('DO_THE_MODEL_IO [THE_MAIN_LOOP]',myThid)
248    
249     C-- Save state for restarts
250     C Note: (jmc: is it still the case after ckp35 ?)
251     C =====
252     C Because of the ordering of the timestepping code and
253     C tendency term code at end of loop model arrays hold
254     C U,V,T,S at "time-level" N but gu, gv, gs, gt, guNM1,...
255     C at "time-level" N+1/2 (guNM1 at "time-level" N+1/2 is
256     C gu at "time-level" N-1/2) and etaN at "time-level" N+1/2.
257     C where N = I+timeLevBase-1
258     C Thus a checkpoint contains U.0000000000, GU.0000000001 and
259     C etaN.0000000001 in the indexing scheme used for the model
260     C "state" files. This example is referred to as a checkpoint
261     C at time level 1
262     CALL TIMER_START('WRITE_CHECKPOINT [THE_MAIN_LOOP]',myThid)
263     CALL WRITE_CHECKPOINT(
264     & .FALSE., myTime, myIter, myThid )
265     CALL TIMER_STOP ('WRITE_CHECKPOINT [THE_MAIN_LOOP]',myThid)
266    
267     #endif /* ALLOW_AUTODIFF_TAMC */
268    
269     END

  ViewVC Help
Powered by ViewVC 1.1.22