/[MITgcm]/MITgcm/pkg/mom_fluxform/mom_fluxform.F
ViewVC logotype

Diff of /MITgcm/pkg/mom_fluxform/mom_fluxform.F

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

revision 1.8 by jmc, Sun Jan 26 21:18:50 2003 UTC revision 1.15 by jmc, Sat Oct 11 16:37:55 2003 UTC
# Line 25  C where ${\bf v}=(u,v,w)$ and $\tau$, th Line 25  C where ${\bf v}=(u,v,w)$ and $\tau$, th
25  C stresses as well as internal viscous stresses.  C stresses as well as internal viscous stresses.
26  CEOI  CEOI
27    
28  #include "CPP_OPTIONS.h"  #include "MOM_FLUXFORM_OPTIONS.h"
29    
30  CBOP  CBOP
31  C !ROUTINE: MOM_FLUXFORM  C !ROUTINE: MOM_FLUXFORM
# Line 33  C !ROUTINE: MOM_FLUXFORM Line 33  C !ROUTINE: MOM_FLUXFORM
33  C !INTERFACE: ==========================================================  C !INTERFACE: ==========================================================
34        SUBROUTINE MOM_FLUXFORM(        SUBROUTINE MOM_FLUXFORM(
35       I        bi,bj,iMin,iMax,jMin,jMax,k,kUp,kDown,       I        bi,bj,iMin,iMax,jMin,jMax,k,kUp,kDown,
36       I        phi_hyd,KappaRU,KappaRV,       I        dPhihydX,dPhiHydY,KappaRU,KappaRV,
37       U        fVerU, fVerV,       U        fVerU, fVerV,
38       I        myTime,myIter,myThid)       I        myTime,myIter,myThid)
39    
# Line 58  C  iMin,iMax,jMin,jMAx  :: loop ranges Line 58  C  iMin,iMax,jMin,jMAx  :: loop ranges
58  C  k                    :: vertical level  C  k                    :: vertical level
59  C  kUp                  :: =1 or 2 for consecutive k  C  kUp                  :: =1 or 2 for consecutive k
60  C  kDown                :: =2 or 1 for consecutive k  C  kDown                :: =2 or 1 for consecutive k
61  C  phi_hyd              :: hydrostatic pressure (perturbation)  C  dPhiHydX,Y           :: Gradient (X & Y dir.) of Hydrostatic Potential
62  C  KappaRU              :: vertical viscosity  C  KappaRU              :: vertical viscosity
63  C  KappaRV              :: vertical viscosity  C  KappaRV              :: vertical viscosity
64  C  fVerU                :: vertical flux of U, 2 1/2 dim for pipe-lining  C  fVerU                :: vertical flux of U, 2 1/2 dim for pipe-lining
# Line 68  C  myIter               :: current time- Line 68  C  myIter               :: current time-
68  C  myThid               :: thread number  C  myThid               :: thread number
69        INTEGER bi,bj,iMin,iMax,jMin,jMax        INTEGER bi,bj,iMin,iMax,jMin,jMax
70        INTEGER k,kUp,kDown        INTEGER k,kUp,kDown
71        _RL phi_hyd(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL dPhiHydX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72          _RL dPhiHydY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73        _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL KappaRU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
74        _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL KappaRV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
75        _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fVerU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
# Line 178  C     Initialise intermediate terms Line 179  C     Initialise intermediate terms
179          fMer(i,j) = 0.          fMer(i,j) = 0.
180          rTransU(i,j) = 0.          rTransU(i,j) = 0.
181          rTransV(i,j) = 0.          rTransV(i,j) = 0.
182    #ifdef ALLOW_AUTODIFF_TAMC
183    C- jmc: this is wrong, but at least with #ifdef/endif TAMC, it does not break
184    C       the forward code ; (same thing in mom_vectinv)
185            fVerU(i,j,1) = 0. _d 0
186            fVerU(i,j,2) = 0. _d 0
187            fVerV(i,j,1) = 0. _d 0
188            fVerV(i,j,2) = 0. _d 0
189    #endif
190         ENDDO         ENDDO
191        ENDDO        ENDDO
192    
# Line 292  C---  Calculate vertical transports (at Line 301  C---  Calculate vertical transports (at
301  C---- Zonal momentum equation starts here  C---- Zonal momentum equation starts here
302    
303  C     Bi-harmonic term del^2 U -> v4F  C     Bi-harmonic term del^2 U -> v4F
304        IF (momViscosity)        IF (momViscosity .AND. viscA4.NE.0. )
305       & CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)       & CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)
306    
307  C---  Calculate mean and eddy fluxes between cells for zonal flow.  C---  Calculate mean and eddy fluxes between cells for zonal flow.
# Line 325  C     Laplacian and bi-harmonic term Line 334  C     Laplacian and bi-harmonic term
334       & CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,vF,myThid)       & CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,vF,myThid)
335    
336  C     Combine fluxes -> fMer  C     Combine fluxes -> fMer
337        DO j=jMin,jMax        DO j=jMin,jMax+1
338         DO i=iMin,iMax         DO i=iMin,iMax
339          fMer(i,j) = vDudyFac*aF(i,j) + AhDudyFac*vF(i,j)          fMer(i,j) = vDudyFac*aF(i,j) + AhDudyFac*vF(i,j)
340         ENDDO         ENDDO
# Line 348  C     Combine fluxes Line 357  C     Combine fluxes
357         ENDDO         ENDDO
358        ENDDO        ENDDO
359    
 C---  Hydrostatic term ( -1/rhoConst . dphi/dx )  
       IF (momPressureForcing) THEN  
        DO j=jMin,jMax  
         DO i=iMin,iMax  
          pf(i,j) = - _recip_dxC(i,j,bi,bj)  
      &    *(phi_hyd(i,j,k)-phi_hyd(i-1,j,k))  
         ENDDO  
        ENDDO  
       ENDIF  
   
