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

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

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


Revision 1.36 - (hide annotations) (download)
Mon Jan 2 21:17:02 2006 UTC (18 years, 5 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint58
Changes since 1.35: +1 -30 lines
o Fix I/O inconsistency in pkg/rbcs: replace precFloat32 by readBinaryPrec
o Remove 3-dim. relaxation code from pkg/exf (now use only pkg/rbcs)
o Thanks to Tom Haine for testing!

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

  ViewVC Help
Powered by ViewVC 1.1.22