model/src/momentum_correction_step.F

a462ae7a37 Jean*0001 #include "PACKAGES_CONFIG.h"
                0002 #include "CPP_OPTIONS.h"
                0003
                0004 CBOP
                0005 C     !ROUTINE: MOMENTUM_CORRECTION_STEP
                0006 C     !INTERFACE:
bfdbc242ac Jean*0007       SUBROUTINE MOMENTUM_CORRECTION_STEP( myTime, myIter, myThid )
a462ae7a37 Jean*0008 C     !DESCRIPTION: \bv
                0009 C     *==========================================================*
f3ec658fc9 Jean*0010 C     | SUBROUTINE MOMENTUM_CORRECTION_STEP
a462ae7a37 Jean*0011 C     *==========================================================*
                0012 C     |1rst Part : Update U,V.
                0013 C     |
                0014 C     | The arrays used for time stepping are cycled.
                0015 C     | Momentum:
                0016 C     |           V(n) = Gv(n) - dt * grad Eta
                0017 C     |
                0018 C     |part1: update U,V
f549345f7f Jean*0019 C     |  U*,V* (contained in gU,gV) have the surface
                0020 C     |     pressure gradient term added and the result stored
                0021 C     |     in U,V (contained in uVel, vVel)
                0022 C     |
a462ae7a37 Jean*0023 C     |part2: Adjustments
f549345f7f Jean*0024 C     |   o Filter  U,V (Shapiro Filter, Zonal_Filter)
a462ae7a37 Jean*0025 C     *==========================================================*
                0026 C     \ev
                0027
                0028 C     !USES:
                0029       IMPLICIT NONE
                0030 C     == Global variables ===
                0031 #include "SIZE.h"
                0032 #include "EEPARAMS.h"
                0033 #include "PARAMS.h"
                0034 #include "DYNVARS.h"
                0035
                0036 #ifdef ALLOW_SHAP_FILT
                0037 #include "SHAP_FILT.h"
                0038 #endif
f3ec658fc9 Jean*0039 #ifdef ALLOW_ZONAL_FILT
a462ae7a37 Jean*0040 #include "ZONAL_FILT.h"
                0041 #endif
                0042
                0043 C     !INPUT/OUTPUT PARAMETERS:
f549345f7f Jean*0044 C     myTime :: Current time in simulation
                0045 C     myIter :: Current iteration number in simulation
bfdbc242ac Jean*0046 C     myThid :: my Thread Id. number
a462ae7a37 Jean*0047       _RL myTime
                0048       INTEGER myIter
                0049       INTEGER myThid
                0050
                0051 C     !LOCAL VARIABLES:
                0052       _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0053       _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
6c16fd0209 Jean*0054       INTEGER iMin, iMax
                0055       INTEGER jMin, jMax
                0056       INTEGER bi, bj
                0057       INTEGER i, j
a462ae7a37 Jean*0058 CEOP
                0059
                0060       DO bj=myByLo(myThid),myByHi(myThid)
                0061        DO bi=myBxLo(myThid),myBxHi(myThid)
                0062
6c16fd0209 Jean*0063         IF ( momStepping ) THEN
a462ae7a37 Jean*0064 C--     Set up work arrays that need valid initial values
6c16fd0209 Jean*0065           DO j=1-OLy,sNy+OLy
                0066            DO i=1-OLx,sNx+OLx
                0067             phiSurfX(i,j) = 0.
                0068             phiSurfY(i,j) = 0.
                0069            ENDDO
                0070           ENDDO
a462ae7a37 Jean*0071
                0072 C       Loop range: Gradients of Eta are evaluated so valid
                0073 C       range is all but first row and column in overlaps.
6c16fd0209 Jean*0074           iMin = 1-OLx+1
                0075           iMax = sNx+OLx
                0076           jMin = 1-OLy+1
                0077           jMax = sNy+OLy
a462ae7a37 Jean*0078
                0079 C-      Calculate gradient of surface Potentiel
6c16fd0209 Jean*0080           CALL CALC_GRAD_PHI_SURF(
                0081      I         bi, bj, iMin, iMax, jMin, jMax,
                0082      I         etaN,
                0083      O         phiSurfX, phiSurfY,
                0084      I         myThid )
                0085
                0086 C-      Update velocity fields:  V(n) = V** - dt * grad Eta
                0087           CALL CORRECTION_STEP(
                0088      I         bi, bj, iMin, iMax, jMin, jMax,
                0089      I         phiSurfX, phiSurfY,
                0090      I         myTime, myIter, myThid )
f549345f7f Jean*0091         ENDIF
a462ae7a37 Jean*0092
f549345f7f Jean*0093 #ifdef ALLOW_OBCS
6c16fd0209 Jean*0094         IF ( useOBCS ) THEN
f549345f7f Jean*0095           CALL OBCS_APPLY_UV( bi, bj, 0, uVel, vVel, myThid )
                0096         ENDIF
                0097 #endif /* ALLOW_OBCS */
a462ae7a37 Jean*0098
                0099 C--    End of 1rst bi,bj loop
                0100        ENDDO
                0101       ENDDO
                0102
                0103 C--- 2nd Part : Adjustment.
                0104
                0105 C--   Filter (and exchange)
                0106 #ifdef ALLOW_SHAP_FILT
6c16fd0209 Jean*0107       IF ( useSHAP_FILT ) THEN
f3ec658fc9 Jean*0108        IF ( .NOT.shap_filt_uvStar ) THEN
                0109         CALL TIMER_START('SHAP_FILT_UV       [MOM_CORR_STEP]',myThid)
                0110         CALL SHAP_FILT_APPLY_UV( uVel, vVel, myTime, myIter, myThid )
                0111         CALL TIMER_STOP ('SHAP_FILT_UV       [MOM_CORR_STEP]',myThid)
                0112        ENDIF
a462ae7a37 Jean*0113       ENDIF
f3ec658fc9 Jean*0114 #endif
a462ae7a37 Jean*0115 #ifdef ALLOW_ZONAL_FILT
6c16fd0209 Jean*0116       IF ( useZONAL_FILT ) THEN
f3ec658fc9 Jean*0117        IF ( .NOT.zonal_filt_uvStar ) THEN
                0118         CALL TIMER_START('ZONAL_FILT_UV      [MOM_CORR_STEP]',myThid)
                0119         CALL ZONAL_FILT_APPLY_UV( uVel, vVel, myThid )
                0120         CALL TIMER_STOP ('ZONAL_FILT_UV      [MOM_CORR_STEP]',myThid)
                0121        ENDIF
a462ae7a37 Jean*0122       ENDIF
f3ec658fc9 Jean*0123 #endif
a462ae7a37 Jean*0124
35c76859f0 Jean*0125 C--   Try to fix restart (Pb at machine trucation level accumulating in wVel):
                0126 C     Apply EXCH to horiz velocity before integrating continuity (-> wVel)
                0127       IF ( applyExchUV_early )
                0128      &  CALL EXCH_UV_3D_RL( uVel, vVel, .TRUE., Nr, myThid )
                0129
a462ae7a37 Jean*0130       RETURN
                0131       END