360  C--   Tendency is minus divergence of the fluxes + coriolis + pressure term  C--   Tendency is minus divergence of the fluxes + coriolis + pressure term
361        DO j=jMin,jMax        DO j=jMin,jMax
362         DO i=iMin,iMax         DO i=iMin,iMax
# Line 373  C--   Tendency is minus divergence of th Line 372  C--   Tendency is minus divergence of th
372       &   +fMer(i,j+1)          - fMer(i  ,j)       &   +fMer(i,j+1)          - fMer(i  ,j)
373       &   +fVerU(i,j,kUp)*rkFac - fVerU(i,j,kDown)*rkFac       &   +fVerU(i,j,kUp)*rkFac - fVerU(i,j,kDown)*rkFac
374       &   )       &   )
375       & _PHM( +phxFac * pf(i,j) )       &  - phxFac*dPhiHydX(i,j)
376         ENDDO         ENDDO
377        ENDDO        ENDDO
378    
# Line 418  C-    No-slip BCs impose a drag at botto Line 417  C-    No-slip BCs impose a drag at botto
417         ENDDO         ENDDO
418        ENDIF        ENDIF
419    
420  C--   Forcing term  C--   Forcing term (moved to timestep.F)
421        IF (momForcing)  c     IF (momForcing)
422       &  CALL EXTERNAL_FORCING_U(  c    &  CALL EXTERNAL_FORCING_U(
423       I     iMin,iMax,jMin,jMax,bi,bj,k,  c    I     iMin,iMax,jMin,jMax,bi,bj,k,
424       I     myTime,myThid)  c    I     myTime,myThid)
425    
426  C--   Metric terms for curvilinear grid systems  C--   Metric terms for curvilinear grid systems
427        IF (useNHMTerms) THEN        IF (useNHMTerms) THEN
# Line 454  C--   Set du/dt on boundaries to zero Line 453  C--   Set du/dt on boundaries to zero
453  C---- Meridional momentum equation starts here  C---- Meridional momentum equation starts here
454    
455  C     Bi-harmonic term del^2 V -> v4F  C     Bi-harmonic term del^2 V -> v4F
456        IF (momViscosity)        IF (momViscosity .AND. viscA4.NE.0. )
457       & CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)       & CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)
458    
459  C---  Calculate mean and eddy fluxes between cells for meridional flow.  C---  Calculate mean and eddy fluxes between cells for meridional flow.
# Line 471  C     Laplacian and bi-harmonic terms -> Line 470  C     Laplacian and bi-harmonic terms ->
470    
471  C     Combine fluxes -> fZon  C     Combine fluxes -> fZon
472        DO j=jMin,jMax        DO j=jMin,jMax
473         DO i=iMin,iMax         DO i=iMin,iMax+1
474          fZon(i,j) = uDvdxFac*aF(i,j) + AhDvdxFac*vF(i,j)          fZon(i,j) = uDvdxFac*aF(i,j) + AhDvdxFac*vF(i,j)
475         ENDDO         ENDDO
476        ENDDO        ENDDO
# Line 510  C     Combine fluxes -> fVerV Line 509  C     Combine fluxes -> fVerV
509         ENDDO         ENDDO
510        ENDDO        ENDDO
511    
 C---  Hydorstatic term (-1/rhoConst . dphi/dy )  
       IF (momPressureForcing) THEN  
        DO j=jMin,jMax  
         DO i=iMin,iMax  
          pF(i,j) = -_recip_dyC(i,j,bi,bj)  
      &    *(phi_hyd(i,j,k)-phi_hyd(i,j-1,k))  
         ENDDO  
        ENDDO  
       ENDIF  
   
