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

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

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

revision 1.36 by jmc, Sat Jul 22 03:53:13 2006 UTC revision 1.41 by jmc, Thu Aug 23 22:25:41 2007 UTC
# Line 30  C     === Global variables === Line 30  C     === Global variables ===
30  #include "SURFACE.h"  #include "SURFACE.h"
31  #ifdef ALLOW_SEAICE  #ifdef ALLOW_SEAICE
32  #include "SEAICE.h"  #include "SEAICE.h"
33    #include "SEAICE_PARAMS.h"
34  #endif /* ALLOW_SEAICE */  #endif /* ALLOW_SEAICE */
35  #ifdef ALLOW_SHELFICE  #ifdef ALLOW_SHELFICE
36  #include "SHELFICE.h"  #include "SHELFICE.h"
# Line 58  CEOP Line 59  CEOP
59  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
60        _RL tmpFac        _RL tmpFac
61  #endif /* ALLOW_DIAGNOSTICS */  #endif /* ALLOW_DIAGNOSTICS */
62    #ifdef ALLOW_SALT_PLUME
63          _RL saltPlume
64    #endif /* ALLOW_SALT_PLUME */
65    
66        IF ( usingPCoords ) THEN        IF ( usingPCoords ) THEN
67         ks        = Nr         ks        = Nr
# Line 70  C---+----1----+----2----+----3----+----4 Line 74  C---+----1----+----2----+----3----+----4
74        IF ( doThetaClimRelax .OR. doSaltClimRelax ) THEN        IF ( doThetaClimRelax .OR. doSaltClimRelax ) THEN
75  C--   Start with surface restoring term :  C--   Start with surface restoring term :
76    
        DO j = jMin, jMax  
         DO i = iMin, iMax  
