Changeset 4688 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90

Timestamp:

2014-06-25T01:39:59+02:00 (10 years ago)

Author:

clem

Message:

new version of LIM3 with perfect conservation of heat, see ticket #1352

File:

• trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90 (modified) (6 diffs)

trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90

@@ -53 +53 @@
       !
       INTEGER  ::   ji, jk     ! dummy loop indices 
-      REAL(wp) ::   zsold, iflush, iaccrbo, igravdr, isnowic, i_ice_switch,  ztmelts   ! local scalars
-      REAL(wp) ::   zaaa, zbbb, zccc, zdiscrim   ! local scalars
-      REAL(wp), POINTER, DIMENSION(:) ::   ze_init, zhiold, zsiold
+      REAL(wp) ::   iflush, igravdr   ! local scalars
       !!---------------------------------------------------------------------
 
-      CALL wrk_alloc( jpij, ze_init, zhiold, zsiold )
-
+      !---------------------------------------------------------
+      !  0) Update ice salinity from snow-ice and bottom growth
+      !---------------------------------------------------------
+      DO ji = kideb, kiut
+         sm_i_b(ji) = sm_i_b(ji) + dsm_i_se_1d(ji) + dsm_i_si_1d(ji)
+      END DO
+ 
       !------------------------------------------------------------------------------|
       ! 1) Constant salinity, constant in time                                       |
@@ -74 +77 @@
       !  Module 2 : Constant salinity varying in time                                |
       !------------------------------------------------------------------------------|
-
       IF(  num_sal == 2  ) THEN
-
-         !---------------------------------
-         ! Thickness at previous time step
-         !---------------------------------
-         DO ji = kideb, kiut
-            zhiold(ji) = ht_i_b(ji) - dh_i_bott(ji) - dh_snowice(ji) - dh_i_surf(ji)
-            zsiold(ji) = sm_i_b(ji)
-         END DO
-
-         !---------------------
-         ! Global heat content
-         !---------------------
-         ze_init(:)  =  0._wp
-         DO jk = 1, nlay_i
-            DO ji = kideb, kiut
-               ze_init(ji) = ze_init(ji) + q_i_b(ji,jk) * ht_i_b(ji) / REAL (nlay_i )
-            END DO
-         END DO
 
          DO ji = kideb, kiut
@@ -99 +83 @@
             ! Switches 
             !----------
-            iflush       =         MAX( 0._wp , SIGN( 1.0 , t_su_b(ji) - rtt )        )    ! =1 if summer 
-            igravdr      =         MAX( 0._wp , SIGN( 1.0 , t_bo_b(ji) - t_su_b(ji) ) )    ! =1 if t_su < t_bo
-            iaccrbo      =         MAX( 0._wp , SIGN( 1.0 , dh_i_bott(ji) )           )    ! =1 if bottom accretion
-            i_ice_switch = 1._wp - MAX ( 0._wp , SIGN( 1._wp , - ht_i_b(ji) + 1.e-2 ) )
-            isnowic      = 1._wp - MAX ( 0._wp , SIGN( 1._wp , - dh_snowice(ji) ) ) * i_ice_switch   ! =1 if snow ice formation
+            iflush  = MAX( 0._wp , SIGN( 1._wp , t_su_b(ji) - rtt )        )    ! =1 if summer 
+            igravdr = MAX( 0._wp , SIGN( 1._wp , t_bo_b(ji) - t_su_b(ji) ) )    ! =1 if t_su < t_bo
 
             !---------------------
             ! Salinity tendencies
             !---------------------
-            !                                   ! drainage by gravity drainage
+            ! drainage by gravity drainage
             dsm_i_gd_1d(ji) = - igravdr * MAX( sm_i_b(ji) - sal_G , 0._wp ) / time_G * rdt_ice 
-            !                                   ! drainage by flushing  
+            ! drainage by flushing  
             dsm_i_fl_1d(ji) = - iflush  * MAX( sm_i_b(ji) - sal_F , 0._wp ) / time_F * rdt_ice
 
@@ -120 +101 @@
             sm_i_b(ji) = sm_i_b(ji) + dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji)
 
-            ! if no ice, salinity = 0.1
-            i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) )
-            sm_i_b(ji)   = i_ice_switch * sm_i_b(ji) + s_i_min * ( 1._wp - i_ice_switch )
-
-            !----------------------------
-            ! Heat flux - brine drainage
-            !----------------------------
-            fhbri_1d(ji) = 0._wp
-
             !----------------------------
             ! Salt flux - brine drainage
             !----------------------------
-            sfx_bri_1d(ji) = sfx_bri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) * ( sm_i_b(ji) - zsiold(ji) ) * r1_rdtice
+            sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoic * a_i_b(ji) * ht_i_b(ji) * ( dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji) ) * r1_rdtice
 
          END DO
@@ -138 +110 @@
          ! Salinity profile
          CALL lim_var_salprof1d( kideb, kiut )
-
-
-         ! Only necessary for conservation check since salinity is modified
-         !--------------------
-         ! Temperature update
-         !--------------------
-         DO jk = 1, nlay_i
-            DO ji = kideb, kiut
-               ztmelts    =  -tmut*s_i_b(ji,jk) + rtt
-               !Conversion q(S,T) -> T (second order equation)
-               zaaa         =  cpic
-               zbbb         =  ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk) / rhoic - lfus
-               zccc         =  lfus * ( ztmelts - rtt )
-               zdiscrim     =  SQRT(  MAX( zbbb*zbbb - 4.0*zaaa*zccc, 0._wp )  )
-               t_i_b(ji,jk) =  rtt - ( zbbb + zdiscrim ) / ( 2.0 *zaaa )
-            END DO
-         END DO
          !
       ENDIF 
@@ -161 +116 @@
       !  Module 3 : Profile of salinity, constant in time                            |
       !------------------------------------------------------------------------------|
-
       IF(  num_sal == 3  )   CALL lim_var_salprof1d( kideb, kiut )
 
-      !
-      CALL wrk_dealloc( jpij, ze_init, zhiold, zsiold )
       !
    END SUBROUTINE lim_thd_sal