/[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.23 by jmc, Sat Jul 30 22:07:00 2005 UTC revision 1.31 by heimbach, Thu Dec 8 15:44:34 2005 UTC
# Line 121  C     uDudxFac, AhDudxFac, etc ... indiv Line 121  C     uDudxFac, AhDudxFac, etc ... indiv
121        _RL  rTransU(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL  rTransU(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
122        _RL  rTransV(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL  rTransV(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
123        _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
124  c     _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
125  c     _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
126  c     _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
127  c     _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
128  c     _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
129  c     _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
130        _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
131        _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
132        _RL  uDudxFac        _RL  uDudxFac
# Line 145  c     _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+O Line 145  c     _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+O
145        _RL  ArDvdrFac        _RL  ArDvdrFac
146        _RL  fvFac        _RL  fvFac
147        _RL  mtFacV        _RL  mtFacV
148        LOGICAL bottomDragTerms        _RL  sideMaskFac
149          LOGICAL bottomDragTerms,harmonic,biharmonic,useVariableViscosity
150  CEOP  CEOP
151    
152  C     Initialise intermediate terms  C     Initialise intermediate terms
# Line 165  C     Initialise intermediate terms Line 166  C     Initialise intermediate terms
166          tension(i,j)= 0.          tension(i,j)= 0.
167          guDiss(i,j) = 0.          guDiss(i,j) = 0.
168          gvDiss(i,j) = 0.          gvDiss(i,j) = 0.
169    #ifdef ALLOW_AUTODIFF_TAMC
170            vort3(i,j)   = 0. _d 0
171            strain(i,j)  = 0. _d 0
172            tension(i,j) = 0. _d 0
173    #endif
174         ENDDO         ENDDO
175        ENDDO        ENDDO
176    
# Line 193  C     o V momentum equation Line 199  C     o V momentum equation
199          ArDvdrFac  = 0.          ArDvdrFac  = 0.
200        ENDIF        ENDIF
201    
202    C note: using standard stencil (no mask) results in under-estimating
203    C       vorticity at a no-slip boundary by a factor of 2 = sideDragFactor
204          IF ( no_slip_sides ) THEN
205            sideMaskFac = sideDragFactor
206          ELSE
207            sideMaskFac = 0. _d 0
208          ENDIF
209    
210        IF (     no_slip_bottom        IF (     no_slip_bottom
211       &    .OR. bottomDragQuadratic.NE.0.       &    .OR. bottomDragQuadratic.NE.0.
212       &    .OR. bottomDragLinear.NE.0.) THEN       &    .OR. bottomDragLinear.NE.0.) THEN
# Line 231  C     Calculate velocity field "volume t Line 245  C     Calculate velocity field "volume t
245         ENDDO         ENDDO
246        ENDDO        ENDDO
247    
248        IF (bottomDragTerms) THEN        CALL MOM_CALC_KE(bi,bj,k,2,uFld,vFld,KE,myThid)
249          CALL MOM_CALC_KE(bi,bj,k,3,uFld,vFld,KE,myThid)        IF ( momViscosity) THEN
250        ENDIF          CALL MOM_CALC_HDIV(bi,bj,k,2,uFld,vFld,hDiv,myThid)
251            CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
252        IF (viscAstrain.NE.0. .OR. viscAtension.NE.0.) THEN          CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,tension,myThid)
253           CALL MOM_CALC_TENSION(bi,bj,k,uFld,vFld,          CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,strain,myThid)
254       O                         tension,          DO j=1-Oly,sNy+Oly
255       I                         myThid)           DO i=1-Olx,sNx+Olx
256           CALL MOM_CALC_STRAIN(bi,bj,k,uFld,vFld,hFacZ,             IF ( hFacZ(i,j).EQ.0. ) THEN
257       O                        strain,               vort3(i,j)  = sideMaskFac*vort3(i,j)
258       I                        myThid)               strain(i,j) = sideMaskFac*strain(i,j)
259               ENDIF
260             ENDDO
261            ENDDO
262    #ifdef ALLOW_DIAGNOSTICS
263            IF ( useDiagnostics ) THEN
264              CALL DIAGNOSTICS_FILL(hDiv,   'momHDiv ',k,1,2,bi,bj,myThid)
265              CALL DIAGNOSTICS_FILL(vort3,  'momVort3',k,1,2,bi,bj,myThid)
266              CALL DIAGNOSTICS_FILL(tension,'Tension ',k,1,2,bi,bj,myThid)
267              CALL DIAGNOSTICS_FILL(strain, 'Strain  ',k,1,2,bi,bj,myThid)
268            ENDIF
269    #endif
270        ENDIF        ENDIF
271    
272  C---  First call (k=1): compute vertical adv. flux fVerU(kUp) & fVerV(kUp)  C---  First call (k=1): compute vertical adv. flux fVerU(kUp) & fVerV(kUp)
# Line 270  C---  Calculate vertical transports (at Line 295  C---  Calculate vertical transports (at
295       I                        myTime, myIter, myThid)       I                        myTime, myIter, myThid)
296        ENDIF        ENDIF
297    
298  c     IF (momViscosity) THEN        IF (momViscosity) THEN
299  c    &  CALL MOM_CALC_VISCOSITY(bi,bj,k,         CALL MOM_CALC_VISC(
300  c    I                         uFld,vFld,       I        bi,bj,k,
301  c    O                         viscAh_D,viscAh_Z,myThid)       O        viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
302         O        harmonic,biharmonic,useVariableViscosity,
303         I        hDiv,vort3,tension,strain,KE,hFacZ,
304         I        myThid)
305          ENDIF
306    
307  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
308    
# Line 314  C--   Tendency is minus divergence of th Line 343  C--   Tendency is minus divergence of th
343           ENDDO           ENDDO
344          ENDDO          ENDDO
345    
346    #ifdef ALLOW_DIAGNOSTICS
347            IF ( useDiagnostics ) THEN
348              CALL DIAGNOSTICS_FILL(fZon,'ADVx_Um ',k,1,2,bi,bj,myThid)
349              CALL DIAGNOSTICS_FILL(fMer,'ADVy_Um ',k,1,2,bi,bj,myThid)
350              CALL DIAGNOSTICS_FILL(fVerU(1-Olx,1-Oly,kUp),
351         &                               'ADVrE_Um',k,1,2,bi,bj,myThid)
352            ENDIF
353    #endif
354    
355  #ifdef NONLIN_FRSURF  #ifdef NONLIN_FRSURF
356  C-- account for 3.D divergence of the flow in rStar coordinate:  C-- account for 3.D divergence of the flow in rStar coordinate:
357    # ifndef DISABLE_RSTAR_CODE
358          IF ( select_rStar.GT.0 ) THEN          IF ( select_rStar.GT.0 ) THEN
359           DO j=jMin,jMax           DO j=jMin,jMax
360            DO i=iMin,iMax            DO i=iMin,iMax
# Line 333  C-- account for 3.D divergence of the fl Line 372  C-- account for 3.D divergence of the fl
372            ENDDO            ENDDO
373           ENDDO           ENDDO
374          ENDIF          ENDIF
375    # endif /* DISABLE_RSTAR_CODE */
376  #endif /* NONLIN_FRSURF */  #endif /* NONLIN_FRSURF */
377    
378        ELSE        ELSE
# Line 350  C-    endif momAdvection. Line 390  C-    endif momAdvection.
390  C---  Calculate eddy fluxes (dissipation) between cells for zonal flow.  C---  Calculate eddy fluxes (dissipation) between cells for zonal flow.
391    
392  C     Bi-harmonic term del^2 U -> v4F  C     Bi-harmonic term del^2 U -> v4F
393          IF ( viscA4.NE.0. )          IF (biharmonic)
394       &  CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)       &  CALL MOM_U_DEL2U(bi,bj,k,uFld,hFacZ,v4f,myThid)
395    
396  C     Laplacian and bi-harmonic terms, Zonal  Fluxes -> fZon  C     Laplacian and bi-harmonic terms, Zonal  Fluxes -> fZon
397          CALL MOM_U_XVISCFLUX(bi,bj,k,uFld,v4F,fZon,myThid)          CALL MOM_U_XVISCFLUX(bi,bj,k,uFld,v4F,fZon,
398         I    viscAh_D,viscA4_D,myThid)
399    
400  C     Laplacian and bi-harmonic termis, Merid Fluxes -> fMer  C     Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
401          CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,fMer,myThid)          CALL MOM_U_YVISCFLUX(bi,bj,k,uFld,v4F,hFacZ,fMer,
402         I    viscAh_Z,viscA4_Z,myThid)
403    
404  C     Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw  C     Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
405         IF (.NOT.implicitViscosity) THEN         IF (.NOT.implicitViscosity) THEN
# Line 383  C--   Tendency is minus divergence of th Line 425  C--   Tendency is minus divergence of th
425           ENDDO           ENDDO
426          ENDDO          ENDDO
427    
428    #ifdef ALLOW_DIAGNOSTICS
429            IF ( useDiagnostics ) THEN
430              CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Um',k,1,2,bi,bj,myThid)
431              CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Um',k,1,2,bi,bj,myThid)
432              IF (.NOT.implicitViscosity)
433         &    CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Um',k,1,2,bi,bj,myThid)
434            ENDIF
435    #endif
436    
437  C-- No-slip and drag BCs appear as body forces in cell abutting topography  C-- No-slip and drag BCs appear as body forces in cell abutting topography
438          IF (no_slip_sides) THEN          IF (no_slip_sides) THEN
439  C-     No-slip BCs impose a drag at walls...  C-     No-slip BCs impose a drag at walls...
440           CALL MOM_U_SIDEDRAG(bi,bj,k,uFld,v4F,hFacZ,vF,myThid)           CALL MOM_U_SIDEDRAG(
441         I        bi,bj,k,
442         I        uFld, v4f, hFacZ,
443         I        viscAh_Z,viscA4_Z,
444         I        harmonic,biharmonic,useVariableViscosity,
445         O        vF,
446         I        myThid)
447           DO j=jMin,jMax           DO j=jMin,jMax
448            DO i=iMin,iMax            DO i=iMin,iMax
449             gUdiss(i,j) = gUdiss(i,j) + vF(i,j)             gUdiss(i,j) = gUdiss(i,j) + vF(i,j)
# Line 473  C--   Tendency is minus divergence of th Line 530  C--   Tendency is minus divergence of th
530       &      +( fMer(i,  j)     - fMer(i,j-1) )*vDvdyFac       &      +( fMer(i,  j)     - fMer(i,j-1) )*vDvdyFac
531       &      +(fVerV(i,j,kDown) - fVerV(i,j,kUp))*rkSign*rVelDvdrFac       &      +(fVerV(i,j,kDown) - fVerV(i,j,kUp))*rkSign*rVelDvdrFac
532       &     )       &     )
533         ENDDO           ENDDO
534        ENDDO          ENDDO
535    
536    #ifdef ALLOW_DIAGNOSTICS
537            IF ( useDiagnostics ) THEN
538              CALL DIAGNOSTICS_FILL(fZon,'ADVx_Vm ',k,1,2,bi,bj,myThid)
539              CALL DIAGNOSTICS_FILL(fMer,'ADVy_Vm ',k,1,2,bi,bj,myThid)
540              CALL DIAGNOSTICS_FILL(fVerV(1-Olx,1-Oly,kUp),
541         &                               'ADVrE_Vm',k,1,2,bi,bj,myThid)
542            ENDIF
543    #endif
544    
545  #ifdef NONLIN_FRSURF  #ifdef NONLIN_FRSURF
546  C-- account for 3.D divergence of the flow in rStar coordinate:  C-- account for 3.D divergence of the flow in rStar coordinate:
547    # ifndef DISABLE_RSTAR_CODE
548          IF ( select_rStar.GT.0 ) THEN          IF ( select_rStar.GT.0 ) THEN
549           DO j=jMin,jMax           DO j=jMin,jMax
550            DO i=iMin,iMax            DO i=iMin,iMax
# Line 495  C-- account for 3.D divergence of the fl Line 562  C-- account for 3.D divergence of the fl
562            ENDDO            ENDDO
563           ENDDO           ENDDO
564          ENDIF          ENDIF
565    # endif /* DISABLE_RSTAR_CODE */
566  #endif /* NONLIN_FRSURF */  #endif /* NONLIN_FRSURF */
567    
568        ELSE        ELSE
# Line 511  C-    endif momAdvection. Line 579  C-    endif momAdvection.
579        IF (momViscosity) THEN        IF (momViscosity) THEN
580  C---  Calculate eddy fluxes (dissipation) between cells for meridional flow.  C---  Calculate eddy fluxes (dissipation) between cells for meridional flow.
581  C     Bi-harmonic term del^2 V -> v4F  C     Bi-harmonic term del^2 V -> v4F
582          IF ( viscA4.NE.0. )          IF (biharmonic)
583       &  CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)       &  CALL MOM_V_DEL2V(bi,bj,k,vFld,hFacZ,v4f,myThid)
584    
585  C     Laplacian and bi-harmonic terms, Zonal  Fluxes -> fZon  C     Laplacian and bi-harmonic terms, Zonal  Fluxes -> fZon
586          CALL MOM_V_XVISCFLUX(bi,bj,k,vFld,v4f,hFacZ,fZon,myThid)          CALL MOM_V_XVISCFLUX(bi,bj,k,vFld,v4f,hFacZ,fZon,
587         I    viscAh_Z,viscA4_Z,myThid)
588    
589  C     Laplacian and bi-harmonic termis, Merid Fluxes -> fMer  C     Laplacian and bi-harmonic termis, Merid Fluxes -> fMer
590          CALL MOM_V_YVISCFLUX(bi,bj,k,vFld,v4f,fMer,myThid)          CALL MOM_V_YVISCFLUX(bi,bj,k,vFld,v4f,fMer,
591         I    viscAh_D,viscA4_D,myThid)
592    
593  C     Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw  C     Eddy component of vertical flux (interior component only) -> fVrUp & fVrDw
594         IF (.NOT.implicitViscosity) THEN         IF (.NOT.implicitViscosity) THEN
# Line 544  C--   Tendency is minus divergence of th Line 614  C--   Tendency is minus divergence of th
614           ENDDO           ENDDO
615          ENDDO          ENDDO
616    
617    #ifdef ALLOW_DIAGNOSTICS
618            IF ( useDiagnostics ) THEN
619              CALL DIAGNOSTICS_FILL(fZon, 'VISCx_Vm',k,1,2,bi,bj,myThid)
620              CALL DIAGNOSTICS_FILL(fMer, 'VISCy_Vm',k,1,2,bi,bj,myThid)
621              IF (.NOT.implicitViscosity)
622         &    CALL DIAGNOSTICS_FILL(fVrUp,'VISrE_Vm',k,1,2,bi,bj,myThid)
623            ENDIF
624    #endif
625    
626  C-- No-slip and drag BCs appear as body forces in cell abutting topography  C-- No-slip and drag BCs appear as body forces in cell abutting topography
627        IF (no_slip_sides) THEN        IF (no_slip_sides) THEN
628  C-     No-slip BCs impose a drag at walls...  C-     No-slip BCs impose a drag at walls...
629           CALL MOM_V_SIDEDRAG(bi,bj,k,vFld,v4F,hFacZ,vF,myThid)           CALL MOM_V_SIDEDRAG(
630         I        bi,bj,k,
631         I        vFld, v4f, hFacZ,
632         I        viscAh_Z,viscA4_Z,
633         I        harmonic,biharmonic,useVariableViscosity,
634         O        vF,
635         I        myThid)
636           DO j=jMin,jMax           DO j=jMin,jMax
637            DO i=iMin,iMax            DO i=iMin,iMax
638             gvDiss(i,j) = gvDiss(i,j) + vF(i,j)             gvDiss(i,j) = gvDiss(i,j) + vF(i,j)
# Line 614  c     ELSE Line 699  c     ELSE
699            gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)            gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj)+fuFac*cf(i,j)
700           ENDDO           ENDDO
701          ENDDO          ENDDO
702    #ifdef ALLOW_DIAGNOSTICS
703            IF ( useDiagnostics )
704         &    CALL DIAGNOSTICS_FILL(cf,'Um_Cori ',k,1,2,bi,bj,myThid)
705    #endif
706          CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)          CALL MOM_V_CORIOLIS(bi,bj,k,uFld,cf,myThid)
707          DO j=jMin,jMax          DO j=jMin,jMax
708           DO i=iMin,iMax           DO i=iMin,iMax
709            gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)            gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj)+fvFac*cf(i,j)
710           ENDDO           ENDDO
711          ENDDO          ENDDO
712    #ifdef ALLOW_DIAGNOSTICS
713            IF ( useDiagnostics )
714         &    CALL DIAGNOSTICS_FILL(cf,'Vm_Cori ',k,1,2,bi,bj,myThid)
715    #endif
716        ENDIF        ENDIF
717    
718        IF (nonHydrostatic.OR.quasiHydrostatic) THEN        IF (nonHydrostatic.OR.quasiHydrostatic) THEN
# Line 641  C--   Set du/dt & dv/dt on boundaries to Line 734  C--   Set du/dt & dv/dt on boundaries to
734         ENDDO         ENDDO
735        ENDDO        ENDDO
736    
737    #ifdef ALLOW_DIAGNOSTICS
738          IF ( useDiagnostics ) THEN
739            CALL DIAGNOSTICS_FILL(KE,    'momKE   ',k,1,2,bi,bj,myThid)
740            CALL DIAGNOSTICS_FILL(gU(1-Olx,1-Oly,k,bi,bj),
741         &                               'Um_Advec',k,1,2,bi,bj,myThid)
742            CALL DIAGNOSTICS_FILL(gV(1-Olx,1-Oly,k,bi,bj),
743         &                               'Vm_Advec',k,1,2,bi,bj,myThid)
744           IF (momViscosity) THEN
745            CALL DIAGNOSTICS_FILL(guDiss,'Um_Diss ',k,1,2,bi,bj,myThid)
746            CALL DIAGNOSTICS_FILL(gvDiss,'Vm_Diss ',k,1,2,bi,bj,myThid)
747           ENDIF
748          ENDIF
749    #endif /* ALLOW_DIAGNOSTICS */
750    
751        RETURN        RETURN
752        END        END
Legend:
Removed from v.1.23  
changed lines
  Added in v.1.31
  ViewVC Help
Powered by ViewVC 1.1.22