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

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

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

revision 1.17 by mlosch, Wed Sep 25 19:36:50 2002 UTC revision 1.30 by heimbach, Tue Mar 1 18:55:13 2005 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4    #include "PACKAGES_CONFIG.h"
5  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
6    #ifdef ALLOW_OBCS
7    # include "OBCS_OPTIONS.h"
8    #endif
9    
10  CBOP  CBOP
11  C     !ROUTINE: EXTERNAL_FORCING_U  C     !ROUTINE: EXTERNAL_FORCING_U
# Line 48  C     number of surface interface layer Line 52  C     number of surface interface layer
52        INTEGER kSurface        INTEGER kSurface
53  CEOP  CEOP
54    
55        if ( buoyancyRelation .eq. 'OCEANICP' ) then        IF ( fluidIsAir ) THEN
56           kSurface = 0
57          ELSEIF ( usingPCoords ) THEN
58         kSurface = Nr         kSurface = Nr
59        else        ELSE
60         kSurface = 1         kSurface = 1
61        endif        ENDIF
62    
63  C--   Forcing term  C--   Forcing term
64    #ifdef ALLOW_AIM
65          IF ( useAIM ) CALL AIM_TENDENCY_APPLY_U(
66         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
67         &                      myCurrentTime, myThid )
68    #endif /* ALLOW_AIM */
69    C AMM
70    #ifdef ALLOW_FIZHI
71          IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_U(
72         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
73         &                      myCurrentTime, myThid )
74    #endif /* ALLOW_FIZHI */
75    C AMM
76    
77  C     Add windstress momentum impulse into the top-layer  C     Add windstress momentum impulse into the top-layer
78        IF ( kLev .EQ. kSurface ) THEN        IF ( kLev .EQ. kSurface ) THEN
79         DO j=jMin,jMax         DO j=jMin,jMax
80          DO i=iMin,iMax          DO i=iMin,iMax
81           gU(i,j,kLev,bi,bj) = gU(i,j,kLev,bi,bj)           gU(i,j,kLev,bi,bj) = gU(i,j,kLev,bi,bj)
82       &   +foFacMom*surfaceTendencyU(i,j,bi,bj)       &   +foFacMom*surfaceForcingU(i,j,bi,bj)
83       &   *_maskW(i,j,kLev,bi,bj)       &   *recip_drF(kLev)*recip_hFacW(i,j,kLev,bi,bj)
84          ENDDO          ENDDO
85         ENDDO         ENDDO
86        ENDIF        ENDIF
87    
88    #if (defined (ALLOW_TAU_EDDY))
89           CALL TAUEDDY_EXTERNAL_FORCING_U(
90         I           iMin, iMax, jMin, jMax,bi,bj,kLev,
91         I           myCurrentTime,myThid)
92    #endif
93    
94  #if (defined (ALLOW_OBCS) && defined (ALLOW_OBCS_SPONGE))  #if (defined (ALLOW_OBCS) && defined (ALLOW_OBCS_SPONGE))
95        IF (useOBCS) THEN        IF (useOBCS) THEN
96         CALL OBCS_SPONGE_U(         CALL OBCS_SPONGE_U(
# Line 121  C     number of surface interface layer Line 146  C     number of surface interface layer
146        INTEGER kSurface        INTEGER kSurface
147  CEOP  CEOP
148    
149        if ( buoyancyRelation .eq. 'OCEANICP' ) then        IF ( fluidIsAir ) THEN
150           kSurface = 0
151          ELSEIF ( usingPCoords ) THEN
152         kSurface = Nr         kSurface = Nr
153        else        ELSE
154         kSurface = 1         kSurface = 1
155        endif        ENDIF
156    
157  C--   Forcing term  C--   Forcing term
158    #ifdef ALLOW_AIM
159          IF ( useAIM ) CALL AIM_TENDENCY_APPLY_V(
160         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
161         &                      myCurrentTime, myThid )
162    #endif /* ALLOW_AIM */
163    
164    C AMM
165    #ifdef ALLOW_FIZHI
166          IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_V(
167         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
168         &                      myCurrentTime, myThid )
169    #endif /* ALLOW_FIZHI */
170    C AMM
171  C     Add windstress momentum impulse into the top-layer  C     Add windstress momentum impulse into the top-layer
172        IF ( kLev .EQ. kSurface ) THEN        IF ( kLev .EQ. kSurface ) THEN
173         DO j=jMin,jMax         DO j=jMin,jMax
174          DO i=iMin,iMax          DO i=iMin,iMax
175           gV(i,j,kLev,bi,bj) = gV(i,j,kLev,bi,bj)           gV(i,j,kLev,bi,bj) = gV(i,j,kLev,bi,bj)
176       &   +foFacMom*surfaceTendencyV(i,j,bi,bj)       &   +foFacMom*surfaceForcingV(i,j,bi,bj)
177       &   *_maskS(i,j,kLev,bi,bj)       &   *recip_drF(kLev)*recip_hFacS(i,j,kLev,bi,bj)
178          ENDDO          ENDDO
179         ENDDO         ENDDO
180        ENDIF        ENDIF
181    
182    #if (defined (ALLOW_TAU_EDDY))
183           CALL TAUEDDY_EXTERNAL_FORCING_V(
184         I           iMin, iMax, jMin, jMax,bi,bj,kLev,
185         I           myCurrentTime,myThid)
186    #endif
187    
188  #if (defined (ALLOW_OBCS) && defined (ALLOW_OBCS_SPONGE))  #if (defined (ALLOW_OBCS) && defined (ALLOW_OBCS_SPONGE))
189        IF (useOBCS) THEN        IF (useOBCS) THEN
190         CALL OBCS_SPONGE_V(         CALL OBCS_SPONGE_V(
# Line 174  C     == Global data == Line 220  C     == Global data ==
220  #include "GRID.h"  #include "GRID.h"
221  #include "DYNVARS.h"  #include "DYNVARS.h"
222  #include "FFIELDS.h"  #include "FFIELDS.h"
 #ifdef SHORTWAVE_HEATING  
       integer two  
       _RL minusone  
       parameter (two=2,minusone=-1.)  
       _RL swfracb(two)  
 #endif  
223    
224  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
225  C     == Routine arguments ==  C     == Routine arguments ==
# Line 199  C     Loop counters Line 239  C     Loop counters
239        INTEGER I, J        INTEGER I, J
240  C     number of surface interface layer  C     number of surface interface layer
241        INTEGER kSurface        INTEGER kSurface
242    #ifdef SHORTWAVE_HEATING
243          integer two
244          _RL minusone
245          parameter (two=2,minusone=-1.)
246          _RL swfracb(two)
247          INTEGER kp1
248    #endif
249  CEOP  CEOP
250    
251        if ( buoyancyRelation .eq. 'OCEANICP' ) then        IF ( fluidIsAir ) THEN
252           kSurface = 0
253          ELSEIF ( usingPCoords ) THEN
254         kSurface = Nr         kSurface = Nr
255        else        ELSE
256         kSurface = 1         kSurface = 1
257        endif        ENDIF
258    
259  C--   Forcing term  C--   Forcing term
260    #ifdef ALLOW_AIM
261          IF ( useAIM ) CALL AIM_TENDENCY_APPLY_T(
262         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
263         &                      myCurrentTime, myThid )
264    #endif /* ALLOW_AIM */
265    
266    C AMM
267    #ifdef ALLOW_FIZHI
268          IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_T(
269         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
270         &                      myCurrentTime, myThid )
271    #endif /* ALLOW_FIZHI */
272    C AMM
273    
274  C     Add heat in top-layer  C     Add heat in top-layer
275        IF ( kLev .EQ. kSurface ) THEN        IF ( kLev .EQ. kSurface ) THEN
276         DO j=jMin,jMax         DO j=jMin,jMax
277          DO i=iMin,iMax          DO i=iMin,iMax
278           gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)           gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)
279       &     +maskC(i,j,kLev,bi,bj)*surfaceTendencyT(i,j,bi,bj)       &     +surfaceForcingT(i,j,bi,bj)
280         &     *recip_drF(kLev)*recip_hFacC(i,j,kLev,bi,bj)
281          ENDDO          ENDDO
282         ENDDO         ENDDO
283        ENDIF        ENDIF
284    
285  #ifdef SHORTWAVE_HEATING  #ifdef SHORTWAVE_HEATING
286  C Penetrating SW radiation  C Penetrating SW radiation
287          kp1 = klev+1
288        swfracb(1)=abs(rF(klev))        swfracb(1)=abs(rF(klev))
289        swfracb(2)=abs(rF(klev+1))        swfracb(2)=abs(rF(klev+1))
290        call SWFRAC(        CALL SWFRAC(
291       I     two,minusone,       I     two,minusone,
292       I     myCurrentTime,myThid,       I     myCurrentTime,myThid,
293       O     swfracb)       U     swfracb)
294          IF (klev.EQ.Nr) THEN
295            kp1 = klev
296            swfracb(2)=0. _d 0
297          ENDIF
298        DO j=jMin,jMax        DO j=jMin,jMax
299         DO i=iMin,iMax         DO i=iMin,iMax
300          gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)          gT(i,j,klev,bi,bj) = gT(i,j,klev,bi,bj)
301       &   -maskC(i,j,klev,bi,bj)*Qsw(i,j,bi,bj)*(swfracb(1)-swfracb(2))       &   -Qsw(i,j,bi,bj)*(swfracb(1)*maskC(i,j,klev,bi,bj)
302       &    *recip_Cp*recip_rhoConst*recip_drF(klev)       &                   -swfracb(2)*maskC(i,j,kp1, bi,bj))
303         &    *recip_Cp*recip_rhoConst
304         &    *recip_drF(klev)*recip_hFacC(i,j,kLev,bi,bj)
305         ENDDO         ENDDO
306        ENDDO        ENDDO
307  #endif  #endif
# Line 291  C     number of surface interface layer Line 362  C     number of surface interface layer
362        INTEGER kSurface        INTEGER kSurface
363  CEOP  CEOP
364    
365        if ( buoyancyRelation .eq. 'OCEANICP' ) then        IF ( fluidIsAir ) THEN
366           kSurface = 0
367          ELSEIF ( usingPCoords ) THEN
368         kSurface = Nr         kSurface = Nr
369        else        ELSE
370         kSurface = 1         kSurface = 1
371        endif        ENDIF
   
372    
373  C--   Forcing term  C--   Forcing term
374    #ifdef ALLOW_AIM
375          IF ( useAIM ) CALL AIM_TENDENCY_APPLY_S(
376         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
377         &                      myCurrentTime, myThid )
378    #endif /* ALLOW_AIM */
379    
380    C AMM
381    #ifdef ALLOW_FIZHI
382          IF ( useFIZHI ) CALL FIZHI_TENDENCY_APPLY_S(
383         &                      iMin,iMax, jMin,jMax, bi,bj, kLev,
384         &                      myCurrentTime, myThid )
385    #endif /* ALLOW_FIZHI */
386    C AMM
387    
388  C     Add fresh-water in top-layer  C     Add fresh-water in top-layer
389        IF ( kLev .EQ. kSurface ) THEN        IF ( kLev .EQ. kSurface ) THEN
390         DO j=jMin,jMax         DO j=jMin,jMax
391          DO i=iMin,iMax          DO i=iMin,iMax
392           gS(i,j,kLev,bi,bj)=gS(i,j,kLev,bi,bj)           gS(i,j,kLev,bi,bj)=gS(i,j,kLev,bi,bj)
393       &   +maskC(i,j,kLev,bi,bj)*surfaceTendencyS(i,j,bi,bj)       &     +surfaceForcingS(i,j,bi,bj)
394         &     *recip_drF(kLev)*recip_hFacC(i,j,kLev,bi,bj)
395          ENDDO          ENDDO
396         ENDDO         ENDDO
397        ENDIF        ENDIF
Legend:
Removed from v.1.17  
changed lines
  Added in v.1.30
  ViewVC Help
Powered by ViewVC 1.1.22