INTERFACE:
subroutine CanopyFluxes (c, p)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 globals
use clm_varcon, only : sb, cpair, hvap, vkc, grav, denice, denh2o, tfrz, &
csoilc
use QSatMod, only : QSat
use FrictionVelocityMod, only : MoninObukIni, FrictionVelocity
ARGUMENTS:
implicit none
type (column_type), target, intent(inout) :: c !column-level data
type (pft_type), target, intent(inout) :: p !pft-level data
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_co2 and forc_o2. Keeping the same hardwired values as in CLM2 to assure bit-for-bit results in the first comparisons.LOCAL VARIABLES:
local pointers to implicit in scalars
real(r8), pointer :: 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)
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 :: 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_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 !atmospheric reference height (m)
real(r8), pointer :: forc_hgt_u !observational height of wind [m]
real(r8), pointer :: forc_hgt_t !observational height of temperature [m]
real(r8), pointer :: forc_hgt_q !observational height of humidity [m]
real(r8), pointer :: forc_rho !density (kg/m**3)
real(r8), pointer :: forc_co2 !atmospheric CO2 concentration (Pa)
real(r8), pointer :: forc_o2 !atmospheric O2 concentration (Pa)
real(r8), pointer :: forc_lwrad !downward infrared (longwave) radiation (W/m**2)
real(r8), pointer :: smpmax !wilting point potential in mm
real(r8), pointer :: displa !displacement height (m)
real(r8), pointer :: dleaf !characteristic leaf dimension (m)
real(r8), pointer :: parsun !average absorbed PAR for sunlit leaves (W/m**2)
real(r8), pointer :: parsha !average absorbed PAR for shaded leaves (W/m**2)
real(r8), pointer :: qe25 !quantum efficiency at 25C (umol CO2 / umol photon)
real(r8), pointer :: vcmx25 !max rate of carboxylation at 25C (umol CO2/m**2/s)
real(r8), pointer :: mp !slope of conductance-to-photosynthesis relationship
real(r8), pointer :: c3psn !photosynthetic pathway: 0. = c4, 1. = c3
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
integer , pointer :: frac_veg_nosno !frac of veg not covered by snow (0 OR 1 now) [-]
real(r8), pointer :: sabv !solar radiation absorbed by vegetation (W/m**2)
local pointers to implicit inout scalars
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 :: h2ocan !canopy water (mm H2O)
real(r8), pointer :: annpsnpot !annual potential photosynthesis (umol CO2 /m**2)
real(r8), pointer :: annpsn !annual photosynthesis (umol CO2 /m**2)
local pointers to implicit out scalars
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 :: u10 !10-m wind (m/s) (for dust model)
real(r8), pointer :: fv !friction velocity (m/s) (for dust model)
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)
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)
local pointers to implicit in arrays
real(r8), dimension(:), pointer :: watsat !volumetric soil water at saturation (porosity)
real(r8), dimension(:), pointer :: h2osoi_ice !ice lens (kg/m2)
real(r8), dimension(:), pointer :: h2osoi_liq !liquid water (kg/m2)
real(r8), dimension(:), pointer :: dz !layer depth (m)
real(r8), dimension(:), pointer :: t_soisno !soil temperature (Kelvin)
real(r8), dimension(:), pointer :: sucsat !minimum soil suction (mm)
real(r8), dimension(:), pointer :: bsw !Clapp and Hornberger "b"
real(r8), dimension(:), pointer :: rootfr !fraction of roots in each soil layer
local pointers to implicit out arrays
real(r8), dimension(:), pointer :: rootr !effective fraction of roots in each soil layer