77  #ifdef ALLOW_SEAICE  #ifdef ALLOW_SEAICE
78           IF ( useSEAICE ) THEN
79  C     Do not restore under sea-ice  C     Do not restore under sea-ice
80            DO j = jMin, jMax
81             DO i = iMin, iMax
82  C     Heat Flux (restoring term) :  C     Heat Flux (restoring term) :
83            surfaceForcingT(i,j,bi,bj) =            surfaceForcingT(i,j,bi,bj) =
84       &      -lambdaThetaClimRelax(i,j,bi,bj) * (1-AREA(i,j,1,bi,bj))       &      -lambdaThetaClimRelax(i,j,bi,bj)*(1.-AREA(i,j,1,bi,bj))
85       &         *(theta(i,j,ks,bi,bj)-SST(i,j,bi,bj))       &         *(theta(i,j,ks,bi,bj)-SST(i,j,bi,bj))
86       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)
87  C     Salt Flux (restoring term) :  C     Salt Flux (restoring term) :
88            surfaceForcingS(i,j,bi,bj) =            surfaceForcingS(i,j,bi,bj) =
89       &      -lambdaSaltClimRelax(i,j,bi,bj) * (1-AREA(i,j,1,bi,bj))       &      -lambdaSaltClimRelax(i,j,bi,bj) *(1.-AREA(i,j,1,bi,bj))
90       &         *(salt(i,j,ks,bi,bj)-SSS(i,j,bi,bj))       &         *(salt(i,j,ks,bi,bj)-SSS(i,j,bi,bj))
91       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)
92  #else /* ifndef ALLOW_SEAICE */           ENDDO
93            ENDDO
94           ELSE
95    #endif /* ALLOW_SEAICE */
96            DO j = jMin, jMax
97             DO i = iMin, iMax
98  C     Heat Flux (restoring term) :  C     Heat Flux (restoring term) :
99            surfaceForcingT(i,j,bi,bj) =            surfaceForcingT(i,j,bi,bj) =
100       &      -lambdaThetaClimRelax(i,j,bi,bj)       &      -lambdaThetaClimRelax(i,j,bi,bj)
# Line 95  C     Salt Flux (restoring term) : Line 105  C     Salt Flux (restoring term) :
105       &      -lambdaSaltClimRelax(i,j,bi,bj)       &      -lambdaSaltClimRelax(i,j,bi,bj)
106       &         *(salt(i,j,ks,bi,bj)-SSS(i,j,bi,bj))       &         *(salt(i,j,ks,bi,bj)-SSS(i,j,bi,bj))
107       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)       &         *drF(ks)*_hFacC(i,j,ks,bi,bj)
108  #endif /* ALLOW_SEAICE */           ENDDO
109          ENDDO          ENDDO
110         ENDDO  #ifdef ALLOW_SEAICE
111           ENDIF
112    #endif /* ALLOW_SEAICE */
113    
114  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
115  #ifdef NONLIN_FRSURF  #ifdef NONLIN_FRSURF
# Line 134  C     of surf.thickness (hFac): account Line 146  C     of surf.thickness (hFac): account
146         IF ( useDiagnostics ) THEN         IF ( useDiagnostics ) THEN
147    
148  C     tRelax (temperature relaxation [W/m2], positive <-> increasing Theta)  C     tRelax (temperature relaxation [W/m2], positive <-> increasing Theta)
149          tmpFac = HeatCapacity_Cp*recip_horiVertRatio*rhoConst          tmpFac = HeatCapacity_Cp*rUnit2mass
150          CALL DIAGNOSTICS_SCALE_FILL(          CALL DIAGNOSTICS_SCALE_FILL(
151       &           surfaceForcingT(1-OLx,1-OLy,bi,bj),tmpFac,1,       &           surfaceForcingT(1-OLx,1-OLy,bi,bj),tmpFac,1,
152       &                             'TRELAX  ',0, 1,2,bi,bj,myThid)       &                             'TRELAX  ',0, 1,2,bi,bj,myThid)
153    
154  C     sRelax (salt relaxation [g/m2/s], positive <-> increasing Salt)  C     sRelax (salt relaxation [g/m2/s], positive <-> increasing Salt)
155          tmpFac = recip_horiVertRatio*rhoConst          tmpFac = rUnit2mass
156          CALL DIAGNOSTICS_SCALE_FILL(          CALL DIAGNOSTICS_SCALE_FILL(
157       &           surfaceForcingS(1-OLx,1-OLy,bi,bj),tmpFac,1,       &           surfaceForcingS(1-OLx,1-OLy,bi,bj),tmpFac,1,
158       &                             'SRELAX  ',0, 1,2,bi,bj,myThid)       &                             'SRELAX  ',0, 1,2,bi,bj,myThid)
# Line 170  C--   Surface Fluxes : Line 182  C--   Surface Fluxes :
182    
183  C     Zonal wind stress fu:  C     Zonal wind stress fu:
184            surfaceForcingU(i,j,bi,bj) =            surfaceForcingU(i,j,bi,bj) =
185       &      fu(i,j,bi,bj)*horiVertRatio*recip_rhoConst       &      fu(i,j,bi,bj)*mass2rUnit
186  C     Meridional wind stress fv:  C     Meridional wind stress fv:
187            surfaceForcingV(i,j,bi,bj) =            surfaceForcingV(i,j,bi,bj) =
188       &      fv(i,j,bi,bj)*horiVertRatio*recip_rhoConst       &      fv(i,j,bi,bj)*mass2rUnit
189  C     Net heat flux Qnet:  C     Net heat flux Qnet:
190            surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)            surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
191       &       - ( Qnet(i,j,bi,bj)       &       - ( Qnet(i,j,bi,bj)
192  #ifdef SHORTWAVE_HEATING  #ifdef SHORTWAVE_HEATING
193       &          -Qsw(i,j,bi,bj)       &          -Qsw(i,j,bi,bj)
194  #endif  #endif
195       &         ) *recip_Cp*horiVertRatio*recip_rhoConst       &         ) *recip_Cp*mass2rUnit
196  C     Net Salt Flux :  C     Net Salt Flux :
197            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
198       &      -saltFlux(i,j,bi,bj)*horiVertRatio*recip_rhoConst       &      -saltFlux(i,j,bi,bj)*mass2rUnit
199    #ifdef ALLOW_SALT_PLUME
200    C saltPlume is the amount of salt rejected by ice while freezing;
201    C it is here subtracted from surfaceForcingS and will be redistributed
202    C to multiple vertical levels later on as per Duffy et al. (GRL 1999)
203              saltPlume = 0.
204    #ifdef ALLOW_SEAICE
205              IF ( saltFlux(i,j,bi,bj) .GT. 0. .AND.
206         &         salt(i,j,ks,bi,bj)  .GT. SEAICE_salinity ) THEN
207               saltPlume = (salt(i,j,ks,bi,bj)-SEAICE_salinity) *
208         &          saltFlux(i,j,bi,bj) / salt(i,j,ks,bi,bj)
209              ENDIF
210    #endif /* ALLOW_SEAICE */
211              surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
212         &      -saltPlume*mass2rUnit
213    #endif /* ALLOW_SALT_PLUME */
214    
215           ENDDO           ENDDO
216        ENDDO        ENDDO
# Line 202  C      (needed for SSH forcing, whether Line 229  C      (needed for SSH forcing, whether
229        ENDIF        ENDIF
230    
231  #ifdef EXACT_CONSERV  #ifdef EXACT_CONSERV
232  c NB: synchronous time step: PmEpR lag 1 time step behind EmPmR  C NB: synchronous time step: PmEpR lag 1 time step behind EmPmR
233  c     to stay consitent with volume change (=d/dt etaH).  C     to stay consitent with volume change (=d/dt etaH).
234        IF ( staggerTimeStep ) THEN        IF ( staggerTimeStep ) THEN
235          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
236           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
# Line 215  c     to stay consitent with volume chan Line 242  c     to stay consitent with volume chan
242        IF ( (nonlinFreeSurf.GT.0 .OR. usingPCoords)        IF ( (nonlinFreeSurf.GT.0 .OR. usingPCoords)
243       &     .AND. useRealFreshWaterFlux ) THEN       &     .AND. useRealFreshWaterFlux ) THEN
244    
245  c-  NonLin_FrSurf and RealFreshWaterFlux : PmEpR effectively changes  C--   NonLin_FrSurf and RealFreshWaterFlux : PmEpR effectively changes
246  c   the water column height ; temp., salt, (tracer) flux associated  C     the water column height ; temp., salt, (tracer) flux associated
247  c   with this input/output of water is added here to the surface tendency.  C     with this input/output of water is added here to the surface tendency.
 c  
