New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 11804 for NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90 – NEMO

Ignore:
Timestamp:
2019-10-25T17:25:19+02:00 (4 years ago)
Author:
laurent
Message:

Use of "TURB_FLUXES@…" inside sbcblk.F90, positive latent HF & positive cool-skin temperature now allowed!

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90

    r11772 r11804  
    520520      !! Estimate of potential temperature at z=rn_zqt, based on adiabatic lapse-rate 
    521521      !!    (see Josey, Gulev & Yu, 2013) / doi=10.1016/B978-0-12-391851-2.00005-2 
    522       !!    (since reanalysis products provide T at z, not theta !) 
     522      !!    (since reanalysis products provide absolute temperature "T" at z, not theta !) 
    523523      ztpot = sf(jp_tair)%fnow(:,:,1) + gamma_moist( sf(jp_tair)%fnow(:,:,1), zqair(:,:) ) * rn_zqt 
    524524 
     
    594594      END IF ! IF ( ln_skin_cs .OR. ln_skin_wl ) 
    595595 
    596  
    597  
    598  
    599       !                          ! Compute true air density : 
    600       IF( ABS(rn_zu - rn_zqt) > 0.01 ) THEN     ! At zu: (probably useless to remove rhoa*grav*rn_zu from SLP...) 
    601          rhoa(:,:) = rho_air( t_zu(:,:)              , q_zu(:,:)              , sf(jp_slp)%fnow(:,:,1) ) 
    602       ELSE                                      ! At zt: 
    603          rhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), zqair(:,:), sf(jp_slp)%fnow(:,:,1) ) 
    604       END IF 
    605  
    606596      !!      CALL iom_put( "Cd_oce", Cd_atm)  ! output value of pure ocean-atm. transfer coef. 
    607597      !!      CALL iom_put( "Ch_oce", Ch_atm)  ! output value of pure ocean-atm. transfer coef. 
    608  
    609       DO jj = 1, jpj             ! tau module, i and j component 
     598       
     599      IF( ABS(rn_zu - rn_zqt) < 0.1_wp ) THEN 
     600         !! If zu == zt, then ensuring once for all that: 
     601         t_zu(:,:) = ztpot(:,:) 
     602         q_zu(:,:) = zqair(:,:) 
     603      END IF 
     604       
     605 
     606      !  Turbulent fluxes over ocean  => TURB_FLUXES @ sbcblk_phy.F90 
     607      ! ------------------------------------------------------------- 
     608       
     609      CALL TURB_FLUXES( rn_zu, zst(:,:), zsq(:,:), t_zu(:,:), q_zu(:,:), Cd_atm(:,:), Ch_atm(:,:), Ce_atm(:,:), & 
     610         &                 wndm(:,:), zU_zu(:,:), sf(jp_slp)%fnow(:,:,1), & 
     611         &                 taum(:,:), zqsb(:,:), zqla(:,:),               & 
     612         &                 pEvap=zevap(:,:), prhoa=rhoa(:,:) ) 
     613       
     614      zqla(:,:)  =  zqla(:,:) * tmask(:,:,1) 
     615      zqsb(:,:)  =  zqsb(:,:) * tmask(:,:,1) 
     616      taum(:,:)  =  taum(:,:) * tmask(:,:,1) 
     617      zevap(:,:) = zevap(:,:) * tmask(:,:,1) 
     618       
     619      DO jj = 1, jpj             ! tau i and j component on T-grid points 
    610620         DO ji = 1, jpi 
    611             zztmp = rhoa(ji,jj)  * zU_zu(ji,jj) * Cd_atm(ji,jj)   ! using bulk wind speed 
    612             taum  (ji,jj) = zztmp * wndm  (ji,jj) 
     621            zztmp = taum(ji,jj) / wndm(ji,jj) 
    613622            zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) 
    614623            zwnd_j(ji,jj) = zztmp * zwnd_j(ji,jj) 
    615624         END DO 
    616625      END DO 
    617  
    618626      !                          ! add the HF tau contribution to the wind stress module 
    619627      IF( lhftau )   taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1) 
     
    634642      CALL lbc_lnk_multi( 'sbcblk', utau, 'U', -1., vtau, 'V', -1. ) 
    635643 
    636       !  Turbulent fluxes over ocean 
    637       ! ----------------------------- 
    638  
    639       ! zqla used as temporary array, for rho*U (common term of bulk formulae): 
    640       zqla(:,:) = rhoa(:,:) * zU_zu(:,:) * tmask(:,:,1) 
    641  
    642       IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 
    643          !! q_air and t_air are given at 10m (wind reference height) 
    644          zevap(:,:) = rn_efac*MAX( 0._wp,  zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - zqair(:,:)) ) ! Evaporation, using bulk wind speed 
    645          zqsb (:,:) = cp_air(zqair(:,:)) * zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - ztpot(:,:)  ) ! Sensible Heat, using bulk wind speed 
    646       ELSE 
    647          !! q_air and t_air are not given at 10m (wind reference height) 
    648          ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 
    649          zevap(:,:) = rn_efac*MAX( 0._wp,  zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - q_zu(:,:)) )  ! Evaporation, using bulk wind speed 
    650          zqsb (:,:) = cp_air(zqair(:,:)) * zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - t_zu(:,:) )   ! Sensible Heat, using bulk wind speed 
    651       ENDIF 
    652  
    653       zqla(:,:) = L_vap(zst(:,:)) * zevap(:,:)     ! Latent Heat flux 
    654  
    655644 
    656645      ! ----------------------------------------------------------------------------- ! 
     
    658647      ! ----------------------------------------------------------------------------- ! 
    659648 
    660       !! LB: now after Turbulent fluxes because must use the skin temperature rather that the SST ! (zst is skin temperature if ln_skin_cs==.TRUE. .OR. ln_skin_wl==.TRUE.) 
     649      !! LB: now moved after Turbulent fluxes because must use the skin temperature rather that the SST ! (zst is skin temperature if ln_skin_cs==.TRUE. .OR. ln_skin_wl==.TRUE.) 
    661650      zqlw(:,:) = emiss_w * ( sf(jp_qlw)%fnow(:,:,1) - stefan*zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:) ) * tmask(:,:,1)   ! Net radiative longwave flux 
    662651 
     
    681670         &         - sf(jp_prec)%fnow(:,:,1) * rn_pfac  ) * tmask(:,:,1) 
    682671      ! 
    683       qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:)                                &   ! Downward Non Solar 
     672      qns(:,:) = zqlw(:,:) + zqsb(:,:) + zqla(:,:)                                &   ! Downward Non Solar 
    684673         &     - sf(jp_snow)%fnow(:,:,1) * rn_pfac * rLfus                        &   ! remove latent melting heat for solid precip 
    685674         &     - zevap(:,:) * pst(:,:) * rcp                                      &   ! remove evap heat content at SST !LB??? pst is Celsius !? 
     
    691680      ! 
    692681#if defined key_si3 
    693       qns_oce(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:)                                ! non solar without emp (only needed by SI3) 
     682      qns_oce(:,:) = zqlw(:,:) + zqsb(:,:) + zqla(:,:)                                ! non solar without emp (only needed by SI3) 
    694683      qsr_oce(:,:) = qsr(:,:) 
    695684#endif 
     
    698687      CALL iom_put( "rho_air"  ,   rhoa )                 ! output air density (kg/m^3) !#LB 
    699688      CALL iom_put( "qlw_oce"  ,   zqlw )                 ! output downward longwave heat over the ocean 
    700       CALL iom_put( "qsb_oce"  , - zqsb )                 ! output downward sensible heat over the ocean 
    701       CALL iom_put( "qla_oce"  , - zqla )                 ! output downward latent   heat over the ocean 
     689      CALL iom_put( "qsb_oce"  ,   zqsb )                 ! output downward sensible heat over the ocean 
     690      CALL iom_put( "qla_oce"  ,   zqla )                 ! output downward latent   heat over the ocean 
    702691      CALL iom_put( "evap_oce" ,  zevap )                 ! evaporation 
    703       CALL iom_put( "qemp_oce" ,   qns-zqlw+zqsb+zqla )   ! output downward heat content of E-P over the ocean 
     692      CALL iom_put( "qemp_oce" ,   qns-zqlw-zqsb-zqla )   ! output downward heat content of E-P over the ocean 
    704693      CALL iom_put( "qns_oce"  ,   qns  )                 ! output downward non solar heat over the ocean 
    705694      CALL iom_put( "qsr_oce"  ,   qsr  )                 ! output downward solar heat over the ocean 
Note: See TracChangeset for help on using the changeset viewer.