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

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

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


Revision 1.15 - (show annotations) (download)
Mon Mar 22 15:54:04 1999 UTC (25 years, 2 months ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint20
Changes since 1.14: +37 -2 lines
Modifications for non-hydrostatic ability + updates for open-boundaries.

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_masks_etc.F,v 1.14 1998/12/10 00:16:16 adcroft Exp $
2
3 #include "CPP_OPTIONS.h"
4
5 CStartOfInterface
6 SUBROUTINE INI_MASKS_ETC( myThid )
7 C /==========================================================\
8 C | SUBROUTINE INI_MASKS_ETC |
9 C | o Initialise masks and topography factors |
10 C |==========================================================|
11 C | These arrays are used throughout the code and describe |
12 C | the topography of the domain through masks (0s and 1s) |
13 C | and fractional height factors (0<hFac<1). The latter |
14 C | distinguish between the lopped-cell and full-step |
15 C | topographic representations. |
16 C \==========================================================/
17 IMPLICIT NONE
18
19 C === Global variables ===
20 #include "SIZE.h"
21 #include "EEPARAMS.h"
22 #include "PARAMS.h"
23 #include "GRID.h"
24
25 C == Routine arguments ==
26 C myThid - Number of this instance of INI_MASKS_ETC
27 INTEGER myThid
28 CEndOfInterface
29
30 C == Local variables ==
31 C bi,bj - Loop counters
32 C I,J,K
33 INTEGER bi, bj
34 INTEGER I, J, K
35 #ifdef ALLOW_NONHYDROSTATIC
36 INTEGER Km1
37 _RL hFacUpper,hFacLower
38 #endif
39
40 C Calculate quantities derived from XY depth map
41 DO bj = myByLo(myThid), myByHi(myThid)
42 DO bi = myBxLo(myThid), myBxHi(myThid)
43 DO J=1,sNy
44 DO I=1,sNx
45 C Inverse of depth
46 IF ( h(i,j,bi,bj) .EQ. 0. _d 0 ) THEN
47 recip_H(i,j,bi,bj) = 0. _d 0
48 ELSE
49 recip_H(i,j,bi,bj) = 1. _d 0 / abs( H(i,j,bi,bj) )
50 ENDIF
51 depthInK(i,j,bi,bj) = 0.
52 ENDDO
53 ENDDO
54 ENDDO
55 ENDDO
56 _EXCH_XY_R4( recip_H, myThid )
57
58 C Calculate lopping factor hFacC
59 DO bj=myByLo(myThid), myByHi(myThid)
60 DO bi=myBxLo(myThid), myBxHi(myThid)
61 DO K=1, Nr
62 DO J=1,sNy
63 DO I=1,sNx
64 C Round depths within a small fraction of layer depth to that
65 C layer depth.
66 IF ( ABS(H(I,J,bi,bj)-rF(K)) .LT.
67 & 1. _d -6*ABS(rF(K)) .AND.
68 & ABS(H(I,J,bi,bj)-rF(K)) .LT.
69 & 1. _d -6*ABS(H(I,J,bi,bj)) )THEN
70 H(I,J,bi,bj) = rF(K)
71 ENDIF
72 IF ( H(I,J,bi,bj)*rkFac .GE. rF(K)*rkFac ) THEN
73 C Top of cell is below base of domain
74 hFacC(I,J,K,bi,bj) = 0.
75 ELSEIF ( H(I,J,bi,bj)*rkFac .LE. rF(K+1)*rkFac ) THEN
76 C Base of domain is below bottom of this cell
77 hFacC(I,J,K,bi,bj) = 1.
78 ELSE
79 C Base of domain is in this cell
80 C Set hFac to the fraction of the cell that is open.
81 hFacC(I,J,K,bi,bj) =
82 & (rF(K)*rkFac-H(I,J,bi,bj)*rkFac)*recip_drF(K)
83 ENDIF
84 C Impose minimum fraction
85 IF (hFacC(I,J,K,bi,bj).LT.hFacMin) THEN
86 IF (hFacC(I,J,K,bi,bj).LT.hFacMin*0.5) THEN
87 hFacC(I,J,K,bi,bj)=0.
88 ELSE
89 hFacC(I,J,K,bi,bj)=hFacMin
90 ENDIF
91 ENDIF
92 C Impose minimum size (dimensional)
93 IF (drF(k)*hFacC(I,J,K,bi,bj).LT.hFacMinDr) THEN
94 IF (drF(k)*hFacC(I,J,K,bi,bj).LT.hFacMinDr*0.5) THEN
95 hFacC(I,J,K,bi,bj)=0.
96 ELSE
97 hFacC(I,J,K,bi,bj)=hFacMinDr*recip_drF(k)
98 ENDIF
99 ENDIF
100 depthInK(i,j,bi,bj) = depthInK(i,j,bi,bj)
101 & +hFacC(i,j,k,bi,bj)
102 ENDDO
103 ENDDO
104 ENDDO
105 ENDDO
106 ENDDO
107 _EXCH_XYZ_R4(hFacC , myThid )
108 _EXCH_XY_R4( depthInK, myThid )
109
110 CALL PLOT_FIELD_XYRS( depthInK,
111 & 'Model Depths K Index' , 1, myThid )
112
113 C hFacW and hFacS (at U and V points)
114 DO bj=myByLo(myThid), myByHi(myThid)
115 DO bi=myBxLo(myThid), myBxHi(myThid)
116 DO K=1, Nr
117 DO J=1,sNy
118 DO I=1,sNx
119 hFacW(I,J,K,bi,bj)=
120 & MIN(hFacC(I,J,K,bi,bj),hFacC(I-1,J,K,bi,bj))
121 hFacS(I,J,K,bi,bj)=
122 & MIN(hFacC(I,J,K,bi,bj),hFacC(I,J-1,K,bi,bj))
123 ENDDO
124 ENDDO
125 ENDDO
126 ENDDO
127 ENDDO
128 _EXCH_XYZ_R4(hFacW , myThid )
129 _EXCH_XYZ_R4(hFacS , myThid )
130
131 C Masks and reciprocals of hFac[CWS]
132 DO bj = myByLo(myThid), myByHi(myThid)
133 DO bi = myBxLo(myThid), myBxHi(myThid)
134 DO K=1,Nr
135 DO J=1,sNy
136 DO I=1,sNx
137 IF (HFacC(I,J,K,bi,bj) .NE. 0. _d 0 ) THEN
138 recip_HFacC(I,J,K,bi,bj) = 1. _d 0 / HFacC(I,J,K,bi,bj)
139 ELSE
140 recip_HFacC(I,J,K,bi,bj) = 0. _d 0
141 ENDIF
142 IF (HFacW(I,J,K,bi,bj) .NE. 0. _d 0 ) THEN
143 recip_HFacW(I,J,K,bi,bj) = 1. _d 0 / HFacW(I,J,K,bi,bj)
144 maskW(I,J,K,bi,bj) = 1. _d 0
145 ELSE
146 recip_HFacW(I,J,K,bi,bj) = 0. _d 0
147 maskW(I,J,K,bi,bj) = 0.0 _d 0
148 ENDIF
149 IF (HFacS(I,J,K,bi,bj) .NE. 0. _d 0 ) THEN
150 recip_HFacS(I,J,K,bi,bj) = 1. _d 0 / HFacS(I,J,K,bi,bj)
151 maskS(I,J,K,bi,bj) = 1. _d 0
152 ELSE
153 recip_HFacS(I,J,K,bi,bj) = 0. _d 0
154 maskS(I,J,K,bi,bj) = 0. _d 0
155 ENDIF
156 ENDDO
157 ENDDO
158 ENDDO
159 ENDDO
160 ENDDO
161 _EXCH_XYZ_R4(recip_HFacC , myThid )
162 _EXCH_XYZ_R4(recip_HFacW , myThid )
163 _EXCH_XYZ_R4(recip_HFacS , myThid )
164 _EXCH_XYZ_R4(maskW , myThid )
165 _EXCH_XYZ_R4(maskS , myThid )
166
167 C Calculate recipricols grid lengths
168 DO bj = myByLo(myThid), myByHi(myThid)
169 DO bi = myBxLo(myThid), myBxHi(myThid)
170 DO J=1,sNy
171 DO I=1,sNx
172 recip_dxG(I,J,bi,bj)=1.d0/dxG(I,J,bi,bj)
173 recip_dyG(I,J,bi,bj)=1.d0/dyG(I,J,bi,bj)
174 recip_dxC(I,J,bi,bj)=1.d0/dxC(I,J,bi,bj)
175 recip_dyC(I,J,bi,bj)=1.d0/dyC(I,J,bi,bj)
176 recip_dxF(I,J,bi,bj)=1.d0/dxF(I,J,bi,bj)
177 recip_dyF(I,J,bi,bj)=1.d0/dyF(I,J,bi,bj)
178 recip_dxV(I,J,bi,bj)=1.d0/dxV(I,J,bi,bj)
179 recip_dyU(I,J,bi,bj)=1.d0/dyU(I,J,bi,bj)
180 ENDDO
181 ENDDO
182 ENDDO
183 ENDDO
184 _EXCH_XY_R4(recip_dxG, myThid )
185 _EXCH_XY_R4(recip_dyG, myThid )
186 _EXCH_XY_R4(recip_dxC, myThid )
187 _EXCH_XY_R4(recip_dyC, myThid )
188 _EXCH_XY_R4(recip_dxF, myThid )
189 _EXCH_XY_R4(recip_dyF, myThid )
190 _EXCH_XY_R4(recip_dxV, myThid )
191 _EXCH_XY_R4(recip_dyU, myThid )
192
193 #ifdef ALLOW_NONHYDROSTATIC
194 C-- Calculate the reciprocal hfac distance/volume for W cells
195 DO bj = myByLo(myThid), myByHi(myThid)
196 DO bi = myBxLo(myThid), myBxHi(myThid)
197 DO K=1,Nr
198 Km1=max(K-1,1)
199 hFacUpper=drF(Km1)/(drF(Km1)+drF(K))
200 IF (Km1.EQ.K) hFacUpper=0.
201 hFacLower=drF(K)/(drF(Km1)+drF(K))
202 DO J=1,sNy
203 DO I=1,sNx
204 IF (hFacC(I,J,K,bi,bj).NE.0.) THEN
205 IF (hFacC(I,J,K,bi,bj).LE.0.5) THEN
206 recip_hFacU(I,J,K,bi,bj)=
207 & hFacUpper+hFacLower*hFacC(I,J,K,bi,bj)
208 ELSE
209 recip_hFacU(I,J,K,bi,bj)=1.
210 ENDIF
211 ELSE
212 recip_hFacU(I,J,K,bi,bj)=0.
213 ENDIF
214 IF (recip_hFacU(I,J,K,bi,bj).NE.0.)
215 & recip_hFacU(I,J,K,bi,bj)=1./recip_hFacU(I,J,K,bi,bj)
216 ENDDO
217 ENDDO
218 ENDDO
219 ENDDO
220 ENDDO
221 _EXCH_XY_R4(recip_hFacU, myThid )
222 #endif
223 C
224 RETURN
225 END

  ViewVC Help
Powered by ViewVC 1.1.22