/[MITgcm]/MITgcm_contrib/PRM/multi_comp_setup/cg/code/set_ddtvars.F
ViewVC logotype

Contents of /MITgcm_contrib/PRM/multi_comp_setup/cg/code/set_ddtvars.F

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


Revision 1.7 - (show annotations) (download)
Sun Jan 4 17:41:54 2009 UTC (16 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint63m, checkpoint63n, HEAD
Changes since 1.6: +19 -8 lines
add option to only apply gU-fg tendency (and not gV-fg)

1 C $Header: /u/gcmpack/MITgcm_contrib/PRM/multi_comp_setup/cg/code/set_ddtvars.F,v 1.6 2008/05/04 20:39:32 jmc Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #ifdef ALLOW_MYPACKAGE
6 # include "MYPACKAGE_OPTIONS.h"
7 #else
8 # include "CPP_OPTIONS.h"
9 #endif
10
11 CBOP
12 C !ROUTINE: SET_DDTVARS
13 C !INTERFACE:
14 SUBROUTINE SET_DDTVARS( gUVelC, gVVelC, gThetaC, gSaltC,
15 I cnx, cny, cnr,
16 I myTime, myIter, myThid )
17 C !DESCRIPTION: \bv
18 C *==========================================================*
19 C | S/R SET_DDTVARS
20 C | o Set tendencies to be used as additional forcing
21 C *==========================================================*
22 C *==========================================================*
23 C \ev
24
25 C !USES:
26 IMPLICIT NONE
27
28 C == Global data ==
29 #include "SIZE.h"
30 #include "EEPARAMS.h"
31 #include "PARAMS.h"
32 #include "GRID.h"
33 #include "ORIENTATION.h"
34 #include "EXT_MOM_TEND.h"
35 #ifdef ALLOW_MYPACKAGE
36 # include "MYPACKAGE.h"
37 #endif
38
39 C !INPUT/OUTPUT PARAMETERS:
40 C == Routine arguments ==
41 C myTime :: Current time in simulation
42 C myIter :: Current iteration number
43 C myThid :: Thread Id number
44 INTEGER cnx, cny, cnr
45 REAL*8 gUVelC( cnx, cny, cnr)
46 REAL*8 gVVelC( cnx, cny, cnr)
47 REAL*8 gThetaC(cnx, cny, cnr)
48 REAL*8 gSaltC( cnx, cny, cnr)
49 _RL myTime
50 INTEGER myIter
51 INTEGER myThid
52
53 C !FUNCTIONS:
54 LOGICAL DIFFERENT_MULTIPLE
55 EXTERNAL DIFFERENT_MULTIPLE
56
57 C !LOCAL VARIABLES:
58 C == Local variables ==
59 C i,j,k :: Loop counters
60 INTEGER i,j,k
61 INTEGER bi,bj
62 C Variables used for I/O
63 CHARACTER*(10) suff
64 C Variables used for smoothing momentum tendencies:
65 _RL viscAhDivDt, viscAhVortDt
66 _RL locDiv (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
67 _RL locVort(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
68 _RL locFx (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
69 _RL locFy (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
70 CEOP
71
72 viscAhDivDt = myPa_param1 * deltaTmom
73 viscAhVortDt = myPa_param2 * deltaTmom
74
75 IF ( cnx .NE. sNx ) THEN
76 STOP 'SET_DTTVARS cnx NE sNx'
77 ENDIF
78 IF ( cny .NE. sNy ) THEN
79 STOP 'SET_DTTVARS cny NE sNy'
80 ENDIF
81 IF ( cnr .NE. Nr ) THEN
82 STOP 'SET_DTTVARS cnr NE Nr'
83 ENDIF
84
85 DO bj = myByLo(myThid), myByHi(myThid)
86 DO bi = myBxLo(myThid), myBxHi(myThid)
87
88 C--- Temp. & Salt Tendency from CG : just do a copy
89 DO k=1,Nr
90 DO j=1,sNy
91 DO i=1,sNx
92 myPa_TendScal1(i,j,k,bi,bj) = gThetaC(i,j,k)
93 myPa_TendScal2(i,j,k,bi,bj) = gsaltC (i,j,k)
94 ENDDO
95 ENDDO
96 ENDDO
97
98 C--- Momentum Tendency from CG
99 C-- 1) rotate to get u,v aline with CG axes
100 DO k=1,Nr
101 DO j=1-Oly,sNy+Oly
102 DO i=1-Olx,sNx+Olx
103 ext_gu(i,j,k,bi,bj) = 0. _d 0
104 ext_gv(i,j,k,bi,bj) = 0. _d 0
105 ENDDO
106 ENDDO
107 IF ( myPa_doSwitch1 ) THEN
108 DO j=1,sNy
109 DO i=1,sNx
110 IF ( velDvRMS(i,j,bi,bj).GT.0. _d 0 ) THEN
111 ext_gu(i,j,k,bi,bj) = gUVelC(i,j,k)*cAngleFG(i,j,bi,bj)
112 ext_gv(i,j,k,bi,bj) = gUVelC(i,j,k)*sAngleFG(i,j,bi,bj)
113 ENDIF
114 ENDDO
115 ENDDO
116 ELSE
117 DO j=1,sNy
118 DO i=1,sNx
119 IF ( velDvRMS(i,j,bi,bj).GT.0. _d 0 ) THEN
120 ext_gu(i,j,k,bi,bj) = gUVelC(i,j,k)*cAngleFG(i,j,bi,bj)
121 & - gVVelC(i,j,k)*sAngleFG(i,j,bi,bj)
122 ext_gv(i,j,k,bi,bj) = gUVelC(i,j,k)*sAngleFG(i,j,bi,bj)
123 & + gVVelC(i,j,k)*cAngleFG(i,j,bi,bj)
124 ENDIF
125 ENDDO
126 ENDDO
127 ENDIF
128 ENDDO
129
130 C- end bi,bj loops
131 ENDDO
132 ENDDO
133
134 C-- 2) Fill the overlap : ext_gu & ext_gv are both on A-grid
135 CALL EXCH_UV_AGRID_3D_RL( ext_gu, ext_gv,
136 & .TRUE., Nr, myThid )
137
138 DO bj = myByLo(myThid), myByHi(myThid)
139 DO bi = myBxLo(myThid), myBxHi(myThid)
140 DO k=1,Nr
141 C-- 3) average to C-grid
142 DO j=1-OLy,sNy+OLy
143 DO i=2-OLx,sNx+OLx
144 myPa_TendVelU(i,j,k,bi,bj) = maskW(i,j,k,bi,bj)
145 & *( ext_gu(i-1,j,k,bi,bj)
146 & +ext_gu( i ,j,k,bi,bj) )*0.5 _d 0
147 ENDDO
148 ENDDO
149 DO j=2-OLy,sNy+OLy
150 DO i=1-OLx,sNx+OLx
151 myPa_TendVelV(i,j,k,bi,bj) = maskS(i,j,k,bi,bj)
152 & *( ext_gv(i,j-1,k,bi,bj)
153 & +ext_gv(i, j ,k,bi,bj) )*0.5 _d 0
154 ENDDO
155 ENDDO
156 C-- 4) Smooth Momentum tendency :
157 C apply horizontal viscosity (*time-step = viscAhDivDt) on Divergence
158 C tendency and on vorticity tend. (viscosity*time-step = viscAhVortDt)
159 IF ( viscAhDivDt.GT.0. .OR. viscAhVortDt.GT.0. ) THEN
160 C Fluxes
161 DO j=2-OLy,sNy+OLy
162 DO i=2-OLx,sNx+OLx
163 locFx(i,j) = myPa_TendVelU(i,j,k,bi,bj)
164 & *dyG(i,j,bi,bj)*hFacW(i,j,k,bi,bj)
165 locFy(i,j) = myPa_TendVelV(i,j,k,bi,bj)
166 & *dxG(i,j,bi,bj)*hFacS(i,j,k,bi,bj)
167 ENDDO
168 ENDDO
169 C Divergence
170 DO j=2-OLy,sNy+OLy-1
171 DO i=2-OLx,sNx+OLx-1
172 locDiv(i,j) =
173 & ( ( locFx(i+1,j)-locFx(i,j) )
174 & +( locFy(i,j+1)-locFy(i,j) )
175 & )*recip_rA(i,j,bi,bj)*recip_hFacC(i,j,k,bi,bj)
176 ENDDO
177 ENDDO
178 C Vorticity
179 DO j=2-OLy,sNy+OLy
180 DO i=2-OLx,sNx+OLx
181 locVort(i,j) =
182 & ( ( myPa_TendVelV(i,j,k,bi,bj)*dyC(i,j,bi,bj)
183 & -myPa_TendVelV(i-1,j,k,bi,bj)*dyC(i-1,j,bi,bj) )
184 & -( myPa_TendVelU(i,j,k,bi,bj)*dxC(i,j,bi,bj)
185 & -myPa_TendVelU(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj) )
186 & )*recip_rAz(I,J,bi,bj)
187 ENDDO
188 ENDDO
189 C Apply to C-grid tendencies:
190 DO j=3-OLy,sNy+OLy-1
191 DO i=3-OLx,sNx+OLx-1
192 myPa_TendVelU(i,j,k,bi,bj) = myPa_TendVelU(i,j,k,bi,bj)
193 & + ( viscAhDivDt * ( locDiv(i,j) - locDiv(i-1,j) )
194 & * recip_dxC(i,j,bi,bj)
195 & - viscAhVortDt* ( locVort(i,j+1)-locVort(i,j) )
196 & * recip_dyG(i,j,bi,bj)
197 & )*maskW(i,j,k,bi,bj)
198 myPa_TendVelV(i,j,k,bi,bj) = myPa_TendVelV(i,j,k,bi,bj)
199 & + ( viscAhDivDt * ( locDiv(i,j) - locDiv(i,j-1) )
200 & * recip_dyC(i,j,bi,bj)
201 & + viscAhVortDt* ( locVort(i+1,j)-locVort(i,j) )
202 & * recip_dxG(i,j,bi,bj)
203 & )*maskS(i,j,k,bi,bj)
204 ENDDO
205 ENDDO
206 C- end smoothing & end of k loop
207 ENDIF
208 ENDDO
209
210 C- end bi,bj loops
211 ENDDO
212 ENDDO
213
214 C-- Diagnostics
215 #ifdef ALLOW_DIAGNOSTICS
216 IF ( useDiagnostics ) THEN
217 CALL DIAGNOSTICS_FILL ( myPa_TendScal1,'MYPadTdt',
218 & 0, Nr, 0, 1, 1, myThid )
219 CALL DIAGNOSTICS_FILL ( myPa_TendScal2,'MYPadSdt',
220 & 0, Nr, 0, 1, 1, myThid )
221 CALL DIAGNOSTICS_FILL ( myPa_TendVelU, 'MYPadUdt',
222 & 0, Nr, 0, 1, 1, myThid )
223 CALL DIAGNOSTICS_FILL ( myPa_TendVelV, 'MYPadVdt',
224 & 0, Nr, 0, 1, 1, myThid )
225 C- This is a hack (to avoid changing mypackage_diagnostics_init.F)
226 CALL DIAGNOSTICS_FILL ( ext_gu, 'MYPaStaU',
227 & 0, Nr, 0, 1, 1, myThid )
228 CALL DIAGNOSTICS_FILL ( ext_gv, 'MYPaStaV',
229 & 0, Nr, 0, 1, 1, myThid )
230 ENDIF
231 #endif /* ALLOW_DIAGNOSTICS */
232
233 C-- Write snap-shot output:
234 IF ( DIFFERENT_MULTIPLE(dumpFreq,myTime,deltaTClock)
235 & .OR. dumpInitAndLast.AND.( myTime.EQ.endTime .OR.
236 & myTime.EQ.startTime )
237 & ) THEN
238 WRITE(suff,'(I10.10)') myIter
239 CALL WRITE_FLD_XYZ_RL( 'ext_gT.', suff, myPa_TendScal1,
240 & myIter, myThid )
241 CALL WRITE_FLD_XYZ_RL( 'ext_gS.', suff, myPa_TendScal2,
242 & myIter, myThid )
243 CALL WRITE_FLD_XYZ_RL( 'ext_gU.', suff, ext_gu,
244 & myIter, myThid )
245 CALL WRITE_FLD_XYZ_RL( 'ext_gV.', suff, ext_gv,
246 & myIter, myThid )
247 ENDIF
248
249 RETURN
250 END

  ViewVC Help
Powered by ViewVC 1.1.22