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

Contents of /MITgcm/model/src/calc_div_ghat.F

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


Revision 1.14 - (show annotations) (download)
Mon Jun 25 20:38:15 2001 UTC (22 years, 11 months ago) by ljmc
Branch: MAIN
Changes since 1.13: +3 -1 lines
the default is now to call the filter after solve_for_pressure

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/calc_div_ghat.F,v 1.13 2001/03/06 16:51:02 jmc Exp $
2 C $Name: $
3
4 #include "CPP_OPTIONS.h"
5 #undef ALLOW_ZONAL_FILT
6 #undef ALLOW_SHAP_FILT
7
8 C /==========================================================\
9 C | S/R CALC_DIV_GHAT |
10 C | o Form the right hand-side of the surface pressure eqn. |
11 C |==========================================================|
12 C \==========================================================/
13 SUBROUTINE CALC_DIV_GHAT(
14 I bi,bj,iMin,iMax,jMin,jMax,K,
15 I xA,yA,
16 U cg2d_b,
17 I myThid)
18
19 IMPLICIT NONE
20
21 C == Global variables ==
22 #include "SIZE.h"
23 #include "DYNVARS.h"
24 #include "FFIELDS.h"
25 #include "EEPARAMS.h"
26 #include "PARAMS.h"
27 #include "GRID.h"
28 #ifdef ALLOW_NONHYDROSTATIC
29 #include "CG3D.h"
30 #endif
31
32 C == Routine arguments ==
33 C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
34 C results will be set.
35 C k - Index of layer.
36 C xA, yA - Cell face areas
37 C cg2d_b - Conjugate Gradient 2-D solver : Right-hand side vector
38 C myThid - Instance number for this call of CALC_DIV_GHAT
39 INTEGER bi,bj,iMin,iMax,jMin,jMax
40 INTEGER K
41 _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
42 _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
43 _RL cg2d_b(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
44 INTEGER myThid
45
46 C == Local variables ==
47 INTEGER i,j
48 _RL pf (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49
50 C-- Pressure equation source term
51 C Term is the vertical integral of the divergence of the
52 C time tendency terms along with a relaxation term that
53 C pulls div(U) + dh/dt back toward zero.
54
55 IF (implicDiv2Dflow.EQ.1.) then
56 C Fully Implicit treatment of the Barotropic Flow Divergence
57 DO j=1,sNy
58 DO i=1,sNx+1
59 pf(i,j) = xA(i,j)*gUNm1(i,j,k,bi,bj) / deltaTmom
60 ENDDO
61 ENDDO
62 ELSE
63 C Explicit+Implicit part of the Barotropic Flow Divergence
64 C => Filtering of uVel,vVel is necessary
65 #ifdef ALLOW_ZONAL_FILT
66 IF (zonal_filt_lat.LT.90.) THEN
67 CALL ZONAL_FILTER(
68 & uVel, hFacW, 1-1, sNy+1, k, k, bi, bj, 1, myThid)
69 CALL ZONAL_FILTER(
70 & vVel, hFacS, 1-1, sNy+1, k, k, bi, bj, 2, myThid)
71 ENDIF
72 #endif
73 DO j=1,sNy
74 DO i=1,sNx+1
75 pf(i,j) = ( implicDiv2Dflow * gUNm1(i,j,k,bi,bj)
76 & + (1.-implicDiv2Dflow) * uVel(i,j,k,bi,bj)
77 & ) * xA(i,j) / deltaTmom
78 ENDDO
79 ENDDO
80 ENDIF
81 DO j=1,sNy
82 DO i=1,sNx
83 cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj) +
84 & pf(i+1,j)-pf(i,j)
85 ENDDO
86 ENDDO
87
88 #ifdef ALLOW_NONHYDROSTATIC
89 IF (nonHydrostatic) THEN
90 DO j=1,sNy
91 DO i=1,sNx
92 cg3d_b(i,j,k,bi,bj) = pf(i+1,j)-pf(i,j)
93 ENDDO
94 ENDDO
95 ENDIF
96 #endif
97
98 IF (implicDiv2Dflow.EQ.1.) then
99 C Fully Implicit treatment of the Barotropic Flow Divergence
100 DO j=1,sNy+1
101 DO i=1,sNx
102 pf(i,j) = yA(i,j)*gVNm1(i,j,k,bi,bj) / deltatmom
103 ENDDO
104 ENDDO
105 ELSE
106 C Explicit+Implicit part of the Barotropic Flow Divergence
107 DO j=1,sNy+1
108 DO i=1,sNx
109 pf(i,j) = ( implicDiv2Dflow * gVNm1(i,j,k,bi,bj)
110 & + (1.-implicDiv2Dflow) * vVel(i,j,k,bi,bj)
111 & ) * yA(i,j) / deltaTmom
112 ENDDO
113 ENDDO
114 ENDIF
115 DO j=1,sNy
116 DO i=1,sNx
117 cg2d_b(i,j,bi,bj) = cg2d_b(i,j,bi,bj) +
118 & pf(i,j+1)-pf(i,j)
119 ENDDO
120 ENDDO
121
122 #ifdef ALLOW_NONHYDROSTATIC
123 IF (nonHydrostatic) THEN
124 DO j=1,sNy
125 DO i=1,sNx
126 cg3d_b(i,j,k,bi,bj) = cg3d_b(i,j,k,bi,bj) +
127 & pf(i,j+1)-pf(i,j)
128 ENDDO
129 ENDDO
130 ENDIF
131 #endif
132
133 RETURN
134 END

  ViewVC Help
Powered by ViewVC 1.1.22