248    
249         IF (temp_EvPrRn.NE.UNSET_RL) THEN         IF (temp_EvPrRn.NE.UNSET_RL) THEN
250          DO j = jMin, jMax          DO j = jMin, jMax
251           DO i = iMin, iMax           DO i = iMin, iMax
252            surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)            surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
253       &      + PmEpR(i,j,bi,bj)       &      + PmEpR(i,j,bi,bj)
254       &         *( temp_EvPrRn - theta(i,j,ks,bi,bj) )       &          *( temp_EvPrRn - theta(i,j,ks,bi,bj) )
255       &         *convertEmP2rUnit       &          *convertEmP2rUnit
256           ENDDO           ENDDO
257          ENDDO          ENDDO
258         ENDIF         ENDIF
# Line 236  c Line 262  c
262           DO i = iMin, iMax           DO i = iMin, iMax
263            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
264       &      + PmEpR(i,j,bi,bj)       &      + PmEpR(i,j,bi,bj)
265       &         *( salt_EvPrRn - salt(i,j,ks,bi,bj) )       &          *( salt_EvPrRn - salt(i,j,ks,bi,bj) )
266       &         *convertEmP2rUnit       &          *convertEmP2rUnit
267           ENDDO           ENDDO
268          ENDDO          ENDDO
269         ENDIF         ENDIF
# Line 248  C---+----1----+----2----+----3----+----4 Line 274  C---+----1----+----2----+----3----+----4
274        IF (.TRUE.) THEN        IF (.TRUE.) THEN
275  #endif /* EXACT_CONSERV */  #endif /* EXACT_CONSERV */
276    
277  c- EmPmR does not really affect the water column height (for tracer budget)  C--   EmPmR does not really affect the water column height (for tracer budget)
278  c   and is converted to a salt tendency.  C     and is converted to a salt tendency.
279    
280         IF (convertFW2Salt .EQ. -1.) THEN         IF (convertFW2Salt .EQ. -1.) THEN
281  c- converts EmPmR to salinity tendency using surface local salinity  C-    use local surface tracer field to calculate forcing term:
282          DO j = jMin, jMax  
283           DO i = iMin, iMax          IF (temp_EvPrRn.NE.UNSET_RL) THEN
284            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)  C     account for Rain/Evap heat content (temp_EvPrRn) using local SST
285       &      + EmPmR(i,j,bi,bj)*salt(i,j,ks,bi,bj)           DO j = jMin, jMax
286       &         *convertEmP2rUnit            DO i = iMin, iMax
287               surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
288         &       + EmPmR(i,j,bi,bj)
289         &           *( theta(i,j,ks,bi,bj) - temp_EvPrRn )
290         &           *convertEmP2rUnit
291              ENDDO
292           ENDDO           ENDDO
293          ENDDO          ENDIF
294            IF (salt_EvPrRn.NE.UNSET_RL) THEN
295    C     converts EmPmR to salinity tendency using surface local salinity
296             DO j = jMin, jMax
297              DO i = iMin, iMax
298               surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
299         &       + EmPmR(i,j,bi,bj)
300         &           *( salt(i,j,ks,bi,bj) - salt_EvPrRn )
301         &           *convertEmP2rUnit
302              ENDDO
303             ENDDO
304            ENDIF
305    
306         ELSE         ELSE
307  c- converts EmPmR to virtual salt flux using uniform salinity (default=35)  C-    use uniform tracer value to calculate forcing term:
308          DO j = jMin, jMax  
309           DO i = iMin, iMax          IF (temp_EvPrRn.NE.UNSET_RL) THEN
310            surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)  C     account for Rain/Evap heat content (temp_EvPrRn) assuming uniform SST (=tRef)
311       &      + EmPmR(i,j,bi,bj)*convertFW2Salt           DO j = jMin, jMax
312       &         *convertEmP2rUnit            DO i = iMin, iMax
313               surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
314         &       + EmPmR(i,j,bi,bj)
315         &           *( tRef(ks) - temp_EvPrRn )
316         &           *convertEmP2rUnit
317              ENDDO
318           ENDDO           ENDDO
319          ENDDO          ENDIF
320            IF (salt_EvPrRn.NE.UNSET_RL) THEN
321    C     converts EmPmR to virtual salt flux using uniform salinity (default=35)
322             DO j = jMin, jMax
323              DO i = iMin, iMax
324               surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
325         &       + EmPmR(i,j,bi,bj)
326         &           *( convertFW2Salt - salt_EvPrRn )
327         &           *convertEmP2rUnit
328              ENDDO
329             ENDDO
330            ENDIF
331    
332    C-    end local-surface-tracer / uniform-value distinction
333         ENDIF         ENDIF
334    
335        ENDIF        ENDIF
Legend:
Removed from v.1.36  
changed lines
  Added in v.1.41
  ViewVC Help
Powered by ViewVC 1.1.22