512  C--   Tendency is minus divergence of the fluxes + coriolis + pressure term  C--   Tendency is minus divergence of the fluxes + coriolis + pressure term
513        DO j=jMin,jMax        DO j=jMin,jMax
514         DO i=iMin,iMax         DO i=iMin,iMax
# Line 535  C--   Tendency is minus divergence of th Line 524  C--   Tendency is minus divergence of th
524       &   +fMer(i,j  )          - fMer(i,j-1)       &   +fMer(i,j  )          - fMer(i,j-1)
525       &   +fVerV(i,j,kUp)*rkFac - fVerV(i,j,kDown)*rkFac       &   +fVerV(i,j,kUp)*rkFac - fVerV(i,j,kDown)*rkFac
526       &   )       &   )
527       & _PHM( +phyFac*pf(i,j) )       &  - phyFac*dPhiHydY(i,j)
528         ENDDO         ENDDO
529        ENDDO        ENDDO
530    
# Line 580  C-    No-slip BCs impose a drag at botto Line 569  C-    No-slip BCs impose a drag at botto
569         ENDDO         ENDDO
570        ENDIF        ENDIF
571    
572  C--   Forcing term  C--   Forcing term (moved to timestep.F)
573        IF (momForcing)  c     IF (momForcing)
574       & CALL EXTERNAL_FORCING_V(  c    & CALL EXTERNAL_FORCING_V(
575       I     iMin,iMax,jMin,jMax,bi,bj,k,  c    I     iMin,iMax,jMin,jMax,bi,bj,k,
576       I     myTime,myThid)  c    I     myTime,myThid)
577    
578  C--   Metric terms for curvilinear grid systems  C--   Metric terms for curvilinear grid systems
579        IF (useNHMTerms) THEN        IF (useNHMTerms) THEN
# Line 614  C--   Set dv/dt on boundaries to zero Line 603  C--   Set dv/dt on boundaries to zero
603    
604  C--   Coriolis term  C--   Coriolis term
605  C     Note. As coded here, coriolis will not work with "thin walls"  C     Note. As coded here, coriolis will not work with "thin walls"
606  #ifdef INCLUDE_CD_CODE  c     IF (useCDscheme) THEN
607        CALL MOM_CDSCHEME(bi,bj,k,phi_hyd,myThid)  c       CALL MOM_CDSCHEME(bi,bj,k,dPhiHydX,dPhiHydY,myThid)
608  #else  c     ELSE
609        CALL MOM_U_CORIOLIS(bi,bj,k,vFld,cf,myThid)        IF (.NOT.useCDscheme) THEN
610        DO j=jMin,jMax          CALL MOM_U_CORIOLIS(bi,bj,k,vFld,cf,myThid)
611         DO i=iMin,iMax          DO j=jMin,jMax
612          gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)           DO i=iMin,iMax
613         ENDDO            gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
614        ENDDO           ENDDO
615        CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)          ENDDO
616        DO j=jMin,jMax          CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)
617         DO i=iMin,iMax          DO j=jMin,jMax
618          gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)           DO i=iMin,iMax
619         ENDDO            gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)
620        ENDDO           ENDDO
621  #endif /* INCLUDE_CD_CODE */          ENDDO
622          ENDIF
623    
624        IF (nonHydrostatic.OR.quasiHydrostatic) THEN        IF (nonHydrostatic.OR.quasiHydrostatic) THEN
625         CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,cf,myThid)         CALL MOM_U_CORIOLIS_NH(bi,bj,k,wVel,cf,myThid)
626         DO j=jMin,jMax         DO j=jMin,jMax
Legend:
Removed from v.1.8  
changed lines
  Added in v.1.15
  ViewVC Help
Powered by ViewVC 1.1.22