Changeset 5357 for branches/2015/dev_r5218_CNRS17_coupling/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90

Timestamp:

2015-06-04T20:39:20+02:00 (9 years ago)

Author:

clem

Message:

LIM3: change the interface between the ice and atm for both coupled and forced modes. Some work still needs to be done to deal with sublimation in coupled mode.

File:

• branches/2015/dev_r5218_CNRS17_coupling/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90 (modified) (9 diffs)

branches/2015/dev_r5218_CNRS17_coupling/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90

@@ -29 +29 @@
    PRIVATE
 
-   PUBLIC   lim_thd_dh   ! called by lim_thd
+   PUBLIC   lim_thd_dh      ! called by lim_thd
+   PUBLIC   lim_thd_snwblow ! called in sbcblk/sbccpl and here
+
+   INTERFACE lim_thd_snwblow
+      MODULE PROCEDURE lim_thd_snwblow_1d, lim_thd_snwblow_2d
+   END INTERFACE
 
    !!----------------------------------------------------------------------
@@ -71 +76 @@
       REAL(wp) ::   zfdum       
       REAL(wp) ::   zfracs       ! fractionation coefficient for bottom salt entrapment
-      REAL(wp) ::   zcoeff       ! dummy argument for snowfall partitioning over ice and leads
-      REAL(wp) ::   zs_snic  ! snow-ice salinity
+      REAL(wp) ::   zs_snic      ! snow-ice salinity
       REAL(wp) ::   zswi1        ! switch for computation of bottom salinity
       REAL(wp) ::   zswi12       ! switch for computation of bottom salinity
@@ -103 +107 @@
       REAL(wp), POINTER, DIMENSION(:) ::   zqh_s       ! total snow heat content (J.m-2)
       REAL(wp), POINTER, DIMENSION(:) ::   zq_s        ! total snow enthalpy     (J.m-3)
+      REAL(wp), POINTER, DIMENSION(:) ::   zsnw        ! distribution of snow after wind blowing
 
       REAL(wp) :: zswitch_sal
@@ -118 +123 @@
 
       CALL wrk_alloc( jpij, zqprec, zq_su, zq_bo, zf_tt, zq_rema )
-      CALL wrk_alloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s )
+      CALL wrk_alloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s, zsnw )
       CALL wrk_alloc( jpij, nlay_i, zdeltah, zh_i )
       CALL wrk_alloc( jpij, nlay_i, icount )
@@ -219 +224 @@
 
       zdeltah(:,:) = 0._wp
+      CALL lim_thd_snwblow( 1. - at_i_1d, zsnw ) ! snow distribution over ice after wind blowing
       DO ji = kideb, kiut
          !-----------
@@ -224 +230 @@
          !-----------
          ! thickness change
-         zcoeff = ( 1._wp - ( 1._wp - at_i_1d(ji) )**rn_betas ) / at_i_1d(ji) 
-         zdh_s_pre(ji) = zcoeff * sprecip_1d(ji) * rdt_ice * r1_rhosn
-         ! enthalpy of the precip (>0, J.m-3) (tatm_ice is now in K)
-         zqprec   (ji) = rhosn * ( cpic * ( rt0 - MIN( tatm_ice_1d(ji), rt0_snow) ) + lfus )   
+         zdh_s_pre(ji) = zsnw(ji) * sprecip_1d(ji) * rdt_ice * r1_rhosn / at_i_1d(ji)
+         ! enthalpy of the precip (>0, J.m-3)
+         zqprec   (ji) = - qprec_ice_1d(ji)   
          IF( sprecip_1d(ji) == 0._wp ) zqprec(ji) = 0._wp
          ! heat flux from snow precip (>0, W.m-2)
@@ -280 +285 @@
       ! clem comment: ice should also sublimate
       zdeltah(:,:) = 0._wp
