scale_land_sfc_thick-slab.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
module scale_land_sfc_thick_slab
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_stdio
  use scale_grid_index
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: land_sfc_thick_slab_setup
  public :: land_sfc_thick_slab

  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine land_sfc_thick_slab_setup( LAND_TYPE )
    use scale_process, only: &
       prc_mpistop
    implicit none

    character(len=*), intent(in) :: land_type

    logical :: dummy

    namelist / param_land_sfc_thick_slab / &
       dummy

    integer :: ierr
    !---------------------------------------------------------------------------

    if( io_l ) write(io_fid_log,*)
    if( io_l ) write(io_fid_log,*) '++++++ Module[THICK-SLAB] / Categ[LAND SFC] / Origin[SCALElib]'

    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_land_sfc_thick_slab,iostat=ierr)
    if( ierr < 0 ) then !--- missing
       if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
       write(*,*) 'xxx Not appropriate names in namelist PARAM_LAND_SFC_THICK_SLAB. Check!'
       call prc_mpistop
    endif
    if( io_nml ) write(io_fid_nml,nml=param_land_sfc_thick_slab)

    return
  end subroutine land_sfc_thick_slab_setup

  !-----------------------------------------------------------------------------
  subroutine land_sfc_thick_slab( &
        LST_t,      &
        ZMFLX,      &
        XMFLX,      &
        YMFLX,      &
        SHFLX,      &
        LHFLX,      &
        GHFLX,      &
        U10,        &
        V10,        &
        T2,         &
        Q2,         &
        TMPA,       &
        PRSA,       &
        WA,         &
        UA,         &
        VA,         &
        RHOA,       &
        QVA,        &
        Z1,         &
        PBL,        &
        PRSS,       &
        LWD,        &
        SWD,        &
        TG,         &
        LST,        &
        QVEF,       &
        ALB_LW,     &
        ALB_SW,     &
        DZG,        &
        TCS,        &
        Z0M,        &
        Z0H,        &
        Z0E,        &
        dt          )
    use scale_process, only: &
      prc_mpistop
    use scale_const, only: &
      pre00 => const_pre00, &
      rdry  => const_rdry,  &
      cpdry => const_cpdry, &
      stb   => const_stb
    use scale_landuse, only: &
      landuse_fact_land
    use scale_atmos_hydrometeor, only: &
      hydrometeor_lhv => atmos_hydrometeor_lhv
    use scale_atmos_saturation, only: &
      qsat => atmos_saturation_pres2qsat_all
    use scale_bulkflux, only: &
      bulkflux
    implicit none

    ! arguments
    real(RP), intent(out) :: lst_t(ia,ja) ! tendency of LST
    real(RP), intent(out) :: zmflx(ia,ja) ! z-momentum flux at the surface [kg/m/s2]
    real(RP), intent(out) :: xmflx(ia,ja) ! x-momentum flux at the surface [kg/m/s2]
    real(RP), intent(out) :: ymflx(ia,ja) ! y-momentum flux at the surface [kg/m/s2]
    real(RP), intent(out) :: shflx(ia,ja) ! sensible heat flux at the surface [J/m2/s]
    real(RP), intent(out) :: lhflx(ia,ja) ! latent heat flux at the surface [J/m2/s]
    real(RP), intent(out) :: ghflx(ia,ja) ! ground heat flux at the surface [J/m2/s]
    real(RP), intent(out) :: u10  (ia,ja) ! velocity u at 10m [m/s]
    real(RP), intent(out) :: v10  (ia,ja) ! velocity v at 10m [m/s]
    real(RP), intent(out) :: t2   (ia,ja) ! temperature at 2m [K]
    real(RP), intent(out) :: q2   (ia,ja) ! water vapor at 2m [kg/kg]

    real(RP), intent(in) :: tmpa(ia,ja) ! temperature at the lowest atmospheric layer [K]
    real(RP), intent(in) :: prsa(ia,ja) ! pressure at the lowest atmospheric layer [Pa]
    real(RP), intent(in) :: wa  (ia,ja) ! velocity w at the lowest atmospheric layer [m/s]
    real(RP), intent(in) :: ua  (ia,ja) ! velocity u at the lowest atmospheric layer [m/s]
    real(RP), intent(in) :: va  (ia,ja) ! velocity v at the lowest atmospheric layer [m/s]
    real(RP), intent(in) :: rhoa(ia,ja) ! density at the lowest atmospheric layer [kg/m3]
    real(RP), intent(in) :: qva (ia,ja) ! ratio of water vapor mass to total mass at the lowest atmospheric layer [kg/kg]
    real(RP), intent(in) :: z1  (ia,ja) ! cell center height at the lowest atmospheric layer [m]
    real(RP), intent(in) :: pbl (ia,ja) ! the top of atmospheric mixing layer [m]
    real(RP), intent(in) :: prss(ia,ja) ! pressure at the surface [Pa]
    real(RP), intent(in) :: lwd (ia,ja) ! downward long-wave radiation flux at the surface [J/m2/s]
    real(RP), intent(in) :: swd (ia,ja) ! downward short-wave radiation flux at the surface [J/m2/s]

    real(RP), intent(in) :: tg    (ia,ja) ! soil temperature [K]
    real(RP), intent(in) :: lst   (ia,ja) ! land surface temperature [K]
    real(RP), intent(in) :: qvef  (ia,ja) ! efficiency of evaporation (0-1)
    real(RP), intent(in) :: alb_lw(ia,ja) ! surface albedo for LW (0-1)
    real(RP), intent(in) :: alb_sw(ia,ja) ! surface albedo for SW (0-1)
    real(RP), intent(in) :: dzg   (ia,ja) ! soil depth [m]
    real(RP), intent(in) :: tcs   (ia,ja) ! thermal conductivity for soil [J/m/K/s]
    real(RP), intent(in) :: z0m   (ia,ja) ! roughness length for momemtum [m]
    real(RP), intent(in) :: z0h   (ia,ja) ! roughness length for heat [m]
    real(RP), intent(in) :: z0e   (ia,ja) ! roughness length for vapor [m]
    real(DP), intent(in) :: dt            ! delta time

    ! works
    real(RP) :: lst1(ia,ja)

    real(RP) :: res    ! residual

    real(RP) :: ustar  ! friction velocity [m]
    real(RP) :: tstar  ! friction potential temperature [K]
    real(RP) :: qstar  ! friction water vapor mass ratio [kg/kg]
    real(RP) :: uabs   ! modified absolute velocity [m/s]
    real(RP) :: qvsat  ! saturation water vapor mixing ratio at surface [kg/kg]
    real(RP) :: qvs    ! water vapor mixing ratio at surface [kg/kg]

    real(RP) :: fracu10 ! calculation parameter for U10 [-]
    real(RP) :: fract2  ! calculation parameter for T2 [-]
    real(RP) :: fracq2  ! calculation parameter for Q2 [-]

    real(RP) :: lhv(ia,ja)    ! latent heat of vaporization [J/kg]

    integer  :: i, j
    !---------------------------------------------------------------------------

    if( io_l ) write(io_fid_log,*) '*** Land surface step: Thick-Slab'

    call hydrometeor_lhv( lhv(:,:), tmpa(:,:) )

    ! update surface temperature
    do j = js, je
    do i = is, ie
      lst1(i,j) = tg(i,j)
      lst_t(i,j) = ( lst1(i,j) - lst(i,j) ) / dt
    end do
    end do

    ! calculate surface flux
    do j = js, je
    do i = is, ie

      if( landuse_fact_land(i,j) > 0.0_rp ) then

        call qsat( qvsat,     & ! [OUT]
                   lst1(i,j), & ! [IN]
                   prss(i,j)  ) ! [IN]

        qvs  = ( 1.0_rp - qvef(i,j) ) * qva(i,j) + qvef(i,j) * qvsat

        call bulkflux( &
            ustar,     & ! [OUT]
            tstar,     & ! [OUT]
            qstar,     & ! [OUT]
            uabs,      & ! [OUT]
            fracu10,   & ! [OUT]
            fract2,    & ! [OUT]
            fracq2,    & ! [OUT]
            tmpa(i,j), & ! [IN]
            lst1(i,j), & ! [IN]
            prsa(i,j), & ! [IN]
            prss(i,j), & ! [IN]
            qva(i,j), & ! [IN]
            qvs,       & ! [IN]
            ua(i,j), & ! [IN]
            va(i,j), & ! [IN]
            z1(i,j), & ! [IN]
            pbl(i,j), & ! [IN]
            z0m(i,j), & ! [IN]
            z0h(i,j), & ! [IN]
            z0e(i,j)  ) ! [IN]

        zmflx(i,j) = -rhoa(i,j) * ustar**2 / uabs * wa(i,j)
        xmflx(i,j) = -rhoa(i,j) * ustar**2 / uabs * ua(i,j)
        ymflx(i,j) = -rhoa(i,j) * ustar**2 / uabs * va(i,j)

        shflx(i,j) = - rhoa(i,j) * ustar * tstar &
                   * cpdry * ( prss(i,j) / pre00 )**( rdry/cpdry )
        lhflx(i,j) = - rhoa(i,j) * ustar * qstar * lhv(i,j)

        ghflx(i,j) = -2.0_rp * tcs(i,j) * ( lst1(i,j) - tg(i,j) ) / dzg(i,j)

        ! calculation for residual
        res = ( 1.0_rp - alb_sw(i,j) ) * swd(i,j) &
            + ( 1.0_rp - alb_lw(i,j) ) * ( lwd(i,j) - stb * lst1(i,j)**4 ) &
            - shflx(i,j) - lhflx(i,j) + ghflx(i,j)

        ! put residual in ground heat flux
        ghflx(i,j) = ghflx(i,j) - res

        ! diagnostic variables considering unstable/stable state
        !U10(i,j) = FracU10 * UA(i,j)
        !V10(i,j) = FracU10 * VA(i,j)
        !T2 (i,j) = ( 1.0_RP - FracT2 ) * LST1(i,j) + FracT2 * TMPA(i,j)
        !Q2 (i,j) = ( 1.0_RP - FracQ2 ) * QVS       + FracQ2 * QVA (i,j)

        ! diagnostic variables for neutral state
        u10(i,j) = ua(i,j) * log( 10.0_rp / z0m(i,j) ) / log( z1(i,j) / z0m(i,j) )
        v10(i,j) = va(i,j) * log( 10.0_rp / z0m(i,j) ) / log( z1(i,j) / z0m(i,j) )
        t2(i,j) = lst1(i,j) + ( tmpa(i,j) - lst1(i,j) ) * ( log(  2.0_rp / z0m(i,j) ) * log(  2.0_rp / z0h(i,j) ) ) &
                                                         / ( log( z1(i,j) / z0m(i,j) ) * log( z1(i,j) / z0h(i,j) ) )
        q2(i,j) = qvs       + (  qva(i,j) - qvs       ) * ( log(  2.0_rp / z0m(i,j) ) * log(  2.0_rp / z0e(i,j) ) ) &
                                                         / ( log( z1(i,j) / z0m(i,j) ) * log( z1(i,j) / z0e(i,j) ) )

      else

        ! not calculate surface flux
        zmflx(i,j) = 0.0_rp
        xmflx(i,j) = 0.0_rp
        ymflx(i,j) = 0.0_rp
        shflx(i,j) = 0.0_rp
        lhflx(i,j) = 0.0_rp
        ghflx(i,j) = 0.0_rp
        u10(i,j) = 0.0_rp
        v10(i,j) = 0.0_rp
        t2(i,j) = 0.0_rp
        q2(i,j) = 0.0_rp

      end if

    end do
    end do

    return
  end subroutine land_sfc_thick_slab

end module scale_land_sfc_thick_slab