/[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.3 - (hide annotations) (download)
Tue May 12 19:44:58 2009 UTC (15 years ago) by jmc
Branch: MAIN
CVS Tags: checkpoint61o, checkpoint61r, checkpoint61p, checkpoint61q
Changes since 1.2: +4 -1 lines
new header files "W2_EXCH2_SIZE.h" (taken out of W2_EXCH2_TOPOLOGY.h)
             and "W2_EXCH2_BUFFER.h" (taken out of W2_EXCH2_PARAMS.h)

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

  ViewVC Help
Powered by ViewVC 1.1.22