/[MITgcm]/MITgcm/model/src/external_forcing.F
ViewVC logotype

Contents of /MITgcm/model/src/external_forcing.F

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


Revision 1.35 - (show annotations) (download)
Mon Dec 19 19:09:35 2005 UTC (18 years, 5 months ago) by stephd
Branch: MAIN
CVS Tags: checkpoint57z_post
Changes since 1.34: +15 -1 lines
o add RBCS call, to relax T/S to 3-D field

1 C $Header: /u/gcmpack/MITgcm/model/src/external_forcing.F,v 1.34 2005/12/15 17:47:54 jmc Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6 #ifdef ALLOW_EXF
7 # include "EXF_OPTIONS.h"
8 #endif
9
10 CBOP
11 C !ROUTINE: EXTERNAL_FORCING_U
12 C !INTERFACE:
13 SUBROUTINE EXTERNAL_FORCING_U(
14 I iMin,iMax, jMin,jMax, bi,bj, kLev,
15 I myTime, myThid )
16 C !DESCRIPTION: \bv
17 C *==========================================================*
18 C | S/R EXTERNAL_FORCING_U
19 C | o Contains problem specific forcing for zonal velocity.
20 C *==========================================================*
21 C | Adds terms to gU for forcing by external sources
22 C | e.g. wind stress, bottom friction etc ...
23 C *==========================================================*
24 C \ev
25
26 C !USES:
27 IMPLICIT NONE
28 C == Global data ==
29 #include "SIZE.h"
30 #include "EEPARAMS.h"
31 #include "PARAMS.h"
32 #include "GRID.h"
33 #include "DYNVARS.h"
34 #include "FFIELDS.h"
35
36 C !INPUT/OUTPUT PARAMETERS:
37 C == Routine arguments ==
38 C iMin,iMax :: Working range of x-index for applying forcing.
39 C jMin,jMax :: Working range of y-index for applying forcing.
40 C bi,bj :: Current tile indices
41 C kLev :: Current vertical level index
42 C myTime :: Current time in simulation
43 C myThid :: Thread Id number
44 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
45 _RL myTime
46 INTEGER myThid
47
48 C !LOCAL VARIABLES:
49 C == Local variables ==
50 C i,j :: Loop counters
51 C kSurface :: index of surface layer
52 INTEGER i, j
53 INTEGER kSurface
54 CEOP
55
56 IF ( fluidIsAir ) THEN
57 kSurface = 0
58 ELSEIF ( usingPCoords ) THEN
59 kSurface = Nr
60 ELSE
61 kSurface = 1
62 ENDIF
63
64 C-- Forcing term
65 #ifdef ALLOW_AIM
66 IF ( useAIM ) CALL AIM_TENDENCY_APPLY_U(
67 & iMin,iMax, jMin,jMax, bi,bj, kLev,
68 & myTime, myThid )
69 #endif /* ALLOW_AIM */
70
71 #ifdef ALLOW_FIZHI
72 IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_U(
73 & iMin,iMax, jMin,jMax, bi,bj, kLev,
74 & myTime, myThid )
75 #endif /* ALLOW_FIZHI */
76
77 C Add windstress momentum impulse into the top-layer
78 IF ( kLev .EQ. kSurface ) THEN
79 c DO j=1,sNy
80 C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNy+1]
81 DO j=0,sNy+1
82 DO i=1,sNx+1
83 gU(i,j,kLev,bi,bj) = gU(i,j,kLev,bi,bj)
84 & +foFacMom*surfaceForcingU(i,j,bi,bj)
85 & *recip_drF(kLev)*recip_hFacW(i,j,kLev,bi,bj)
86 ENDDO
87 ENDDO
88 ENDIF
89
90 #if (defined (ALLOW_TAU_EDDY))
91 CALL TAUEDDY_EXTERNAL_FORCING_U(
92 I iMin,iMax, jMin,jMax, bi,bj, kLev,
93 I myTime, myThid )
94 #endif
95
96 #ifdef ALLOW_OBCS
97 IF (useOBCS) THEN
98 CALL OBCS_SPONGE_U(
99 I iMin,iMax, jMin,jMax, bi,bj, kLev,
100 I myTime, myThid )
101 ENDIF
102 #endif
103
104 RETURN
105 END
106
107 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
108 CBOP
109 C !ROUTINE: EXTERNAL_FORCING_V
110 C !INTERFACE:
111 SUBROUTINE EXTERNAL_FORCING_V(
112 I iMin,iMax, jMin,jMax, bi,bj, kLev,
113 I myTime, myThid )
114 C !DESCRIPTION: \bv
115 C *==========================================================*
116 C | S/R EXTERNAL_FORCING_V
117 C | o Contains problem specific forcing for merid velocity.
118 C *==========================================================*
119 C | Adds terms to gV for forcing by external sources
120 C | e.g. wind stress, bottom friction etc ...
121 C *==========================================================*
122 C \ev
123
124 C !USES:
125 IMPLICIT NONE
126 C == Global data ==
127 #include "SIZE.h"
128 #include "EEPARAMS.h"
129 #include "PARAMS.h"
130 #include "GRID.h"
131 #include "DYNVARS.h"
132 #include "FFIELDS.h"
133
134 C !INPUT/OUTPUT PARAMETERS:
135 C == Routine arguments ==
136 C iMin,iMax :: Working range of x-index for applying forcing.
137 C jMin,jMax :: Working range of y-index for applying forcing.
138 C bi,bj :: Current tile indices
139 C kLev :: Current vertical level index
140 C myTime :: Current time in simulation
141 C myThid :: Thread Id number
142 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
143 _RL myTime
144 INTEGER myThid
145
146 C !LOCAL VARIABLES:
147 C == Local variables ==
148 C i,j :: Loop counters
149 C kSurface :: index of surface layer
150 INTEGER i, j
151 INTEGER kSurface
152 CEOP
153
154 IF ( fluidIsAir ) THEN
155 kSurface = 0
156 ELSEIF ( usingPCoords ) THEN
157 kSurface = Nr
158 ELSE
159 kSurface = 1
160 ENDIF
161
162 C-- Forcing term
163 #ifdef ALLOW_AIM
164 IF ( useAIM ) CALL AIM_TENDENCY_APPLY_V(
165 & iMin,iMax, jMin,jMax, bi,bj, kLev,
166 & myTime, myThid )
167 #endif /* ALLOW_AIM */
168
169 #ifdef ALLOW_FIZHI
170 IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_V(
171 & iMin,iMax, jMin,jMax, bi,bj, kLev,
172 & myTime, myThid )
173 #endif /* ALLOW_FIZHI */
174
175 C Add windstress momentum impulse into the top-layer
176 IF ( kLev .EQ. kSurface ) THEN
177 DO j=1,sNy+1
178 c DO i=1,sNx
179 C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNx+1]
180 DO i=0,sNx+1
181 gV(i,j,kLev,bi,bj) = gV(i,j,kLev,bi,bj)
182 & +foFacMom*surfaceForcingV(i,j,bi,bj)
183 & *recip_drF(kLev)*recip_hFacS(i,j,kLev,bi,bj)
184 ENDDO
185 ENDDO
186 ENDIF
187
188 #if (defined (ALLOW_TAU_EDDY))
189 CALL TAUEDDY_EXTERNAL_FORCING_V(
190 I iMin,iMax, jMin,jMax, bi,bj, kLev,
191 I myTime, myThid )
192 #endif
193
194 #ifdef ALLOW_OBCS
195 IF (useOBCS) THEN
196 CALL OBCS_SPONGE_V(
197 I iMin,iMax, jMin,jMax, bi,bj, kLev,
198 I myTime, myThid )
199 ENDIF
200 #endif
201
202 RETURN
203 END
204
205 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
206 CBOP
207 C !ROUTINE: EXTERNAL_FORCING_T
208 C !INTERFACE:
209 SUBROUTINE EXTERNAL_FORCING_T(
210 I iMin,iMax, jMin,jMax, bi,bj, kLev,
211 I myTime, myThid )
212 C !DESCRIPTION: \bv
213 C *==========================================================*
214 C | S/R EXTERNAL_FORCING_T
215 C | o Contains problem specific forcing for temperature.
216 C *==========================================================*
217 C | Adds terms to gT for forcing by external sources
218 C | e.g. heat flux, climatalogical relaxation, etc ...
219 C *==========================================================*
220 C \ev
221
222 C !USES:
223 IMPLICIT NONE
224 C == Global data ==
225 #include "SIZE.h"
226 #include "EEPARAMS.h"
227 #include "PARAMS.h"
228 #include "GRID.h"
229 #include "DYNVARS.h"
230 #include "FFIELDS.h"
231
232 C !INPUT/OUTPUT PARAMETERS:
233 C == Routine arguments ==
234 C iMin,iMax :: Working range of x-index for applying forcing.
235 C jMin,jMax :: Working range of y-index for applying forcing.
236 C bi,bj :: Current tile indices
237 C kLev :: Current vertical level index
238 C myTime :: Current time in simulation
239 C myThid :: Thread Id number
240 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
241 _RL myTime
242 INTEGER myThid
243
244 C !LOCAL VARIABLES:
245 C == Local variables ==
246 C i,j :: Loop counters
247 C kSurface :: index of surface layer
248 INTEGER i, j
249 INTEGER kSurface
250 CEOP
251 #ifdef SHORTWAVE_HEATING
252 integer two
253 _RL minusone
254 parameter (two=2,minusone=-1.)
255 _RL swfracb(two)
256 INTEGER kp1
257 #endif
258
259 IF ( fluidIsAir ) THEN
260 kSurface = 0
261 ELSEIF ( usingPCoords ) THEN
262 kSurface = Nr
263 ELSE
264 kSurface = 1
265 ENDIF
266
267 C-- Forcing term
268 #ifdef ALLOW_AIM
269 IF ( useAIM ) CALL AIM_TENDENCY_APPLY_T(
270 & iMin,iMax, jMin,jMax, bi,bj, kLev,
271 & myTime, myThid )
272 #endif /* ALLOW_AIM */
273
274 #ifdef ALLOW_FIZHI
275 IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_T(
276 & iMin,iMax, jMin,jMax, bi,bj, kLev,
277 & myTime, myThid )
278 #endif /* ALLOW_FIZHI */
279
280 C Add heat in top-layer
281 IF ( kLev .EQ. kSurface ) THEN
282 DO j=1,sNy
283 DO i=1,sNx
284 gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)
285 & +surfaceForcingT(i,j,bi,bj)
286 & *recip_drF(kLev)*recip_hFacC(i,j,kLev,bi,bj)
287 ENDDO
288 ENDDO
289 ENDIF
290
291 #ifdef SHORTWAVE_HEATING
292 C Penetrating SW radiation
293 c IF ( usePenetratingSW ) THEN
294 swfracb(1)=abs(rF(klev))
295 swfracb(2)=abs(rF(klev+1))
296 CALL SWFRAC(
297 I two,minusone,
298 I myTime,myThid,
299 U swfracb)
300 kp1 = klev+1
301 IF (klev.EQ.Nr) THEN
302 kp1 = klev
303 swfracb(2)=0. _d 0
304 ENDIF
305 DO j=1,sNy
306 DO i=1,sNx
307 gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)
308 & -Qsw(i,j,bi,bj)*(swfracb(1)*maskC(i,j,klev,bi,bj)
309 & -swfracb(2)*maskC(i,j,kp1, bi,bj))
310 & *recip_Cp*recip_rhoConst
311 & *recip_drF(klev)*recip_hFacC(i,j,kLev,bi,bj)
312 ENDDO
313 ENDDO
314 c ENDIF
315 #endif
316
317 #ifdef ALLOW_RBCS
318 if (useRBCS) then
319 call RBCS_ADD_TENDENCY(bi,bj,klev, 1,
320 & myTime, myThid )
321 endif
322 #endif
323
324 #ifdef ALLOW_CLIMTEMP_RELAXATION
325 IF ( tauThetaClimRelax3Dim .NE. 0. ) THEN
326 DO j=1,sNy
327 DO i=1,sNx
328 gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)
329 & -1./tauThetaClimRelax3Dim
330 & *(theta(i,j,klev,bi,bj)-thetaStar(i,j,klev,bi,bj))
331 & *hFacC(i,j,klev,bi,bj)*recip_hFacC(i,j,kLev,bi,bj)
332 ENDDO
333 ENDDO
334 ENDIF
335 #endif
336
337 #ifdef ALLOW_OBCS
338 IF (useOBCS) THEN
339 CALL OBCS_SPONGE_T(
340 I iMin,iMax, jMin,jMax, bi,bj, kLev,
341 I myTime, myThid )
342 ENDIF
343 #endif
344
345 RETURN
346 END
347
348 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
349 CBOP
350 C !ROUTINE: EXTERNAL_FORCING_S
351 C !INTERFACE:
352 SUBROUTINE EXTERNAL_FORCING_S(
353 I iMin,iMax, jMin,jMax, bi,bj, kLev,
354 I myTime, myThid )
355
356 C !DESCRIPTION: \bv
357 C *==========================================================*
358 C | S/R EXTERNAL_FORCING_S
359 C | o Contains problem specific forcing for merid velocity.
360 C *==========================================================*
361 C | Adds terms to gS for forcing by external sources
362 C | e.g. fresh-water flux, climatalogical relaxation, etc ...
363 C *==========================================================*
364 C \ev
365
366 C !USES:
367 IMPLICIT NONE
368 C == Global data ==
369 #include "SIZE.h"
370 #include "EEPARAMS.h"
371 #include "PARAMS.h"
372 #include "GRID.h"
373 #include "DYNVARS.h"
374 #include "FFIELDS.h"
375
376 C !INPUT/OUTPUT PARAMETERS:
377 C == Routine arguments ==
378 C iMin,iMax :: Working range of x-index for applying forcing.
379 C jMin,jMax :: Working range of y-index for applying forcing.
380 C bi,bj :: Current tile indices
381 C kLev :: Current vertical level index
382 C myTime :: Current time in simulation
383 C myThid :: Thread Id number
384 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
385 _RL myTime
386 INTEGER myThid
387
388 C !LOCAL VARIABLES:
389 C == Local variables ==
390 C i,j :: Loop counters
391 C kSurface :: index of surface layer
392 INTEGER i, j
393 INTEGER kSurface
394 CEOP
395
396 IF ( fluidIsAir ) THEN
397 kSurface = 0
398 ELSEIF ( usingPCoords ) THEN
399 kSurface = Nr
400 ELSE
401 kSurface = 1
402 ENDIF
403
404 C-- Forcing term
405 #ifdef ALLOW_AIM
406 IF ( useAIM ) CALL AIM_TENDENCY_APPLY_S(
407 & iMin,iMax, jMin,jMax, bi,bj, kLev,
408 & myTime, myThid )
409 #endif /* ALLOW_AIM */
410
411 #ifdef ALLOW_FIZHI
412 IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_S(
413 & iMin,iMax, jMin,jMax, bi,bj, kLev,
414 & myTime, myThid )
415 #endif /* ALLOW_FIZHI */
416
417 C Add fresh-water in top-layer
418 IF ( kLev .EQ. kSurface ) THEN
419 DO j=1,sNy
420 DO i=1,sNx
421 gS(i,j,kLev,bi,bj)=gS(i,j,kLev,bi,bj)
422 & +surfaceForcingS(i,j,bi,bj)
423 & *recip_drF(kLev)*recip_hFacC(i,j,kLev,bi,bj)
424 ENDDO
425 ENDDO
426 ENDIF
427
428 #ifdef ALLOW_RBCS
429 if (useRBCS) then
430 call RBCS_ADD_TENDENCY(bi,bj,klev, 2,
431 & myTime, myThid )
432 endif
433 #endif
434
435 #ifdef ALLOW_CLIMSALT_RELAXATION
436 IF ( tauSaltClimRelax3Dim .NE. 0. ) THEN
437 DO j=1,sNy
438 DO i=1,sNx
439 gS(i,j,klev,bi,bj) = gS(i,j,klev,bi,bj)
440 & -1./tauSaltClimRelax3Dim
441 & *(salt(i,j,klev,bi,bj)-saltStar(i,j,klev,bi,bj))
442 & *hFacC(i,j,klev,bi,bj)*recip_hFacC(i,j,kLev,bi,bj)
443 ENDDO
444 ENDDO
445 ENDIF
446 #endif
447
448 #ifdef ALLOW_OBCS
449 IF (useOBCS) THEN
450 CALL OBCS_SPONGE_S(
451 I iMin,iMax, jMin,jMax, bi,bj, kLev,
452 I myTime, myThid )
453 ENDIF
454 #endif
455
456 RETURN
457 END

  ViewVC Help
Powered by ViewVC 1.1.22