Changeset 3962 for branches/2013/dev_r3406_CNRS_LIM3

Timestamp:

2013-07-09T14:22:30+02:00 (11 years ago)

Author:

gm

Message:

dev_r3406_CNRS_LIM3: fix a bug in the freshwater budget correction, see ticket #1116

Location:

branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO

Files:

• LIM_SRC_3/limsbc.F90 (modified) (5 diffs)
• OPA_SRC/SBC/sbcfwb.F90 (modified) (6 diffs)
• OPA_SRC/oce.F90 (modified) (1 diff)

branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

@@ -9 +9 @@
    !!            3.3  ! 2010-05 (G. Madec) decrease ocean & ice reference salinities in the Baltic sea
    !!                 !                  + simplification of the ice-ocean stress calculation
-   !!            4.0  ! 2011-02 (G. Madec) dynamical allocation
+   !!            3.4  ! 2011-02 (G. Madec) dynamical allocation
    !!             -   ! 2012    (D. Iovino) salt flux change
    !!             -   ! 2012-05 (C. Rousset) add penetration solar flux
+   !!            3.5  ! 2012-10 (A. Coward, G. Madec) salt fluxes ; ice+snow mass
    !!----------------------------------------------------------------------
 #if defined key_lim3
@@ -37 +38 @@
    USE cpl_oasis3, ONLY : lk_cpl
    USE traqsr           ! clem: add penetration of solar flux into the calculation of heat budget
-   USE lib_fortran      ! to use key_nosignedzero
+   USE oce,        ONLY : sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass, sshu_b, sshv_b, sshu_n, sshv_n, sshf_n
+   USE dom_ice,    ONLY : tms
+   USE lib_fortran      ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)  
 
    IMPLICIT NONE
@@ -173 +176 @@
             ! qdtcn Energy from the turbulent oceanic heat flux heat flux coming in the lead
 
-            IF ( num_sal == 2 ) zfcm2(ji,jj) = zfcm2(ji,jj) + &
-               fhbri(ji,jj) ! new contribution due to brine drainage 
+            IF ( num_sal == 2 ) zfcm2(ji,jj) = zfcm2(ji,jj) + fhbri(ji,jj) ! new contribution due to brine drainage 
 
             ! bottom radiative component is sent to the computation of the
@@ -281 +283 @@
       ENDIF
 
+      !-----------------------------------------------!
+      !   mass of snow and ice per unit area          !
+      !-----------------------------------------------!
+      IF( nn_ice_embd /= 0 ) THEN                               ! embedded sea-ice (mass required)
+         snwice_mass_b(:,:) = snwice_mass(:,:)                  ! save mass from the previous ice time step
+         !                                                      ! new mass per unit area
+         snwice_mass  (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  ) 
+         !                                                      ! time evolution of snow+ice mass
+         snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_rdtice
+      ENDIF
 
       !-----------------------------------------------!
@@ -420 +432 @@
       oatte(:,:) = 1._wp
       !
+      !                                      ! sea ice  with mass exchange
+         snwice_mass  (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  )
+         snwice_mass_b(:,:) = snwice_mass(:,:)
+      !
+      !
    END SUBROUTINE lim_sbc_init
 

branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90

@@ -59 +59 @@
       !!                =3 global mean of emp set to zero at each nn_fsbc time step
       !!                   & spread out over erp area depending its sign
+      !! Note: if mass exchanges between ice and ocean, it is taken into account 
+      !!       when computing the budget 
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
@@ -87 +89 @@
          !
          area = glob_sum( e1e2t(:,:) )           ! interior global domain surface
+         !
+!!#if ! defined key_lim2 &&  ! defined key_lim3 && ! defined key_cice 
+#if ! defined key_lim3 
+         snwice_mass_b(:,:) = 0._wp              ! no sea-ice model is being used : no snow+ice mass
+         snwice_mass  (:,:) = 0._wp
+         snwice_fmass (:,:) = 0._wp 
+#endif
+         !
       ENDIF
       
@@ -95 +105 @@
          !
          IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN
-            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) ) ) / area   ! sum over the global domain
+            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) - snwice_fmass(:,:) ) ) / area   ! sum over the global domain
             emp (:,:) = emp (:,:) - z_fwf 
             erp (:,:) = erp (:,:) - z_fwf 
@@ -120 +130 @@
          ikty = 365 * 86400 / rdttra(1)    !!bug  use of 365 days leap year or 360d year !!!!!!!
          IF( MOD( kt, ikty ) == 0 ) THEN
-            a_fwb_b = a_fwb
-            a_fwb   = glob_sum( e1e2t(:,:) * sshn(:,:) )   ! sum over the global domain
+            a_fwb_b = a_fwb                           ! mean sea level taking into account the ice+snow
+                                                      ! sum over the global domain
+            a_fwb   = glob_sum( e1e2t(:,:) * ( sshn(:,:) + snwice_mass(:,:) * r1_rau0 ) )
             a_fwb   = a_fwb * 1.e+3 / ( area * 86400. * 365. )     ! convert in Kg/m3/s = mm/s
 !!gm        !                                                      !!bug 365d year 
@@ -151 +162 @@
             zsurf_neg = glob_sum( e1e2t(:,:)*ztmsk_neg(:,:) )   ! Area filled by <0 and >0 erp 
             zsurf_pos = glob_sum( e1e2t(:,:)*ztmsk_pos(:,:) )
-            !                                                  ! fwf global mean 
-            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) ) ) / area
+            !                                                  ! fwf global mean (excluding ocean to ice/snow exchanges) 
+            z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) - snwice_fmass(:,:) ) ) / area
             !            
             IF( z_fwf < 0._wp ) THEN         ! spread out over >0 erp area to increase evaporation
@@ -166 +177 @@
             z_fwf_nsrf =  zsum_fwf / ( zsurf_tospread + rsmall )
             !                                                  ! weight to respect erp field 2D structure 
-            zsum_erp = glob_sum( ztmsk_tospread(:,:) * erp(:,:) * e1e2t(:,:) )
+            zsum_erp   = glob_sum( ztmsk_tospread(:,:) * erp(:,:) * e1e2t(:,:) )
             z_wgt(:,:) = ztmsk_tospread(:,:) * erp(:,:) / ( zsum_erp + rsmall )
             !                                                  ! final correction term to apply

branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/oce.F90

@@ -47 +47 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   gru , grv    !: horizontal gradient of rd at bottom u-point
 
+   !! arrays relating to embedding ice in the ocean. These arrays need to be declared 
+   !! even if no ice model is required. In the no ice model or traditional levitating 
+   !! ice cases they contain only zeros
+   !! ---------------------
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   snwice_mass        !: mass of snow and ice at current  ice time step   [Kg/m2]
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   snwice_mass_b      !: mass of snow and ice at previous ice time step   [Kg/m2]
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   snwice_fmass       !: time evolution of mass of snow+ice               [Kg/m2/s]
+
    !!----------------------------------------------------------------------
    !! NEMO/OPA 4.0 , NEMO Consortium (2011)