/[MITgcm]/MITgcm/pkg/exch2/exch2_uv_cgrid_3d_rx.template
ViewVC logotype

Annotation of /MITgcm/pkg/exch2/exch2_uv_cgrid_3d_rx.template

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


Revision 1.4 - (hide annotations) (download)
Sun Jun 28 00:57:51 2009 UTC (14 years, 11 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61s, checkpoint61z, checkpoint61x, checkpoint61y
Changes since 1.3: +4 -22 lines
-always call exch2_*_cube, not exch-1 anymore, if useCubedSphereExchange=F
-add bj in exch2 arrays and S/R.

1 jmc 1.4 C $Header: /u/gcmpack/MITgcm/pkg/exch2/exch2_uv_cgrid_3d_rx.template,v 1.3 2009/05/12 19:44:58 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "CPP_EEOPTIONS.h"
5     #include "W2_OPTIONS.h"
6    
7     CBOP
8     C !ROUTINE: EXCH2_UV_CGRID_3D_RX
9    
10     C !INTERFACE:
11     SUBROUTINE EXCH2_UV_CGRID_3D_RX(
12     U uPhi, vPhi,
13     I withSigns, myNz, myThid )
14    
15     C !DESCRIPTION:
16     C*=====================================================================*
17     C Purpose: SUBROUTINE EXCH2_UV_CGRID_3D_RX
18     C handle exchanges for a 3D vector field on a C-grid.
19     C
20     C Input:
21     C uPhi(lon,lat,levs,bi,bj) :: first component of vector
22     C vPhi(lon,lat,levs,bi,bj) :: second component of vector
23     C withSigns (logical) :: true to use sign of components
24     C myNz :: 3rd dimension of input arrays uPhi,vPhi
25     C myThid :: my Thread Id number
26     C
27     C Output: uPhi and vPhi are updated (halo regions filled)
28     C
29     C Calls: exch_RX (exch2_RX1_cube) - for each component
30     C
31     C*=====================================================================*
32    
33     C !USES:
34     IMPLICIT NONE
35    
36     #include "SIZE.h"
37     #include "EEPARAMS.h"
38 jmc 1.3 #include "W2_EXCH2_SIZE.h"
39 jmc 1.1 #include "W2_EXCH2_TOPOLOGY.h"
40 jmc 1.3 #ifdef W2_FILL_NULL_REGIONS
41 jmc 1.1 #include "W2_EXCH2_PARAMS.h"
42 jmc 1.3 #endif
43 jmc 1.1
44     C !INPUT/OUTPUT PARAMETERS:
45     C == Argument list variables ==
46     INTEGER myNz
47     _RX uPhi(1-OLx:sNx+OLx,1-OLy:sNy+OLy,myNz,nSx,nSy)
48     _RX vPhi(1-OLx:sNx+OLx,1-OLy:sNy+OLy,myNz,nSx,nSy)
49     LOGICAL withSigns
50     INTEGER myThid
51    
52     C !LOCAL VARIABLES:
53     C == Local variables ==
54     C i,j,k,bi,bj :: loop indices.
55     C OL[wens] :: Overlap extents in west, east, north, south.
56     C exchWidth[XY] :: Extent of regions that will be exchanged.
57     C uLoc,vLoc :: local copy of the vector components with haloes filled.
58    
59     INTEGER i,j,k,bi,bj
60     INTEGER OLw, OLe, OLn, OLs, exchWidthX, exchWidthY
61     _RX uLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
62     _RX vLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
63     _RX negOne
64     INTEGER myTile, myFace
65     CEOP
66    
67     OLw = OLx
68     OLe = OLx
69     OLn = OLy
70     OLs = OLy
71     exchWidthX = OLx
72     exchWidthY = OLy
73     negOne = 1.
74     IF (withSigns) negOne = -1.
75    
76     C-- First call the exchanges for the two components
77    
78     CALL EXCH2_RX1_CUBE( uPhi, 'T ',
79     I OLw, OLe, OLs, OLn, myNz,
80     I exchWidthX, exchWidthY,
81     I FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
82     CALL EXCH2_RX1_CUBE( uPhi, 'T ',
83     I OLw, OLe, OLs, OLn, myNz,
84     I exchWidthX, exchWidthY,
85     I FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
86    
87     CALL EXCH2_RX1_CUBE( vPhi, 'T ',
88     I OLw, OLe, OLs, OLn, myNz,
89     I exchWidthX, exchWidthY,
90     I FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
91     CALL EXCH2_RX1_CUBE( vPhi, 'T ',
92     I OLw, OLe, OLs, OLn, myNz,
93     I exchWidthX, exchWidthY,
94     I FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
95    
96     C- note: can substitute the low-level S/R calls above with:
97     c CALL EXCH2_3D_RX( uPhi, myNz, myThid )
98     c CALL EXCH2_3D_RX( vPhi, myNz, myThid )
99    
100 jmc 1.4 IF ( useCubedSphereExchange ) THEN
101    
102 jmc 1.1 C-- Then, depending on which tile we are, we may need
103     C 1) to switch u and v components and also to switch the signs
104     C 2) to shift the index along the face edge.
105     C 3) ensure that near-corner halo regions is filled in the correct order
106     C (i.e. with velocity component already available after 1 exch)
107     C- note: because of index shift, the order really matter:
108     C odd faces, do North 1rst and then West;
109     C even faces, do East 1rst and then South.
110    
111     C-- Loops on tile indices:
112     DO bj = myByLo(myThid), myByHi(myThid)
113     DO bi = myBxLo(myThid), myBxHi(myThid)
114    
115     C- Now choose what to do at each edge of the halo based on which face
116     C (we assume that bj is always=1)
117 jmc 1.4 myTile = W2_myTileList(bi,bj)
118 jmc 1.1 myFace = exch2_myFace(myTile)
119    
120     C-- Loops on level index:
121     DO k = 1,myNz
122    
123     C- First we copy the 2 components info into local dummy arrays uLoc,vLoc
124     DO j = 1-OLy,sNy+OLy
125     DO i = 1-OLx,sNx+OLx
126     uLoc(i,j) = uPhi(i,j,k,bi,bj)
127     vLoc(i,j) = vPhi(i,j,k,bi,bj)
128     ENDDO
129     ENDDO
130    
131     C- odd faces share disposition of all sections of the halo
132     IF ( MOD(myFace,2).EQ.1 ) THEN
133     C- North:
134     IF (exch2_isNedge(myTile).EQ.1) THEN
135     C switch u <- v , reverse the sign & shift i+1 <- i
136     DO j = 1,exchWidthY
137     DO i = 1-OLx,sNx+OLx-1
138     uPhi(i+1,sNy+j,k,bi,bj) = vLoc(i,sNy+j)*negOne
139     ENDDO
140     ENDDO
141     C switch v <- u , keep the sign
142     DO j = 1,exchWidthY
143     DO i = 1-OLx,sNx+OLx
144     vPhi(i,sNy+j,k,bi,bj) = uLoc(i,sNy+j)
145     ENDDO
146     ENDDO
147     ENDIF
148     C- South (nothing to change)
149     c IF (exch2_isSedge(myTile).EQ.1) THEN
150     c DO j = 1,exchWidthY
151     c DO i = 1-OLx,sNx+OLx
152     c uPhi(i,1-j,k,bi,bj) = uLoc(i,1-j)
153     c vPhi(i,1-j,k,bi,bj) = vLoc(i,1-j)
154     c ENDDO
155     c ENDDO
156     c ENDIF
157     C- East (nothing to change)
158     c IF (exch2_isEedge(myTile).EQ.1) THEN
159     c DO j = 1-OLy,sNy+OLy
160     c DO i = 1,exchWidthX
161     c uPhi(sNx+i,j,k,bi,bj) = uLoc(sNx+i,j)
162     c vPhi(sNx+i,j,k,bi,bj) = vLoc(sNx+i,j)
163     c ENDDO
164     c ENDDO
165     c ENDIF
166     C- West:
167     IF (exch2_isWedge(myTile).EQ.1) THEN
168     C switch u <- v , keep the sign
169     DO j = 1-OLy,sNy+OLy
170     DO i = 1,exchWidthX
171     uPhi(1-i,j,k,bi,bj) = vLoc(1-i,j)
172     ENDDO
173     ENDDO
174     C switch v <- u , reverse the sign & shift j+1 <- j
175     DO j = 1-OLy,sNy+OLy-1
176     DO i = 1,exchWidthX
177     vPhi(1-i,j+1,k,bi,bj) = uLoc(1-i,j)*negOne
178     ENDDO
179     ENDDO
180     ENDIF
181    
182     ELSE
183     C- Now the even faces (share disposition of all sections of the halo)
184    
185     C- East:
186     IF (exch2_isEedge(myTile).EQ.1) THEN
187     C switch u <- v , keep the sign
188     DO j = 1-OLy,sNy+OLy
189     DO i = 1,exchWidthX
190     uPhi(sNx+i,j,k,bi,bj) = vLoc(sNx+i,j)
191     ENDDO
192     ENDDO
193     C switch v <- u , reverse the sign & shift j+1 <- j
194     DO j = 1-OLy,sNy+OLy-1
195     DO i = 1,exchWidthX
196     vPhi(sNx+i,j+1,k,bi,bj) = uLoc(sNx+i,j)*negOne
197     ENDDO
198     ENDDO
199     ENDIF
200     C- West (nothing to change)
201     c IF (exch2_isWedge(myTile).EQ.1) THEN
202     c DO j = 1-OLy,sNy+OLy
203     c DO i = 1,exchWidthX
204     c uPhi(1-i,j,k,bi,bj) = uLoc(1-i,j)
205     c vPhi(1-i,j,k,bi,bj) = vLoc(1-i,j)
206     c ENDDO
207     c ENDDO
208     c ENDIF
209     C- North (nothing to change)
210     c IF (exch2_isNedge(myTile).EQ.1) THEN
211     c DO j = 1,exchWidthY
212     c DO i = 1-OLx,sNx+OLx
213     c uPhi(i,sNy+j,k,bi,bj) = uLoc(i,sNy+j)
214     c vPhi(i,sNy+j,k,bi,bj) = vLoc(i,sNy+j)
215     c ENDDO
216     c ENDDO
217     c ENDIF
218     C- South:
219     IF (exch2_isSedge(myTile).EQ.1) THEN
220     C switch u <- v , reverse the sign & shift i+1 <- i
221     DO j = 1,exchWidthY
222     DO i = 1-OLx,sNx+OLx-1
223     uPhi(i+1,1-j,k,bi,bj) = vLoc(i,1-j)*negOne
224     ENDDO
225     ENDDO
226     C switch v <- u , keep the sign
227     DO j = 1,exchWidthY
228     DO i = 1-OLx,sNx+OLx
229     vPhi(i,1-j,k,bi,bj) = uLoc(i,1-j)
230     ENDDO
231     ENDDO
232     ENDIF
233    
234     C- end odd / even faces
235     ENDIF
236    
237     C-- end of Loops on level index k.
238     ENDDO
239    
240     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
241     C-- Now fix edges near cube-corner
242    
243     IF ( exch2_isWedge(myTile) .EQ. 1 .AND.
244     & exch2_isSedge(myTile) .EQ. 1 ) THEN
245     IF ( MOD(myFace,2).EQ.1 ) THEN
246     DO k=1,myNz
247     DO i=1,OLx
248     vPhi(1-i,1,k,bi,bj) = uPhi(1,1-i,k,bi,bj)*negOne
249     ENDDO
250     ENDDO
251     ELSE
252     DO k=1,myNz
253     DO i=1,OLx
254     uPhi(1,1-i,k,bi,bj) = vPhi(1-i,1,k,bi,bj)*negOne
255     ENDDO
256     ENDDO
257     ENDIF
258     ENDIF
259    
260     IF ( exch2_isEedge(myTile) .EQ. 1 .AND.
261     & exch2_isSedge(myTile) .EQ. 1 ) THEN
262     IF ( MOD(myFace,2).EQ.1 ) THEN
263     DO k=1,myNz
264     DO i=1,OLx
265     uPhi(sNx+1,1-i,k,bi,bj) = vPhi(sNx+i,1,k,bi,bj)
266     ENDDO
267     ENDDO
268     ELSE
269     DO k=1,myNz
270     DO i=1,OLx
271     vPhi(sNx+i,1,k,bi,bj) = uPhi(sNx+1,1-i,k,bi,bj)
272     ENDDO
273     ENDDO
274     ENDIF
275     ENDIF
276    
277     IF ( exch2_isEedge(myTile) .EQ. 1 .AND.
278     & exch2_isNedge(myTile) .EQ. 1 ) THEN
279     IF ( MOD(myFace,2).EQ.1 ) THEN
280     DO k=1,myNz
281     DO i=1,OLx
282     vPhi(sNx+i,sNy+1,k,bi,bj)=uPhi(sNx+1,sNy+i,k,bi,bj)*negOne
283     ENDDO
284     ENDDO
285     ELSE
286     DO k=1,myNz
287     DO i=1,OLx
288     uPhi(sNx+1,sNy+i,k,bi,bj)=vPhi(sNx+i,sNy+1,k,bi,bj)*negOne
289     ENDDO
290     ENDDO
291     ENDIF
292     ENDIF
293    
294     IF ( exch2_isWedge(myTile) .EQ. 1 .AND.
295     & exch2_isNedge(myTile) .EQ. 1 ) THEN
296     IF ( MOD(myFace,2).EQ.1 ) THEN
297     DO k=1,myNz
298     DO i=1,OLx
299     uPhi(1,sNy+i,k,bi,bj) = vPhi(1-i,sNy+1,k,bi,bj)
300     ENDDO
301     ENDDO
302     ELSE
303     DO k=1,myNz
304     DO i=1,OLx
305     vPhi(1-i,sNy+1,k,bi,bj) = uPhi(1,sNy+i,k,bi,bj)
306     ENDDO
307     ENDDO
308     ENDIF
309     ENDIF
310    
311     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
312    
313     C-- Now zero out the null areas that should not be used in the numerics
314     C Also add one valid u,v value next to the corner, that allows
315     C to compute vorticity on a wider stencil (e.g., vort3(0,1) & (1,0))
316    
317     IF ( exch2_isWedge(myTile) .EQ. 1 .AND.
318     & exch2_isSedge(myTile) .EQ. 1 ) THEN
319     C Zero SW corner points
320     DO k=1,myNz
321     #ifdef W2_FILL_NULL_REGIONS
322     DO j=1-OLy,0
323     DO i=1-OLx,0
324     uPhi(i,j,k,bi,bj)=e2FillValue_RX
325     ENDDO
326     ENDDO
327     DO j=1-OLy,0
328     DO i=1-OLx,0
329     vPhi(i,j,k,bi,bj)=e2FillValue_RX
330     ENDDO
331     ENDDO
332     #endif
333     uPhi(0,0,k,bi,bj)=vPhi(1,0,k,bi,bj)
334     vPhi(0,0,k,bi,bj)=uPhi(0,1,k,bi,bj)
335     ENDDO
336     ENDIF
337    
338     IF ( exch2_isWedge(myTile) .EQ. 1 .AND.
339     & exch2_isNedge(myTile) .EQ. 1 ) THEN
340     C Zero NW corner points
341     DO k=1,myNz
342     #ifdef W2_FILL_NULL_REGIONS
343     DO j=sNy+1,sNy+OLy
344     DO i=1-OLx,0
345     uPhi(i,j,k,bi,bj)=e2FillValue_RX
346     ENDDO
347     ENDDO
348     DO j=sNy+2,sNy+OLy
349     DO i=1-OLx,0
350     vPhi(i,j,k,bi,bj)=e2FillValue_RX
351     ENDDO
352     ENDDO
353     #endif
354     uPhi(0,sNy+1,k,bi,bj)= vPhi(1,sNy+2,k,bi,bj)*negOne
355     vPhi(0,sNy+2,k,bi,bj)= uPhi(0,sNy,k,bi,bj)*negOne
356     ENDDO
357     ENDIF
358    
359     IF ( exch2_isEedge(myTile) .EQ. 1 .AND.
360     & exch2_isSedge(myTile) .EQ. 1 ) THEN
361     C Zero SE corner points
362     DO k=1,myNz
363     #ifdef W2_FILL_NULL_REGIONS
364     DO j=1-OLy,0
365     DO i=sNx+2,sNx+OLx
366     uPhi(i,j,k,bi,bj)=e2FillValue_RX
367     ENDDO
368     ENDDO
369     DO j=1-OLy,0
370     DO i=sNx+1,sNx+OLx
371     vPhi(i,j,k,bi,bj)=e2FillValue_RX
372     ENDDO
373     ENDDO
374     #endif
375     uPhi(sNx+2,0,k,bi,bj)= vPhi(sNx,0,k,bi,bj)*negOne
376     vPhi(sNx+1,0,k,bi,bj)= uPhi(sNx+2,1,k,bi,bj)*negOne
377     ENDDO
378     ENDIF
379    
380     IF ( exch2_isEedge(myTile) .EQ. 1 .AND.
381     & exch2_isNedge(myTile) .EQ. 1 ) THEN
382     C Zero NE corner points
383     DO k=1,myNz
384     #ifdef W2_FILL_NULL_REGIONS
385     DO j=sNy+1,sNy+OLy
386     DO i=sNx+2,sNx+OLx
387     uPhi(i,j,k,bi,bj)=e2FillValue_RX
388     ENDDO
389     ENDDO
390     DO j=sNy+2,sNy+OLy
391     DO i=sNx+1,sNx+OLx
392     vPhi(i,j,k,bi,bj)=e2FillValue_RX
393     ENDDO
394     ENDDO
395     #endif
396     uPhi(sNx+2,sNy+1,k,bi,bj)=vPhi(sNx,sNy+2,k,bi,bj)
397     vPhi(sNx+1,sNy+2,k,bi,bj)=uPhi(sNx+2,sNy,k,bi,bj)
398     ENDDO
399     ENDIF
400     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
401    
402     C-- end of Loops on tile indices (bi,bj).
403     ENDDO
404     ENDDO
405    
406     C--- using or not using CubedSphereExchange: end
407     ENDIF
408    
409     RETURN
410     END
411    
412     CEH3 ;;; Local Variables: ***
413     CEH3 ;;; mode:fortran ***
414     CEH3 ;;; End: ***

  ViewVC Help
Powered by ViewVC 1.1.22