Changeset 12064 for NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk

Timestamp:

2019-12-05T11:59:14+01:00 (4 years ago)

Author:

laurent

Message:

Added sanity/consistence check on air humidity.

File:

• NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90 (modified) (9 diffs)

NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90

@@ -107 +107 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   cdn_oce, chn_oce, cen_oce ! needed by Lupkes 2015 bulk scheme
 
-   !LB:
    LOGICAL  ::   ln_skin_cs     ! use the cool-skin (only available in ECMWF and COARE algorithms) !LB
    LOGICAL  ::   ln_skin_wl     ! use the warm-layer parameterization (only available in ECMWF and COARE algorithms) !LB
@@ -119 +118 @@
    INTEGER, PARAMETER :: np_humi_dpt = 2
    INTEGER, PARAMETER :: np_humi_rlh = 3
-   !LB.
 
    INTEGER  ::   nblk           ! choice of the bulk algorithm
@@ -216 +214 @@
             & CALL ctl_stop( 'sbc_blk_init: Please set "nn_fsbc" to 1 when using cool-skin/warm-layer param.')
       END IF
-      
+
       IF( ln_skin_wl ) THEN
          !! Check if the frequency of downwelling solar flux input makes sense and if ln_dm2dc=T if it is daily!
@@ -224 +222 @@
             & CALL ctl_stop( 'sbc_blk_init: Please set ln_dm2dc=T for warm-layer param. (ln_skin_wl) to work properly' )
       END IF
-      
+
       ioptio = 0
       IF( ln_humi_sph ) THEN
@@ -369 +367 @@
       INTEGER, INTENT(in) ::   kt   ! ocean time step
       !!---------------------------------------------------------------------
+      REAL(wp) :: ztmp
+      !!---------------------------------------------------------------------
       !
       CALL fld_read( kt, nn_fsbc, sf )             ! input fields provided at the current time-step
-      !
+
+      ! Sanity/consistence test on humidity at first time step to detect potential screw-up:
+      IF( kt == nit000 ) THEN
+         WRITE(numout,*) ''
+         ztmp = SUM(tmask(:,:,1)) ! number of ocean points on local proc domain
+         IF( ztmp > 8._wp ) THEN ! test only on proc domains with at least 8 ocean points!
+            ztmp = SUM(sf(jp_humi)%fnow(:,:,1)*tmask(:,:,1))/ztmp ! mean humidity over ocean on proc
+            SELECT CASE( nhumi )
+            CASE( np_humi_sph ) ! specific humidity => expect: 0. <= something < 0.065 [kg/kg] (0.061 is saturation at 45degC !!!)
+               IF(  (ztmp < 0._wp) .OR. (ztmp > 0.065)  ) ztmp = -1._wp
+            CASE( np_humi_dpt ) ! dew-point temperature => expect: 110. <= something < 320. [K]
+               IF( (ztmp < 110._wp).OR.(ztmp > 320._wp) ) ztmp = -1._wp
+            CASE( np_humi_rlh ) ! relative humidity => expect: 0. <= something < 100. [%]
+               IF(  (ztmp < 0._wp) .OR.(ztmp > 100._wp) ) ztmp = -1._wp
+            END SELECT
+            IF(ztmp < 0._wp) THEN
+               WRITE(numout,'("   Mean humidity value found on proc #",i5.5," is: ",f)'), narea, ztmp
+               CALL ctl_stop( 'STOP', 'Something is wrong with air humidity!!!', &
+                  &   ' ==> check the unit in your input files'       , &
+                  &   ' ==> check consistence of namelist choice: specific? relative? dew-point?', &
+                  &   ' ==> ln_humi_sph -> [kg/kg] | ln_humi_rlh -> [%] | ln_humi_dpt -> [K] !!!' )
+            END IF
+         END IF
+         WRITE(numout,*) ' === Sanity/consistence test on air humidity sucessfuly passed! ==='
+         WRITE(numout,*) ''
+      END IF !IF( kt == nit000 )
+
       !                                            ! compute the surface ocean fluxes using bulk formulea
       IF( MOD( kt - 1, nn_fsbc ) == 0 )   CALL blk_oce( kt, sf, sst_m, ssu_m, ssv_m )
@@ -522 +548 @@
       END DO
 
-      
+
       !! Time to call the user-selected bulk parameterization for
       !!  ==  transfer coefficients  ==!   Cd, Ch, Ce at T-point, and more...
@@ -592 +618 @@
       zwnd_i(:,:) = Cd_atm(:,:) * zwnd_i(:,:)
       zwnd_j(:,:) = Cd_atm(:,:) * zwnd_j(:,:)
-      
+
       !                          ! add the HF tau contribution to the wind stress module
       IF( lhftau )   taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1)
@@ -664 +690 @@
       CALL iom_put( 'precip' , tprecip )                 ! Total precipitation
       !
-      IF ( nn_ice == 0 ) THEN
+      IF( nn_ice == 0 ) THEN
          CALL iom_put( "qemp_oce" ,   qns-zqlw-zqsb-zqla )   ! output downward heat content of E-P over the ocean
          CALL iom_put( "qns_oce"  ,   qns  )                 ! output downward non solar heat over the ocean
@@ -735 +761 @@
       !
       ! Make ice-atm. drag dependent on ice concentration
-      IF    ( ln_Cd_L12 ) THEN   ! calculate new drag from Lupkes(2012) equations
+      IF( ln_Cd_L12 ) THEN   ! calculate new drag from Lupkes(2012) equations
          CALL Cdn10_Lupkes2012( Cd_atm )
          Ch_atm(:,:) = Cd_atm(:,:)       ! momentum and heat transfer coef. are considered identical