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

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

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


Revision 1.1 - (show annotations) (download)
Fri Apr 28 17:22:44 2017 UTC (8 years, 2 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, HEAD
Error occurred while calculating annotation data.
add sea ice momentum advection code

1 C $Header: $
2 C $Name: $
3
4 #include "SEAICE_OPTIONS.h"
5 #ifdef ALLOW_AUTODIFF
6 # include "AUTODIFF_OPTIONS.h"
7 #endif
8
9 CBOP
10 C !ROUTINE: SEAICE_MOM_ADVECTION
11
12 C !INTERFACE: ==========================================================
13 SUBROUTINE SEAICE_MOM_ADVECTION(
14 I bi,bj,iMin,iMax,jMin,jMax,
15 I uIceLoc, vIceLoc,
16 O gU, gV,
17 I myTime, myIter, myThid )
18 C *==========================================================*
19 C | S/R SEAICE_MOM_ADVECTION |
20 C | o Form the advection of sea ice momentum to be added to |
21 C | the right hand-side of the momentum equation. |
22 C *==========================================================*
23 C | Most of the code is take from S/R MOM_VECINV and reuses |
24 C | code from mom_vecinv and mom_common |
25 C *==========================================================*
26 IMPLICIT NONE
27
28 C == Global variables ==
29 #include "SIZE.h"
30 #include "EEPARAMS.h"
31 #include "PARAMS.h"
32 #include "GRID.h"
33 #include "SEAICE_SIZE.h"
34 #include "SEAICE_PARAMS.h"
35 #ifdef ALLOW_AUTODIFF_TAMC
36 # include "tamc.h"
37 # include "tamc_keys.h"
38 #endif
39
40 C == Routine arguments ==
41 C bi,bj :: current tile indices
42 C iMin,iMax,jMin,jMax :: loop ranges
43 C uIceLoc ::
44 C vIceLoc ::
45
46 C gU :: advection tendency (all explicit terms), u component
47 C gV :: advection tendency (all explicit terms), v component
48 C myTime :: current time
49 C myIter :: current time-step number
50 C myThid :: my Thread Id number
51 INTEGER bi,bj
52 INTEGER iMin,iMax,jMin,jMax
53 _RL uIceLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
54 _RL vIceLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
55 _RL gU(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
56 _RL gV(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
57 _RL myTime
58 INTEGER myIter
59 INTEGER myThid
60 CEOP
61
62 #ifdef SEAICE_ALLOW_MOM_ADVECTION
63
64 C == Functions ==
65 LOGICAL DIFFERENT_MULTIPLE
66 EXTERNAL DIFFERENT_MULTIPLE
67
68 C == Local variables ==
69 _RL uCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
70 _RL vCf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
71 _RS hFacZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72 _RS r_hFacZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73 _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74 _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 _RL KE (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76 _RL vort3 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
77 C i,j :: Loop counters
78 C k :: surface level index
79 INTEGER i,j,k
80 C later these will be run time parameters
81 CML LOGICAL SEAICEhighOrderVorticity, SEAICEupwindVorticity
82 CML LOGICAL SEAICEuseAbsVorticity,
83 LOGICAL vorticityFlag
84 #ifdef ALLOW_AUTODIFF_TAMC
85 INTEGER imomkey
86 #endif
87
88 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
89
90 #ifdef ALLOW_AUTODIFF_TAMC
91 act0 = k - 1
92 max0 = Nr
93 act1 = bi - myBxLo(myThid)
94 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
95 act2 = bj - myByLo(myThid)
96 max2 = myByHi(myThid) - myByLo(myThid) + 1
97 act3 = myThid - 1
98 max3 = nTx*nTy
99 act4 = ikey_dynamics - 1
100 imomkey = (act0 + 1)
101 & + act1*max0
102 & + act2*max0*max1
103 & + act3*max0*max1*max2
104 & + act4*max0*max1*max2*max3
105 #endif /* ALLOW_AUTODIFF_TAMC */
106
107 CML SEAICEselectKEscheme = selectKEscheme
108 CML SEAICEselectVortScheme = selectVortScheme
109 CML SEAICEhighOrderVorticity = highOrderVorticity
110 CML SEAICEupwindVorticity = upwindVorticity
111 CML SEAICEuseAbsVorticity = useAbsVorticity
112 CML SEAICEuseJamartMomAdv = useJamartMomAdv
113
114 C-- Initialise intermediate terms
115 DO j=1-OLy,sNy+OLy
116 DO i=1-OLx,sNx+OLx
117 uCf(i,j) = 0.
118 vCf(i,j) = 0.
119 gU(i,j) = 0.
120 gV(i,j) = 0.
121 vort3(i,j) = 0.
122 KE(i,j) = 0.
123 #ifdef ALLOW_AUTODIFF
124 hFacZ(i,j) = 0. _d 0
125 #endif
126 ENDDO
127 ENDDO
128
129 k = 1
130 C-- Calculate open water fraction at vorticity points
131 CALL MOM_CALC_HFACZ(bi,bj,k,hFacZ,r_hFacZ,myThid)
132
133 C Make local copies of horizontal flow field
134 DO j=1-OLy,sNy+OLy
135 DO i=1-OLx,sNx+OLx
136 uFld(i,j) = uIceLoc(i,j,bi,bj)
137 vFld(i,j) = vIceLoc(i,j,bi,bj)
138 ENDDO
139 ENDDO
140
141 #ifdef ALLOW_AUTODIFF_TAMC
142 CADJ STORE ufld(:,:) =
143 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
144 CADJ STORE vfld(:,:) =
145 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
146 CADJ STORE hFacZ(:,:) =
147 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
148 CADJ STORE r_hFacZ(:,:) =
149 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
150 #endif
151
152 CALL MOM_CALC_KE(bi,bj,k,SEAICEselectKEscheme,uFld,vFld,KE,myThid)
153
154 CALL MOM_CALC_RELVORT3(bi,bj,k,uFld,vFld,hFacZ,vort3,myThid)
155
156 CMLC- calculate absolute vorticity
157 CML IF (useAbsVorticity) THEN
158 CML DO j=1-Oly,sNy+Oly
159 CML DO i=1-Olx,sNx+Olx
160 CML vort3(i,j) = vort3(i,j) + fCoriG(i,j,bi,bj)
161 CML ENDDO
162 CML ENDDO
163 CML ENDIF
164
165 #ifdef ALLOW_AUTODIFF_TAMC
166 CADJ STORE ke(:,:) =
167 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
168 CADJ STORE vort3(:,:) =
169 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
170 #endif
171
172 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
173
174 #ifdef ALLOW_AUTODIFF_TAMC
175 CADJ STORE ucf(:,:) =
176 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
177 CADJ STORE vcf(:,:) =
178 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
179 #endif
180
181 C-- Horizontal advection of relative (or absolute) vorticity
182 vorticityFlag = SEAICEhighOrderVorticity.OR.SEAICEupwindVorticity
183 IF ( vorticityFlag ) THEN
184 CALL MOM_VI_U_CORIOLIS_C4(bi,bj,k,SEAICEselectVortScheme,
185 & SEAICEhighOrderVorticity,
186 & SEAICEupwindVorticity,
187 & vFld,vort3,r_hFacZ,
188 & uCf,myThid)
189 ELSE
190 CALL MOM_VI_U_CORIOLIS(bi,bj,k,SEAICEselectVortScheme,
191 & SEAICEuseJamartMomAdv,
192 & vFld,vort3,hFacZ,r_hFacZ,
193 & uCf,myThid)
194 ENDIF
195 DO j=jMin,jMax
196 DO i=iMin,iMax
197 gU(i,j) = gU(i,j)+uCf(i,j)
198 ENDDO
199 ENDDO
200 IF ( vorticityFlag ) THEN
201 CALL MOM_VI_V_CORIOLIS_C4(bi,bj,k,SEAICEselectVortScheme,
202 & SEAICEhighOrderVorticity,
203 & SEAICEupwindVorticity,
204 & uFld,vort3,r_hFacZ,
205 & vCf,myThid)
206 ELSE
207 CALL MOM_VI_V_CORIOLIS(bi,bj,k,SEAICEselectVortScheme,
208 & SEAICEuseJamartMomAdv,
209 & uFld,vort3,hFacZ,r_hFacZ,
210 & vCf,myThid)
211 ENDIF
212 DO j=jMin,jMax
213 DO i=iMin,iMax
214 gV(i,j) = gV(i,j)+vCf(i,j)
215 ENDDO
216 ENDDO
217
218 #ifdef ALLOW_AUTODIFF_TAMC
219 CADJ STORE ucf(:,:) =
220 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
221 CADJ STORE vcf(:,:) =
222 CADJ & comlev1_bibj_lsr, key = imomkey, byte = isbyte
223 #endif
224
225 #ifdef ALLOW_DIAGNOSTICS
226 IF ( useDiagnostics ) THEN
227 CALL DIAGNOSTICS_FILL(uCf,'SIuAdvZ3',k,1,2,bi,bj,myThid)
228 CALL DIAGNOSTICS_FILL(vCf,'SIvAdvZ3',k,1,2,bi,bj,myThid)
229 ENDIF
230 #endif /* ALLOW_DIAGNOSTICS */
231
232 C-- Bernoulli term
233 CALL MOM_VI_U_GRAD_KE(bi,bj,k,KE,uCf,myThid)
234 DO j=jMin,jMax
235 DO i=iMin,iMax
236 gU(i,j) = gU(i,j)+uCf(i,j)
237 ENDDO
238 ENDDO
239 CALL MOM_VI_V_GRAD_KE(bi,bj,k,KE,vCf,myThid)
240 DO j=jMin,jMax
241 DO i=iMin,iMax
242 gV(i,j) = gV(i,j)+vCf(i,j)
243 ENDDO
244 ENDDO
245 #ifdef ALLOW_DIAGNOSTICS
246 IF ( useDiagnostics ) THEN
247 CALL DIAGNOSTICS_FILL(uCf,'SIKEx ',k,1,2,bi,bj,myThid)
248 CALL DIAGNOSTICS_FILL(vCf,'SIKEy ',k,1,2,bi,bj,myThid)
249 ENDIF
250 #endif /* ALLOW_DIAGNOSTICS */
251
252 C-- Set du/dt & dv/dt on boundaries to zero
253 C apply masks for interior (important when we have open boundaries)
254 DO j=jMin,jMax
255 DO i=iMin,iMax
256 gU(i,j) = gU(i,j)*maskInW(i,j,bi,bj)
257 gV(i,j) = gV(i,j)*maskInS(i,j,bi,bj)
258 ENDDO
259 ENDDO
260
261 #ifdef ALLOW_DIAGNOSTICS
262 IF ( useDiagnostics ) THEN
263 CALL DIAGNOSTICS_FILL(KE, 'SImomKE ',k,1,2,bi,bj,myThid)
264 CALL DIAGNOSTICS_FILL(gU, 'SIuMmAdv',k,1,2,bi,bj,myThid)
265 CALL DIAGNOSTICS_FILL(gV, 'SIvMmAdv',k,1,2,bi,bj,myThid)
266 ENDIF
267 #endif /* ALLOW_DIAGNOSTICS */
268
269 #endif /* SEAICE_ALLOW_MOM_ADVECTION */
270
271 RETURN
272 END

  ViewVC Help
Powered by ViewVC 1.1.22