/[MITgcm]/MITgcm/pkg/mom_common/mom_calc_3d_strain.F
ViewVC logotype

Annotation of /MITgcm/pkg/mom_common/mom_calc_3d_strain.F

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

Revision 1.2 - (hide annotations) (download)
Wed Nov 6 00:37:11 2013 UTC (10 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint65, HEAD
Changes since 1.1: +35 -36 lines 
empty routine when ALLOW_SMAG_3D is not defined
1 jmc 1.2 C $Header: /u/gcmpack/MITgcm/pkg/mom_common/mom_calc_3d_strain.F,v 1.1 2013/11/05 13:31:50 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "MOM_COMMON_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7     CBOP
8     C !ROUTINE: MOM_CALC_3D_STRAIN
9    
10     C !INTERFACE:
11     SUBROUTINE MOM_CALC_3D_STRAIN(
12     O str11, str22, str33, str12, str13, str23,
13     I bi, bj, myThid )
14    
15     C !DESCRIPTION:
16     C Calculates the strain tensor of the 3-D flow field
17    
18     C !USES:
19     IMPLICIT NONE
20    
21     C == Global variables ==
22     #include "SIZE.h"
23     #include "EEPARAMS.h"
24     #include "PARAMS.h"
25     #include "GRID.h"
26     #include "DYNVARS.h"
27    
28     C !INPUT PARAMETERS:
29     C bi, bj :: tile indices
30     C myThid :: my Thread Id number
31     INTEGER bi, bj
32     INTEGER myThid
33    
34     C !OUTPUT PARAMETERS:
35     C str11 :: strain component Vxx @ grid-cell center
36     C str22 :: strain component Vyy @ grid-cell center
37     C str33 :: strain component Vzz @ grid-cell center
38     C str12 :: strain component Vxy @ grid-cell corner
39     C str13 :: strain component Vxz @ above uVel
40     C str23 :: strain component Vyz @ above vVel
41     _RL str11(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
42     _RL str22(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
43     _RL str33(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
44     _RL str12(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
45     _RL str13(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
46     _RL str23(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1)
47    
48 jmc 1.2 #ifdef ALLOW_SMAG_3D
49 jmc 1.1 C !LOCAL VARIABLES:
50     C i, j, k :: loop indices
51     INTEGER i, j, k
52     INTEGER kp1
53     _RL maskp1
54     LOGICAL freeSlip3d
55     CEOP
56     freeSlip3d = .NOT.( no_slip_sides .AND. no_slip_bottom )
57    
58     DO k=1,Nr
59     kp1 = MIN(k+1,Nr)
60     maskp1 = oneRL
61     IF ( k.EQ.Nr ) maskp1 = zeroRL
62    
63     C- Fills up array edges:
64     i = sNx+OLx
65     DO j=1-OLy,sNy+OLy
66     str11(i,j,k) = 0. _d 0
67     ENDDO
68     j = sNy+OLy
69     DO i=1-OLx,sNx+OLx
70     str22(i,j,k) = 0. _d 0
71     ENDDO
72     i = 1-OLx
73     DO j=1-OLy,sNy+OLy
74     str12(i,j,k) = 0. _d 0
75     str13(i,j,k) = 0. _d 0
76     ENDDO
77     j = 1-OLy
78     DO i=1-OLx,sNx+OLx
79     str12(i,j,k) = 0. _d 0
80     str23(i,j,k) = 0. _d 0
81     ENDDO
82    
83     C str11 = u_x
84     DO j=1-OLy,sNy+OLy
85     DO i=1-OLx,sNx+OLx-1
86 jmc 1.2 str11(i,j,k) = recip_dxF(i,j,bi,bj)
87     & *( uVel(i+1, j , k ,bi,bj)-uVel( i , j , k ,bi,bj) )
88 jmc 1.1 ENDDO
89     ENDDO
90    
91     C str22 = v_y
92     DO j=1-OLy,sNy+OLy-1
93     DO i=1-OLx,sNx+OLx
94 jmc 1.2 str22(i,j,k) = recip_dyF(i,j,bi,bj)
95     & *( vVel( i ,j+1, k ,bi,bj)-vVel( i , j , k ,bi,bj) )
96 jmc 1.1 ENDDO
97     ENDDO
98    
99     C str33 = w_z
100     DO j=1-OLy,sNy+OLy
101     DO i=1-OLx,sNx+OLx
102 jmc 1.2 str33(i,j,k) = recip_drF(k)*rkSign
103     & *( maskp1*wVel( i , j ,kp1,bi,bj)-wVel( i , j , k ,bi,bj) )
104 jmc 1.1 ENDDO
105     ENDDO
106    
107     C str12 = ( u_y + v_x )/2
108     DO j=2-OLy,sNy+OLy
109     DO i=2-OLx,sNx+OLx
110     str12(i,j,k) = halfRL*(
111 jmc 1.2 & recip_dyU(i,j,bi,bj)
112     & *( uVel( i , j , k ,bi,bj)-uVel( i ,j-1, k ,bi,bj) )
113     & +recip_dxV(i,j,bi,bj)
114     & *( vVel( i , j , k ,bi,bj)-vVel(i-1, j , k ,bi,bj) )
115 jmc 1.1 & )
116     ENDDO
117     ENDDO
118    
119     C str13 & str23 special case: k=1
120     IF ( k.EQ.1 .AND. freeSlip3d ) THEN
121     DO j=1-OLy,sNy+OLy
122     DO i=1-OLx,sNx+OLx
123     str13(i,j,k) = 0. _d 0
124     str23(i,j,k) = 0. _d 0
125     ENDDO
126     ENDDO
127     ELSEIF ( k.EQ.1 ) THEN
128     C-- should put surface wind-stress if z-coords; but right in p-coords:
129     DO j=1-OLy,sNy+OLy
130     DO i=2-OLx,sNx+OLx
131     str13(i,j,k) = halfRL*(
132 jmc 1.2 & recip_drC(k)*rkSign
133     & *( uVel( i , j , k ,bi,bj)*twoRL )
134     & +recip_dxC(i,j,bi,bj)
135     & *( wVel( i , j , k ,bi,bj)-wVel(i-1, j , k ,bi,bj) )
136 jmc 1.1 & )
137     ENDDO
138     ENDDO
139     DO j=2-OLy,sNy+OLy
140     DO i=1-OLx,sNx+OLx
141     str23(i,j,k) = halfRL*(
142 jmc 1.2 & recip_drC(k)*rkSign
143     & *( vVel( i , j , k ,bi,bj)*twoRL )
144     & +recip_dyC(i,j,bi,bj)
145     & *( wVel( i , j , k ,bi,bj)-wVel( i ,j-1, k ,bi,bj) )
146 jmc 1.1 & )
147     ENDDO
148     ENDDO
149     ELSE
150     C str13 = ( u_z + w_x )/2
151     DO j=1-OLy,sNy+OLy
152     DO i=2-OLx,sNx+OLx
153     str13(i,j,k) = halfRL*(
154 jmc 1.2 & recip_drC(k)*rkSign
155     & *( uVel( i , j , k ,bi,bj)-uVel( i , j ,k-1 ,bi,bj) )
156     & +recip_dxC(i,j,bi,bj)
157     & *( wVel( i , j , k ,bi,bj)-wVel(i-1, j , k ,bi,bj) )
158 jmc 1.1 & )
159     ENDDO
160     ENDDO
161     C str23 = ( v_z + w_y )/2
162     DO j=2-OLy,sNy+OLy
163     DO i=1-OLx,sNx+OLx
164     str23(i,j,k) = halfRL*(
165 jmc 1.2 & recip_drC(k)*rkSign
166     & *( vVel( i , j , k ,bi,bj)-vVel( i , j ,k-1,bi,bj) )
167     & +recip_dyC(i,j,bi,bj)
168     & *( wVel( i , j , k ,bi,bj)-wVel( i ,j-1, k ,bi,bj) )
169 jmc 1.1 & )
170     ENDDO
171     ENDDO
172     ENDIF
173    
174     IF ( freeSlip3d ) THEN
175     DO j=2-OLy,sNy+OLy
176     DO i=2-OLx,sNx+OLx
177     str12(i,j,k) = str12(i,j,k)
178     & *maskW(i,j-1,k,bi,bj)*maskW(i,j,k,bi,bj)
179     ENDDO
180     ENDDO
181     IF ( k.GE.2 ) THEN
182     DO j=1-OLy,sNy+OLy
183     DO i=2-OLx,sNx+OLx
184     str13(i,j,k) = str13(i,j,k)
185     & *maskW(i,j,k-1,bi,bj)*maskW(i,j,k,bi,bj)
186     ENDDO
187     ENDDO
188     DO j=2-OLy,sNy+OLy
189     DO i=1-OLx,sNx+OLx
190     str23(i,j,k) = str23(i,j,k)
191     & *maskS(i,j,k-1,bi,bj)*maskS(i,j,k,bi,bj)
192     ENDDO
193     ENDDO
194     ENDIF
195     ENDIF
196    
197     C-- end k loop
198     ENDDO
199    
200     C-- fill-up strain tensor component at the very bottom (k=Nr+1)
201     k = Nr+1
202    
203     DO j=1-OLy,sNy+OLy
204     DO i=1-OLx,sNx+OLx
205     str13(i,j,k) = 0. _d 0
206     str23(i,j,k) = 0. _d 0
207     ENDDO
208     ENDDO
209    
210     IF ( .NOT.freeSlip3d ) THEN
211    
212     C str13 = ( u_z + w_x )/2
213     DO j=1-OLy,sNy+OLy
214     DO i=2-OLx,sNx+OLx
215     str13(i,j,k) =
216 jmc 1.2 & recip_drF(Nr)*rkSign
217     c & recip_drC(k)*rkSign
218     & *( 0. _d 0 - uVel( i , j ,k-1 ,bi,bj) )
219 jmc 1.1 ENDDO
220     ENDDO
221    
222     C str23 = ( v_z + w_y )/2
223     DO j=2-OLy,sNy+OLy
224     DO i=1-OLx,sNx+OLx
225     str23(i,j,k) =
226 jmc 1.2 & recip_drF(Nr)*rkSign
227     c & recip_drC(k)*rkSign
228     & *( 0. _d 0 - vVel( i , j ,k-1,bi,bj) )
229 jmc 1.1 ENDDO
230     ENDDO
231    
232     ENDIF
233    
234     C Special stuff for Cubed Sphere
235     c IF (useCubedSphereExchange) THEN
236     c STOP 'S/R MOM_CALC_3D_STRAIN: should not be used on the cube!'
237     c ENDIF
238    
239 jmc 1.2 #endif /* ALLOW_SMAG_3D */
240 jmc 1.1 RETURN
241     END
  ViewVC Help
Powered by ViewVC 1.1.22