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 11772 for NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk_phy.F90 – NEMO

Ignore:
Timestamp:
2019-10-23T16:04:12+02:00 (5 years ago)
Author:
laurent
Message:

LB: solid updates+improvements of cool-skin/warm-layer capabilty of COARE and ECMWF bulk algorithms!

File:
1 edited

Legend:

Unmodified
Added
Removed

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

    r11638 r11772  
    4545   END INTERFACE cp_air 
    4646 
    47       INTERFACE alpha_sw 
     47   INTERFACE alpha_sw 
    4848      MODULE PROCEDURE alpha_sw_vctr, alpha_sw_sclr 
    4949   END INTERFACE alpha_sw 
     50 
     51   INTERFACE turb_fluxes 
     52      MODULE PROCEDURE turb_fluxes_vctr, turb_fluxes_sclr 
     53   END INTERFACE turb_fluxes 
    5054 
    5155 
    … …  
    6468   PUBLIC update_qnsol_tau 
    6569   PUBLIC alpha_sw 
     70   PUBLIC turb_fluxes 
    6671 
    6772   !!---------------------------------------------------------------------- 
    … …  
    128133      rho_air_sclr = MAX( pslp / (R_dry*ptak * ( 1._wp + rctv0*pqa )) , 0.8_wp ) 
    129134   END FUNCTION rho_air_sclr 
    130     
     135 
    131136 
    132137 
    … …  
    464469 
    465470 
    466    SUBROUTINE UPDATE_QNSOL_TAU( pTs, pqs, pTa, pqa, pust, ptst, pqst, pUb, pslp, prlw, & 
     471   SUBROUTINE UPDATE_QNSOL_TAU( pzu, pTs, pqs, pTa, pqa, pust, ptst, pqst, pwnd, pUb, pslp, prlw, & 
    467472      &                         pQns, pTau,  & 
    468473      &                         Qlat) 
    … …  
    472477      !! ** Author: L. Brodeau, Sept. 2019 / AeroBulk (https://github.com/brodeau/aerobulk/) 
    473478      !!---------------------------------------------------------------------------------- 
     479      REAL(wp),                     INTENT(in)  :: pzu  ! height above the sea-level where all this takes place (normally 10m) 
    474480      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pTs  ! water temperature at the air-sea interface [K] 
    475481      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqs  ! satur. spec. hum. at T=pTs   [kg/kg] 
    476       REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pTa  ! air temperature at z=zu [K] 
    477       REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqa  ! specific humidity at z=zu [kg/kg] 
     482      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pTa  ! potential air temperature at z=pzu [K] 
     483      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqa  ! specific humidity at z=pzu [kg/kg] 
    478484      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pust ! u* 
    479485      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: ptst ! t* 
    480486      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqst ! q* 
    481       REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pUb  ! bulk wind speed at z=zu [m/s] 
     487      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pwnd ! wind speed module at z=pzu [m/s] 
     488      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pUb  ! bulk wind speed at z=pzu (inc. pot. effect of gustiness etc) [m/s] 
    482489      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pslp ! sea-level atmospheric pressure [Pa] 
    483490      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: prlw ! downwelling longwave radiative flux [W/m^2] 
    … …  
    491498         &        zQlat, zQsen, zQlw 
    492499      INTEGER  ::   ji, jj     ! dummy loop indices 
    493       !!----------------------------------------------------------------------------------      
     500      !!---------------------------------------------------------------------------------- 
    494501      DO jj = 1, jpj 
    495502         DO ji = 1, jpi 
    … …  
    502509            zCe = zz0*pqst(ji,jj)/zdq 
    503510 
    504             zUrho = pUb(ji,jj)*MAX(rho_air(pTa(ji,jj), pqa(ji,jj), pslp(ji,jj)), 1._wp)     ! rho*U10 
     511            !zUrho = pUb(ji,jj)*MAX(rho_air(pTa(ji,jj), pqa(ji,jj), pslp(ji,jj)), 1._wp)     ! rho*U10 
    505512            zTs2  = pTs(ji,jj)*pTs(ji,jj) 
    506513 
     514            CALL TURB_FLUXES( pzu, pTs(ji,jj), pqs(ji,jj), pTa(ji,jj), pqa(ji,jj), zCd, zCh, zCe, & 
     515               &              pwnd(ji,jj), pUb(ji,jj), pslp(ji,jj), & 
     516               &              pTau(ji,jj), zQsen, zQlat ) 
     517 
     518 
    507519            ! Wind stress module: 
    508             pTau(ji,jj) = zCd*zUrho*pUb(ji,jj) ! lolo? 
     520            !pTau(ji,jj) = zCd*zUrho*pUb(ji,jj) ! lolo? 
    509521 
    510522            ! Non-Solar heat flux to the ocean: 
    511             zQlat = MIN ( zUrho*zCe*L_vap( pTs(ji,jj)) * zdq , 0._wp )  ! we do not want a Qlat > 0 ! 
    512             zQsen = zUrho*zCh*cp_air(pqa(ji,jj)) * zdt 
     523            !zQlat = MIN ( zUrho*zCe*L_vap( pTs(ji,jj)) * zdq , 0._wp )  ! we do not want a Qlat > 0 ! 
     524            !zQsen = zUrho*zCh*cp_air(pqa(ji,jj)) * zdt 
    513525            zQlw  = emiss_w*(prlw(ji,jj) - stefan*zTs2*zTs2) ! Net longwave flux 
    514526 
    515527            pQns(ji,jj) = zQlat + zQsen + zQlw 
    516              
    517             IF ( PRESENT(Qlat) ) Qlat(ji,jj) = zQlat             
     528 
     529            IF ( PRESENT(Qlat) ) Qlat(ji,jj) = zQlat 
    518530         END DO 
    519531      END DO 
    520532   END SUBROUTINE UPDATE_QNSOL_TAU 
     533 
     534 
     535 
     536 
     537 
     538   SUBROUTINE TURB_FLUXES_VCTR( pzu, pTs, pqs, pTa, pqa, pCd, pCh, pCe, pwnd, pUb, pslp, & 
     539      &                                 pTau, pQsen, pQlat,  pEvap ) 
     540      !!---------------------------------------------------------------------------------- 
     541      REAL(wp),                     INTENT(in)  :: pzu  ! height above the sea-level where all this takes place (normally 10m) 
     542      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pTs  ! water temperature at the air-sea interface [K] 
     543      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqs  ! satur. spec. hum. at T=pTs   [kg/kg] 
     544      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pTa  ! potential air temperature at z=pzu [K] 
     545      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pqa  ! specific humidity at z=pzu [kg/kg] 
     546      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pCd 
     547      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pCh 
     548      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pCe 
     549      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pwnd ! wind speed module at z=pzu [m/s] 
     550      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pUb  ! bulk wind speed at z=pzu (inc. pot. effect of gustiness etc) [m/s] 
     551      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: pslp ! sea-level atmospheric pressure [Pa] 
     552      !! 
     553      REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pTau  ! module of the wind stress [N/m^2] 
     554      REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pQsen !  [W/m^2] 
     555      REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pQlat !  [W/m^2] 
     556      !! 
     557      REAL(wp), DIMENSION(jpi,jpj), INTENT(out), OPTIONAL :: pEvap !  [kg/m^2/s] 
     558      !! 
     559      REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap 
     560      INTEGER  :: ji, jj, jq     ! dummy loop indices 
     561      !!---------------------------------------------------------------------------------- 
     562      DO jj = 1, jpj 
     563         DO ji = 1, jpi 
     564 
     565            !! Need ztaa, absolute temperature at pzu (formula to estimate rho_air needs absolute temperature, not the potential temperature "pTa") 
     566            ztaa = pTa(ji,jj) ! first guess... 
     567            DO jq = 1, 4 
     568               zgamma = gamma_moist( 0.5*(ztaa+pTs(ji,jj)) , pqa(ji,jj) ) 
     569               ztaa = pTa(ji,jj) - zgamma*pzu   ! Absolute temp. is slightly colder... 
     570            END DO 
     571            zrho = rho_air(ztaa, pqa(ji,jj), pslp(ji,jj)) 
     572            zrho = rho_air(ztaa, pqa(ji,jj), pslp(ji,jj)-zrho*grav*pzu) ! taking into account that we are pzu m above the sea level where SLP is given! 
     573 
     574            zUrho = pUb(ji,jj)*MAX(zrho, 1._wp)     ! rho*U10 
     575 
     576            pTau(ji,jj) = zUrho * pCd(ji,jj) * pwnd(ji,jj) ! Wind stress module 
     577 
     578            zevap        = MIN( zUrho * pCe(ji,jj) * (pqa(ji,jj) - pqs(ji,jj)) , 0._wp )   ! we do not want condensation & Qlat > 0 ! 
     579            pQsen(ji,jj) =      zUrho * pCh(ji,jj) * (pTa(ji,jj) - pTs(ji,jj)) * cp_air(pqa(ji,jj)) 
     580            pQlat(ji,jj) =  L_vap(pTs(ji,jj)) * zevap 
     581 
     582            IF ( PRESENT(pEvap) ) pEvap(ji,jj) = - zevap 
     583 
     584         END DO 
     585      END DO 
     586   END SUBROUTINE TURB_FLUXES_VCTR 
     587 
     588 
     589   SUBROUTINE TURB_FLUXES_SCLR( pzu, pTs, pqs, pTa, pqa, pCd, pCh, pCe, pwnd, pUb, pslp, & 
     590      &                                 pTau, pQsen, pQlat,  pEvap ) 
     591      !!---------------------------------------------------------------------------------- 
     592      REAL(wp),                     INTENT(in)  :: pzu  ! height above the sea-level where all this takes place (normally 10m) 
     593      REAL(wp), INTENT(in)  :: pTs  ! water temperature at the air-sea interface [K] 
     594      REAL(wp), INTENT(in)  :: pqs  ! satur. spec. hum. at T=pTs   [kg/kg] 
     595      REAL(wp), INTENT(in)  :: pTa  ! potential air temperature at z=pzu [K] 
     596      REAL(wp), INTENT(in)  :: pqa  ! specific humidity at z=pzu [kg/kg] 
     597      REAL(wp), INTENT(in)  :: pCd 
     598      REAL(wp), INTENT(in)  :: pCh 
     599      REAL(wp), INTENT(in)  :: pCe 
     600      REAL(wp), INTENT(in)  :: pwnd ! wind speed module at z=pzu [m/s] 
     601      REAL(wp), INTENT(in)  :: pUb  ! bulk wind speed at z=pzu (inc. pot. effect of gustiness etc) [m/s] 
     602      REAL(wp), INTENT(in)  :: pslp ! sea-level atmospheric pressure [Pa] 
     603      !! 
     604      REAL(wp), INTENT(out) :: pTau  ! module of the wind stress [N/m^2] 
     605      REAL(wp), INTENT(out) :: pQsen !  [W/m^2] 
     606      REAL(wp), INTENT(out) :: pQlat !  [W/m^2] 
     607      !! 
     608      REAL(wp), INTENT(out), OPTIONAL :: pEvap !  [kg/m^2/s] 
     609      !! 
     610      REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap 
     611      INTEGER  :: jq 
     612      !!---------------------------------------------------------------------------------- 
     613 
     614      !! Need ztaa, absolute temperature at pzu (formula to estimate rho_air needs absolute temperature, not the potential temperature "pTa") 
     615      ztaa = pTa ! first guess... 
     616      DO jq = 1, 4 
     617         zgamma = gamma_moist( 0.5*(ztaa+pTs) , pqa ) 
     618         ztaa = pTa - zgamma*pzu   ! Absolute temp. is slightly colder... 
     619      END DO 
     620      zrho = rho_air(ztaa, pqa, pslp) 
     621      zrho = rho_air(ztaa, pqa, pslp-zrho*grav*pzu) ! taking into account that we are pzu m above the sea level where SLP is given! 
     622 
     623      zUrho = pUb*MAX(zrho, 1._wp)     ! rho*U10 
     624 
     625      pTau = zUrho * pCd * pwnd ! Wind stress module 
     626 
     627      zevap        = MIN( zUrho * pCe * (pqa - pqs) , 0._wp )   ! we do not want condensation & Qlat > 0 ! 
     628      pQsen =      zUrho * pCh * (pTa - pTs) * cp_air(pqa) 
     629      pQlat =  L_vap(pTs) * zevap 
     630 
     631      IF ( PRESENT(pEvap) ) pEvap = - zevap 
     632 
     633   END SUBROUTINE TURB_FLUXES_SCLR 
     634 
     635 
     636 
    521637 
    522638 
    … …  
    529645      !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://github.com/brodeau/aerobulk/) 
    530646      !!---------------------------------------------------------------------------------- 
    531       REAL(wp), DIMENSION(jpi,jpj)             ::   alpha_sw_vctr   ! latent heat of vaporization   [J/kg] 
     647      REAL(wp), DIMENSION(jpi,jpj)             ::   alpha_sw_vctr   ! thermal expansion coefficient of sea-water [1/K] 
    532648      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   psst   ! water temperature                [K] 
    533649      !!---------------------------------------------------------------------------------- 
    … …  
    543659      !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://github.com/brodeau/aerobulk/) 
    544660      !!---------------------------------------------------------------------------------- 
    545       REAL(wp)             ::   alpha_sw_sclr   ! latent heat of vaporization   [J/kg] 
     661      REAL(wp)             ::   alpha_sw_sclr   ! thermal expansion coefficient of sea-water [1/K] 
    546662      REAL(wp), INTENT(in) ::   psst   ! sea-water temperature                   [K] 
    547663      !!---------------------------------------------------------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.