--- MITgcm/model/src/dynamics.F 2001/02/20 15:06:21 1.63 +++ MITgcm/model/src/dynamics.F 2001/07/20 19:16:28 1.73 @@ -1,4 +1,4 @@ -C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/model/src/dynamics.F,v 1.63 2001/02/20 15:06:21 jmc Exp $ +C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/model/src/dynamics.F,v 1.73 2001/07/20 19:16:28 adcroft Exp $ C $Name: $ #include "CPP_OPTIONS.h" @@ -26,22 +26,25 @@ C == Global variables === #include "SIZE.h" #include "EEPARAMS.h" -#include "CG2D.h" #include "PARAMS.h" #include "DYNVARS.h" #include "GRID.h" +#include "TR1.h" #ifdef ALLOW_AUTODIFF_TAMC # include "tamc.h" # include "tamc_keys.h" +# include "FFIELDS.h" +# ifdef ALLOW_KPP +# include "KPP.h" +# endif +# ifdef ALLOW_GMREDI +# include "GMREDI.h" +# endif #endif /* ALLOW_AUTODIFF_TAMC */ -#ifdef ALLOW_KPP -# include "KPP.h" -#endif - -#ifdef INCLUDE_DIAGNOSTICS_INTERFACE_CODE -#include "AVER.h" +#ifdef ALLOW_TIMEAVE +#include "TIMEAVE_STATV.h" #endif C == Routine arguments == @@ -59,8 +62,7 @@ C o uTrans: Zonal transport C o vTrans: Meridional transport C o rTrans: Vertical transport -C maskC,maskUp o maskC: land/water mask for tracer cells -C o maskUp: land/water mask for W points +C maskUp o maskUp: land/water mask for W points C fVer[STUV] o fVer: Vertical flux term - note fVer C is "pipelined" in the vertical C so we need an fVer for each @@ -68,10 +70,9 @@ C rhoK, rhoKM1 - Density at current level, and level above C phiHyd - Hydrostatic part of the potential phiHydi. C In z coords phiHydiHyd is the hydrostatic -C pressure anomaly +C Potential (=pressure/rho0) anomaly C In p coords phiHydiHyd is the geopotential -C surface height -C anomaly. +C surface height anomaly. C phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean) C phiSurfY or geopotentiel (atmos) in X and Y direction C KappaRT, - Total diffusion in vertical for T and S. @@ -82,15 +83,16 @@ C k, kup, - Index for layer above and below. kup and kDown C kDown, km1 are switched with layer to be the appropriate C index into fVerTerm. +C tauAB - Adams-Bashforth timestepping weight: 0=forward ; 1/2=Adams-Bashf. _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RL uTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RL vTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RL rTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy) - _RS maskC (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RS maskUp (1-OLx:sNx+OLx,1-OLy:sNy+OLy) _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2) _RL fVerS (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2) + _RL fVerTr1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2) _RL fVerU (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2) _RL fVerV (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2) _RL phiHyd (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) @@ -105,6 +107,7 @@ _RL sigmaX (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) _RL sigmaY (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) + _RL tauAB C This is currently used by IVDC and Diagnostics _RL ConvectCount (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) @@ -121,16 +124,6 @@ c EXTERNAL DIFFERENT_MULTIPLE Cjmc(end) -#ifdef ALLOW_AUTODIFF_TAMC - INTEGER isbyte - PARAMETER( isbyte = 4 ) - - INTEGER act1, act2, act3, act4 - INTEGER max1, max2, max3 - INTEGER iikey, kkey - INTEGER maximpl -#endif /* ALLOW_AUTODIFF_TAMC */ - C--- The algorithm... C C "Correction Step" @@ -201,7 +194,6 @@ ENDDO rhoKM1 (i,j) = 0. _d 0 rhok (i,j) = 0. _d 0 - maskC (i,j) = 0. _d 0 phiSurfX(i,j) = 0. _d 0 phiSurfY(i,j) = 0. _d 0 ENDDO @@ -218,7 +210,7 @@ #ifdef ALLOW_AUTODIFF_TAMC C-- HPF directive to help TAMC CHPF$ INDEPENDENT, NEW (rTrans,fVerT,fVerS,fVerU,fVerV -CHPF$& ,phiHyd,utrans,vtrans,maskc,xA,yA +CHPF$& ,phiHyd,utrans,vtrans,xA,yA CHPF$& ,KappaRT,KappaRS,KappaRU,KappaRV CHPF$& ) #endif /* ALLOW_AUTODIFF_TAMC */ @@ -245,15 +237,17 @@ C-- Set up work arrays that need valid initial values DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx - rTrans(i,j) = 0. _d 0 - fVerT (i,j,1) = 0. _d 0 - fVerT (i,j,2) = 0. _d 0 - fVerS (i,j,1) = 0. _d 0 - fVerS (i,j,2) = 0. _d 0 - fVerU (i,j,1) = 0. _d 0 - fVerU (i,j,2) = 0. _d 0 - fVerV (i,j,1) = 0. _d 0 - fVerV (i,j,2) = 0. _d 0 + rTrans (i,j) = 0. _d 0 + fVerT (i,j,1) = 0. _d 0 + fVerT (i,j,2) = 0. _d 0 + fVerS (i,j,1) = 0. _d 0 + fVerS (i,j,2) = 0. _d 0 + fVerTr1(i,j,1) = 0. _d 0 + fVerTr1(i,j,2) = 0. _d 0 + fVerU (i,j,1) = 0. _d 0 + fVerU (i,j,2) = 0. _d 0 + fVerV (i,j,1) = 0. _d 0 + fVerV (i,j,2) = 0. _d 0 ENDDO ENDDO @@ -274,13 +268,25 @@ jMax = sNy+OLy +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE theta(:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +CADJ STORE tr1 (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + C-- Start of diagnostic loop DO k=Nr,1,-1 #ifdef ALLOW_AUTODIFF_TAMC C? Patrick, is this formula correct now that we change the loop range? C? Do we still need this? - kkey = (ikey-1)*(Nr-2+1) + (k-2) + 1 +cph kkey formula corrected. +cph Needed for rhok, rhokm1, in the case useGMREDI. + kkey = (ikey-1)*Nr + k +CADJ STORE rhokm1(:,:) = comlev1_bibj_k , key=kkey, byte=isbyte +CADJ STORE rhok (:,:) = comlev1_bibj_k , key=kkey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ C-- Integrate continuity vertically for vertical velocity @@ -302,16 +308,26 @@ C slope terms (e.g. GM/Redi tensor or IVDC diffusivity) c IF ( k.GT.1 .AND. (useGMRedi.OR.ivdc_kappa.NE.0.) ) THEN IF ( useGMRedi .OR. (k.GT.1 .AND. ivdc_kappa.NE.0.) ) THEN +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE theta(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte +CADJ STORE salt (:,:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ CALL FIND_RHO( I bi, bj, iMin, iMax, jMin, jMax, k, k, eosType, I theta, salt, O rhoK, I myThid ) - IF (k.GT.1) CALL FIND_RHO( + IF (k.GT.1) THEN +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE theta(:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte +CADJ STORE salt (:,:,k-1,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + CALL FIND_RHO( I bi, bj, iMin, iMax, jMin, jMax, k-1, k, eosType, I theta, salt, O rhoKm1, I myThid ) + ENDIF CALL GRAD_SIGMA( I bi, bj, iMin, iMax, jMin, jMax, k, I rhoK, rhoKm1, rhoK, @@ -332,6 +348,11 @@ C-- end of diagnostic k loop (Nr:1) ENDDO +#ifdef ALLOW_AUTODIFF_TAMC +cph avoids recomputation of integrate_for_w +CADJ STORE wvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + #ifdef ALLOW_OBCS C-- Calculate future values on open boundaries IF (useOBCS) THEN @@ -346,8 +367,25 @@ CALL EXTERNAL_FORCING_SURF( I bi, bj, iMin, iMax, jMin, jMax, I myThid ) +#ifdef ALLOW_AUTODIFF_TAMC +cph needed for KPP +CADJ STORE surfacetendencyU(:,:,bi,bj) +CADJ & = comlev1_bibj, key=ikey, byte=isbyte +CADJ STORE surfacetendencyV(:,:,bi,bj) +CADJ & = comlev1_bibj, key=ikey, byte=isbyte +CADJ STORE surfacetendencyS(:,:,bi,bj) +CADJ & = comlev1_bibj, key=ikey, byte=isbyte +CADJ STORE surfacetendencyT(:,:,bi,bj) +CADJ & = comlev1_bibj, key=ikey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ #ifdef ALLOW_GMREDI + +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE sigmaX(:,:,:) = comlev1, key=ikey, byte=isbyte +CADJ STORE sigmaY(:,:,:) = comlev1, key=ikey, byte=isbyte +CADJ STORE sigmaR(:,:,:) = comlev1, key=ikey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ C-- Calculate iso-neutral slopes for the GM/Redi parameterisation IF (useGMRedi) THEN DO k=1,Nr @@ -366,6 +404,13 @@ ENDDO #endif /* ALLOW_AUTODIFF_TAMC */ ENDIF + +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE Kwx(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte +CADJ STORE Kwy(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte +CADJ STORE Kwz(:,:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + #endif /* ALLOW_GMREDI */ #ifdef ALLOW_KPP @@ -373,7 +418,22 @@ IF (useKPP) THEN CALL KPP_CALC( I bi, bj, myTime, myThid ) +#ifdef ALLOW_AUTODIFF_TAMC + ELSE + CALL KPP_CALC_DUMMY( + I bi, bj, myTime, myThid ) +#endif /* ALLOW_AUTODIFF_TAMC */ ENDIF + +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE KPPghat (:,:,:,bi,bj) +CADJ & , KPPviscAz (:,:,:,bi,bj) +CADJ & , KPPdiffKzT(:,:,:,bi,bj) +CADJ & , KPPdiffKzS(:,:,:,bi,bj) +CADJ & , KPPfrac (:,: ,bi,bj) +CADJ & = comlev1_bibj, key=ikey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + #endif /* ALLOW_KPP */ #ifdef ALLOW_AUTODIFF_TAMC @@ -383,6 +443,7 @@ CADJ STORE salt (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte CADJ STORE uvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte CADJ STORE vvel (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte +CADJ STORE tr1 (:,:,:,bi,bj) = comlev1_bibj, key = ikey, byte = isbyte #endif /* ALLOW_AUTODIFF_TAMC */ #ifdef ALLOW_AIM @@ -390,7 +451,7 @@ C note(jmc) : phiHyd=0 at this point but is not really used in Molteni Physics IF ( useAIM ) THEN CALL TIMER_START('AIM_DO_ATMOS_PHYS [DYNAMICS]', myThid) - CALL AIM_DO_ATMOS_PHYSICS( phiHyd, myTime, myThid ) + CALL AIM_DO_ATMOS_PHYSICS( phiHyd, bi, bj, myTime, myThid ) CALL TIMER_STOP ('AIM_DO_ATMOS_PHYS [DYNAMICS]', myThid) ENDIF #endif /* ALLOW_AIM */ @@ -398,6 +459,12 @@ C-- Start of thermodynamics loop DO k=Nr,1,-1 +#ifdef ALLOW_AUTODIFF_TAMC +C? Patrick Is this formula correct? +cph Yes, but I rewrote it. +cph Also, the KappaR? need the index and subscript k! + kkey = (ikey-1)*Nr + k +#endif /* ALLOW_AUTODIFF_TAMC */ C-- km1 Points to level above k (=k-1) C-- kup Cycles through 1,2 to point to layer above @@ -412,25 +479,22 @@ jMin = 1-OLy+2 jMax = sNy+OLy-1 -#ifdef ALLOW_AUTODIFF_TAMC -CPatrick Is this formula correct? - kkey = (ikey-1)*(Nr-1+1) + (k-1) + 1 -CADJ STORE rTrans(:,:) = comlev1_bibj_k, key = kkey, byte = isbyte -CADJ STORE KappaRT(:,:,:) = comlev1_bibj_k, key = kkey, byte = isbyte -CADJ STORE KappaRS(:,:,:) = comlev1_bibj_k, key = kkey, byte = isbyte -#endif /* ALLOW_AUTODIFF_TAMC */ - C-- Get temporary terms used by tendency routines CALL CALC_COMMON_FACTORS ( - I bi,bj,iMin,iMax,jMin,jMax,k,km1,kup,kDown, - O xA,yA,uTrans,vTrans,rTrans,maskC,maskUp, + I bi,bj,iMin,iMax,jMin,jMax,k, + O xA,yA,uTrans,vTrans,rTrans,maskUp, I myThid) +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE KappaRT(:,:,k) = comlev1_bibj_k, key=kkey, byte=isbyte +CADJ STORE KappaRS(:,:,k) = comlev1_bibj_k, key=kkey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + #ifdef INCLUDE_CALC_DIFFUSIVITY_CALL C-- Calculate the total vertical diffusivity CALL CALC_DIFFUSIVITY( I bi,bj,iMin,iMax,jMin,jMax,k, - I maskC,maskup, + I maskUp, O KappaRT,KappaRS,KappaRU,KappaRV, I myThid) #endif @@ -440,12 +504,13 @@ IF ( tempStepping ) THEN CALL CALC_GT( I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown, - I xA,yA,uTrans,vTrans,rTrans,maskUp,maskC, + I xA,yA,uTrans,vTrans,rTrans,maskUp, I KappaRT, U fVerT, I myTime, myThid) + tauAB = 0.5d0 + abEps CALL TIMESTEP_TRACER( - I bi,bj,iMin,iMax,jMin,jMax,k, + I bi,bj,iMin,iMax,jMin,jMax,k,tauAB, I theta, gT, U gTnm1, I myIter, myThid) @@ -453,16 +518,31 @@ IF ( saltStepping ) THEN CALL CALC_GS( I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown, - I xA,yA,uTrans,vTrans,rTrans,maskUp,maskC, + I xA,yA,uTrans,vTrans,rTrans,maskUp, I KappaRS, U fVerS, I myTime, myThid) + tauAB = 0.5d0 + abEps CALL TIMESTEP_TRACER( - I bi,bj,iMin,iMax,jMin,jMax,k, + I bi,bj,iMin,iMax,jMin,jMax,k,tauAB, I salt, gS, U gSnm1, I myIter, myThid) ENDIF + IF ( tr1Stepping ) THEN + CALL CALC_GTR1( + I bi,bj,iMin,iMax,jMin,jMax, k,km1,kup,kDown, + I xA,yA,uTrans,vTrans,rTrans,maskUp, + I KappaRT, + U fVerTr1, + I myTime, myThid) + tauAB = 0.5d0 + abEps + CALL TIMESTEP_TRACER( + I bi,bj,iMin,iMax,jMin,jMax,k,tauAB, + I Tr1, gTr1, + U gTr1NM1, + I myIter, myThid) + ENDIF #ifdef ALLOW_OBCS C-- Apply open boundary conditions @@ -485,9 +565,13 @@ #ifdef ALLOW_AUTODIFF_TAMC -CPatrick? What about this one? - maximpl = 6 - iikey = (ikey-1)*maximpl +C? Patrick? What about this one? +cph Keys iikey and idkey don't seem to be needed +cph since storing occurs on different tape for each +cph impldiff call anyways. +cph Thus, common block comlev1_impl isn't needed either. +cph Storing below needed in the case useGMREDI. + iikey = (ikey-1)*maximpl #endif /* ALLOW_AUTODIFF_TAMC */ C-- Implicit diffusion @@ -496,6 +580,7 @@ IF (tempStepping) THEN #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 1 +CADJ STORE gTNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -507,6 +592,7 @@ IF (saltStepping) THEN #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 2 +CADJ STORE gSNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -515,6 +601,17 @@ I myThid ) ENDIF + IF (tr1Stepping) THEN +#ifdef ALLOW_AUTODIFF_TAMC +CADJ STORE gTr1Nm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte +#endif /* ALLOW_AUTODIFF_TAMC */ + CALL IMPLDIFF( + I bi, bj, iMin, iMax, jMin, jMax, + I deltaTtracer, KappaRT, recip_HFacC, + U gTr1Nm1, + I myThid ) + ENDIF + #ifdef ALLOW_OBCS C-- Apply open boundary conditions IF (useOBCS) THEN @@ -533,17 +630,14 @@ jMin = 1-OLy+2 jMax = sNy+OLy-1 -C-- Explicit part of the Surface Pressure Gradient (add in TIMESTEP) +C-- Explicit part of the Surface Potentiel Gradient (add in TIMESTEP) C (note: this loop will be replaced by CALL CALC_GRAD_ETA) IF (implicSurfPress.NE.1.) THEN - DO j=jMin,jMax - DO i=iMin,iMax - phiSurfX(i,j) = _recip_dxC(i,j,bi,bj)*gBaro - & *(cg2d_x(i,j,bi,bj)-cg2d_x(i-1,j,bi,bj)) - phiSurfY(i,j) = _recip_dyC(i,j,bi,bj)*gBaro - & *(cg2d_x(i,j,bi,bj)-cg2d_x(i,j-1,bi,bj)) - ENDDO - ENDDO + CALL CALC_GRAD_PHI_SURF( + I bi,bj,iMin,iMax,jMin,jMax, + I etaN, + O phiSurfX,phiSurfY, + I myThid ) ENDIF C-- Start of dynamics loop @@ -617,6 +711,7 @@ IF (implicitViscosity.AND.momStepping) THEN #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 3 +CADJ STORE gUNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -625,6 +720,7 @@ I myThid ) #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 4 +CADJ STORE gVNm1(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -644,6 +740,7 @@ #ifdef INCLUDE_CD_CODE #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 5 +CADJ STORE vVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -652,6 +749,7 @@ I myThid ) #ifdef ALLOW_AUTODIFF_TAMC idkey = iikey + 6 +CADJ STORE uVelD(:,:,:,bi,bj) = comlev1_bibj , key=ikey, byte=isbyte #endif /* ALLOW_AUTODIFF_TAMC */ CALL IMPLDIFF( I bi, bj, iMin, iMax, jMin, jMax, @@ -670,21 +768,38 @@ c ENDIF Cjmc(end) -#ifdef INCLUDE_DIAGNOSTICS_INTERFACE_CODE +#ifdef ALLOW_TIMEAVE IF (taveFreq.GT.0.) THEN - DO K=1,Nr - CALL TIMEAVER_1FLD_XYZ(phiHyd, phiHydtave, - I deltaTclock, bi, bj, K, myThid) + CALL TIMEAVE_CUMUL_1T(phiHydtave, phiHyd, Nr, + I deltaTclock, bi, bj, myThid) IF (ivdc_kappa.NE.0.) THEN - CALL TIMEAVER_1FLD_XYZ(ConvectCount, ConvectCountTave, - I deltaTclock, bi, bj, K, myThid) + CALL TIMEAVE_CUMULATE(ConvectCountTave, ConvectCount, Nr, + I deltaTclock, bi, bj, myThid) ENDIF - ENDDO ENDIF -#endif /* INCLUDE_DIAGNOSTICS_INTERFACE_CODE */ +#endif /* ALLOW_TIMEAVE */ ENDDO ENDDO +#ifndef EXCLUDE_DEBUGMODE + If (debugMode) THEN + CALL DEBUG_STATS_RL(1,EtaN,'EtaN (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,uVel,'Uvel (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,vVel,'Vvel (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,wVel,'Wvel (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,theta,'Theta (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,salt,'Salt (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,Gu,'Gu (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,Gv,'Gv (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,Gt,'Gt (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,Gs,'Gs (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,GuNm1,'GuNm1 (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,GvNm1,'GvNm1 (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,GtNm1,'GtNm1 (DYNAMICS)',myThid) + CALL DEBUG_STATS_RL(Nr,GsNm1,'GsNm1 (DYNAMICS)',myThid) + ENDIF +#endif + RETURN END