INTERFACE:
subroutine CanopyFluxes(lbg, ubg, lbc, ubc, lbp, ubp, &
num_nolakep, filter_nolakep)
DESCRIPTION:
1. Calculates the leaf temperature: 2. Calculates the leaf fluxes, transpiration, photosynthesis and updates the dew accumulation due to evaporation.
Method: Use the Newton-Raphson iteration to solve for the foliage temperature that balances the surface energy budget:
f(t_veg) = Net radiation - Sensible - Latent = 0 f(t_veg) + d(f)/d(t_veg) * dt_veg = 0 (*)
Note: (1) In solving for t_veg, t_grnd is given from the previous timestep. (2) The partial derivatives of aerodynamical resistances, which cannot be determined analytically, are ignored for d(H)/dT and d(LE)/dT (3) The weighted stomatal resistance of sunlit and shaded foliage is used (4) Canopy air temperature and humidity are derived from => Hc + Hg = Ha => Ec + Eg = Ea (5) Energy loss is due to: numerical truncation of energy budget equation (*); and "ecidif" (see the code) which is dropped into the sensible heat (6) The convergence criteria: the difference, del = t_veg(n+1)-t_veg(n) and del2 = t_veg(n)-t_veg(n-1) less than 0.01 K, and the difference of water flux from the leaf between the iteration step (n+1) and (n) less than 0.1 W/m2; or the iterative steps over 40.
USES:
use shr_kind_mod , only : r8 => shr_kind_r8
use clmtype
use clm_atmlnd , only : clm_a2l
use clm_time_manager , only : get_step_size, get_prev_date
use clm_varpar , only : nlevgrnd, nlevsno
use clm_varcon , only : sb, cpair, hvap, vkc, grav, denice, &
denh2o, tfrz, csoilc, tlsai_crit, alpha_aero, &
isecspday, degpsec
use shr_const_mod , only : SHR_CONST_TKFRZ
use pftvarcon , only : nirrig
use QSatMod , only : QSat
use FrictionVelocityMod, only : FrictionVelocity, MoninObukIni
use spmdMod , only : masterproc
ARGUMENTS:
implicit none
integer, intent(in) :: lbg, ubg ! gridcell bounds
integer, intent(in) :: lbc, ubc ! column bounds
integer, intent(in) :: lbp, ubp ! pft bounds
integer, intent(in) :: num_nolakep ! number of column non-lake points in pft filter
integer, intent(in) :: filter_nolakep(ubp-lbp+1) ! pft filter for non-lake points
CALLED FROM:
subroutine Biogeophysics1 in module Biogeophysics1ModREVISION HISTORY:
15 September 1999: Yongjiu Dai; Initial code 15 December 1999: Paul Houser and Jon Radakovich; F90 Revision 12/19/01, Peter Thornton Changed tg to t_grnd for consistency with other routines 1/29/02, Peter Thornton Migrate to new data structures, new calling protocol. For now co2 and o2 partial pressures are hardwired, but they should be coming in from forc_pco2 and forc_po2. Keeping the same hardwired values as in CLM2 to assure bit-for-bit results in the first comparisons. 27 February 2008: Keith Oleson; Sparse/dense aerodynamic parameters from X. Zeng 6 March 2009: Peter Thornton; Daylength control on Vcmax, from Bill BauerleLOCAL VARIABLES:
local pointers to implicit in variables
integer , pointer :: frac_veg_nosno(:) ! frac of veg not covered by snow (0 OR 1 now) [-]
integer , pointer :: ivt(:) ! pft vegetation type
integer , pointer :: pcolumn(:) ! pft's column index
integer , pointer :: plandunit(:) ! pft's landunit index
integer , pointer :: pgridcell(:) ! pft's gridcell index
real(r8), pointer :: forc_th(:) ! atmospheric potential temperature (Kelvin)
real(r8), pointer :: t_grnd(:) ! ground surface temperature [K]
real(r8), pointer :: thm(:) ! intermediate variable (forc_t+0.0098*forc_hgt_t_pft)
real(r8), pointer :: qg(:) ! specific humidity at ground surface [kg/kg]
real(r8), pointer :: thv(:) ! virtual potential temperature (kelvin)
real(r8), pointer :: z0mv(:) ! roughness length over vegetation, momentum [m]
real(r8), pointer :: z0hv(:) ! roughness length over vegetation, sensible heat [m]
real(r8), pointer :: z0qv(:) ! roughness length over vegetation, latent heat [m]
real(r8), pointer :: z0mg(:) ! roughness length of ground, momentum [m]
real(r8), pointer :: dqgdT(:) ! temperature derivative of "qg"
real(r8), pointer :: htvp(:) ! latent heat of evaporation (/sublimation) [J/kg]
real(r8), pointer :: emv(:) ! ground emissivity
real(r8), pointer :: emg(:) ! vegetation emissivity
real(r8), pointer :: forc_pbot(:) ! atmospheric pressure (Pa)
real(r8), pointer :: forc_pco2(:) ! partial pressure co2 (Pa)
#if (defined C13)
! 4/14/05: PET
! Adding isotope code
real(r8), pointer :: forc_pc13o2(:) ! partial pressure c13o2 (Pa)
#endif
real(r8), pointer :: forc_po2(:) ! partial pressure o2 (Pa)
real(r8), pointer :: forc_q(:) ! atmospheric specific humidity (kg/kg)
real(r8), pointer :: forc_u(:) ! atmospheric wind speed in east direction (m/s)
real(r8), pointer :: forc_v(:) ! atmospheric wind speed in north direction (m/s)
real(r8), pointer :: forc_hgt_u_pft(:) !observational height of wind at pft level [m]
real(r8), pointer :: forc_rho(:) ! density (kg/m**3)
real(r8), pointer :: forc_lwrad(:) ! downward infrared (longwave) radiation (W/m**2)
real(r8), pointer :: displa(:) ! displacement height (m)
real(r8), pointer :: elai(:) ! one-sided leaf area index with burying by snow
real(r8), pointer :: esai(:) ! one-sided stem area index with burying by snow
real(r8), pointer :: fdry(:) ! fraction of foliage that is green and dry [-]
real(r8), pointer :: fwet(:) ! fraction of canopy that is wet (0 to 1)
real(r8), pointer :: laisun(:) ! sunlit leaf area
real(r8), pointer :: laisha(:) ! shaded leaf area
real(r8), pointer :: sabv(:) ! solar radiation absorbed by vegetation (W/m**2)
real(r8), pointer :: watsat(:,:) ! volumetric soil water at saturation (porosity)
real(r8), pointer :: watdry(:,:) ! btran parameter for btran=0
real(r8), pointer :: watopt(:,:) ! btran parameter for btran = 1
real(r8), pointer :: h2osoi_ice(:,:)! ice lens (kg/m2)
real(r8), pointer :: h2osoi_liq(:,:)! liquid water (kg/m2)
real(r8), pointer :: dz(:,:) ! layer depth (m)
real(r8), pointer :: t_soisno(:,:) ! soil temperature (Kelvin)
real(r8), pointer :: sucsat(:,:) ! minimum soil suction (mm)
real(r8), pointer :: bsw(:,:) ! Clapp and Hornberger "b"
real(r8), pointer :: rootfr(:,:) ! fraction of roots in each soil layer
real(r8), pointer :: dleaf(:) ! characteristic leaf dimension (m)
real(r8), pointer :: smpso(:) ! soil water potential at full stomatal opening (mm)
real(r8), pointer :: smpsc(:) ! soil water potential at full stomatal closure (mm)
real(r8), pointer :: frac_sno(:) ! fraction of ground covered by snow (0 to 1)
real(r8), pointer :: htop(:) ! canopy top(m)
real(r8), pointer :: snowdp(:) ! snow height (m)
real(r8), pointer :: soilbeta(:) ! soil wetness relative to field capacity
real(r8), pointer :: lat(:) ! latitude (radians)
real(r8), pointer :: decl(:) ! declination angle (radians)
real(r8), pointer :: max_dayl(:) !maximum daylength for this column (s)
real(r8), pointer :: londeg(:) ! longitude (degrees) (for calculation of local time)
local pointers to implicit inout arguments
real(r8), pointer :: cgrnds(:) ! deriv. of soil sensible heat flux wrt soil temp [w/m2/k]
real(r8), pointer :: cgrndl(:) ! deriv. of soil latent heat flux wrt soil temp [w/m**2/k]
real(r8), pointer :: t_veg(:) ! vegetation temperature (Kelvin)
real(r8), pointer :: t_ref2m(:) ! 2 m height surface air temperature (Kelvin)
real(r8), pointer :: q_ref2m(:) ! 2 m height surface specific humidity (kg/kg)
real(r8), pointer :: t_ref2m_r(:) ! Rural 2 m height surface air temperature (Kelvin)
real(r8), pointer :: rh_ref2m(:) ! 2 m height surface relative humidity (%)
real(r8), pointer :: rh_ref2m_r(:) ! Rural 2 m height surface relative humidity (%)
real(r8), pointer :: h2ocan(:) ! canopy water (mm H2O)
real(r8), pointer :: cisun(:) !sunlit intracellular CO2 (Pa)
real(r8), pointer :: cisha(:) !shaded intracellular CO2 (Pa)
local pointers to implicit out arguments
real(r8), pointer :: rb1(:) ! boundary layer resistance (s/m)
real(r8), pointer :: cgrnd(:) ! deriv. of soil energy flux wrt to soil temp [w/m2/k]
real(r8), pointer :: dlrad(:) ! downward longwave radiation below the canopy [W/m2]
real(r8), pointer :: ulrad(:) ! upward longwave radiation above the canopy [W/m2]
real(r8), pointer :: ram1(:) ! aerodynamical resistance (s/m)
real(r8), pointer :: btran(:) ! transpiration wetness factor (0 to 1)
real(r8), pointer :: rssun(:) ! sunlit stomatal resistance (s/m)
real(r8), pointer :: rssha(:) ! shaded stomatal resistance (s/m)
real(r8), pointer :: psnsun(:) ! sunlit leaf photosynthesis (umol CO2 /m**2/ s)
real(r8), pointer :: psnsha(:) ! shaded leaf photosynthesis (umol CO2 /m**2/ s)
#if (defined C13)
! 4/14/05: PET
! Adding isotope code
real(r8), pointer :: c13_psnsun(:) ! sunlit leaf photosynthesis (umol 13CO2 /m**2/ s)
real(r8), pointer :: c13_psnsha(:) ! shaded leaf photosynthesis (umol 13CO2 /m**2/ s)
! 4/21/05: PET
! Adding isotope code
real(r8), pointer :: rc13_canair(:) !C13O2/C12O2 in canopy air
real(r8), pointer :: rc13_psnsun(:) !C13O2/C12O2 in sunlit canopy psn flux
real(r8), pointer :: rc13_psnsha(:) !C13O2/C12O2 in shaded canopy psn flux
real(r8), pointer :: alphapsnsun(:) !fractionation factor in sunlit canopy psn flux
real(r8), pointer :: alphapsnsha(:) !fractionation factor in shaded canopy psn flux
#endif
real(r8), pointer :: qflx_tran_veg(:) ! vegetation transpiration (mm H2O/s) (+ = to atm)
real(r8), pointer :: dt_veg(:) ! change in t_veg, last iteration (Kelvin)
real(r8), pointer :: qflx_evap_veg(:) ! vegetation evaporation (mm H2O/s) (+ = to atm)
real(r8), pointer :: eflx_sh_veg(:) ! sensible heat flux from leaves (W/m**2) [+ to atm]
real(r8), pointer :: taux(:) ! wind (shear) stress: e-w (kg/m/s**2)
real(r8), pointer :: tauy(:) ! wind (shear) stress: n-s (kg/m/s**2)
real(r8), pointer :: eflx_sh_grnd(:) ! sensible heat flux from ground (W/m**2) [+ to atm]
real(r8), pointer :: qflx_evap_soi(:) ! soil evaporation (mm H2O/s) (+ = to atm)
real(r8), pointer :: fpsn(:) ! photosynthesis (umol CO2 /m**2 /s)
real(r8), pointer :: rootr(:,:) ! effective fraction of roots in each soil layer
real(r8), pointer :: rresis(:,:) ! root resistance by layer (0-1) (nlevgrnd)
real(r8), pointer :: irrig_rate(:) ! current irrigation rate [mm/s]
integer, pointer :: n_irrig_steps_left(:) ! # of time steps for which we still need to irrigate today
!OTHER LOCAL VARIABLES: