Changeset 4688 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.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/limcons.F90 (modified) (4 diffs)

trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90

@@ -7 +7 @@
    !!            3.0  ! 2004-06  (M. Vancoppenolle)   Energy Conservation 
    !!            4.0  ! 2011-02  (G. Madec)  add mpp considerations
+   !!             -   ! 2014-05  (C. Rousset) add lim_cons_hsm
    !!----------------------------------------------------------------------
 #if defined key_lim3
@@ -14 +15 @@
    !!    lim_cons     :   checks whether energy, mass and salt are conserved 
    !!----------------------------------------------------------------------
+   USE phycst         ! physical constants
    USE par_ice        ! LIM-3 parameter
    USE ice            ! LIM-3 variables
@@ -28 +30 @@
    PUBLIC   lim_column_sum_energy
    PUBLIC   lim_cons_check
+   PUBLIC   lim_cons_hsm
 
    !!----------------------------------------------------------------------
@@ -151 +154 @@
    END SUBROUTINE lim_cons_check
 
+
+   SUBROUTINE lim_cons_hsm( icount, cd_routine, zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b )
+      !!-------------------------------------------------------------------
+      !!               ***  ROUTINE lim_cons_hsm ***
+      !!
+      !! ** Purpose : Test the conservation of heat, salt and mass for each routine
+      !!
+      !! ** Method  :
+      !!---------------------------------------------------------------------
+      INTEGER         , INTENT(in)    :: icount      ! determine wether this is the beggining of the routine (0) or the end (1)
+      CHARACTER(len=*), INTENT(in)    :: cd_routine  ! name of the routine
+      REAL(wp)        , INTENT(inout) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 
+      REAL(wp)                        :: zvi,   zsmv,   zei,   zfs,   zfw,   zft
+      REAL(wp)                        :: zvmin, zamin, zamax 
+
+      IF( icount == 0 ) THEN
+
+         zvi_b  = glob_sum( SUM(   v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) )
+         zsmv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) )
+         zei_b  = glob_sum( SUM(   e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) )
+         zfw_b  = glob_sum( - ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) +  &
+            &                   wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) + wfx_sub(:,:) + wfx_spr(:,:) ) * area(:,:) * tms(:,:) )
+         zfs_b  = glob_sum(   ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) +  &
+            &                   sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) )
+         zft_b  = glob_sum(   ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:)  & 
+            &                 - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:) ) * area(:,:) / unit_fac * tms(:,:) )
+
+      ELSEIF( icount == 1 ) THEN
+
+         zfs  = glob_sum(   ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) +  &
+            &                 sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zfs_b
+         zfw  = glob_sum( - ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) +  &
+            &                 wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) + wfx_sub(:,:) + wfx_spr(:,:) ) * area(:,:) * tms(:,:) ) - zfw_b
+         zft  = glob_sum(   ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:)  & 
+            &               - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zft_b
+ 
+         zvi  = ( glob_sum( SUM(   v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zvi_b ) * r1_rdtice - zfw 
+         zsmv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zsmv_b ) * r1_rdtice + ( zfs / rhoic )
+         zei  =   glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zei_b * r1_rdtice + zft
+
+         zvmin = glob_min(v_i)
+         zamax = glob_max(SUM(a_i,dim=3))
+         zamin = glob_min(a_i)
+       
+         IF(lwp) THEN
+            IF ( ABS( zvi    ) >  1.e-4 ) WRITE(numout,*) 'violation volume [kg/day]     (',cd_routine,') = ',(zvi * rday)
+            IF ( ABS( zsmv   ) >  1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (',cd_routine,') = ',(zsmv * rday)
+            IF ( ABS( zei    ) >  1.    ) WRITE(numout,*) 'violation enthalpy [1e9 J]    (',cd_routine,') = ',(zei)
+            IF ( zvmin <  0.            ) WRITE(numout,*) 'violation v_i<0  [m]          (',cd_routine,') = ',(zvmin)
+            IF( cd_routine /= 'limtrp' .AND. cd_routine /= 'limitd_me' .AND. zamax > amax+1.e-10 ) THEN
+                                          WRITE(numout,*) 'violation a_i>amax            (',cd_routine,') = ',zamax
+            ENDIF
+            IF ( zamin <  0.            ) WRITE(numout,*) 'violation a_i<0               (',cd_routine,') = ',zamin
+         ENDIF
+
+      ENDIF
+
+   END SUBROUTINE lim_cons_hsm
+
 #else
    !!----------------------------------------------------------------------