/[MITgcm]/MITgcm/model/src/solve_for_pressure.F
ViewVC logotype

Diff of /MITgcm/model/src/solve_for_pressure.F

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

revision 1.50 by jmc, Mon Dec 19 21:10:34 2005 UTC revision 1.51 by jmc, Tue Dec 20 20:31:28 2005 UTC
# Line 64  C     == Local variables == Line 64  C     == Local variables ==
64        INTEGER numIters        INTEGER numIters
65        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
66  #ifdef ALLOW_NONHYDROSTATIC  #ifdef ALLOW_NONHYDROSTATIC
67        INTEGER ks        INTEGER ks, kp1
68          _RL     maskKp1
69        LOGICAL zeroPsNH        LOGICAL zeroPsNH
70  #endif  #endif
71  CEOP  CEOP
# Line 166  C--   Add source term arising from w=d/d Line 167  C--   Add source term arising from w=d/d
167          IF ( nonHydrostatic .AND. zeroPsNH ) THEN          IF ( nonHydrostatic .AND. zeroPsNH ) THEN
168           DO j=1,sNy           DO j=1,sNy
169            DO i=1,sNx            DO i=1,sNx
170             cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj)             ks = ksurfC(i,j,bi,bj)
171               IF ( ks.LE.Nr ) THEN
172                cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj)
173       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf
174       &         * etaH(i,j,bi,bj)       &         * etaH(i,j,bi,bj)
175             cg3d_b(i,j,1,bi,bj) = cg3d_b(i,j,1,bi,bj)              cg3d_b(i,j,ks,bi,bj) = cg3d_b(i,j,ks,bi,bj)
176       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf
177       &         * etaH(i,j,bi,bj)       &         * etaH(i,j,bi,bj)
178               ENDIF
179            ENDDO            ENDDO
180           ENDDO           ENDDO
181          ELSEIF ( nonHydrostatic ) THEN          ELSEIF ( nonHydrostatic ) THEN
182           DO j=1,sNy           DO j=1,sNy
183            DO i=1,sNx            DO i=1,sNx
184             cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj)             ks = ksurfC(i,j,bi,bj)
185               IF ( ks.LE.Nr ) THEN
186                cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj)
187       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf
188       &         *( etaN(i,j,bi,bj)       &         *( etaN(i,j,bi,bj)
189       &           +phi_nh(i,j,1,bi,bj)*horiVertRatio/gravity )       &           +phi_nh(i,j,ks,bi,bj)*horiVertRatio/gravity )
190             cg3d_b(i,j,1,bi,bj) = cg3d_b(i,j,1,bi,bj)              cg3d_b(i,j,ks,bi,bj) = cg3d_b(i,j,ks,bi,bj)
191       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf       &       -freeSurfFac*_rA(i,j,bi,bj)/deltaTMom/deltaTfreesurf
192       &         *( etaN(i,j,bi,bj)       &         *( etaN(i,j,bi,bj)
193       &           +phi_nh(i,j,1,bi,bj)*horiVertRatio/gravity )       &           +phi_nh(i,j,ks,bi,bj)*horiVertRatio/gravity )
194               ENDIF
195            ENDDO            ENDDO
196           ENDDO           ENDDO
197          ELSEIF ( exactConserv ) THEN          ELSEIF ( exactConserv ) THEN
# Line 338  C Western boundary Line 345  C Western boundary
345           ENDIF           ENDIF
346  #endif /* ALLOW_OBCS */  #endif /* ALLOW_OBCS */
347    
348             IF ( usingZCoords ) THEN
349    C-       Z coordinate: assume surface @ level k=1
350               tmpFac = freeSurfFac
351             ELSE
352    C-       Other than Z coordinate: no assumption on surface level index
353               tmpFac = 0.
354               DO j=1,sNy
355                DO i=1,sNx
356                  ks = ksurfC(i,j,bi,bj)
357                  IF ( ks.LE.Nr ) THEN
358                   cg3d_b(i,j,ks,bi,bj) = cg3d_b(i,j,ks,bi,bj)
359         &              +freeSurfFac*etaN(i,j,bi,bj)/deltaTfreesurf
360         &                          *_rA(i,j,bi,bj)/deltaTmom
361                  ENDIF
362                ENDDO
363               ENDDO
364             ENDIF
365           K=1           K=1
366             kp1 = MIN(k+1,Nr)
367             maskKp1 = 1.
368             IF (k.GE.Nr) maskKp1 = 0.
369           DO j=1,sNy           DO j=1,sNy
370            DO i=1,sNx            DO i=1,sNx
371             cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj)              cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj)
372       &       +drF(K)*dYG(i+1,j,bi,bj)*hFacW(i+1,j,k,bi,bj)*uf(i+1,j)       &       +drF(K)*dyG(i+1,j,bi,bj)*hFacW(i+1,j,k,bi,bj)*uf(i+1,j)
373       &       -drF(K)*dYG( i ,j,bi,bj)*hFacW( i ,j,k,bi,bj)*uf( i ,j)       &       -drF(K)*dyG( i ,j,bi,bj)*hFacW( i ,j,k,bi,bj)*uf( i ,j)
374       &       +drF(K)*dXG(i,j+1,bi,bj)*hFacS(i,j+1,k,bi,bj)*vf(i,j+1)       &       +drF(K)*dxG(i,j+1,bi,bj)*hFacS(i,j+1,k,bi,bj)*vf(i,j+1)
375       &       -drF(K)*dXG(i, j ,bi,bj)*hFacS(i, j ,k,bi,bj)*vf(i, j )       &       -drF(K)*dxG(i, j ,bi,bj)*hFacS(i, j ,k,bi,bj)*vf(i, j )
376       &       +( freeSurfFac*etaN(i,j,bi,bj)/deltaTMom       &       +( tmpFac*etaN(i,j,bi,bj)/deltaTfreesurf
377       &          -wVel(i,j,k+1,bi,bj)       &         -wVel(i,j,kp1,bi,bj)*maskKp1
378       &        )*_rA(i,j,bi,bj)/deltaTmom       &        )*_rA(i,j,bi,bj)/deltaTmom
379            ENDDO            ENDDO
380           ENDDO           ENDDO
381           DO K=2,Nr-1           DO K=2,Nr
382              kp1 = MIN(k+1,Nr)
383              maskKp1 = 1.
384              IF (k.GE.Nr) maskKp1 = 0.
385            DO j=1,sNy            DO j=1,sNy
386             DO i=1,sNx             DO i=1,sNx
387              cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj)              cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj)
388       &       +drF(K)*dYG(i+1,j,bi,bj)*hFacW(i+1,j,k,bi,bj)*uf(i+1,j)       &       +drF(K)*dyG(i+1,j,bi,bj)*hFacW(i+1,j,k,bi,bj)*uf(i+1,j)
389       &       -drF(K)*dYG( i ,j,bi,bj)*hFacW( i ,j,k,bi,bj)*uf( i ,j)       &       -drF(K)*dyG( i ,j,bi,bj)*hFacW( i ,j,k,bi,bj)*uf( i ,j)
390       &       +drF(K)*dXG(i,j+1,bi,bj)*hFacS(i,j+1,k,bi,bj)*vf(i,j+1)       &       +drF(K)*dxG(i,j+1,bi,bj)*hFacS(i,j+1,k,bi,bj)*vf(i,j+1)
391       &       -drF(K)*dXG(i, j ,bi,bj)*hFacS(i, j ,k,bi,bj)*vf(i, j )       &       -drF(K)*dxG(i, j ,bi,bj)*hFacS(i, j ,k,bi,bj)*vf(i, j )
392       &       +( wVel(i,j,k  ,bi,bj)       &       +( wVel(i,j,k  ,bi,bj)*maskC(i,j,k-1,bi,bj)
393       &         -wVel(i,j,k+1,bi,bj)       &         -wVel(i,j,kp1,bi,bj)*maskKp1
394       &        )*_rA(i,j,bi,bj)/deltaTmom       &        )*_rA(i,j,bi,bj)/deltaTmom
395    
396             ENDDO             ENDDO
397            ENDDO            ENDDO
398           ENDDO           ENDDO
          K=Nr  
          DO j=1,sNy  
           DO i=1,sNx  
             cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj)  
      &       +drF(K)*dYG(i+1,j,bi,bj)*hFacW(i+1,j,k,bi,bj)*uf(i+1,j)  
      &       -drF(K)*dYG( i ,j,bi,bj)*hFacW( i ,j,k,bi,bj)*uf( i ,j)  
      &       +drF(K)*dXG(i,j+1,bi,bj)*hFacS(i,j+1,k,bi,bj)*vf(i,j+1)  
      &       -drF(K)*dXG(i, j ,bi,bj)*hFacS(i, j ,k,bi,bj)*vf(i, j )  
      &       +( wVel(i,j,k  ,bi,bj)  
      &        )*_rA(i,j,bi,bj)/deltaTmom  
   
           ENDDO  
          ENDDO  
399    
400  #ifdef ALLOW_OBCS  #ifdef ALLOW_OBCS
401           IF (useOBCS) THEN           IF (useOBCS) THEN
Legend:
Removed from v.1.50  
changed lines
  Added in v.1.51
  ViewVC Help
Powered by ViewVC 1.1.22