/[MITgcm]/MITgcm/pkg/seaice/adi.F
ViewVC logotype

Contents of /MITgcm/pkg/seaice/adi.F

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


Revision 1.12 - (show annotations) (download)
Mon May 3 06:09:39 2004 UTC (20 years, 5 months ago) by dimitri
Branch: MAIN
CVS Tags: HEAD
Changes since 1.11: +1 -1 lines
FILE REMOVED
o removed ADI dynamic solver from pkg/seaice

1 C $Header: /u/gcmpack/MITgcm/pkg/seaice/adi.F,v 1.11 2004/04/28 12:00:53 mlosch Exp $
2 C $Name: $
3
4 #include "SEAICE_OPTIONS.h"
5
6 CStartOfInterface
7 SUBROUTINE adi( myThid )
8 C /==========================================================\
9 C | SUBROUTINE adi |
10 C | o Solve ice momentum equation with an ADI dynamics solver|
11 C | (see Zhang and Rothrock, JGR, 105, 3325-3338, 2000) |
12 C | written by Jinlun Zhang, PSC/UW, Feb-2001 |
13 C | zhang@apl.washington.edu |
14 C |==========================================================|
15 C \==========================================================/
16 IMPLICIT NONE
17
18 C === Global variables ===
19 #include "SIZE.h"
20 #include "EEPARAMS.h"
21 #include "PARAMS.h"
22 #include "SEAICE.h"
23 #include "SEAICE_PARAMS.h"
24 #include "SEAICE_GRID.h"
25
26 C === Routine arguments ===
27 C myThid - Thread no. that called this routine.
28 INTEGER myThid
29 CEndOfInterface
30
31 #ifdef ALLOW_SEAICE
32 #ifdef SEAICE_ALLOW_DYNAMICS
33
34 C === Local variables ===
35 C i,j,bi,bj - Loop counters
36
37 INTEGER i, j, bi, bj, j1, j2, im, jm
38 _RL AA1, AA2, AA3, AA4, AA5, AA6, AA9, RADIUS, RADIUS2
39
40 _RL AU (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
41 _RL BU (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
42 _RL CU (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
43 _RL AV (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
44 _RL BV (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45 _RL CV (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
46 _RL DELXY (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47 _RL DELXR (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48 _RL DELYR (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49 _RL DELX2 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50 _RL DELY2 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51 _RL ETAMEAN (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
52 _RL ZETAMEAN(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
53 _RL FXY (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
54 _RL FXY1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
55
56 _RL URT(1-OLx:sNx+OLx), VRT(1-OLy:sNy+OLy)
57 _RL CUU(1-OLx:sNx+OLx), CVV(1-OLy:sNy+OLy)
58
59 C SET SOME VALUES
60 RADIUS=6370. _d 3
61 RADIUS2=ONE/(RADIUS*RADIUS)
62
63 C SOLVE FOR UICE
64 C FIRST HALF FIRST
65
66 c DO bj=myByLo(myThid),myByHi(myThid)
67 c DO bi=myBxLo(myThid),myBxHi(myThid)
68 c DO j=1,sNy
69 c DO i=1,sNx
70 c FORCEX(I,J,bi,bj)=FORCEX(I,J,bi,bj)*UVM(I,J,bi,bj)
71 c FORCEY(I,J,bi,bj)=FORCEY(I,J,bi,bj)*UVM(I,J,bi,bj)
72 c ENDDO
73 c ENDDO
74 c ENDDO
75 c ENDDO
76
77 C-- Update overlap regions
78 CALL EXCH_UV_XY_RL(FORCEX,FORCEY,.TRUE.,myThid)
79 _EXCH_XY_R8(DRAGS, myThid)
80 _EXCH_XY_R8(DRAGA, myThid)
81 _EXCH_XY_R8(AMASS, myThid)
82
83 c$taf loop = parallel
84 DO bj=myByLo(myThid),myByHi(myThid)
85 c$taf loop = parallel
86 DO bi=myBxLo(myThid),myBxHi(myThid)
87
88 DO J=1-OLy+1,sNy+OLy-1
89 DO I=1-OLx+1,sNx+OLx-1
90 DELXY(I,J)=HALF/(DXUICE(I,J,bi,bj)*DYUICE(I,J,bi,bj))
91 DELXR(I,J)=HALF/(DXUICE(I,J,bi,bj)*RADIUS)
92 DELX2(I,J)=HALF/(DXUICE(I,J,bi,bj)*DXUICE(I,J,bi,bj))
93 ETAMEAN(I,J)=QUART*(ETA(I,J-1,bi,bj)+ETA(I-1,J-1,bi,bj)
94 & +ETA(I,J,bi,bj)+ETA(I-1,J,bi,bj))
95 ZETAMEAN(I,J)=QUART*(ZETA(I,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
96 & +ZETA(I,J,bi,bj)+ZETA(I-1,J,bi,bj))
97
98 FXY(I,J)=DRAGA(I,J,bi,bj)*VICEC(I,J,bi,bj)+FORCEX(I,J,bi,bj)
99 3 +HALF*(ZETA(I,J,bi,bj)*(VICEC(I+1,J+1,bi,bj)
100 3 +VICEC(I,J+1,bi,bj)
101 3 -VICEC(I+1,J,bi,bj)-VICEC(I,J,bi,bj))
102 3 +ZETA(I,J-1,bi,bj)*(VICEC(I+1,J,bi,bj)
103 3 +VICEC(I,J,bi,bj)-VICEC(I+1,J-1,bi,bj)
104 3 -VICEC(I,J-1,bi,bj))+ZETA(I-1,J,bi,bj)
105 3 *(VICEC(I,J,bi,bj)+VICEC(I-1,J,bi,bj)
106 3 -VICEC(I,J+1,bi,bj)-VICEC(I-1,J+1,bi,bj))
107 3 +ZETA(I-1,J-1,bi,bj)*(VICEC(I,J-1,bi,bj)
108 3 +VICEC(I-1,J-1,bi,bj)-VICEC(I,J,bi,bj)
109 3 -VICEC(I-1,J,bi,bj)))*DELXY(I,J)
110 3 *RECIP_CSUICE(I,J,bi,bj)
111 3
112 4 -HALF*(ETA(I,J,bi,bj)*(VICEC(I+1,J+1,bi,bj)
113 4 +VICEC(I,J+1,bi,bj)
114 4 -VICEC(I+1,J,bi,bj)-VICEC(I,J,bi,bj))
115 4 +ETA(I,J-1,bi,bj)*(VICEC(I+1,J,bi,bj)
116 4 +VICEC(I,J,bi,bj)-VICEC(I+1,J-1,bi,bj)
117 4 -VICEC(I,J-1,bi,bj))+ETA(I-1,J,bi,bj)
118 4 *(VICEC(I,J,bi,bj)+VICEC(I-1,J,bi,bj)
119 4 -VICEC(I,J+1,bi,bj)-VICEC(I-1,J+1,bi,bj))
120 4 +ETA(I-1,J-1,bi,bj)*(VICEC(I,J-1,bi,bj)
121 4 +VICEC(I-1,J-1,bi,bj)-VICEC(I,J,bi,bj)
122 4 -VICEC(I-1,J,bi,bj)))*DELXY(I,J)
123 4 *RECIP_CSUICE(I,J,bi,bj)
124 4
125 5 +HALF*(VICEC(I+1,J,bi,bj)-VICEC(I-1,J,bi,bj))
126 5 *(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
127 5 -ETA(I-1,J-1,bi,bj)-ETA(I,J-1,bi,bj))*DELXY(I,J)
128 5 *RECIP_CSUICE(I,J,bi,bj)+HALF*ETAMEAN(I,J)
129 5 *((VICEC(I+1,J+1,bi,bj)
130 5 -VICEC(I-1,J+1,bi,bj))*RECIP_CSUICE(I,J+1,bi,bj)
131 5 -(VICEC(I+1,J-1,bi,bj)-VICEC(I-1,J-1,bi,bj))
132 5 *RECIP_CSUICE(I,J-1,bi,bj))*DELXY(I,J)
133 5
134 6 -((ZETA(I,J,bi,bj)+ZETA(I,J-1,bi,bj)
135 6 -ZETA(I-1,J-1,bi,bj)-ZETA(I-1,J,bi,bj))
136 6 +(ETA(I,J,bi,bj)+ETA(I,J-1,bi,bj)
137 6 -ETA(I-1,J-1,bi,bj)-ETA(I-1,J,bi,bj)))
138 6 *TNGICE(I,J,bi,bj)*VICEC(I,J,bi,bj)
139 & *DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
140 6 -(ETAMEAN(I,J)+ZETAMEAN(I,J))*TNGICE(I,J,bi,bj)
141 6 *(VICEC(I+1,J,bi,bj)-VICEC(I-1,J,bi,bj))
142 6 *DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
143 6
144 7 -ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)
145 7 *(VICEC(I+1,J,bi,bj)-VICEC(I-1,J,bi,bj))
146 7 *DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
147 END DO
148 END DO
149
150 DO J=1-OLy+1,sNy+OLy-1
151 DO I=1-OLx+1,sNx+OLx-1
152 DELY2(I,J)=HALF/(DYUICE(I,J,bi,bj)*DYUICE(I,J,bi,bj))
153 DELYR(I,J)=HALF/(DYUICE(I,J,bi,bj)*RADIUS)
154 AA1=((ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj))
155 & *RECIP_CSUICE(I,J,bi,bj)
156 & +(ETA(I,J,bi,bj)+ZETA(I,J,bi,bj))
157 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
158 AA2=((ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj))
159 & *RECIP_CSUICE(I,J,bi,bj)
160 & +(ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj))
161 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
162 AA3=ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
163 AA4=ETA(I-1,J-1,bi,bj)+ETA(I,J-1,bi,bj)
164 AA5=-(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)-ETA(I-1,J-1,bi,bj)
165 & -ETA(I,J-1,bi,bj))*TNGICE(I,J,bi,bj)
166 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
167 AU(I,J)=-AA2*DELX2(I,J)*UVM(I,J,bi,bj)
168 BU(I,J)=((AA1+AA2)*DELX2(I,J)+AA6*RADIUS2
169 & +AMASS(I,J,bi,bj)/SEAICE_DT*TWO
170 & +DRAGS(I,J,bi,bj))*UVM(I,J,bi,bj)
171 & +(ONE-UVM(I,J,bi,bj))
172 CU(I,J)=-AA1*DELX2(I,J)*UVM(I,J,bi,bj)
173 END DO
174 END DO
175
176 DO J=1-OLy+1,sNy+OLy-1
177 AU(1-OLx+1,J)=ZERO
178 CU(sNx+OLx-1,J)=ZERO
179 CU(1-OLx+1,J)=CU(1-OLx+1,J)/BU(1-OLx+1,J)
180 END DO
181
182 c$taf loop = parallel
183 DO 1200 J=1-OLy+1,sNy+OLy-1
184 DO I=1-OLx+1,sNx+OLx-1
185
186 AA1=((ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj))
187 & *RECIP_CSUICE(I,J,bi,bj)
188 & +(ETA(I,J,bi,bj)+ZETA(I,J,bi,bj))
189 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
190 AA2=((ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj))
191 & *RECIP_CSUICE(I,J,bi,bj)
192 & +(ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj))
193 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
194 AA3=ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
195 AA4=ETA(I-1,J-1,bi,bj)+ETA(I,J-1,bi,bj)
196 AA5=-(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)-ETA(I-1,J-1,bi,bj)
197 & -ETA(I,J-1,bi,bj))*TNGICE(I,J,bi,bj)
198 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
199
200 IF(I.EQ.1-OLx+1) THEN
201 AA9=AA2*DELX2(I,J)*UICEC(I-1,J,bi,bj)*UVM(I,J,bi,bj)
202 ELSE IF(I.EQ.sNx+OLx-1) THEN
203 AA9=AA1*DELX2(I,J)*UICEC(I+1,J,bi,bj)*UVM(I,J,bi,bj)
204 ELSE
205 AA9=ZERO
206 END IF
207
208 URT(I)=AA9+FXY(I,J)-AA5*DELYR(I,J)*UICE(I,J,2,bi,bj)
209 1 -(AA3+AA4)*DELY2(I,J)*UICE(I,J,2,bi,bj)
210 1 +(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj))
211 1 *UICE(I,J+1,2,bi,bj)*DELY2(I,J)
212 2 +(ETA(I-1,J-1,bi,bj)+ETA(I,J-1,bi,bj))
213 2 *UICE(I,J-1,2,bi,bj)*DELY2(I,J)
214 3 +ETAMEAN(I,J)*DELYR(I,J)*(UICE(I,J+1,2,bi,bj)
215 3 *TNGICE(I,J+1,bi,bj)
216 3 -UICE(I,J-1,2,bi,bj)*TNGICE(I,J-1,bi,bj))
217 4 -ETAMEAN(I,J)*DELYR(I,J)*TWO*TNGICE(I,J,bi,bj)
218 4 *(UICE(I,J+1,2,bi,bj)
219 4 -UICE(I,J-1,2,bi,bj))
220 URT(I)=(URT(I)+AMASS(I,J,bi,bj)/SEAICE_DT
221 & *UICE(I,J,2,bi,bj)*TWO)*UVM(I,J,bi,bj)
222 END DO
223
224 DO I=1-OLx+1,sNx+OLx-1
225 CUU(I)=CU(I,J)
226 END DO
227 URT(1-OLx+1)=URT(1-OLx+1)/BU(1-OLx+1,J)
228 DO I=1-OLx+2,sNx+OLx-1
229 IM=I-1
230 CUU(I)=CUU(I)/(BU(I,J)-AU(I,J)*CUU(IM))
231 URT(I)=(URT(I)-AU(I,J)*URT(IM))/(BU(I,J)-AU(I,J)*CUU(IM))
232 END DO
233 DO I=1-OLx+1,sNx+OLx-1-1
234 J1=sNx+OLx-1-I
235 J2=J1+1
236 URT(J1)=URT(J1)-CUU(J1)*URT(J2)
237 END DO
238 DO I=1-OLx+1,sNx+OLx-1
239 UICE(I,J,1,bi,bj)=URT(I)
240 END DO
241 1200 CONTINUE
242
243 c DO J=1,sNy
244 c DO I=1,sNx
245 c UICE(I,J,3,bi,bj)=UICE(I,J,1,bi,bj)
246 c END DO
247 c END DO
248
249 ENDDO
250 ENDDO
251
252 CALL SEAICE_EXCH( UICE, myThid )
253
254 C NOW SECOND HALF
255
256 c$taf loop = parallel
257 DO bj=myByLo(myThid),myByHi(myThid)
258 c$taf loop = parallel
259 DO bi=myBxLo(myThid),myBxHi(myThid)
260
261 DO I=1-OLx+1,sNx+OLx-1
262 DO J=1-OLy+1,sNy+OLy-1
263
264 AA1=ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
265 AA2=ETA(I-1,J-1,bi,bj)+ETA(I,J-1,bi,bj)
266 AA5=-(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
267 & -ETA(I-1,J-1,bi,bj)-ETA(I,J-1,bi,bj))*TNGICE(I,J,bi,bj)
268 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
269
270 AV(I,J)=(-AA2*DELY2(I,J)+ETAMEAN(I,J)*DELYR(I,J)
271 & *(TNGICE(I,J-1,bi,bj)
272 & -TWO*TNGICE(I,J,bi,bj)))*UVM(I,J,bi,bj)
273 BV(I,J)=((AA1+AA2)*DELY2(I,J)+AA5*DELYR(I,J)+AA6*RADIUS2
274 & +AMASS(I,J,bi,bj)/SEAICE_DT*TWO
275 & +DRAGS(I,J,bi,bj))*UVM(I,J,bi,bj)
276 & +(ONE-UVM(I,J,bi,bj))
277 CV(I,J)=(-AA1*DELY2(I,J)-ETAMEAN(I,J)*DELYR(I,J)
278 & *(TNGICE(I,J+1,bi,bj)
279 & -TWO*TNGICE(I,J,bi,bj)))*UVM(I,J,bi,bj)
280 END DO
281 END DO
282
283 DO I=1-OLx+1,sNx+OLx-1
284 AV(I,1-OLy+1)=ZERO
285 CV(I,sNy+OLy-1)=ZERO
286 CV(I,1-OLy+1)=CV(I,1-OLy+1)/BV(I,1-OLy+1)
287 END DO
288
289 DO I=1-OLx+1,sNx+OLx-1
290 DO J=1-OLy+1,sNy+OLy-1
291
292 AA1=((ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj))
293 & *RECIP_CSUICE(I,J,bi,bj)
294 & +(ETA(I,J,bi,bj)+ZETA(I,J,bi,bj))
295 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
296 AA2=((ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj))
297 & *RECIP_CSUICE(I,J,bi,bj)
298 & +(ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj))
299 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
300
301 AA3=ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
302 AA4=ETA(I-1,J-1,bi,bj)+ETA(I,J-1,bi,bj)
303
304 IF(J.EQ.1-OLy+1) THEN
305 AA9=( AA4*DELY2(I,J)*UICEC(I,J-1,bi,bj)
306 & -ETAMEAN(I,J)*DELYR(I,J)*(TNGICE(I,J-1,bi,bj)
307 & -TWO*TNGICE(I,J,bi,bj))
308 & *UICEC(I,J-1,bi,bj) )*UVM(I,J,bi,bj)
309 ELSE IF(J.EQ.sNy+OLy-1) THEN
310 AA9=( AA3*DELY2(I,J)*UICEC(I,J+1,bi,bj)
311 & +ETAMEAN(I,J)*DELYR(I,J)*(TNGICE(I,J+1,bi,bj)
312 & -TWO*TNGICE(I,J,bi,bj))
313 & *UICEC(I,J+1,bi,bj) )*UVM(I,J,bi,bj)
314 ELSE
315 AA9=ZERO
316 END IF
317
318 FXY1(I,J)=AA9+AMASS(I,J,bi,bj)/SEAICE_DT*UICE(I,J,1,bi,bj)*TWO
319 5 -(AA1+AA2)*DELX2(I,J)*UICE(I,J,1,bi,bj)
320 6 +((ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj)
321 6 +ETA(I,J,bi,bj)+ZETA(I,J,bi,bj))
322 6 *UICE(I+1,J,1,bi,bj)
323 6 +(ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
324 6 + ETA(I-1,J ,bi,bj)+ZETA(I-1,J ,bi,bj))
325 6 *UICE(I-1,J,1,bi,bj))*DELX2(I,J)
326 6 *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
327
328 END DO
329 END DO
330
331 DO 1300 I=1-OLx+1,sNx+OLx-1
332 DO J=1-OLy+1,sNy+OLy-1
333 VRT(J)=FXY(I,J)+FXY1(I,J)
334 VRT(J)=VRT(J)*UVM(I,J,bi,bj)
335 END DO
336
337 DO J=1-OLy+1,sNy+OLy-1
338 CVV(J)=CV(I,J)
339 END DO
340 VRT(1-OLy+1)=VRT(1-OLy+1)/BV(I,1-OLy+1)
341 DO J=1-OLy+2,sNy+OLy-1
342 JM=J-1
343 CVV(J)=CVV(J)/(BV(I,J)-AV(I,J)*CVV(JM))
344 VRT(J)=(VRT(J)-AV(I,J)*VRT(JM))/(BV(I,J)-AV(I,J)*CVV(JM))
345 END DO
346 DO J=1-OLy+1,sNy+OLy-1-1
347 J1=sNy+OLy-1-J
348 J2=J1+1
349 VRT(J1)=VRT(J1)-CVV(J1)*VRT(J2)
350 END DO
351 DO J=1-OLy+1,sNy+OLy-1
352 UICE(I,J,1,bi,bj)=VRT(J)
353 END DO
354 1300 CONTINUE
355
356 ENDDO
357 ENDDO
358
359 C SOLVE FOR VICE
360 C FIRST HALF FIRST
361
362 c$taf loop = parallel
363 DO bj=myByLo(myThid),myByHi(myThid)
364 c$taf loop = parallel
365 DO bi=myBxLo(myThid),myBxHi(myThid)
366
367 DO I=1-OLx+1,sNx+OLx-1
368 DO J=1-OLy+1,sNy+OLy-1
369
370 FXY(I,J)=-DRAGA(I,J,bi,bj)*UICEC(I,J,bi,bj)+FORCEY(I,J,bi,bj)
371 3 +(HALF*(UICEC(I+1,J,bi,bj)-UICEC(I-1,J,bi,bj))
372 3 *(ZETA(I-1,J,bi,bj)+ZETA(I,J,bi,bj)
373 3 -ZETA(I-1,J-1,bi,bj)-ZETA(I,J-1,bi,bj))*DELXY(I,J)
374 3 *RECIP_CSUICE(I,J,bi,bj)+HALF*ZETAMEAN(I,J)
375 3 *((UICEC(I+1,J+1,bi,bj)
376 3 -UICEC(I-1,J+1,bi,bj))*RECIP_CSUICE(I,J+1,bi,bj)
377 3 -(UICEC(I+1,J-1,bi,bj)-UICEC(I-1,J-1,bi,bj))
378 3 *RECIP_CSUICE(I,J-1,bi,bj))*DELXY(I,J))
379 3
380 4 -(HALF*(UICEC(I+1,J,bi,bj)-UICEC(I-1,J,bi,bj))
381 4 *(ETA(I-1,J,bi,bj)+ETA(I,J,bi,bj)
382 4 -ETA(I-1,J-1,bi,bj)-ETA(I,J-1,bi,bj))*DELXY(I,J)
383 4 *RECIP_CSUICE(I,J,bi,bj)+HALF*ETAMEAN(I,J)
384 4 *((UICEC(I+1,J+1,bi,bj)
385 4 -UICEC(I-1,J+1,bi,bj))*RECIP_CSUICE(I,J+1,bi,bj)
386 4 -(UICEC(I+1,J-1,bi,bj)-UICEC(I-1,J-1,bi,bj))
387 4 *RECIP_CSUICE(I,J-1,bi,bj))*DELXY(I,J))
388 4
389 5 +HALF*(ETA(I,J,bi,bj)*(UICEC(I+1,J+1,bi,bj)
390 5 +UICEC(I,J+1,bi,bj)
391 5 -UICEC(I+1,J,bi,bj)-UICEC(I,J,bi,bj))+ETA(I,J-1,bi,bj)
392 5 *(UICEC(I+1,J,bi,bj)
393 5 +UICEC(I,J,bi,bj)-UICEC(I+1,J-1,bi,bj)-UICEC(I,J-1,bi,bj))
394 5 +ETA(I-1,J,bi,bj)
395 5 *(UICEC(I,J,bi,bj)+UICEC(I-1,J,bi,bj)-UICEC(I,J+1,bi,bj)
396 5 -UICEC(I-1,J+1,bi,bj))
397 5 +ETA(I-1,J-1,bi,bj)*(UICEC(I,J-1,bi,bj)
398 5 +UICEC(I-1,J-1,bi,bj)
399 5 -UICEC(I,J,bi,bj)
400 5 -UICEC(I-1,J,bi,bj)))*DELXY(I,J)*RECIP_CSUICE(I,J,bi,bj)
401 5
402 6 +(ETA(I,J,bi,bj)+ETA(I,J-1,bi,bj)
403 6 -ETA(I-1,J-1,bi,bj)-ETA(I-1,J,bi,bj))
404 6 *TNGICE(I,J,bi,bj)*UICEC(I,J,bi,bj)
405 6 *DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
406 6 +ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*(UICEC(I+1,J,bi,bj)
407 6 -UICEC(I-1,J,bi,bj))*DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
408 6
409 7 +ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)*(UICEC(I+1,J,bi,bj)
410 7 -UICEC(I-1,J,bi,bj))*DELXR(I,J)*RECIP_CSUICE(I,J,bi,bj)
411
412 END DO
413 END DO
414
415 DO I=1-OLx+1,sNx+OLx-1
416 DO J=1-OLy+1,sNy+OLy-1
417 AA1=ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj)
418 & +ETA(I,J,bi,bj)+ZETA(I,J,bi,bj)
419 AA2=ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
420 & +ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj)
421 AA3=(ETA(I,J-1,bi,bj)*RECIP_CSUICE(I,J,bi,bj)+ETA(I,J,bi,bj)
422 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
423 AA4=(ETA(I-1,J-1,bi,bj)+ETA(I-1,J,bi,bj))
424 & *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
425 AA5=((ZETA(I-1,J,bi,bj)-ETA(I-1,J,bi,bj))
426 & +(ZETA(I,J,bi,bj)-ETA(I,J,bi,bj))
427 & -(ZETA(I-1,J-1,bi,bj)-ETA(I-1,J-1,bi,bj))
428 & -(ZETA(I,J-1,bi,bj)
429 & -ETA(I,J-1,bi,bj)))*TNGICE(I,J,bi,bj)
430 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
431
432 AV(I,J)=(-AA2*DELY2(I,J)-(ZETAMEAN(I,J)-ETAMEAN(I,J))
433 & *TNGICE(I,J-1,bi,bj)
434 & *DELYR(I,J)-ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)
435 & *DELYR(I,J))*UVM(I,J,bi,bj)
436 BV(I,J)=((AA1+AA2)*DELY2(I,J)+AA5*DELYR(I,J)+AA6*RADIUS2
437 & +AMASS(I,J,bi,bj)/SEAICE_DT*TWO
438 & +DRAGS(I,J,bi,bj))*UVM(I,J,bi,bj)
439 & +(ONE-UVM(I,J,bi,bj))
440 CV(I,J)=(-AA1*DELY2(I,J)+(ZETAMEAN(I,J)-ETAMEAN(I,J))
441 & *TNGICE(I,J+1,bi,bj)
442 & *DELYR(I,J)+ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)
443 & *DELYR(I,J))*UVM(I,J,bi,bj)
444 END DO
445 END DO
446 DO I=1-OLx+1,sNx+OLx-1
447 AV(I,1-OLy+1)=ZERO
448 CV(I,sNy+OLy-1)=ZERO
449 CV(I,1-OLy+1)=CV(I,1-OLy+1)/BV(I,1-OLy+1)
450 END DO
451
452 c$taf loop = parallel
453 DO 1301 I=1-OLx+1,sNx+OLx-1
454 DO J=1-OLy+1,sNy+OLy-1
455
456 AA1=ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj)
457 & +ETA(I,J,bi,bj)+ZETA(I,J,bi,bj)
458 AA2=ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
459 & +ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj)
460 AA3=(ETA(I,J-1,bi,bj)*RECIP_CSUICE(I,J,bi,bj)+ETA(I,J,bi,bj)
461 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
462 AA4=(ETA(I-1,J-1,bi,bj)+ETA(I-1,J,bi,bj))
463 & *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
464 AA5=((ZETA(I-1,J,bi,bj)-ETA(I-1,J,bi,bj))
465 & +(ZETA(I,J,bi,bj)-ETA(I,J,bi,bj))
466 & -(ZETA(I-1,J-1,bi,bj)-ETA(I-1,J-1,bi,bj))
467 & -(ZETA(I,J-1,bi,bj)-ETA(I,J-1,bi,bj)))*TNGICE(I,J,bi,bj)
468 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
469
470 IF(J.EQ.1-OLy+1) THEN
471 AA9=(AA2*DELY2(I,J)+(ZETAMEAN(I,J)-ETAMEAN(I,J))
472 & *TNGICE(I,J-1,bi,bj)*DELYR(I,J)
473 & +ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)*DELYR(I,J))
474 & *VICEC(I,J-1,bi,bj)*UVM(I,J,bi,bj)
475 ELSE IF(J.EQ.sNy+OLy-1) THEN
476 AA9=(AA1*DELY2(I,J)-(ZETAMEAN(I,J)-ETAMEAN(I,J))
477 & *TNGICE(I,J+1,bi,bj)*DELYR(I,J)
478 & -ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)*DELYR(I,J))
479 & *VICEC(I,J+1,bi,bj)*UVM(I,J,bi,bj)
480 ELSE
481 AA9=ZERO
482 END IF
483
484 VRT(J)=AA9+FXY(I,J)-(AA3+AA4)*DELX2(I,J)*VICE(I,J,2,bi,bj)
485 6 +((ETA(I,J-1,bi,bj)*RECIP_CSUICE(I,J,bi,bj)+ETA(I,J,bi,bj)
486 6 *RECIP_CSUICE(I,J,bi,bj))*VICE(I+1,J,2,bi,bj)*DELX2(I,J)
487 7 +(ETA(I-1,J-1,bi,bj)+ETA(I-1,J,bi,bj))
488 7 *VICE(I-1,J,2,bi,bj)*DELX2(I,J))
489 7 *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
490 VRT(J)=(VRT(J)+AMASS(I,J,bi,bj)/SEAICE_DT
491 & *VICE(I,J,2,bi,bj)*TWO)
492 & *UVM(I,J,bi,bj)
493 END DO
494
495 DO J=1-OLy+1,sNy+OLy-1
496 CVV(J)=CV(I,J)
497 END DO
498 VRT(1-OLy+1)=VRT(1-OLy+1)/BV(I,1-OLy+1)
499 DO J=1-OLy+2,sNy+OLy-1
500 JM=J-1
501 CVV(J)=CVV(J)/(BV(I,J)-AV(I,J)*CVV(JM))
502 VRT(J)=(VRT(J)-AV(I,J)*VRT(JM))/(BV(I,J)-AV(I,J)*CVV(JM))
503 END DO
504 DO J=1-OLy+1,sNy+OLy-1-1
505 J1=sNy+OLy-1-J
506 J2=J1+1
507 VRT(J1)=VRT(J1)-CVV(J1)*VRT(J2)
508 END DO
509 DO J=1-OLy+1,sNy+OLy-1
510 VICE(I,J,1,bi,bj)=VRT(J)
511 END DO
512 1301 CONTINUE
513
514 c DO J=1,sNy
515 c DO I=1,sNx
516 c VICE(I,J,3,bi,bj)=VICE(I,J,1,bi,bj)
517 c END DO
518 c END DO
519
520 ENDDO
521 ENDDO
522
523 CALL SEAICE_EXCH( VICE, myThid )
524
525 C NOW SECOND HALF
526
527 c$taf loop = parallel
528 DO bj=myByLo(myThid),myByHi(myThid)
529 c$taf loop = parallel
530 DO bi=myBxLo(myThid),myBxHi(myThid)
531
532 DO J=1-OLy+1,sNy+OLy-1
533 DO I=1-OLx+1,sNx+OLx-1
534
535 AA1=ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj)
536 & +ETA(I,J,bi,bj)+ZETA(I,J,bi,bj)
537 AA2=ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
538 & +ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj)
539 AA3=(ETA(I,J-1,bi,bj)*RECIP_CSUICE(I,J,bi,bj)+ETA(I,J,bi,bj)
540 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
541 AA4=(ETA(I-1,J-1,bi,bj)+ETA(I-1,J,bi,bj))
542 & *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
543 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
544
545 AU(I,J)=-AA4*DELX2(I,J)*UVM(I,J,bi,bj)
546 BU(I,J)=((AA3+AA4)*DELX2(I,J)+AA6*RADIUS2
547 & +AMASS(I,J,bi,bj)/SEAICE_DT*TWO
548 & +DRAGS(I,J,bi,bj))*UVM(I,J,bi,bj)
549 & +(ONE-UVM(I,J,bi,bj))
550 CU(I,J)=-AA3*DELX2(I,J)*UVM(I,J,bi,bj)
551 END DO
552 END DO
553 DO J=1-OLy+1,sNy+OLy-1
554 AU(1-OLx+1,J)=ZERO
555 CU(sNx+OLx-1,J)=ZERO
556 CU(1-OLx+1,J)=CU(1-OLx+1,J)/BU(1-OLx+1,J)
557 END DO
558
559 DO J=1-OLy+1,sNy+OLy-1
560 DO I=1-OLx+1,sNx+OLx-1
561
562 AA1=ETA(I-1,J,bi,bj)+ZETA(I-1,J,bi,bj)
563 & +ETA(I,J,bi,bj)+ZETA(I,J,bi,bj)
564 AA2=ETA(I-1,J-1,bi,bj)+ZETA(I-1,J-1,bi,bj)
565 & +ETA(I,J-1,bi,bj)+ZETA(I,J-1,bi,bj)
566 AA3=(ETA(I,J-1,bi,bj)*RECIP_CSUICE(I,J,bi,bj)+ETA(I,J,bi,bj)
567 & *RECIP_CSUICE(I,J,bi,bj))*RECIP_CSUICE(I,J,bi,bj)
568 AA4=(ETA(I-1,J-1,bi,bj)+ETA(I-1,J,bi,bj))
569 & *RECIP_CSUICE(I,J,bi,bj)*RECIP_CSUICE(I,J,bi,bj)
570 AA5=((ZETA(I-1,J,bi,bj)-ETA(I-1,J,bi,bj))
571 & +(ZETA(I,J,bi,bj)-ETA(I,J,bi,bj))
572 & -(ZETA(I-1,J-1,bi,bj)-ETA(I-1,J-1,bi,bj))
573 & -(ZETA(I,J-1,bi,bj)
574 & -ETA(I,J-1,bi,bj)))*TNGICE(I,J,bi,bj)
575 AA6=TWO*ETAMEAN(I,J)*TNGICE(I,J,bi,bj)*TNGICE(I,J,bi,bj)
576
577 IF(I.EQ.1-OLx+1) THEN
578 AA9=AA4*DELX2(I,J)*VICEC(I-1,J,bi,bj)*UVM(I,J,bi,bj)
579 ELSE IF(I.EQ.sNx+OLx-1) THEN
580 AA9=AA3*DELX2(I,J)*VICEC(I+1,J,bi,bj)*UVM(I,J,bi,bj)
581 ELSE
582 AA9=ZERO
583 END IF
584
585 FXY1(I,J)=AA9+AMASS(I,J,bi,bj)/SEAICE_DT
586 1 *VICE(I,J,1,bi,bj)*TWO
587 1 -AA5*DELYR(I,J)*VICE(I,J,1,bi,bj)
588 1 -(AA1+AA2)*DELY2(I,J)*VICE(I,J,1,bi,bj)
589 1 +AA1*DELY2(I,J)*VICE(I,J+1,1,bi,bj)
590 1 -((ZETAMEAN(I,J)-ETAMEAN(I,J))
591 1 *TNGICE(I,J+1,bi,bj)*DELYR(I,J)
592 1 +ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)*DELYR(I,J))
593 1 *VICE(I,J+1,1,bi,bj)
594 2 +AA2*DELY2(I,J)*VICE(I,J-1,1,bi,bj)
595 2 +((ZETAMEAN(I,J)-ETAMEAN(I,J))
596 2 *TNGICE(I,J-1,bi,bj)*DELYR(I,J)
597 2 +ETAMEAN(I,J)*TWO*TNGICE(I,J,bi,bj)*DELYR(I,J))
598 2 *VICE(I,J-1,1,bi,bj)
599 END DO
600 END DO
601
602 DO 1201 J=1-OLy+1,sNy+OLy-1
603 DO I=1-OLx+1,sNx+OLx-1
604 URT(I)=FXY(I,J)+FXY1(I,J)
605 URT(I)=URT(I)*UVM(I,J,bi,bj)
606 END DO
607
608 DO I=1-OLx+1,sNx+OLx-1
609 CUU(I)=CU(I,J)
610 END DO
611 URT(1-OLx+1)=URT(1-OLx+1)/BU(1-OLx+1,J)
612 DO I=1-OLx+2,sNx+OLx-1
613 IM=I-1
614 CUU(I)=CUU(I)/(BU(I,J)-AU(I,J)*CUU(IM))
615 URT(I)=(URT(I)-AU(I,J)*URT(IM))/(BU(I,J)-AU(I,J)*CUU(IM))
616 END DO
617 DO I=1-OLx+1,sNx+OLx-1-1
618 J1=sNx+OLx-1-I
619 J2=J1+1
620 URT(J1)=URT(J1)-CUU(J1)*URT(J2)
621 END DO
622 DO I=1-OLx+1,sNx+OLx-1
623 VICE(I,J,1,bi,bj)=URT(I)
624 END DO
625 1201 CONTINUE
626
627 ENDDO
628 ENDDO
629
630 DO bj=myByLo(myThid),myByHi(myThid)
631 DO bi=myBxLo(myThid),myBxHi(myThid)
632 DO J=1,sNy
633 DO I=1,sNx
634 UICEC(I,J,bi,bj)=UICE(I,J,1,bi,bj)*UVM(I,J,bi,bj)
635 VICEC(I,J,bi,bj)=VICE(I,J,1,bi,bj)*UVM(I,J,bi,bj)
636 END DO
637 END DO
638 ENDDO
639 ENDDO
640
641 C-- Update overlap regions
642 CALL EXCH_UV_XY_RL(UICEC,VICEC,.TRUE.,myThid)
643
644 DO bj=myByLo(myThid),myByHi(myThid)
645 DO bi=myBxLo(myThid),myBxHi(myThid)
646 DO j=1-OLy,sNy+OLy
647 DO i=1-OLx,sNx+OLx
648 UICE(I,J,1,bi,bj)=UICEC(I,J,bi,bj)
649 VICE(I,J,1,bi,bj)=VICEC(I,J,bi,bj)
650 ENDDO
651 ENDDO
652 ENDDO
653 ENDDO
654
655 #endif /* SEAICE_ALLOW_DYNAMICS */
656 #endif /* ALLOW_SEAICE */
657
658 RETURN
659 END

  ViewVC Help
Powered by ViewVC 1.1.22