-      IF( ln_cpl ) THEN
-         ! coupled mode: sublimation already included in emp_ice (to do in limsbc_ice)
-         zdh_s_sub(:)      =  0._wp 
-      ELSE
-         ! forced  mode: snow thickness change due to sublimation
-         DO ji = kideb, kiut
-            zdh_s_sub(ji)  =  MAX( - ht_s_1d(ji) , - qla_ice_1d(ji) / ( rhosn * lsub ) * rdt_ice )
-            ! Heat flux by sublimation [W.m-2], < 0
-            !      sublimate first snow that had fallen, then pre-existing snow
-            zdeltah(ji,1)  = MAX( zdh_s_sub(ji), - zdh_s_pre(ji) )
-            hfx_sub_1d(ji) = hfx_sub_1d(ji) + ( zdeltah(ji,1) * zqprec(ji) + ( zdh_s_sub(ji) - zdeltah(ji,1) ) * q_s_1d(ji,1)  &
-               &                              ) * a_i_1d(ji) * r1_rdtice
-            ! Mass flux by sublimation
-            wfx_sub_1d(ji) =  wfx_sub_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice
-            ! new snow thickness
-            ht_s_1d(ji)    =  MAX( 0._wp , ht_s_1d(ji) + zdh_s_sub(ji) )
-            ! update precipitations after sublimation and correct sublimation
-            zdh_s_pre(ji) = zdh_s_pre(ji) + zdeltah(ji,1)
-            zdh_s_sub(ji) = zdh_s_sub(ji) - zdeltah(ji,1)
-         END DO
-      ENDIF
-
+      ! coupled mode: sublimation is set to 0 (evap_ice = 0) until further notice
+      ! forced  mode: snow thickness change due to sublimation
+      DO ji = kideb, kiut
+         zdh_s_sub(ji)  =  MAX( - ht_s_1d(ji) , - evap_ice_1d(ji) * r1_rhosn * rdt_ice )
+         ! Heat flux by sublimation [W.m-2], < 0
+         !      sublimate first snow that had fallen, then pre-existing snow
+         zdeltah(ji,1)  = MAX( zdh_s_sub(ji), - zdh_s_pre(ji) )
+         hfx_sub_1d(ji) = hfx_sub_1d(ji) + ( zdeltah(ji,1) * zqprec(ji) + ( zdh_s_sub(ji) - zdeltah(ji,1) ) * q_s_1d(ji,1)  &
+            &                              ) * a_i_1d(ji) * r1_rdtice
+         ! Mass flux by sublimation
+         wfx_sub_1d(ji) =  wfx_sub_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice
+         ! new snow thickness
+         ht_s_1d(ji)    =  MAX( 0._wp , ht_s_1d(ji) + zdh_s_sub(ji) )
+         ! update precipitations after sublimation and correct sublimation
+         zdh_s_pre(ji) = zdh_s_pre(ji) + zdeltah(ji,1)
+         zdh_s_sub(ji) = zdh_s_sub(ji) - zdeltah(ji,1)
+      END DO
+      
       ! --- Update snow diags --- !
       DO ji = kideb, kiut
@@ -689 +690 @@
       
       CALL wrk_dealloc( jpij, zqprec, zq_su, zq_bo, zf_tt, zq_rema )
-      CALL wrk_dealloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s )
+      CALL wrk_dealloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zqh_i, zqh_s, zq_s, zsnw )
       CALL wrk_dealloc( jpij, nlay_i, zdeltah, zh_i )
       CALL wrk_dealloc( jpij, nlay_i, icount )
@@ -695 +696 @@
       !
    END SUBROUTINE lim_thd_dh
+
+
+   !!--------------------------------------------------------------------------
+   !! INTERFACE lim_thd_snwblow
+   !! ** Purpose :   Compute distribution of precip over the ice
+   !!--------------------------------------------------------------------------
+   SUBROUTINE lim_thd_snwblow_2d( pin, pout )
+      REAL(wp), DIMENSION(:,:), INTENT(in)  :: pin   ! previous fraction lead ( pfrld or (1. - a_i_b) )
+      REAL(wp), DIMENSION(:,:), INTENT(out) :: pout
+      pout = ( 1._wp - ( pin )**rn_betas )
+   END SUBROUTINE lim_thd_snwblow_2d
+
+   SUBROUTINE lim_thd_snwblow_1d( pin, pout )
+      REAL(wp), DIMENSION(:), INTENT(in)  :: pin
+      REAL(wp), DIMENSION(:), INTENT(out) :: pout
+      pout = ( 1._wp - ( pin )**rn_betas )
+   END SUBROUTINE lim_thd_snwblow_1d
+
    
 #else