/[MITgcm]/MITgcm/pkg/thsice/thsice_advdiff.F
ViewVC logotype

Diff of /MITgcm/pkg/thsice/thsice_advdiff.F

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

revision 1.3 by jmc, Sun Apr 8 18:53:16 2007 UTC revision 1.4 by jmc, Mon Apr 9 15:51:22 2007 UTC
# Line 69  C     vTrIce    :: sea-ice volume transp Line 69  C     vTrIce    :: sea-ice volume transp
69  C     afx       :: horizontal advective flux, x direction  C     afx       :: horizontal advective flux, x direction
70  C     afy       :: horizontal advective flux, y direction  C     afy       :: horizontal advective flux, y direction
71  C     iceFrc    :: (new) sea-ice fraction  C     iceFrc    :: (new) sea-ice fraction
 C     iceFld    :: (new) effective sea-ice thickness  
72  C     iceVol    :: temporary array used in advection S/R  C     iceVol    :: temporary array used in advection S/R
73  C     oldVol    :: (old) sea-ice volume  C     oldVol    :: (old) sea-ice volume
74  C     msgBuf    :: Informational/error meesage buffer  C     msgBuf    :: Informational/error meesage buffer
# Line 85  C     msgBuf    :: Informational/error m Line 84  C     msgBuf    :: Informational/error m
84        _RL afy       (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL afy       (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85        _RS maskOce   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RS maskOce   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86        _RL iceFrc    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL iceFrc    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL iceFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)  
87        _RL iceVol    (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL iceVol    (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
88        _RL oldVol    (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL oldVol    (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
89        _RL minIcHeff, minIcArea, r_minArea        _RL r_minArea
90        _RL meanCellArea, areaEpsil, vol_Epsil        _RL meanCellArea, areaEpsil, vol_Epsil
91  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
92        CHARACTER*8 diagName        CHARACTER*8 diagName
# Line 120  C     but what matter here is just to ha Line 118  C     but what matter here is just to ha
118        areaEpsil = 1. _d -10 * meanCellArea        areaEpsil = 1. _d -10 * meanCellArea
119        vol_Epsil = 1. _d -15 * meanCellArea        vol_Epsil = 1. _d -15 * meanCellArea
120    
       minIcArea = iceMaskMin  
       minIcHeff = iceMaskMin*himin  
121        r_minArea = 0. _d 0        r_minArea = 0. _d 0
122        IF ( minIcArea.GT.0. _d 0 ) r_minArea = 1. _d 0 / minIcArea        IF ( iceMaskMin.GT.0. _d 0 ) r_minArea = 1. _d 0 / iceMaskMin
123    
124        thSIce_multiDimAdv = .TRUE.        thSIce_multiDimAdv = .TRUE.
125  #ifdef ALLOW_GENERIC_ADVDIFF  #ifdef ALLOW_GENERIC_ADVDIFF
# Line 366  c           WRITE(6,'(A,1P4E14.6)') 'ICE Line 362  c           WRITE(6,'(A,1P4E14.6)') 'ICE
362          ENDIF          ENDIF
363  #endif  #endif
364    
365  C--   Update Ice Fraction, Ice thickness and snow thickness:  C--   Update Ice Fraction: ensure that fraction is > iceMaskMin & < 1
366  C      and adjust sea-ice state if not enough ice.  C      and adjust Ice thickness and snow thickness accordingly
367          DO j=1,sNy          DO j=1,sNy
368           DO i=1,sNx           DO i=1,sNx
369  C-    store new effective ice-thickness            IF ( iceFrc(i,j) .GT. 1. _d 0 ) THEN
           iceFld(i,j) = iceHeight(i,j,bi,bj)*iceFrc(i,j)  
           IF ( iceFld(i,j) .GE. minIcHeff ) THEN  
 C-    where there is enough ice, ensure that Ice fraction is > minIcArea & < 1  
            IF ( iceFrc(i,j) .GT. 1. _d 0 ) THEN  
370              iceMask(i,j,bi,bj)    = 1. _d 0              iceMask(i,j,bi,bj)    = 1. _d 0
371              iceHeight(i,j,bi,bj)  = iceFld(i,j)              iceHeight(i,j,bi,bj)  = iceHeight(i,j,bi,bj) *iceFrc(i,j)
372              snowHeight(i,j,bi,bj) = snowHeight(i,j,bi,bj)*iceFrc(i,j)              snowHeight(i,j,bi,bj) = snowHeight(i,j,bi,bj)*iceFrc(i,j)
373             ELSEIF ( iceFrc(i,j) .LT. minIcArea ) THEN            ELSEIF ( iceFrc(i,j) .LT. iceMaskMin ) THEN
374              iceMask(i,j,bi,bj)    = minIcArea              iceMask(i,j,bi,bj)    = iceMaskMin
375              iceHeight(i,j,bi,bj)  = iceFld(i,j)*r_minArea              iceHeight(i,j,bi,bj)  = iceHeight(i,j,bi,bj)
376         &                             *iceFrc(i,j)*r_minArea
377              snowHeight(i,j,bi,bj) = snowHeight(i,j,bi,bj)              snowHeight(i,j,bi,bj) = snowHeight(i,j,bi,bj)
378       &                             *iceFrc(i,j)*r_minArea       &                             *iceFrc(i,j)*r_minArea
            ELSE  
             iceMask(i,j,bi,bj)    = iceFrc(i,j)  
            ENDIF  
379            ELSE            ELSE
380                iceMask(i,j,bi,bj)    = iceFrc(i,j)
381              ENDIF
382             ENDDO
383            ENDDO
384    C-     adjust sea-ice state if not enough ice.
385            DO j=1,sNy
386             DO i=1,sNx
387              IF ( iceHeight(i,j,bi,bj).LT.himin ) THEN
388  C-    Not enough ice, melt the tiny amount of snow & ice:  C-    Not enough ice, melt the tiny amount of snow & ice:
389  C     and return frsh-water, salt & energy to the ocean (flx > 0 = into ocean)  C     and return fresh-water, salt & energy to the ocean (flx > 0 = into ocean)
390  C- -  Note: using 1rst.Order Upwind, I can get the same results as when  C- -  Note: using 1rst.Order Upwind, I can get the same results as when
391  C     using seaice_advdiff (with SEAICEadvScheme=1) providing I comment  C     using seaice_advdiff (with SEAICEadvScheme=1) providing I comment
392  C     out the following lines (and then loose conservation).  C     out the following lines (and then loose conservation).
393  C- -  C- -
394              oceFWfx(i,j,bi,bj) = ( rhos*snowHeight(i,j,bi,bj)              oceFWfx(i,j,bi,bj) =  ( rhos*snowHeight(i,j,bi,bj)
395       &                            +rhoi*iceHeight(i,j,bi,bj) )       &                             +rhoi*iceHeight(i,j,bi,bj)
396       &                                      *iceFrc(i,j)/thSIce_deltaT       &                            )*iceMask(i,j,bi,bj)/thSIce_deltaT
397              oceSflx(i,j,bi,bj) =saltice*rhoi*iceFld(i,j)/thSIce_deltaT              oceSflx(i,j,bi,bj) =    rhoi*iceHeight(i,j,bi,bj)*saltice
398              oceQnet(i,j,bi,bj) = -qsnow*rhos*snowHeight(i,j,bi,bj)       &                             *iceMask(i,j,bi,bj)/thSIce_deltaT
399       &                                      *iceFrc(i,j)/thSIce_deltaT              oceQnet(i,j,bi,bj) = -( rhos*snowHeight(i,j,bi,bj)*qsnow
400       &                         -( Qice1(i,j,bi,bj)       &                             +rhoi*iceHeight(i,j,bi,bj)
401       &                           +Qice2(i,j,bi,bj) )*0.5 _d 0       &                                  *( Qice1(i,j,bi,bj)
402       &                                 *rhoi*iceFld(i,j)/thSIce_deltaT       &                                    +Qice2(i,j,bi,bj) )*0.5 _d 0
403         &                            )*iceMask(i,j,bi,bj)/thSIce_deltaT
404  C- -  C- -
405  c           flx2oc (i,j) = flx2oc (i,j) +  c           flx2oc (i,j) = flx2oc (i,j) +
406  c           frw2oc (i,j) = frw2oc (i,j) +  c           frw2oc (i,j) = frw2oc (i,j) +

Legend:
Removed from v.1.3  
changed lines
  Added in v.1.4

  ViewVC Help
Powered by ViewVC 1.1.22