diff --git a/build/source/engine/computFlux.f90 b/build/source/engine/computFlux.f90
index 73ed736c9b979dd1fe5fdbdb6d8088c878d2ea01..0b7034aee1689650c6b340e3ad1fd868288824f7 100755
--- a/build/source/engine/computFlux.f90
+++ b/build/source/engine/computFlux.f90
@@ -26,9 +26,9 @@ USE nrtype
! provide access to the derived types to define the data structures
USE data_types,only:&
var_i, & ! data vector (i4b)
- var_d, & ! data vector (dp)
+ var_d, & ! data vector (rkind)
var_ilength, & ! data vector with variable length dimension (i4b)
- var_dlength, & ! data vector with variable length dimension (dp)
+ var_dlength, & ! data vector with variable length dimension (rkind)
model_options ! defines the model decisions
! indices that define elements of the data structures
@@ -71,27 +71,27 @@ USE multiconst,only:&
! look-up values for the choice of groundwater representation (local-column, or single-basin)
USE mDecisions_module,only: &
- localColumn, & ! separate groundwater representation in each local soil column
- singleBasin ! single groundwater store over the entire basin
+ localColumn, & ! separate groundwater representation in each local soil column
+ singleBasin ! single groundwater store over the entire basin
! look-up values for the choice of groundwater parameterization
USE mDecisions_module,only: &
- qbaseTopmodel, & ! TOPMODEL-ish baseflow parameterization
- bigBucket, & ! a big bucket (lumped aquifer model)
- noExplicit ! no explicit groundwater parameterization
+ qbaseTopmodel, & ! TOPMODEL-ish baseflow parameterization
+ bigBucket, & ! a big bucket (lumped aquifer model)
+ noExplicit ! no explicit groundwater parameterization
! look-up values for the form of Richards' equation
USE mDecisions_module,only: &
- moisture, & ! moisture-based form of Richards' equation
- mixdform ! mixed form of Richards' equation
+ moisture, & ! moisture-based form of Richards' equation
+ mixdform ! mixed form of Richards' equation
! look-up values for the choice of boundary conditions for hydrology
USE mDecisions_module,only: &
- prescribedHead, & ! prescribed head (volumetric liquid water content for mixed form of Richards' eqn)
- funcBottomHead, & ! function of matric head in the lower-most layer
- freeDrainage, & ! free drainage
- liquidFlux, & ! liquid water flux
- zeroFlux ! zero flux
+ prescribedHead, & ! prescribed head (volumetric liquid water content for mixed form of Richards' eqn)
+ funcBottomHead, & ! function of matric head in the lower-most layer
+ freeDrainage, & ! free drainage
+ liquidFlux, & ! liquid water flux
+ zeroFlux ! zero flux
implicit none
private
@@ -100,52 +100,53 @@ public::soilCmpres
public::soilCmpresSundials
contains
+
! *********************************************************************************************************
! public subroutine computFlux: compute model fluxes
! *********************************************************************************************************
subroutine computFlux(&
- ! input-output: model control
- nSnow, & ! intent(in): number of snow layers
- nSoil, & ! intent(in): number of soil layers
- nLayers, & ! intent(in): total number of layers
- firstSubStep, & ! intent(in): flag to indicate if we are processing the first sub-step
- firstFluxCall, & ! intent(inout): flag to denote the first flux call
- firstSplitOper, & ! intent(in): flag to indicate if we are processing the first flux call in a splitting operation
- computeVegFlux, & ! intent(in): flag to indicate if we need to compute fluxes over vegetation
- scalarSolution, & ! intent(in): flag to indicate the scalar solution
- insideIDA, & ! intent(in): flag if inside Sundials solver
- drainageMeltPond, & ! intent(in): drainage from the surface melt pond (kg m-2 s-1)
- ! input: state variables
- scalarCanairTempTrial, & ! intent(in): trial value for the temperature of the canopy air space (K)
- scalarCanopyTempTrial, & ! intent(in): trial value for the temperature of the vegetation canopy (K)
- mLayerTempTrial, & ! intent(in): trial value for the temperature of each snow and soil layer (K)
- mLayerMatricHeadLiqTrial, & ! intent(in): trial value for the liquid water matric potential in each soil layer (m)
- mLayerMatricHeadTrial, & ! intent(in) trial vector of total water matric potential (m)
- scalarAquiferStorageTrial,& ! intent(in): trial value of storage of water in the aquifer (m)
- ! input: diagnostic variables defining the liquid water and ice content
- scalarCanopyLiqTrial, & ! intent(in): trial value for the liquid water on the vegetation canopy (kg m-2)
- scalarCanopyIceTrial, & ! intent(in): trial value for the ice on the vegetation canopy (kg m-2)
- mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each snow and soil layer (-)
- mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice in each snow and soil layer (-)
- ! input: data structures
- model_decisions, & ! intent(in): model decisions
- type_data, & ! intent(in): type of vegetation and soil
- attr_data, & ! intent(in): spatial attributes
- mpar_data, & ! intent(in): model parameters
- forc_data, & ! intent(in): model forcing data
- bvar_data, & ! intent(in): average model variables for the entire basin
- prog_data, & ! intent(in): model prognostic variables for a local HRU
- indx_data, & ! intent(in): index data
- ! input-output: data structures
- diag_data, & ! intent(inout): model diagnostic variables for a local HRU
- flux_data, & ! intent(inout): model fluxes for a local HRU
- deriv_data, & ! intent(inout): derivatives in model fluxes w.r.t. relevant state variables
- ! input-output: flux vector and baseflow derivatives
- ixSaturation, & ! intent(inout): index of the lowest saturated layer (NOTE: only computed on the first iteration)
- dBaseflow_dMatric, & ! intent(out): derivative in baseflow w.r.t. matric head (s-1)
- fluxVec, & ! intent(out): flux vector (mixed units)
- ! output: error control
- err,message) ! intent(out): error code and error message
+ ! input-output: model control
+ nSnow, & ! intent(in): number of snow layers
+ nSoil, & ! intent(in): number of soil layers
+ nLayers, & ! intent(in): total number of layers
+ firstSubStep, & ! intent(in): flag to indicate if we are processing the first sub-step
+ firstFluxCall, & ! intent(inout): flag to denote the first flux call
+ firstSplitOper, & ! intent(in): flag to indicate if we are processing the first flux call in a splitting operation
+ computeVegFlux, & ! intent(in): flag to indicate if we need to compute fluxes over vegetation
+ scalarSolution, & ! intent(in): flag to indicate the scalar solution
+ insideIDA, & ! intent(in): flag to indicate inside Sundials Solver (do not require longwave to be balanced)
+ drainageMeltPond, & ! intent(in): drainage from the surface melt pond (kg m-2 s-1)
+ ! input: state variables
+ scalarCanairTempTrial, & ! intent(in): trial value for the temperature of the canopy air space (K)
+ scalarCanopyTempTrial, & ! intent(in): trial value for the temperature of the vegetation canopy (K)
+ mLayerTempTrial, & ! intent(in): trial value for the temperature of each snow and soil layer (K)
+ mLayerMatricHeadLiqTrial, & ! intent(in): trial value for the liquid water matric potential in each soil layer (m)
+ mLayerMatricHeadTrial, & ! intent(in): trial vector of total water matric potential (m)
+ scalarAquiferStorageTrial,& ! intent(in): trial value of storage of water in the aquifer (m)
+ ! input: diagnostic variables defining the liquid water and ice content
+ scalarCanopyLiqTrial, & ! intent(in): trial value for the liquid water on the vegetation canopy (kg m-2)
+ scalarCanopyIceTrial, & ! intent(in): trial value for the ice on the vegetation canopy (kg m-2)
+ mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each snow and soil layer (-)
+ mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice in each snow and soil layer (-)
+ ! input: data structures
+ model_decisions, & ! intent(in): model decisions
+ type_data, & ! intent(in): type of vegetation and soil
+ attr_data, & ! intent(in): spatial attributes
+ mpar_data, & ! intent(in): model parameters
+ forc_data, & ! intent(in): model forcing data
+ bvar_data, & ! intent(in): average model variables for the entire basin
+ prog_data, & ! intent(in): model prognostic variables for a local HRU
+ indx_data, & ! intent(in): index data
+ ! input-output: data structures
+ diag_data, & ! intent(inout): model diagnostic variables for a local HRU
+ flux_data, & ! intent(inout): model fluxes for a local HRU
+ deriv_data, & ! intent(inout): derivatives in model fluxes w.r.t. relevant state variables
+ ! input-output: flux vector and baseflow derivatives
+ ixSaturation, & ! intent(inout): index of the lowest saturated layer (NOTE: only computed on the first iteration)
+ dBaseflow_dMatric, & ! intent(out): derivative in baseflow w.r.t. matric head (s-1)
+ fluxVec, & ! intent(out): flux vector (mixed units)
+ ! output: error control
+ err,message) ! intent(out): error code and error message
! provide access to flux subroutines
USE vegNrgFlux_module,only:vegNrgFlux ! compute energy fluxes over vegetation
USE ssdNrgFlux_module,only:ssdNrgFlux ! compute energy fluxes throughout the snow and soil subdomains
@@ -159,59 +160,60 @@ subroutine computFlux(&
! * dummy variables
! ---------------------------------------------------------------------------------------
! input-output: control
- integer(i4b),intent(in) :: nSnow ! number of snow layers
- integer(i4b),intent(in) :: nSoil ! number of soil layers
- integer(i4b),intent(in) :: nLayers ! total number of layers
- logical(lgt),intent(in) :: firstSubStep ! flag to indicate if we are processing the first sub-step
- logical(lgt),intent(inout) :: firstFluxCall ! flag to indicate if we are processing the first flux call
- logical(lgt),intent(in) :: firstSplitOper ! flag to indicate if we are processing the first flux call in a splitting operation
- logical(lgt),intent(in) :: computeVegFlux ! flag to indicate if computing fluxes over vegetation
- logical(lgt),intent(in) :: scalarSolution ! flag to denote if implementing the scalar solution
- logical(lgt),intent(in) :: insideIDA ! flag if inside Sundials solver
- real(dp),intent(in) :: drainageMeltPond ! drainage from the surface melt pond (kg m-2 s-1)
+ integer(i4b),intent(in) :: nSnow ! number of snow layers
+ integer(i4b),intent(in) :: nSoil ! number of soil layers
+ integer(i4b),intent(in) :: nLayers ! total number of layers
+ logical(lgt),intent(in) :: firstSubStep ! flag to indicate if we are processing the first sub-step
+ logical(lgt),intent(inout) :: firstFluxCall ! flag to indicate if we are processing the first flux call
+ logical(lgt),intent(in) :: firstSplitOper ! flag to indicate if we are processing the first flux call in a splitting operation
+ logical(lgt),intent(in) :: computeVegFlux ! flag to indicate if computing fluxes over vegetation
+ logical(lgt),intent(in) :: scalarSolution ! flag to denote if implementing the scalar solution
+ logical(lgt),intent(in) :: insideIDA ! flag if inside Sundials solver
+ real(rkind),intent(in) :: drainageMeltPond ! drainage from the surface melt pond (kg m-2 s-1)
! input: state variables
- real(dp),intent(in) :: scalarCanairTempTrial ! trial value for temperature of the canopy air space (K)
- real(dp),intent(in) :: scalarCanopyTempTrial ! trial value for temperature of the vegetation canopy (K)
- real(dp),intent(in) :: mLayerTempTrial(:) ! trial value for temperature of each snow/soil layer (K)
- real(dp),intent(in) :: mLayerMatricHeadLiqTrial(:) ! trial value for the liquid water matric potential (m)
- real(dp),intent(in) :: mLayerMatricHeadTrial(:) ! trial value for the total water matric potential (m)
- real(dp),intent(in) :: scalarAquiferStorageTrial ! trial value of aquifer storage (m)
+ real(rkind),intent(in) :: scalarCanairTempTrial ! trial value for temperature of the canopy air space (K)
+ real(rkind),intent(in) :: scalarCanopyTempTrial ! trial value for temperature of the vegetation canopy (K)
+ real(rkind),intent(in) :: mLayerTempTrial(:) ! trial value for temperature of each snow/soil layer (K)
+ real(rkind),intent(in) :: mLayerMatricHeadLiqTrial(:) ! trial value for the liquid water matric potential (m)
+ real(rkind),intent(in) :: mLayerMatricHeadTrial(:) ! trial value for the total water matric potential (m)
+ real(rkind),intent(in) :: scalarAquiferStorageTrial ! trial value of aquifer storage (m)
! input: diagnostic variables
- real(dp),intent(in) :: scalarCanopyLiqTrial ! trial value for mass of liquid water on the vegetation canopy (kg m-2)
- real(dp),intent(in) :: scalarCanopyIceTrial ! trial value for mass of ice on the vegetation canopy (kg m-2)
- real(dp),intent(in) :: mLayerVolFracLiqTrial(:) ! trial value for volumetric fraction of liquid water (-)
- real(dp),intent(in) :: mLayerVolFracIceTrial(:) ! trial value for volumetric fraction of ice (-)
+ real(rkind),intent(in) :: scalarCanopyLiqTrial ! trial value for mass of liquid water on the vegetation canopy (kg m-2)
+ real(rkind),intent(in) :: scalarCanopyIceTrial ! trial value for mass of ice on the vegetation canopy (kg m-2)
+ real(rkind),intent(in) :: mLayerVolFracLiqTrial(:) ! trial value for volumetric fraction of liquid water (-)
+ real(rkind),intent(in) :: mLayerVolFracIceTrial(:) ! trial value for volumetric fraction of ice (-)
! input: data structures
- type(model_options),intent(in) :: model_decisions(:) ! model decisions
- type(var_i), intent(in) :: type_data ! type of vegetation and soil
- type(var_d), intent(in) :: attr_data ! spatial attributes
- type(var_dlength), intent(in) :: mpar_data ! model parameters
- type(var_d), intent(in) :: forc_data ! model forcing data
- type(var_dlength), intent(in) :: bvar_data ! model variables for the local basin
- type(var_dlength), intent(in) :: prog_data ! prognostic variables for a local HRU
- type(var_ilength), intent(in) :: indx_data ! indices defining model states and layers
+ type(model_options),intent(in) :: model_decisions(:) ! model decisions
+ type(var_i), intent(in) :: type_data ! type of vegetation and soil
+ type(var_d), intent(in) :: attr_data ! spatial attributes
+ type(var_dlength), intent(in) :: mpar_data ! model parameters
+ type(var_d), intent(in) :: forc_data ! model forcing data
+ type(var_dlength), intent(in) :: bvar_data ! model variables for the local basin
+ type(var_dlength), intent(in) :: prog_data ! prognostic variables for a local HRU
+ type(var_ilength), intent(in) :: indx_data ! indices defining model states and layers
! input-output: data structures
- type(var_dlength),intent(inout) :: diag_data ! diagnostic variables for a local HRU
- type(var_dlength),intent(inout) :: flux_data ! model fluxes for a local HRU
- type(var_dlength),intent(inout) :: deriv_data ! derivatives in model fluxes w.r.t. relevant state variables
+ type(var_dlength),intent(inout) :: diag_data ! diagnostic variables for a local HRU
+ type(var_dlength),intent(inout) :: flux_data ! model fluxes for a local HRU
+ type(var_dlength),intent(inout) :: deriv_data ! derivatives in model fluxes w.r.t. relevant state variables
! input-output: flux vector and baseflow derivatives
- integer(i4b),intent(inout) :: ixSaturation ! index of the lowest saturated layer (NOTE: only computed on the first iteration)
- real(dp),intent(out) :: dBaseflow_dMatric(:,:) ! derivative in baseflow w.r.t. matric head (s-1)
- real(dp),intent(out) :: fluxVec(:) ! model flux vector (mixed units)
+ integer(i4b),intent(inout) :: ixSaturation ! index of the lowest saturated layer (NOTE: only computed on the first iteration)
+ real(rkind),intent(out) :: dBaseflow_dMatric(:,:) ! derivative in baseflow w.r.t. matric head (s-1)
+ real(rkind),intent(out) :: fluxVec(:) ! model flux vector (mixed units)
! output: error control
- integer(i4b),intent(out) :: err ! error code
- character(*),intent(out) :: message ! error message
+ integer(i4b),intent(out) :: err ! error code
+ character(*),intent(out) :: message ! error message
! ---------------------------------------------------------------------------------------
! * local variables
! ---------------------------------------------------------------------------------------
integer(i4b) :: local_ixGroundwater ! local index for groundwater representation
integer(i4b) :: iLayer ! index of model layers
logical(lgt) :: doVegNrgFlux ! flag to compute the energy flux over vegetation
- real(dp),dimension(nSoil) :: dHydCond_dMatric ! derivative in hydraulic conductivity w.r.t matric head (s-1)
+ real(rkind),dimension(nSoil) :: dHydCond_dMatric ! derivative in hydraulic conductivity w.r.t matric head (s-1)
character(LEN=256) :: cmessage ! error message of downwind routine
- real(dp) :: above_soilLiqFluxDeriv ! derivative in layer above soil (canopy or snow) liquid flux w.r.t. liquid water
- real(dp) :: above_soildLiq_dTk ! derivative of layer above soil (canopy or snow) liquid flux w.r.t. temperature
- real(dp) :: above_soilFracLiq ! fraction of liquid water layer above soil (canopy or snow) (-)
+ real(rkind) :: above_soilLiqFluxDeriv ! derivative in layer above soil (canopy or snow) liquid flux w.r.t. liquid water
+ real(rkind) :: above_soildLiq_dTk ! derivative of layer above soil (canopy or snow) liquid flux w.r.t. temperature
+ real(rkind) :: above_soilFracLiq ! fraction of liquid water layer above soil (canopy or snow) (-)
+
! --------------------------------------------------------------
! initialize error control
err=0; message='computFlux/'
@@ -223,285 +225,283 @@ subroutine computFlux(&
! get the necessary variables for the flux computations
associate(&
- ! model decisions
- ixGroundwater => model_decisions(iLookDECISIONS%groundwatr)%iDecision ,& ! intent(in): [i4b] groundwater parameterization
- ixSpatialGroundwater => model_decisions(iLookDECISIONS%spatial_gw)%iDecision ,& ! intent(in): [i4b] spatial representation of groundwater (local-column or single-basin)
-
- ! domain boundary conditions
- upperBoundTemp => forc_data%var(iLookFORCE%airtemp) ,& ! intent(in): [dp] temperature of the upper boundary of the snow and soil domains (K)
- scalarRainfall => flux_data%var(iLookFLUX%scalarRainfall)%dat(1) ,& ! intent(in): [dp] rainfall rate (kg m-2 s-1)
-
- ! canopy and layer depth
- canopyDepth => diag_data%var(iLookDIAG%scalarCanopyDepth)%dat(1) ,& ! intent(in): [dp ] canopy depth (m)
- mLayerDepth => prog_data%var(iLookPROG%mLayerDepth)%dat ,& ! intent(in): [dp(:)] depth of each layer in the snow-soil sub-domain (m)
-
- ! indices of model state variables for the vegetation subdomain
- ixCasNrg => indx_data%var(iLookINDEX%ixCasNrg)%dat(1) ,& ! intent(in): [i4b] index of canopy air space energy state variable
- ixVegNrg => indx_data%var(iLookINDEX%ixVegNrg)%dat(1) ,& ! intent(in): [i4b] index of canopy energy state variable
- ixVegHyd => indx_data%var(iLookINDEX%ixVegHyd)%dat(1) ,& ! intent(in): [i4b] index of canopy hydrology state variable (mass)
- ixTopNrg => indx_data%var(iLookINDEX%ixTopNrg)%dat(1) ,& ! intent(in): [i4b] index of upper-most energy state in the snow+soil subdomain
- ixTopHyd => indx_data%var(iLookINDEX%ixTopHyd)%dat(1) ,& ! intent(in): [i4b] index of upper-most hydrology state in the snow+soil subdomain
- ixAqWat => indx_data%var(iLookINDEX%ixAqWat)%dat(1) ,& ! intent(in): [i4b] index of water storage in the aquifer
-
- ! indices of model state variables for the snow+soil domain
- ixSnowSoilNrg => indx_data%var(iLookINDEX%ixSnowSoilNrg)%dat ,& ! intent(in): [i4b(:)] indices for energy states in the snow+soil subdomain
- ixSnowSoilHyd => indx_data%var(iLookINDEX%ixSnowSoilHyd)%dat ,& ! intent(in): [i4b(:)] indices for hydrology states in the snow+soil subdomain
- layerType => indx_data%var(iLookINDEX%layerType)%dat ,& ! intent(in): [i4b(:)] type of layer (iname_soil or iname_snow)
-
- ! number of state variables of a specific type
- nSnowSoilNrg => indx_data%var(iLookINDEX%nSnowSoilNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the snow+soil domain
- nSnowOnlyNrg => indx_data%var(iLookINDEX%nSnowOnlyNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the snow domain
- nSoilOnlyNrg => indx_data%var(iLookINDEX%nSoilOnlyNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the soil domain
- nSnowSoilHyd => indx_data%var(iLookINDEX%nSnowSoilHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the snow+soil domain
- nSnowOnlyHyd => indx_data%var(iLookINDEX%nSnowOnlyHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the snow domain
- nSoilOnlyHyd => indx_data%var(iLookINDEX%nSoilOnlyHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the soil domain
-
- ! snow parameters
- snowfrz_scale => mpar_data%var(iLookPARAM%snowfrz_scale)%dat(1) ,& ! intent(in): [dp] scaling parameter for the snow freezing curve (K-1)
-
- ! derivatives
- dPsiLiq_dPsi0 => deriv_data%var(iLookDERIV%dPsiLiq_dPsi0 )%dat ,& ! intent(in): [dp(:)] derivative in liquid water matric pot w.r.t. the total water matric pot (-)
- dPsiLiq_dTemp => deriv_data%var(iLookDERIV%dPsiLiq_dTemp )%dat ,& ! intent(in): [dp(:)] derivative in the liquid water matric potential w.r.t. temperature
- mLayerdTheta_dTk => deriv_data%var(iLookDERIV%mLayerdTheta_dTk)%dat ,& ! intent(in): [dp(:)] derivative of volumetric liquid water content w.r.t. temperature
- dTheta_dTkCanopy => deriv_data%var(iLookDERIV%dTheta_dTkCanopy)%dat(1) ,& ! intent(in): [dp] derivative of volumetric liquid water content w.r.t. temperature
-
- ! number of flux calls
- numFluxCalls => diag_data%var(iLookDIAG%numFluxCalls)%dat(1) ,& ! intent(out): [dp] number of flux calls (-)
-
- ! net fluxes over the vegetation domain
- scalarCanairNetNrgFlux => flux_data%var(iLookFLUX%scalarCanairNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the canopy air space (W m-2)
- scalarCanopyNetNrgFlux => flux_data%var(iLookFLUX%scalarCanopyNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the vegetation canopy (W m-2)
- scalarGroundNetNrgFlux => flux_data%var(iLookFLUX%scalarGroundNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the ground surface (W m-2)
- scalarCanopyNetLiqFlux => flux_data%var(iLookFLUX%scalarCanopyNetLiqFlux)%dat(1) ,& ! intent(out): [dp] net liquid water flux for the vegetation canopy (kg m-2 s-1)
-
- ! net fluxes over the snow+soil domain
- mLayerNrgFlux => flux_data%var(iLookFLUX%mLayerNrgFlux)%dat ,& ! intent(out): [dp] net energy flux for each layer within the snow+soil domain (J m-3 s-1)
- mLayerLiqFluxSnow => flux_data%var(iLookFLUX%mLayerLiqFluxSnow)%dat ,& ! intent(out): [dp] net liquid water flux for each snow layer (s-1)
- mLayerLiqFluxSoil => flux_data%var(iLookFLUX%mLayerLiqFluxSoil)%dat ,& ! intent(out): [dp] net liquid water flux for each soil layer (s-1)
-
- ! evaporative fluxes
- scalarCanopyTranspiration => flux_data%var(iLookFLUX%scalarCanopyTranspiration)%dat(1) ,& ! intent(out): [dp] canopy transpiration (kg m-2 s-1)
- scalarCanopyEvaporation => flux_data%var(iLookFLUX%scalarCanopyEvaporation)%dat(1) ,& ! intent(out): [dp] canopy evaporation/condensation (kg m-2 s-1)
- scalarGroundEvaporation => flux_data%var(iLookFLUX%scalarGroundEvaporation)%dat(1) ,& ! intent(out): [dp] ground evaporation/condensation -- below canopy or non-vegetated (kg m-2 s-1)
- mLayerTranspire => flux_data%var(iLookFLUX%mLayerTranspire)%dat ,& ! intent(out): [dp(:)] transpiration loss from each soil layer (m s-1)
-
- ! fluxes for the snow+soil domain
- iLayerNrgFlux => flux_data%var(iLookFLUX%iLayerNrgFlux)%dat ,& ! intent(out): [dp(0:)] vertical energy flux at the interface of snow and soil layers
- iLayerLiqFluxSnow => flux_data%var(iLookFLUX%iLayerLiqFluxSnow)%dat ,& ! intent(out): [dp(0:)] vertical liquid water flux at snow layer interfaces (-)
- iLayerLiqFluxSoil => flux_data%var(iLookFLUX%iLayerLiqFluxSoil)%dat ,& ! intent(out): [dp(0:)] vertical liquid water flux at soil layer interfaces (-)
- mLayerHydCond => flux_data%var(iLookFLUX%mLayerHydCond)%dat ,& ! intent(out): [dp(:)] hydraulic conductivity in each soil layer (m s-1)
- mLayerBaseflow => flux_data%var(iLookFLUX%mLayerBaseflow)%dat ,& ! intent(out): [dp(:)] baseflow from each soil layer (m s-1)
- scalarSnowDrainage => flux_data%var(iLookFLUX%scalarSnowDrainage)%dat(1) ,& ! intent(out): [dp] drainage from the snow profile (m s-1)
- scalarSoilDrainage => flux_data%var(iLookFLUX%scalarSoilDrainage)%dat(1) ,& ! intent(out): [dp] drainage from the soil profile (m s-1)
- scalarSoilBaseflow => flux_data%var(iLookFLUX%scalarSoilBaseflow)%dat(1) ,& ! intent(out): [dp] total baseflow from the soil profile (m s-1)
-
- ! infiltration
- scalarInfilArea => diag_data%var(iLookDIAG%scalarInfilArea )%dat(1) ,& ! intent(out): [dp] fraction of unfrozen area where water can infiltrate (-)
- scalarFrozenArea => diag_data%var(iLookDIAG%scalarFrozenArea )%dat(1) ,& ! intent(out): [dp] fraction of area that is considered impermeable due to soil ice (-)
- scalarSoilControl => diag_data%var(iLookDIAG%scalarSoilControl )%dat(1) ,& ! intent(out): [dp] soil control on infiltration, zero or one
- scalarMaxInfilRate => flux_data%var(iLookFLUX%scalarMaxInfilRate)%dat(1) ,& ! intent(out): [dp] maximum infiltration rate (m s-1)
- scalarInfiltration => flux_data%var(iLookFLUX%scalarInfiltration)%dat(1) ,& ! intent(out): [dp] infiltration of water into the soil profile (m s-1)
- scalarFracLiqVeg => diag_data%var(iLookDIAG%scalarFracLiqVeg)%dat(1) ,& ! intent(inout): [dp] fraction of liquid water on vegetation (-)
- mLayerFracLiqSnow => diag_data%var(iLookDIAG%mLayerFracLiqSnow)%dat ,& ! intent(inout): [dp(:)] fraction of liquid water in each snow layer (-)
-
- ! boundary fluxes in the soil domain
- scalarThroughfallRain => flux_data%var(iLookFLUX%scalarThroughfallRain)%dat(1) ,& ! intent(out): [dp] rain that reaches the ground without ever touching the canopy (kg m-2 s-1)
- scalarCanopyLiqDrainage => flux_data%var(iLookFLUX%scalarCanopyLiqDrainage)%dat(1) ,& ! intent(out): [dp] drainage of liquid water from the vegetation canopy (kg m-2 s-1)
- scalarRainPlusMelt => flux_data%var(iLookFLUX%scalarRainPlusMelt)%dat(1) ,& ! intent(out): [dp] rain plus melt (m s-1)
- scalarSurfaceRunoff => flux_data%var(iLookFLUX%scalarSurfaceRunoff)%dat(1) ,& ! intent(out): [dp] surface runoff (m s-1)
- scalarExfiltration => flux_data%var(iLookFLUX%scalarExfiltration)%dat(1) ,& ! intent(out): [dp] exfiltration from the soil profile (m s-1)
- mLayerColumnOutflow => flux_data%var(iLookFLUX%mLayerColumnOutflow)%dat ,& ! intent(out): [dp(:)] column outflow from each soil layer (m3 s-1)
-
- ! fluxes for the aquifer
- scalarAquiferTranspire => flux_data%var(iLookFLUX%scalarAquiferTranspire)%dat(1) ,& ! intent(out): [dp] transpiration loss from the aquifer (m s-1
- scalarAquiferRecharge => flux_data%var(iLookFLUX%scalarAquiferRecharge)%dat(1) ,& ! intent(out): [dp] recharge to the aquifer (m s-1)
- scalarAquiferBaseflow => flux_data%var(iLookFLUX%scalarAquiferBaseflow)%dat(1) ,& ! intent(out): [dp] total baseflow from the aquifer (m s-1)
-
- ! total runoff
- scalarTotalRunoff => flux_data%var(iLookFLUX%scalarTotalRunoff)%dat(1) ,& ! intent(out): [dp] total runoff (m s-1)
-
- ! derivatives in net vegetation energy fluxes w.r.t. relevant state variables
- dCanairNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. canopy air temperature
- dCanairNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. canopy temperature
- dCanairNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. ground temperature
- dCanopyNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. canopy air temperature
- dCanopyNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. canopy temperature
- dCanopyNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. ground temperature
- dCanopyNetFlux_dCanWat => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in net canopy fluxes w.r.t. canopy liquid water content
- dGroundNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. canopy air temperature
- dGroundNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. canopy temperature
- dGroundNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. ground temperature
- dGroundNetFlux_dCanWat => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in net ground fluxes w.r.t. canopy liquid water content
-
- ! derivatives in evaporative fluxes w.r.t. relevant state variables
- dCanopyEvaporation_dTCanair => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy air temperature
- dCanopyEvaporation_dTCanopy => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy temperature
- dCanopyEvaporation_dTGround => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. ground temperature
- dCanopyEvaporation_dCanWat => deriv_data%var(iLookDERIV%dCanopyEvaporation_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy liquid water content
- dGroundEvaporation_dTCanair => deriv_data%var(iLookDERIV%dGroundEvaporation_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy air temperature
- dGroundEvaporation_dTCanopy => deriv_data%var(iLookDERIV%dGroundEvaporation_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy temperature
- dGroundEvaporation_dTGround => deriv_data%var(iLookDERIV%dGroundEvaporation_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. ground temperature
- dGroundEvaporation_dCanWat => deriv_data%var(iLookDERIV%dGroundEvaporation_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy liquid water content
-
- ! derivatives in transpiration
- dCanopyTrans_dTCanair => deriv_data%var(iLookDERIV%dCanopyTrans_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTCanopy => deriv_data%var(iLookDERIV%dCanopyTrans_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTGround => deriv_data%var(iLookDERIV%dCanopyTrans_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dCanWat => deriv_data%var(iLookDERIV%dCanopyTrans_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy total water content (s-1)
-
- ! derivatives in canopy water w.r.t canopy temperature
- dCanLiq_dTcanopy => deriv_data%var(iLookDERIV%dCanLiq_dTcanopy )%dat(1) ,& ! intent(out): [dp] derivative of canopy liquid storage w.r.t. temperature
-
- ! derivatives in canopy liquid fluxes w.r.t. canopy water
- scalarCanopyLiqDeriv => deriv_data%var(iLookDERIV%scalarCanopyLiqDeriv )%dat(1) ,& ! intent(out): [dp] derivative in (throughfall + drainage) w.r.t. canopy liquid water
- scalarThroughfallRainDeriv => deriv_data%var(iLookDERIV%scalarThroughfallRainDeriv )%dat(1) ,& ! intent(out): [dp] derivative in throughfall w.r.t. canopy liquid water
- scalarCanopyLiqDrainageDeriv => deriv_data%var(iLookDERIV%scalarCanopyLiqDrainageDeriv)%dat(1) ,& ! intent(out): [dp] derivative in canopy drainage w.r.t. canopy liquid water
-
- ! derivatives in energy fluxes at the interface of snow+soil layers w.r.t. temperature in layers above and below
- dNrgFlux_dTempAbove => deriv_data%var(iLookDERIV%dNrgFlux_dTempAbove )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. temperature in the layer above
- dNrgFlux_dTempBelow => deriv_data%var(iLookDERIV%dNrgFlux_dTempBelow )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. temperature in the layer below
-
- ! derivatives in energy fluxes at the interface of snow+soil layers w.r.t. water state in layers above and below
- dNrgFlux_dWatAbove => deriv_data%var(iLookDERIV%dNrgFlux_dWatAbove )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. water state in the layer above
- dNrgFlux_dWatBelow => deriv_data%var(iLookDERIV%dNrgFlux_dWatBelow )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. water state in the layer below
-
- ! derivative in liquid water fluxes at the interface of snow layers w.r.t. volumetric liquid water content in the layer above
- iLayerLiqFluxSnowDeriv => deriv_data%var(iLookDERIV%iLayerLiqFluxSnowDeriv )%dat ,& ! intent(out): [dp(:)] derivative in vertical liquid water flux at layer interfaces
-
- ! derivative in liquid water fluxes for the soil domain w.r.t hydrology state variables
- dVolTot_dPsi0 => deriv_data%var(iLookDERIV%dVolTot_dPsi0 )%dat ,& ! intent(out): [dp(:)] derivative in total water content w.r.t. total water matric potential
- dq_dHydStateAbove => deriv_data%var(iLookDERIV%dq_dHydStateAbove )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer above
- dq_dHydStateBelow => deriv_data%var(iLookDERIV%dq_dHydStateBelow )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer below
- dq_dHydStateLayerSurfVec => deriv_data%var(iLookDERIV%dq_dHydStateLayerSurfVec )%dat ,& ! intent(out): [dp(:)] change in the flux in soil surface interface w.r.t. state variables in layers
- mLayerdTheta_dPsi => deriv_data%var(iLookDERIV%mLayerdTheta_dPsi )%dat ,& ! intent(out): [dp(:)] derivative in the soil water characteristic w.r.t. psi
- mLayerdPsi_dTheta => deriv_data%var(iLookDERIV%mLayerdPsi_dTheta )%dat ,& ! intent(out): [dp(:)] derivative in the soil water characteristic w.r.t. theta
- dCompress_dPsi => deriv_data%var(iLookDERIV%dCompress_dPsi )%dat ,& ! intent(out): [dp(:)] derivative in compressibility w.r.t matric head
-
- ! derivative in baseflow flux w.r.t. aquifer storage
- dBaseflow_dAquifer => deriv_data%var(iLookDERIV%dBaseflow_dAquifer )%dat(1) ,& ! intent(out): [dp(:)] erivative in baseflow flux w.r.t. aquifer storage (s-1)
-
- ! derivative in liquid water fluxes for the soil domain w.r.t energy state variables
- dq_dNrgStateAbove => deriv_data%var(iLookDERIV%dq_dNrgStateAbove )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer above
- dq_dNrgStateBelow => deriv_data%var(iLookDERIV%dq_dNrgStateBelow )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer below
- dq_dNrgStateLayerSurfVec => deriv_data%var(iLookDERIV%dq_dNrgStateLayerSurfVec )%dat ,& ! intent(out): [dp(:)] change in the flux in soil surface interface w.r.t. state variables in layers
-
- ! derivatives in soil transpiration w.r.t. canopy state variables
- mLayerdTrans_dTCanair => deriv_data%var(iLookDERIV%mLayerdTrans_dTCanair )%dat ,& !intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy air temperature
- mLayerdTrans_dTCanopy => deriv_data%var(iLookDERIV%mLayerdTrans_dTCanopy )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy temperature
- mLayerdTrans_dTGround => deriv_data%var(iLookDERIV%mLayerdTrans_dTGround )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. ground temperature
- mLayerdTrans_dCanWat => deriv_data%var(iLookDERIV%mLayerdTrans_dCanWat )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy total water
-
- ! derivatives in aquifer transpiration w.r.t. canopy state variables
- dAquiferTrans_dTCanair => deriv_data%var(iLookDERIV%dAquiferTrans_dTCanair )%dat(1) ,& !intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy air temperature
- dAquiferTrans_dTCanopy => deriv_data%var(iLookDERIV%dAquiferTrans_dTCanopy )%dat(1) ,& ! intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy temperature
- dAquiferTrans_dTGround => deriv_data%var(iLookDERIV%dAquiferTrans_dTGround )%dat(1) ,& ! intent(out): derivatives in the aquifer transpiration flux w.r.t. ground temperature
- dAquiferTrans_dCanWat => deriv_data%var(iLookDERIV%dAquiferTrans_dCanWat )%dat(1) & ! intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy total water
-
- ) ! association to data in structures
-
- ! *****
- ! * PRELIMINARIES...
- ! ******************
-
- !print*, '***** nSnowSoilNrg, nSnowOnlyNrg, nSoilOnlyNrg, nSnowSoilHyd, nSnowOnlyHyd, nSoilOnlyHyd = ', &
- ! nSnowSoilNrg, nSnowOnlyNrg, nSoilOnlyNrg, nSnowSoilHyd, nSnowOnlyHyd, nSoilOnlyHyd
-
- ! increment the number of flux calls
- numFluxCalls = numFluxCalls+1
-
- ! modify the groundwater representation for this single-column implementation
- select case(ixSpatialGroundwater)
- case(singleBasin); local_ixGroundwater = noExplicit ! force no explicit representation of groundwater at the local scale
- case(localColumn); local_ixGroundwater = ixGroundwater ! go with the specified decision
- case default; err=20; message=trim(message)//'unable to identify spatial representation of groundwater'; return
- end select ! (modify the groundwater representation for this single-column implementation)
-
- ! initialize liquid water fluxes throughout the snow and soil domains
- ! NOTE: used in the energy routines, which is called before the hydrology routines
- if(firstFluxCall)then
- if(nSnow > 0) iLayerLiqFluxSnow(0:nSnow) = 0._dp
- iLayerLiqFluxSoil(0:nSoil) = 0._dp
- end if
+ ! model decisions
+ ixGroundwater => model_decisions(iLookDECISIONS%groundwatr)%iDecision ,& ! intent(in): [i4b] groundwater parameterization
+ ixSpatialGroundwater => model_decisions(iLookDECISIONS%spatial_gw)%iDecision ,& ! intent(in): [i4b] spatial representation of groundwater (local-column or single-basin)
+
+ ! domain boundary conditions
+ upperBoundTemp => forc_data%var(iLookFORCE%airtemp) ,& ! intent(in): [dp] temperature of the upper boundary of the snow and soil domains (K)
+ scalarRainfall => flux_data%var(iLookFLUX%scalarRainfall)%dat(1) ,& ! intent(in): [dp] rainfall rate (kg m-2 s-1)
+
+ ! canopy and layer depth
+ canopyDepth => diag_data%var(iLookDIAG%scalarCanopyDepth)%dat(1) ,& ! intent(in): [dp ] canopy depth (m)
+ mLayerDepth => prog_data%var(iLookPROG%mLayerDepth)%dat ,& ! intent(in): [dp(:)] depth of each layer in the snow-soil sub-domain (m)
+
+ ! indices of model state variables for the vegetation subdomain
+ ixCasNrg => indx_data%var(iLookINDEX%ixCasNrg)%dat(1) ,& ! intent(in): [i4b] index of canopy air space energy state variable
+ ixVegNrg => indx_data%var(iLookINDEX%ixVegNrg)%dat(1) ,& ! intent(in): [i4b] index of canopy energy state variable
+ ixVegHyd => indx_data%var(iLookINDEX%ixVegHyd)%dat(1) ,& ! intent(in): [i4b] index of canopy hydrology state variable (mass)
+ ixTopNrg => indx_data%var(iLookINDEX%ixTopNrg)%dat(1) ,& ! intent(in): [i4b] index of upper-most energy state in the snow+soil subdomain
+ ixTopHyd => indx_data%var(iLookINDEX%ixTopHyd)%dat(1) ,& ! intent(in): [i4b] index of upper-most hydrology state in the snow+soil subdomain
+ ixAqWat => indx_data%var(iLookINDEX%ixAqWat)%dat(1) ,& ! intent(in): [i4b] index of water storage in the aquifer
+
+ ! indices of model state variables for the snow+soil domain
+ ixSnowSoilNrg => indx_data%var(iLookINDEX%ixSnowSoilNrg)%dat ,& ! intent(in): [i4b(:)] indices for energy states in the snow+soil subdomain
+ ixSnowSoilHyd => indx_data%var(iLookINDEX%ixSnowSoilHyd)%dat ,& ! intent(in): [i4b(:)] indices for hydrology states in the snow+soil subdomain
+ layerType => indx_data%var(iLookINDEX%layerType)%dat ,& ! intent(in): [i4b(:)] type of layer (iname_soil or iname_snow)
+
+ ! number of state variables of a specific type
+ nSnowSoilNrg => indx_data%var(iLookINDEX%nSnowSoilNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the snow+soil domain
+ nSnowOnlyNrg => indx_data%var(iLookINDEX%nSnowOnlyNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the snow domain
+ nSoilOnlyNrg => indx_data%var(iLookINDEX%nSoilOnlyNrg )%dat(1) ,& ! intent(in): [i4b] number of energy state variables in the soil domain
+ nSnowSoilHyd => indx_data%var(iLookINDEX%nSnowSoilHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the snow+soil domain
+ nSnowOnlyHyd => indx_data%var(iLookINDEX%nSnowOnlyHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the snow domain
+ nSoilOnlyHyd => indx_data%var(iLookINDEX%nSoilOnlyHyd )%dat(1) ,& ! intent(in): [i4b] number of hydrology variables in the soil domain
+
+ ! snow parameters
+ snowfrz_scale => mpar_data%var(iLookPARAM%snowfrz_scale)%dat(1) ,& ! intent(in): [dp] scaling parameter for the snow freezing curve (K-1)
+
+ ! derivatives
+ dPsiLiq_dPsi0 => deriv_data%var(iLookDERIV%dPsiLiq_dPsi0 )%dat ,& ! intent(in): [dp(:)] derivative in liquid water matric pot w.r.t. the total water matric pot (-)
+ dPsiLiq_dTemp => deriv_data%var(iLookDERIV%dPsiLiq_dTemp )%dat ,& ! intent(in): [dp(:)] derivative in the liquid water matric potential w.r.t. temperature
+ mLayerdTheta_dTk => deriv_data%var(iLookDERIV%mLayerdTheta_dTk)%dat ,& ! intent(in): [dp(:)] derivative of volumetric liquid water content w.r.t. temperature
+ dTheta_dTkCanopy => deriv_data%var(iLookDERIV%dTheta_dTkCanopy)%dat(1) ,& ! intent(in): [dp] derivative of volumetric liquid water content w.r.t. temperature
+
+ ! number of flux calls
+ numFluxCalls => diag_data%var(iLookDIAG%numFluxCalls)%dat(1) ,& ! intent(out): [dp] number of flux calls (-)
+
+ ! net fluxes over the vegetation domain
+ scalarCanairNetNrgFlux => flux_data%var(iLookFLUX%scalarCanairNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the canopy air space (W m-2)
+ scalarCanopyNetNrgFlux => flux_data%var(iLookFLUX%scalarCanopyNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the vegetation canopy (W m-2)
+ scalarGroundNetNrgFlux => flux_data%var(iLookFLUX%scalarGroundNetNrgFlux)%dat(1) ,& ! intent(out): [dp] net energy flux for the ground surface (W m-2)
+ scalarCanopyNetLiqFlux => flux_data%var(iLookFLUX%scalarCanopyNetLiqFlux)%dat(1) ,& ! intent(out): [dp] net liquid water flux for the vegetation canopy (kg m-2 s-1)
+
+ ! net fluxes over the snow+soil domain
+ mLayerNrgFlux => flux_data%var(iLookFLUX%mLayerNrgFlux)%dat ,& ! intent(out): [dp] net energy flux for each layer within the snow+soil domain (J m-3 s-1)
+ mLayerLiqFluxSnow => flux_data%var(iLookFLUX%mLayerLiqFluxSnow)%dat ,& ! intent(out): [dp] net liquid water flux for each snow layer (s-1)
+ mLayerLiqFluxSoil => flux_data%var(iLookFLUX%mLayerLiqFluxSoil)%dat ,& ! intent(out): [dp] net liquid water flux for each soil layer (s-1)
+
+ ! evaporative fluxes
+ scalarCanopyTranspiration => flux_data%var(iLookFLUX%scalarCanopyTranspiration)%dat(1) ,& ! intent(out): [dp] canopy transpiration (kg m-2 s-1)
+ scalarCanopyEvaporation => flux_data%var(iLookFLUX%scalarCanopyEvaporation)%dat(1) ,& ! intent(out): [dp] canopy evaporation/condensation (kg m-2 s-1)
+ scalarGroundEvaporation => flux_data%var(iLookFLUX%scalarGroundEvaporation)%dat(1) ,& ! intent(out): [dp] ground evaporation/condensation -- below canopy or non-vegetated (kg m-2 s-1)
+ mLayerTranspire => flux_data%var(iLookFLUX%mLayerTranspire)%dat ,& ! intent(out): [dp(:)] transpiration loss from each soil layer (m s-1)
+
+ ! fluxes for the snow+soil domain
+ iLayerNrgFlux => flux_data%var(iLookFLUX%iLayerNrgFlux)%dat ,& ! intent(out): [dp(0:)] vertical energy flux at the interface of snow and soil layers
+ iLayerLiqFluxSnow => flux_data%var(iLookFLUX%iLayerLiqFluxSnow)%dat ,& ! intent(out): [dp(0:)] vertical liquid water flux at snow layer interfaces (-)
+ iLayerLiqFluxSoil => flux_data%var(iLookFLUX%iLayerLiqFluxSoil)%dat ,& ! intent(out): [dp(0:)] vertical liquid water flux at soil layer interfaces (-)
+ mLayerHydCond => flux_data%var(iLookFLUX%mLayerHydCond)%dat ,& ! intent(out): [dp(:)] hydraulic conductivity in each soil layer (m s-1)
+ mLayerBaseflow => flux_data%var(iLookFLUX%mLayerBaseflow)%dat ,& ! intent(out): [dp(:)] baseflow from each soil layer (m s-1)
+ scalarSnowDrainage => flux_data%var(iLookFLUX%scalarSnowDrainage)%dat(1) ,& ! intent(out): [dp] drainage from the snow profile (m s-1)
+ scalarSoilDrainage => flux_data%var(iLookFLUX%scalarSoilDrainage)%dat(1) ,& ! intent(out): [dp] drainage from the soil profile (m s-1)
+ scalarSoilBaseflow => flux_data%var(iLookFLUX%scalarSoilBaseflow)%dat(1) ,& ! intent(out): [dp] total baseflow from the soil profile (m s-1)
+
+ ! infiltration
+ scalarInfilArea => diag_data%var(iLookDIAG%scalarInfilArea )%dat(1) ,& ! intent(out): [dp] fraction of unfrozen area where water can infiltrate (-)
+ scalarFrozenArea => diag_data%var(iLookDIAG%scalarFrozenArea )%dat(1) ,& ! intent(out): [dp] fraction of area that is considered impermeable due to soil ice (-)
+ scalarSoilControl => diag_data%var(iLookDIAG%scalarSoilControl )%dat(1) ,& ! intent(out): [dp] soil control on infiltration, zero or one
+ scalarMaxInfilRate => flux_data%var(iLookFLUX%scalarMaxInfilRate)%dat(1) ,& ! intent(out): [dp] maximum infiltration rate (m s-1)
+ scalarInfiltration => flux_data%var(iLookFLUX%scalarInfiltration)%dat(1) ,& ! intent(out): [dp] infiltration of water into the soil profile (m s-1)
+ scalarFracLiqVeg => diag_data%var(iLookDIAG%scalarFracLiqVeg)%dat(1) ,& ! intent(inout): [dp] fraction of liquid water on vegetation (-)
+ mLayerFracLiqSnow => diag_data%var(iLookDIAG%mLayerFracLiqSnow)%dat ,& ! intent(inout): [dp(:)] fraction of liquid water in each snow layer (-)
+
+ ! boundary fluxes in the soil domain
+ scalarThroughfallRain => flux_data%var(iLookFLUX%scalarThroughfallRain)%dat(1) ,& ! intent(out): [dp] rain that reaches the ground without ever touching the canopy (kg m-2 s-1)
+ scalarCanopyLiqDrainage => flux_data%var(iLookFLUX%scalarCanopyLiqDrainage)%dat(1) ,& ! intent(out): [dp] drainage of liquid water from the vegetation canopy (kg m-2 s-1)
+ scalarRainPlusMelt => flux_data%var(iLookFLUX%scalarRainPlusMelt)%dat(1) ,& ! intent(out): [dp] rain plus melt (m s-1)
+ scalarSurfaceRunoff => flux_data%var(iLookFLUX%scalarSurfaceRunoff)%dat(1) ,& ! intent(out): [dp] surface runoff (m s-1)
+ scalarExfiltration => flux_data%var(iLookFLUX%scalarExfiltration)%dat(1) ,& ! intent(out): [dp] exfiltration from the soil profile (m s-1)
+ mLayerColumnOutflow => flux_data%var(iLookFLUX%mLayerColumnOutflow)%dat ,& ! intent(out): [dp(:)] column outflow from each soil layer (m3 s-1)
+
+ ! fluxes for the aquifer
+ scalarAquiferTranspire => flux_data%var(iLookFLUX%scalarAquiferTranspire)%dat(1) ,& ! intent(out): [dp] transpiration loss from the aquifer (m s-1
+ scalarAquiferRecharge => flux_data%var(iLookFLUX%scalarAquiferRecharge)%dat(1) ,& ! intent(out): [dp] recharge to the aquifer (m s-1)
+ scalarAquiferBaseflow => flux_data%var(iLookFLUX%scalarAquiferBaseflow)%dat(1) ,& ! intent(out): [dp] total baseflow from the aquifer (m s-1)
+
+ ! total runoff
+ scalarTotalRunoff => flux_data%var(iLookFLUX%scalarTotalRunoff)%dat(1) ,& ! intent(out): [dp] total runoff (m s-1)
+
+ ! derivatives in net vegetation energy fluxes w.r.t. relevant state variables
+ dCanairNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. canopy air temperature
+ dCanairNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. canopy temperature
+ dCanairNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dCanairNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy air space flux w.r.t. ground temperature
+ dCanopyNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. canopy air temperature
+ dCanopyNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. canopy temperature
+ dCanopyNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dCanopyNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net canopy flux w.r.t. ground temperature
+ dCanopyNetFlux_dCanWat => deriv_data%var(iLookDERIV%dCanopyNetFlux_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in net canopy fluxes w.r.t. canopy total water content
+ dGroundNetFlux_dCanairTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanairTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. canopy air temperature
+ dGroundNetFlux_dCanopyTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanopyTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. canopy temperature
+ dGroundNetFlux_dGroundTemp => deriv_data%var(iLookDERIV%dGroundNetFlux_dGroundTemp )%dat(1) ,& ! intent(out): [dp] derivative in net ground flux w.r.t. ground temperature
+ dGroundNetFlux_dCanWat => deriv_data%var(iLookDERIV%dGroundNetFlux_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in net ground fluxes w.r.t. canopy total water content
+
+ ! derivatives in evaporative fluxes w.r.t. relevant state variables
+ dCanopyEvaporation_dTCanair => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy air temperature
+ dCanopyEvaporation_dTCanopy => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy temperature
+ dCanopyEvaporation_dTGround => deriv_data%var(iLookDERIV%dCanopyEvaporation_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. ground temperature
+ dCanopyEvaporation_dCanWat => deriv_data%var(iLookDERIV%dCanopyEvaporation_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in canopy evaporation w.r.t. canopy total water content
+ dGroundEvaporation_dTCanair => deriv_data%var(iLookDERIV%dGroundEvaporation_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy air temperature
+ dGroundEvaporation_dTCanopy => deriv_data%var(iLookDERIV%dGroundEvaporation_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy temperature
+ dGroundEvaporation_dTGround => deriv_data%var(iLookDERIV%dGroundEvaporation_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. ground temperature
+ dGroundEvaporation_dCanWat => deriv_data%var(iLookDERIV%dGroundEvaporation_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in ground evaporation w.r.t. canopy total water content
+
+ ! derivatives in transpiration
+ dCanopyTrans_dTCanair => deriv_data%var(iLookDERIV%dCanopyTrans_dTCanair )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTCanopy => deriv_data%var(iLookDERIV%dCanopyTrans_dTCanopy )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTGround => deriv_data%var(iLookDERIV%dCanopyTrans_dTGround )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dCanWat => deriv_data%var(iLookDERIV%dCanopyTrans_dCanWat )%dat(1) ,& ! intent(out): [dp] derivative in canopy transpiration w.r.t. canopy total water content (s-1)
+
+ ! derivatives in canopy water w.r.t canopy temperature
+ dCanLiq_dTcanopy => deriv_data%var(iLookDERIV%dCanLiq_dTcanopy )%dat(1) ,& ! intent(out): [dp] derivative of canopy liquid storage w.r.t. temperature
+
+ ! derivatives in canopy liquid fluxes w.r.t. canopy water
+ scalarCanopyLiqDeriv => deriv_data%var(iLookDERIV%scalarCanopyLiqDeriv )%dat(1) ,& ! intent(out): [dp] derivative in (throughfall + drainage) w.r.t. canopy liquid water
+ scalarThroughfallRainDeriv => deriv_data%var(iLookDERIV%scalarThroughfallRainDeriv )%dat(1) ,& ! intent(out): [dp] derivative in throughfall w.r.t. canopy liquid water
+ scalarCanopyLiqDrainageDeriv => deriv_data%var(iLookDERIV%scalarCanopyLiqDrainageDeriv)%dat(1) ,& ! intent(out): [dp] derivative in canopy drainage w.r.t. canopy liquid water
+
+ ! derivatives in energy fluxes at the interface of snow+soil layers w.r.t. temperature in layers above and below
+ dNrgFlux_dTempAbove => deriv_data%var(iLookDERIV%dNrgFlux_dTempAbove )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. temperature in the layer above
+ dNrgFlux_dTempBelow => deriv_data%var(iLookDERIV%dNrgFlux_dTempBelow )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. temperature in the layer below
+
+ ! derivatives in energy fluxes at the interface of snow+soil layers w.r.t. water state in layers above and below
+ dNrgFlux_dWatAbove => deriv_data%var(iLookDERIV%dNrgFlux_dWatAbove )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. water state in the layer above
+ dNrgFlux_dWatBelow => deriv_data%var(iLookDERIV%dNrgFlux_dWatBelow )%dat ,& ! intent(out): [dp(:)] derivatives in the flux w.r.t. water state in the layer below
+
+ ! derivative in liquid water fluxes at the interface of snow layers w.r.t. volumetric liquid water content in the layer above
+ iLayerLiqFluxSnowDeriv => deriv_data%var(iLookDERIV%iLayerLiqFluxSnowDeriv )%dat ,& ! intent(out): [dp(:)] derivative in vertical liquid water flux at layer interfaces
+
+ ! derivative in liquid water fluxes for the soil domain w.r.t hydrology state variables
+ dVolTot_dPsi0 => deriv_data%var(iLookDERIV%dVolTot_dPsi0 )%dat ,& ! intent(out): [dp(:)] derivative in total water content w.r.t. total water matric potential
+ dq_dHydStateAbove => deriv_data%var(iLookDERIV%dq_dHydStateAbove )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer above
+ dq_dHydStateBelow => deriv_data%var(iLookDERIV%dq_dHydStateBelow )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer below
+ dq_dHydStateLayerSurfVec => deriv_data%var(iLookDERIV%dq_dHydStateLayerSurfVec )%dat ,& ! intent(out): [dp(:)] change in the flux in soil surface interface w.r.t. state variables in layers
+ mLayerdTheta_dPsi => deriv_data%var(iLookDERIV%mLayerdTheta_dPsi )%dat ,& ! intent(out): [dp(:)] derivative in the soil water characteristic w.r.t. psi
+ mLayerdPsi_dTheta => deriv_data%var(iLookDERIV%mLayerdPsi_dTheta )%dat ,& ! intent(out): [dp(:)] derivative in the soil water characteristic w.r.t. theta
+ dCompress_dPsi => deriv_data%var(iLookDERIV%dCompress_dPsi )%dat ,& ! intent(out): [dp(:)] derivative in compressibility w.r.t matric head
+
+ ! derivative in baseflow flux w.r.t. aquifer storage
+ dBaseflow_dAquifer => deriv_data%var(iLookDERIV%dBaseflow_dAquifer )%dat(1) ,& ! intent(out): [dp(:)] erivative in baseflow flux w.r.t. aquifer storage (s-1)
+
+ ! derivative in liquid water fluxes for the soil domain w.r.t energy state variables
+ dq_dNrgStateAbove => deriv_data%var(iLookDERIV%dq_dNrgStateAbove )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer above
+ dq_dNrgStateBelow => deriv_data%var(iLookDERIV%dq_dNrgStateBelow )%dat ,& ! intent(out): [dp(:)] change in flux at layer interfaces w.r.t. states in the layer below
+ dq_dNrgStateLayerSurfVec => deriv_data%var(iLookDERIV%dq_dNrgStateLayerSurfVec )%dat ,& ! intent(out): [dp(:)] change in the flux in soil surface interface w.r.t. state variables in layers
+
+ ! derivatives in soil transpiration w.r.t. canopy state variables
+ mLayerdTrans_dTCanair => deriv_data%var(iLookDERIV%mLayerdTrans_dTCanair )%dat ,& !intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy air temperature
+ mLayerdTrans_dTCanopy => deriv_data%var(iLookDERIV%mLayerdTrans_dTCanopy )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy temperature
+ mLayerdTrans_dTGround => deriv_data%var(iLookDERIV%mLayerdTrans_dTGround )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. ground temperature
+ mLayerdTrans_dCanWat => deriv_data%var(iLookDERIV%mLayerdTrans_dCanWat )%dat ,& ! intent(out): derivatives in the soil layer transpiration flux w.r.t. canopy total water
+
+ ! derivatives in aquifer transpiration w.r.t. canopy state variables
+ dAquiferTrans_dTCanair => deriv_data%var(iLookDERIV%dAquiferTrans_dTCanair )%dat(1) ,& !intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy air temperature
+ dAquiferTrans_dTCanopy => deriv_data%var(iLookDERIV%dAquiferTrans_dTCanopy )%dat(1) ,& ! intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy temperature
+ dAquiferTrans_dTGround => deriv_data%var(iLookDERIV%dAquiferTrans_dTGround )%dat(1) ,& ! intent(out): derivatives in the aquifer transpiration flux w.r.t. ground temperature
+ dAquiferTrans_dCanWat => deriv_data%var(iLookDERIV%dAquiferTrans_dCanWat )%dat(1) & ! intent(out): derivatives in the aquifer transpiration flux w.r.t. canopy total water
+
+ ) ! association to data in structures
+
+ ! *****
+ ! * PRELIMINARIES...
+ ! ******************
+
+ ! increment the number of flux calls
+ numFluxCalls = numFluxCalls+1
+
+ ! modify the groundwater representation for this single-column implementation
+ select case(ixSpatialGroundwater)
+ case(singleBasin); local_ixGroundwater = noExplicit ! force no explicit representation of groundwater at the local scale
+ case(localColumn); local_ixGroundwater = ixGroundwater ! go with the specified decision
+ case default; err=20; message=trim(message)//'unable to identify spatial representation of groundwater'; return
+ end select ! (modify the groundwater representation for this single-column implementation)
+
+ ! initialize liquid water fluxes throughout the snow and soil domains
+ ! NOTE: used in the energy routines, which is called before the hydrology routines
+ if(firstFluxCall)then
+ if(nSnow > 0) iLayerLiqFluxSnow(0:nSnow) = 0._rkind
+ end if
+
+ ! *****
+ ! * CALCULATE ENERGY FLUXES OVER VEGETATION...
+ ! *********************************************
+
+ ! identify the need to calculate the energy flux over vegetation
+ doVegNrgFlux = (ixCasNrg/=integerMissing .or. ixVegNrg/=integerMissing .or. ixTopNrg/=integerMissing)
+
+ ! check if there is a need to calculate the energy fluxes over vegetation
+ if(doVegNrgFlux)then
+
+ ! derivative in canopy liquid storage w.r.t. canopy temperature
+ dCanLiq_dTcanopy = dTheta_dTkCanopy*iden_water*canopyDepth ! kg m-2 K-1
+
+ ! calculate the energy fluxes over vegetation
+ call vegNrgFlux(&
+ ! input: model control
+ firstSubStep, & ! intent(in): flag to indicate if we are processing the first sub-step
+ firstFluxCall, & ! intent(in): flag to indicate if we are processing the first flux call
+ computeVegFlux, & ! intent(in): flag to indicate if we need to compute fluxes over vegetation
+ insideIDA, & ! intent(in): flag to indicate inside Sundials Solver (do not require longwave to be balanced)
+ ! input: model state variables
+ upperBoundTemp, & ! intent(in): temperature of the upper boundary (K) --> NOTE: use air temperature
+ scalarCanairTempTrial, & ! intent(in): trial value of the canopy air space temperature (K)
+ scalarCanopyTempTrial, & ! intent(in): trial value of canopy temperature (K)
+ mLayerTempTrial(1), & ! intent(in): trial value of ground temperature (K)
+ scalarCanopyIceTrial, & ! intent(in): trial value of mass of ice on the vegetation canopy (kg m-2)
+ scalarCanopyLiqTrial, & ! intent(in): trial value of mass of liquid water on the vegetation canopy (kg m-2)
+ ! input: model derivatives
+ dCanLiq_dTcanopy, & ! intent(in): derivative in canopy liquid storage w.r.t. canopy temperature (kg m-2 K-1)
+ ! input/output: data structures
+ type_data, & ! intent(in): type of vegetation and soil
+ forc_data, & ! intent(in): model forcing data
+ mpar_data, & ! intent(in): model parameters
+ indx_data, & ! intent(in): index data
+ prog_data, & ! intent(in): model prognostic variables for a local HRU
+ diag_data, & ! intent(inout): model diagnostic variables for a local HRU
+ flux_data, & ! intent(inout): model fluxes for a local HRU
+ bvar_data, & ! intent(in): model variables for the local basin
+ model_decisions, & ! intent(in): model decisions
+ ! output: liquid water fluxes associated with evaporation/transpiration
+ scalarCanopyTranspiration, & ! intent(out): canopy transpiration (kg m-2 s-1)
+ scalarCanopyEvaporation, & ! intent(out): canopy evaporation/condensation (kg m-2 s-1)
+ scalarGroundEvaporation, & ! intent(out): ground evaporation/condensation -- below canopy or non-vegetated (kg m-2 s-1)
+ ! output: fluxes
+ scalarCanairNetNrgFlux, & ! intent(out): net energy flux for the canopy air space (W m-2)
+ scalarCanopyNetNrgFlux, & ! intent(out): net energy flux for the vegetation canopy (W m-2)
+ scalarGroundNetNrgFlux, & ! intent(out): net energy flux for the ground surface (W m-2)
+ ! output: flux derivatives
+ dCanairNetFlux_dCanairTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. canopy air temperature (W m-2 K-1)
+ dCanairNetFlux_dCanopyTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. canopy temperature (W m-2 K-1)
+ dCanairNetFlux_dGroundTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. ground temperature (W m-2 K-1)
+ dCanopyNetFlux_dCanairTemp, & ! intent(out): derivative in net canopy flux w.r.t. canopy air temperature (W m-2 K-1)
+ dCanopyNetFlux_dCanopyTemp, & ! intent(out): derivative in net canopy flux w.r.t. canopy temperature (W m-2 K-1)
+ dCanopyNetFlux_dGroundTemp, & ! intent(out): derivative in net canopy flux w.r.t. ground temperature (W m-2 K-1)
+ dGroundNetFlux_dCanairTemp, & ! intent(out): derivative in net ground flux w.r.t. canopy air temperature (W m-2 K-1)
+ dGroundNetFlux_dCanopyTemp, & ! intent(out): derivative in net ground flux w.r.t. canopy temperature (W m-2 K-1)
+ dGroundNetFlux_dGroundTemp, & ! intent(out): derivative in net ground flux w.r.t. ground temperature (W m-2 K-1)
+ ! output: liquid water flux derivatives (canopy evap)
+ dCanopyEvaporation_dCanWat, & ! intent(out): derivative in canopy evaporation w.r.t. canopy total water content (s-1)
+ dCanopyEvaporation_dTCanair, & ! intent(out): derivative in canopy evaporation w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dCanopyEvaporation_dTCanopy, & ! intent(out): derivative in canopy evaporation w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dCanopyEvaporation_dTGround, & ! intent(out): derivative in canopy evaporation w.r.t. ground temperature (kg m-2 s-1 K-1)
+ ! output: liquid water flux derivatives (ground evap)
+ dGroundEvaporation_dCanWat, & ! intent(out): derivative in ground evaporation w.r.t. canopy total water content (s-1)
+ dGroundEvaporation_dTCanair, & ! intent(out): derivative in ground evaporation w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dGroundEvaporation_dTCanopy, & ! intent(out): derivative in ground evaporation w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dGroundEvaporation_dTGround, & ! intent(out): derivative in ground evaporation w.r.t. ground temperature (kg m-2 s-1 K-1)
+ ! output: transpiration derivatives
+ dCanopyTrans_dCanWat, & ! intent(out): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
+ dCanopyTrans_dTCanair, & ! intent(out): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTCanopy, & ! intent(out): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTGround, & ! intent(out): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
+ ! output: cross derivative terms
+ dCanopyNetFlux_dCanWat, & ! intent(out): derivative in net canopy fluxes w.r.t. canopy total water content (J kg-1 s-1)
+ dGroundNetFlux_dCanWat, & ! intent(out): derivative in net ground fluxes w.r.t. canopy total water content (J kg-1 s-1)
+ ! output: error control
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif ! (check for errors)
- ! *****
- ! * CALCULATE ENERGY FLUXES OVER VEGETATION...
- ! *********************************************
-
- ! identify the need to calculate the energy flux over vegetation
- doVegNrgFlux = (ixCasNrg/=integerMissing .or. ixVegNrg/=integerMissing .or. ixTopNrg/=integerMissing)
-
- ! check if there is a need to calculate the energy fluxes over vegetation
- if(doVegNrgFlux)then
-
- ! derivative in canopy liquid storage w.r.t. canopy temperature
- dCanLiq_dTcanopy = dTheta_dTkCanopy*iden_water*canopyDepth ! kg m-2 K-1
- ! calculate the energy fluxes over vegetation
- call vegNrgFlux(&
- ! input: model control
- firstSubStep, & ! intent(in): flag to indicate if we are processing the first sub-step
- firstFluxCall, & ! intent(in): flag to indicate if we are processing the first flux call
- computeVegFlux, & ! intent(in): flag to indicate if we need to compute fluxes over vegetation
- insideIDA, & ! intent(in): flag to indicate if we need longwave to be balanced
- ! input: model state variables
- upperBoundTemp, & ! intent(in): temperature of the upper boundary (K) --> NOTE: use air temperature
- scalarCanairTempTrial, & ! intent(in): trial value of the canopy air space temperature (K)
- scalarCanopyTempTrial, & ! intent(in): trial value of canopy temperature (K)
- mLayerTempTrial(1), & ! intent(in): trial value of ground temperature (K)
- scalarCanopyIceTrial, & ! intent(in): trial value of mass of ice on the vegetation canopy (kg m-2)
- scalarCanopyLiqTrial, & ! intent(in): trial value of mass of liquid water on the vegetation canopy (kg m-2)
- ! input: model derivatives
- dCanLiq_dTcanopy, & ! intent(in): derivative in canopy liquid storage w.r.t. canopy temperature (kg m-2 K-1)
- ! input/output: data structures
- type_data, & ! intent(in): type of vegetation and soil
- forc_data, & ! intent(in): model forcing data
- mpar_data, & ! intent(in): model parameters
- indx_data, & ! intent(in): index data
- prog_data, & ! intent(in): model prognostic variables for a local HRU
- diag_data, & ! intent(inout): model diagnostic variables for a local HRU
- flux_data, & ! intent(inout): model fluxes for a local HRU
- bvar_data, & ! intent(in): model variables for the local basin
- model_decisions, & ! intent(in): model decisions
- ! output: liquid water fluxes associated with evaporation/transpiration
- scalarCanopyTranspiration, & ! intent(out): canopy transpiration (kg m-2 s-1)
- scalarCanopyEvaporation, & ! intent(out): canopy evaporation/condensation (kg m-2 s-1)
- scalarGroundEvaporation, & ! intent(out): ground evaporation/condensation -- below canopy or non-vegetated (kg m-2 s-1)
- ! output: fluxes
- scalarCanairNetNrgFlux, & ! intent(out): net energy flux for the canopy air space (W m-2)
- scalarCanopyNetNrgFlux, & ! intent(out): net energy flux for the vegetation canopy (W m-2)
- scalarGroundNetNrgFlux, & ! intent(out): net energy flux for the ground surface (W m-2)
- ! output: flux derivatives
- dCanairNetFlux_dCanairTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. canopy air temperature (W m-2 K-1)
- dCanairNetFlux_dCanopyTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. canopy temperature (W m-2 K-1)
- dCanairNetFlux_dGroundTemp, & ! intent(out): derivative in net canopy air space flux w.r.t. ground temperature (W m-2 K-1)
- dCanopyNetFlux_dCanairTemp, & ! intent(out): derivative in net canopy flux w.r.t. canopy air temperature (W m-2 K-1)
- dCanopyNetFlux_dCanopyTemp, & ! intent(out): derivative in net canopy flux w.r.t. canopy temperature (W m-2 K-1)
- dCanopyNetFlux_dGroundTemp, & ! intent(out): derivative in net canopy flux w.r.t. ground temperature (W m-2 K-1)
- dGroundNetFlux_dCanairTemp, & ! intent(out): derivative in net ground flux w.r.t. canopy air temperature (W m-2 K-1)
- dGroundNetFlux_dCanopyTemp, & ! intent(out): derivative in net ground flux w.r.t. canopy temperature (W m-2 K-1)
- dGroundNetFlux_dGroundTemp, & ! intent(out): derivative in net ground flux w.r.t. ground temperature (W m-2 K-1)
- ! output: liquid water flux derivatives (canopy evap)
- dCanopyEvaporation_dCanWat, & ! intent(out): derivative in canopy evaporation w.r.t. canopy total water content (s-1)
- dCanopyEvaporation_dTCanair, & ! intent(out): derivative in canopy evaporation w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dCanopyEvaporation_dTCanopy, & ! intent(out): derivative in canopy evaporation w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dCanopyEvaporation_dTGround, & ! intent(out): derivative in canopy evaporation w.r.t. ground temperature (kg m-2 s-1 K-1)
- ! output: liquid water flux derivatives (ground evap)
- dGroundEvaporation_dCanWat, & ! intent(out): derivative in ground evaporation w.r.t. canopy total water content (s-1)
- dGroundEvaporation_dTCanair, & ! intent(out): derivative in ground evaporation w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dGroundEvaporation_dTCanopy, & ! intent(out): derivative in ground evaporation w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dGroundEvaporation_dTGround, & ! intent(out): derivative in ground evaporation w.r.t. ground temperature (kg m-2 s-1 K-1)
- ! output: transpiration derivatives
- dCanopyTrans_dCanWat, & ! intent(out): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
- dCanopyTrans_dTCanair, & ! intent(out): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTCanopy, & ! intent(out): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTGround, & ! intent(out): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
- ! output: cross derivative terms
- dCanopyNetFlux_dCanWat, & ! intent(out): derivative in net canopy fluxes w.r.t. canopy total water content (J kg-1 s-1)
- dGroundNetFlux_dCanWat, & ! intent(out): derivative in net ground fluxes w.r.t. canopy total water content (J kg-1 s-1)
- ! output: error control
- err,cmessage) ! intent(out): error control
-
- ! check fluxes
- if(globalPrintFlag)then
+ ! check fluxes
+ if(globalPrintFlag)then
print*, '**'
write(*,'(a,1x,10(f30.20))') 'canopyDepth = ', canopyDepth
write(*,'(a,1x,10(f30.20))') 'mLayerDepth(1:2) = ', mLayerDepth(1:2)
@@ -512,444 +512,425 @@ subroutine computFlux(&
write(*,'(a,1x,10(f30.20))') 'scalarCanopyNetNrgFlux = ', scalarCanopyNetNrgFlux
write(*,'(a,1x,10(f30.20))') 'scalarGroundNetNrgFlux = ', scalarGroundNetNrgFlux
write(*,'(a,1x,10(f30.20))') 'dGroundNetFlux_dGroundTemp = ', dGroundNetFlux_dGroundTemp
- endif ! if checking fluxes
-
- endif ! if calculating the energy fluxes over vegetation
-
- ! *****
- ! * CALCULATE ENERGY FLUXES THROUGH THE SNOW-SOIL DOMAIN...
- ! **********************************************************
-
- ! check the need to compute energy fluxes throughout the snow+soil domain
- if(nSnowSoilNrg>0)then
-
- ! calculate energy fluxes at layer interfaces through the snow and soil domain
- call ssdNrgFlux(&
- ! input: model control
- (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
- .true., & ! intent(in): flag indicating if derivatives are desired
- ! input: fluxes and derivatives at the upper boundary
- scalarGroundNetNrgFlux, & ! intent(in): total flux at the ground surface (W m-2)
- dGroundNetFlux_dGroundTemp, & ! intent(in): derivative in total ground surface flux w.r.t. ground temperature (W m-2 K-1)
- ! input: liquid water fluxes throughout the snow and soil domains
- iLayerLiqFluxSnow, & ! intent(in): liquid flux at the interface of each snow layer (m s-1)
- iLayerLiqFluxSoil, & ! intent(in): liquid flux at the interface of each soil layer (m s-1)
- ! input: trial value of model state variabes
- mLayerTempTrial, & ! intent(in): trial temperature at the current iteration (K)
- mLayerMatricHeadTrial, & ! intent(in): trial value for the total water matric potential in each soil layer (m)
- mLayerVolFracLiqTrial, & ! intent(in): trial volumetric fraction of liquid water at the current iteration(-)
- mLayerVolFracIceTrial, & ! intent(in): trial volumetric fraction of ice water at the current iteration(-)
- ! input: pre-computed derivatives
- mLayerdTheta_dTk, & ! intent(in): derivative in volumetric liquid water content w.r.t. temperature (K-1)
- mLayerFracLiqSnow, & ! intent(in): fraction of liquid water (-)
- ! input-output: data structures
- mpar_data, & ! intent(in): model parameters
- indx_data, & ! intent(in): model indices
- prog_data, & ! intent(in): model prognostic variables for a local HRU
- diag_data, & ! intent(in): model diagnostic variables for a local HRU
- flux_data, & ! intent(inout): model fluxes for a local HRU
- ! output: fluxes and derivatives at all layer interfaces
- iLayerNrgFlux, & ! intent(out): energy flux at the layer interfaces (W m-2)
- dNrgFlux_dTempAbove, & ! intent(out): derivatives in the flux w.r.t. temperature in the layer above (W m-2 K-1)
- dNrgFlux_dTempBelow, & ! intent(out): derivatives in the flux w.r.t. temperature in the layer below (W m-2 K-1)
- dNrgFlux_dWatAbove, & ! intent(out): derivatives in the flux w.r.t. water state in the layer above
- dNrgFlux_dWatBelow, & ! intent(out): derivatives in the flux w.r.t. water state in the layer below
- ! output: error control
- err,cmessage) ! intent(out): error control
- if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
-
- ! calculate net energy fluxes for each snow and soil layer (J m-3 s-1)
- do iLayer=1,nLayers
- mLayerNrgFlux(iLayer) = -(iLayerNrgFlux(iLayer) - iLayerNrgFlux(iLayer-1))/mLayerDepth(iLayer)
+ endif ! if checking fluxes
+
+ endif ! if calculating the energy fluxes over vegetation
+
+ ! *****
+ ! * CALCULATE ENERGY FLUXES THROUGH THE SNOW-SOIL DOMAIN...
+ ! **********************************************************
+
+ ! check the need to compute energy fluxes throughout the snow+soil domain
+ if(nSnowSoilNrg>0)then
+
+ ! calculate energy fluxes at layer interfaces through the snow and soil domain
+ call ssdNrgFlux(&
+ ! input: model control
+ (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
+ .true., & ! intent(in): flag indicating if derivatives are desired
+ ! input: fluxes and derivatives at the upper boundary
+ scalarGroundNetNrgFlux, & ! intent(in): total flux at the ground surface (W m-2)
+ dGroundNetFlux_dGroundTemp, & ! intent(in): derivative in total ground surface flux w.r.t. ground temperature (W m-2 K-1)
+ ! input: liquid water fluxes throughout the snow and soil domains
+ iLayerLiqFluxSnow, & ! intent(in): liquid flux at the interface of each snow layer (m s-1)
+ iLayerLiqFluxSoil, & ! intent(in): liquid flux at the interface of each soil layer (m s-1)
+ ! input: trial value of model state variables
+ mLayerTempTrial, & ! intent(in): trial temperature at the current iteration (K)
+ mLayerMatricHeadTrial, & ! intent(in): trial value for the total water matric potential in each soil layer (m)
+ mLayerVolFracLiqTrial, & ! intent(in): trial volumetric fraction of liquid water at the current iteration(-)
+ mLayerVolFracIceTrial, & ! intent(in): trial volumetric fraction of ice water at the current iteration(-)
+ ! input: pre-computed derivatives
+ mLayerdTheta_dTk, & ! intent(in): derivative in volumetric liquid water content w.r.t. temperature (K-1)
+ mLayerFracLiqSnow, & ! intent(in): fraction of liquid water (-)
+ ! input-output: data structures
+ mpar_data, & ! intent(in): model parameters
+ indx_data, & ! intent(in): model indices
+ prog_data, & ! intent(in): model prognostic variables for a local HRU
+ diag_data, & ! intent(in): model diagnostic variables for a local HRU
+ flux_data, & ! intent(inout): model fluxes for a local HRU
+ ! output: fluxes and derivatives at all layer interfaces
+ iLayerNrgFlux, & ! intent(out): energy flux at the layer interfaces (W m-2)
+ dNrgFlux_dTempAbove, & ! intent(out): derivatives in the flux w.r.t. temperature in the layer above (W m-2 K-1)
+ dNrgFlux_dTempBelow, & ! intent(out): derivatives in the flux w.r.t. temperature in the layer below (W m-2 K-1)
+ dNrgFlux_dWatAbove, & ! intent(out): derivatives in the flux w.r.t. water state in the layer above
+ dNrgFlux_dWatBelow, & ! intent(out): derivatives in the flux w.r.t. water state in the layer below
+ ! output: error control
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
+
+ ! calculate net energy fluxes for each snow and soil layer (J m-3 s-1)
+ do iLayer=1,nLayers
+ mLayerNrgFlux(iLayer) = -(iLayerNrgFlux(iLayer) - iLayerNrgFlux(iLayer-1))/mLayerDepth(iLayer)
+ if(globalPrintFlag)then
+ if(iLayer < 10) write(*,'(a,1x,i4,1x,10(f25.15,1x))') 'iLayer, iLayerNrgFlux(iLayer-1:iLayer), mLayerNrgFlux(iLayer) = ', iLayer, iLayerNrgFlux(iLayer-1:iLayer), mLayerNrgFlux(iLayer)
+ endif
+ end do
+
+ endif ! if computing energy fluxes throughout the snow+soil domain
+
+ ! *****
+ ! * CALCULATE THE LIQUID FLUX THROUGH VEGETATION...
+ ! **************************************************
+
+ ! check the need to compute the liquid water fluxes through vegetation
+ if(ixVegHyd/=integerMissing)then
+
+ ! calculate liquid water fluxes through vegetation
+ call vegLiqFlux(&
+ ! input
+ computeVegFlux, & ! intent(in): flag to denote if computing energy flux over vegetation
+ scalarCanopyLiqTrial, & ! intent(in): trial mass of liquid water on the vegetation canopy at the current iteration (kg m-2)
+ scalarRainfall, & ! intent(in): rainfall rate (kg m-2 s-1)
+ ! input-output: data structures
+ mpar_data, & ! intent(in): model parameters
+ diag_data, & ! intent(in): local HRU diagnostic model variables
+ ! output
+ scalarThroughfallRain, & ! intent(out): rain that reaches the ground without ever touching the canopy (kg m-2 s-1)
+ scalarCanopyLiqDrainage, & ! intent(out): drainage of liquid water from the vegetation canopy (kg m-2 s-1)
+ scalarThroughfallRainDeriv, & ! intent(out): derivative in throughfall w.r.t. canopy liquid water (s-1)
+ scalarCanopyLiqDrainageDeriv, & ! intent(out): derivative in canopy drainage w.r.t. canopy liquid water (s-1)
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
+
+ ! calculate the net liquid water flux for the vegetation canopy
+ scalarCanopyNetLiqFlux = scalarRainfall + scalarCanopyEvaporation - scalarThroughfallRain - scalarCanopyLiqDrainage
+
+ ! calculate the total derivative in the downward liquid flux
+ scalarCanopyLiqDeriv = scalarThroughfallRainDeriv + scalarCanopyLiqDrainageDeriv
+
+ ! test
if(globalPrintFlag)then
- if(iLayer < 10) write(*,'(a,1x,i4,1x,10(f25.15,1x))') 'iLayer, iLayerNrgFlux(iLayer-1:iLayer), mLayerNrgFlux(iLayer) = ', iLayer, iLayerNrgFlux(iLayer-1:iLayer), mLayerNrgFlux(iLayer)
+ print*, '**'
+ print*, 'scalarRainfall = ', scalarRainfall
+ print*, 'scalarThroughfallRain = ', scalarThroughfallRain
+ print*, 'scalarCanopyEvaporation = ', scalarCanopyEvaporation
+ print*, 'scalarCanopyLiqDrainage = ', scalarCanopyLiqDrainage
+ print*, 'scalarCanopyNetLiqFlux = ', scalarCanopyNetLiqFlux
+ print*, 'scalarCanopyLiqTrial = ', scalarCanopyLiqTrial
endif
- end do
-
- endif ! if computing energy fluxes throughout the snow+soil domain
-
- ! print*, "After ssdNRGFlux call ", flux_data%var(iLookFLUX%iLayerLiqFluxSoil)%dat
- ! *****
- ! * CALCULATE THE LIQUID FLUX THROUGH VEGETATION...
- ! **************************************************
-
- ! check the need to compute the liquid water fluxes through vegetation
- if(ixVegHyd/=integerMissing)then
-
- ! calculate liquid water fluxes through vegetation
- call vegLiqFlux(&
- ! input
- computeVegFlux, & ! intent(in): flag to denote if computing energy flux over vegetation
- scalarCanopyLiqTrial, & ! intent(in): trial mass of liquid water on the vegetation canopy at the current iteration (kg m-2)
- scalarRainfall, & ! intent(in): rainfall rate (kg m-2 s-1)
- ! input-output: data structures
- mpar_data, & ! intent(in): model parameters
- diag_data, & ! intent(in): local HRU diagnostic model variables
- ! output
- scalarThroughfallRain, & ! intent(out): rain that reaches the ground without ever touching the canopy (kg m-2 s-1)
- scalarCanopyLiqDrainage, & ! intent(out): drainage of liquid water from the vegetation canopy (kg m-2 s-1)
- scalarThroughfallRainDeriv, & ! intent(out): derivative in throughfall w.r.t. canopy liquid water (s-1)
- scalarCanopyLiqDrainageDeriv, & ! intent(out): derivative in canopy drainage w.r.t. canopy liquid water (s-1)
- err,cmessage) ! intent(out): error control
- if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
-
- ! calculate the net liquid water flux for the vegetation canopy
- scalarCanopyNetLiqFlux = scalarRainfall + scalarCanopyEvaporation - scalarThroughfallRain - scalarCanopyLiqDrainage
-
- ! calculate the total derivative in the downward liquid flux
- scalarCanopyLiqDeriv = scalarThroughfallRainDeriv + scalarCanopyLiqDrainageDeriv
-
- ! test
- if(globalPrintFlag)then
- print*, '**'
- print*, 'scalarRainfall = ', scalarRainfall
- print*, 'scalarThroughfallRain = ', scalarThroughfallRain
- print*, 'scalarCanopyEvaporation = ', scalarCanopyEvaporation
- print*, 'scalarCanopyLiqDrainage = ', scalarCanopyLiqDrainage
- print*, 'scalarCanopyNetLiqFlux = ', scalarCanopyNetLiqFlux
- print*, 'scalarCanopyLiqTrial = ', scalarCanopyLiqTrial
- endif
+ endif ! computing the liquid water fluxes through vegetation
+
+ ! *****
+ ! * CALCULATE THE LIQUID FLUX THROUGH SNOW...
+ ! ********************************************
+
+ ! check the need to compute liquid water fluxes through snow
+ if(nSnowOnlyHyd>0)then
+
+ ! compute liquid fluxes through snow
+ call snowLiqFlx(&
+ ! input: model control
+ nSnow, & ! intent(in): number of snow layers
+ firstFluxCall, & ! intent(in): the first flux call (compute variables that are constant over the iterations)
+ (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
+ ! input: forcing for the snow domain
+ scalarThroughfallRain, & ! intent(in): rain that reaches the snow surface without ever touching vegetation (kg m-2 s-1)
+ scalarCanopyLiqDrainage, & ! intent(in): liquid drainage from the vegetation canopy (kg m-2 s-1)
+ ! input: model state vector
+ mLayerVolFracLiqTrial(1:nSnow), & ! intent(in): trial value of volumetric fraction of liquid water at the current iteration (-)
+ ! input-output: data structures
+ indx_data, & ! intent(in): model indices
+ mpar_data, & ! intent(in): model parameters
+ prog_data, & ! intent(in): model prognostic variables for a local HRU
+ diag_data, & ! intent(inout): model diagnostic variables for a local HRU
+ ! output: fluxes and derivatives
+ iLayerLiqFluxSnow(0:nSnow), & ! intent(inout): vertical liquid water flux at layer interfaces (m s-1)
+ iLayerLiqFluxSnowDeriv(0:nSnow), & ! intent(inout): derivative in vertical liquid water flux at layer interfaces (m s-1)
+ ! output: error control
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
- endif ! computing the liquid water fluxes through vegetation
+ ! define forcing for the soil domain
+ scalarRainPlusMelt = iLayerLiqFluxSnow(nSnow) ! drainage from the base of the snowpack
- ! *****
- ! * CALCULATE THE LIQUID FLUX THROUGH SNOW...
- ! ********************************************
-
- ! check the need to compute liquid water fluxes through snow
- if(nSnowOnlyHyd>0)then
-
- ! print*, "scalarThroughfallRain = ", scalarThroughfallRain
- ! print*, "scalarCanopyLiqDrainage = ", scalarCanopyLiqDrainage
- ! print*, "mLayerVolFracLiqTrial(1) =", mLayerVolFracLiqTrial(1)
-
- ! compute liquid fluxes through snow
- call snowLiqFlx(&
- ! input: model control
- nSnow, & ! intent(in): number of snow layers
- firstFluxCall, & ! intent(in): the first flux call (compute variables that are constant over the iterations)
- (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
- ! input: forcing for the snow domain
- scalarThroughfallRain, & ! intent(in): rain that reaches the snow surface without ever touching vegetation (kg m-2 s-1)
- scalarCanopyLiqDrainage, & ! intent(in): liquid drainage from the vegetation canopy (kg m-2 s-1)
- ! input: model state vector
- mLayerVolFracLiqTrial(1:nSnow), & ! intent(in): trial value of volumetric fraction of liquid water at the current iteration (-)
- ! input-output: data structures
- indx_data, & ! intent(in): model indices
- mpar_data, & ! intent(in): model parameters
- prog_data, & ! intent(in): model prognostic variables for a local HRU
- diag_data, & ! intent(inout): model diagnostic variables for a local HRU
- ! output: fluxes and derivatives
- iLayerLiqFluxSnow(0:nSnow), & ! intent(inout): vertical liquid water flux at layer interfaces (m s-1)
- iLayerLiqFluxSnowDeriv(0:nSnow), & ! intent(inout): derivative in vertical liquid water flux at layer interfaces (m s-1)
- ! output: error control
- err,cmessage) ! intent(out): error control
- if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
-
- ! define forcing for the soil domain
- scalarRainPlusMelt = iLayerLiqFluxSnow(nSnow) ! drainage from the base of the snowpack
-
- ! calculate net liquid water fluxes for each soil layer (s-1)
- do iLayer=1,nSnow
- mLayerLiqFluxSnow(iLayer) = -(iLayerLiqFluxSnow(iLayer) - iLayerLiqFluxSnow(iLayer-1))/mLayerDepth(iLayer)
- !write(*,'(a,1x,i4,1x,2(f16.10,1x))') 'iLayer, mLayerLiqFluxSnow(iLayer), iLayerLiqFluxSnow(iLayer-1) = ', &
- ! iLayer, mLayerLiqFluxSnow(iLayer), iLayerLiqFluxSnow(iLayer-1)
- end do
+ ! calculate net liquid water fluxes for each snow layer (s-1)
+ do iLayer=1,nSnow
+ mLayerLiqFluxSnow(iLayer) = -(iLayerLiqFluxSnow(iLayer) - iLayerLiqFluxSnow(iLayer-1))/mLayerDepth(iLayer)
+ end do
- ! compute drainage from the soil zone (needed for mass balance checks)
- scalarSnowDrainage = iLayerLiqFluxSnow(nSnow)
- ! print*, "scalarSnowDrainage = ", scalarSnowDrainage
+ ! compute drainage from the soil zone (needed for mass balance checks)
+ scalarSnowDrainage = iLayerLiqFluxSnow(nSnow)
! save bottom layer of snow derivatives
- above_soilLiqFluxDeriv = iLayerLiqFluxSnowDeriv(nSnow) ! derivative in vertical liquid water flux at bottom snow layer interface
- above_soildLiq_dTk = mLayerdTheta_dTk(nSnow) ! derivative in volumetric liquid water content in bottom snow layer w.r.t. temperature
- above_soilFracLiq = mLayerFracLiqSnow(nSnow) ! fraction of liquid water in bottom snow layer (-)
-
- else
-
- ! define forcing for the soil domain for the case of no snow layers
- ! NOTE: in case where nSnowOnlyHyd==0 AND snow layers exist, then scalarRainPlusMelt is taken from the previous flux evaluation
- if(nSnow==0)then !no snow layers
- scalarRainPlusMelt = (scalarThroughfallRain + scalarCanopyLiqDrainage)/iden_water & ! liquid flux from the canopy (m s-1)
- + drainageMeltPond/iden_water ! melt of the snow without a layer (m s-1)
-
- if(ixVegHyd/=integerMissing)then
- ! save canopy derivatives
- above_soilLiqFluxDeriv = scalarCanopyLiqDeriv/iden_water ! derivative in (throughfall + drainage) w.r.t. canopy liquid water
- above_soildLiq_dTk = dCanLiq_dTcanopy ! derivative of canopy liquid storage w.r.t. temperature
- above_soilFracLiq = scalarFracLiqVeg ! fraction of liquid water in canopy (-)
- else
- above_soilLiqFluxDeriv = 0._dp
- above_soildLiq_dTk = 0._dp
- above_soilFracLiq = 0._dp
- endif
- else ! snow layers, take from previous flux calculation
above_soilLiqFluxDeriv = iLayerLiqFluxSnowDeriv(nSnow) ! derivative in vertical liquid water flux at bottom snow layer interface
above_soildLiq_dTk = mLayerdTheta_dTk(nSnow) ! derivative in volumetric liquid water content in bottom snow layer w.r.t. temperature
above_soilFracLiq = mLayerFracLiqSnow(nSnow) ! fraction of liquid water in bottom snow layer (-)
- endif ! snow layers or not
-
- endif
-
- ! *****
- ! * CALCULATE THE LIQUID FLUX THROUGH SOIL...
- ! ********************************************
-
- ! check the need to calculate the liquid flux through soil
- if(nSoilOnlyHyd>0)then
-
- ! calculate the liquid flux through soil
- call soilLiqFlx(&
- ! input: model control
- nSoil, & ! intent(in): number of soil layers
- firstSplitOper, & ! intent(in): flag indicating first flux call in a splitting operation
- (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
- .true., & ! intent(in): flag indicating if derivatives are desired
- ! input: trial state variables
- mLayerTempTrial(nSnow+1:nLayers), & ! intent(in): trial temperature at the current iteration (K)
- mLayerMatricHeadLiqTrial(1:nSoil), & ! intent(in): liquid water matric potential (m)
- mLayerVolFracLiqTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of liquid water (-)
- mLayerVolFracIceTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of ice (-)
- ! input: pre-computed deriavatives
- mLayerdTheta_dTk(nSnow+1:nLayers), & ! intent(in): derivative in volumetric liquid water content w.r.t. temperature (K-1)
- dPsiLiq_dTemp(1:nSoil), & ! intent(in): derivative in liquid water matric potential w.r.t. temperature (m K-1)
- dCanopyTrans_dCanWat, & ! intent(in): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
- dCanopyTrans_dTCanair, & ! intent(in): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTCanopy, & ! intent(in): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTGround, & ! intent(in): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
- above_soilLiqFluxDeriv, & ! intent(in): derivative in layer above soil (canopy or snow) liquid flux w.r.t. liquid water
- above_soildLiq_dTk, & ! intent(in): derivative of layer above soil (canopy or snow) liquid flux w.r.t. temperature
- above_soilFracLiq, & ! intent(in): fraction of liquid water layer above soil (canopy or snow) (-)
- ! input: fluxes
- scalarCanopyTranspiration, & ! intent(in): canopy transpiration (kg m-2 s-1)
- scalarGroundEvaporation, & ! intent(in): ground evaporation (kg m-2 s-1)
- scalarRainPlusMelt, & ! intent(in): rain plus melt (m s-1)
- ! input-output: data structures
- mpar_data, & ! intent(in): model parameters
- indx_data, & ! intent(in): model indices
- prog_data, & ! intent(inout): model prognostic variables for a local HRU
- diag_data, & ! intent(inout): model diagnostic variables for a local HRU
- flux_data, & ! intent(inout): model fluxes for a local HRU
- ! output: diagnostic variables for surface runoff
- scalarMaxInfilRate, & ! intent(inout): maximum infiltration rate (m s-1)
- scalarInfilArea, & ! intent(inout): fraction of unfrozen area where water can infiltrate (-)
- scalarFrozenArea, & ! intent(inout): fraction of area that is considered impermeable due to soil ice (-)
- scalarSurfaceRunoff, & ! intent(inout): surface runoff (m s-1)
- ! output: diagnostic variables for model layers
- mLayerdTheta_dPsi, & ! intent(inout): derivative in the soil water characteristic w.r.t. psi (m-1)
- mLayerdPsi_dTheta, & ! intent(inout): derivative in the soil water characteristic w.r.t. theta (m)
- dHydCond_dMatric, & ! intent(inout): derivative in hydraulic conductivity w.r.t matric head (s-1)
- ! output: fluxes
- scalarInfiltration, & ! intent(inout): surface infiltration rate (m s-1) -- controls on infiltration only computed for iter==1
- iLayerLiqFluxSoil, & ! intent(inout): liquid fluxes at layer interfaces (m s-1)
- mLayerTranspire, & ! intent(inout): transpiration loss from each soil layer (m s-1)
- mLayerHydCond, & ! intent(inout): hydraulic conductivity in each layer (m s-1)
- ! output: derivatives in fluxes w.r.t. state variables -- matric head or volumetric lquid water -- in the layer above and layer below (m s-1 or s-1)
- dq_dHydStateAbove, & ! intent(inout): derivatives in the flux w.r.t. matric head in the layer above (s-1)
- dq_dHydStateBelow, & ! intent(inout): derivatives in the flux w.r.t. matric head in the layer below (s-1)
- dq_dHydStateLayerSurfVec, & ! intent(inout): derivative in surface infiltration w.r.t. hydrology state in above soil snow or canopy and every soil layer (m s-1 or s-1)
- ! output: derivatives in fluxes w.r.t. energy state variables -- now just temperature -- in the layer above and layer below (m s-1 K-1)
- dq_dNrgStateAbove, & ! intent(inout): derivatives in the flux w.r.t. temperature in the layer above (m s-1 K-1)
- dq_dNrgStateBelow, & ! intent(inout): derivatives in the flux w.r.t. temperature in the layer below (m s-1 K-1)
- dq_dNrgStateLayerSurfVec, & ! intent(inout): derivative in surface infiltration w.r.t. energy state in above soil snow or canopy and every soil layer (m s-1 K-1)
- ! output: derivatives in transpiration w.r.t. canopy state variables
- mLayerdTrans_dTCanair, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy air temperature
- mLayerdTrans_dTCanopy, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy temperature
- mLayerdTrans_dTGround, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. ground temperature
- mLayerdTrans_dCanWat, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy total water
- ! output: error control
- err,cmessage) ! intent(out): error control
- if(err/=0)then
- message=trim(message)//trim(cmessage)
- print*, message
- return
- endif
-
- ! calculate net liquid water fluxes for each soil layer (s-1)
- do iLayer=1,nSoil
- mLayerLiqFluxSoil(iLayer) = -(iLayerLiqFluxSoil(iLayer) - iLayerLiqFluxSoil(iLayer-1))/mLayerDepth(iLayer+nSnow)
- !if(iLayer<8) write(*,'(a,1x,2(i4,1x),3(e20.10),f12.7)') 'iLayerLiqFluxSoil(iLayer-1), iLayerLiqFluxSoil(iLayer), mLayerLiqFluxSoil(iLayer) = ', iLayer-1, iLayer, &
- ! iLayerLiqFluxSoil(iLayer-1), iLayerLiqFluxSoil(iLayer), mLayerLiqFluxSoil(iLayer), mLayerDepth(iLayer+nSnow)
- end do
- ! calculate the soil control on infiltration
- if(nSnow==0) then
- ! * case of infiltration into soil
- if(scalarMaxInfilRate > scalarRainPlusMelt)then ! infiltration is not rate-limited
- scalarSoilControl = (1._dp - scalarFrozenArea)*scalarInfilArea
- else
- scalarSoilControl = 0._dp ! (scalarRainPlusMelt exceeds maximum infiltration rate
- endif
else
- ! * case of infiltration into snow
- scalarSoilControl = 1._dp
- endif
-
- ! compute drainage from the soil zone (needed for mass balance checks)
- scalarSoilDrainage = iLayerLiqFluxSoil(nSoil)
-
- ! expand derivatives to the total water matric potential
- ! NOTE: arrays are offset because computing derivatives in interface fluxes, at the top and bottom of the layer respectively
- if(globalPrintFlag) print*, 'dPsiLiq_dPsi0(1:nSoil) = ', dPsiLiq_dPsi0(1:nSoil)
- dq_dHydStateAbove(1:nSoil) = dq_dHydStateAbove(1:nSoil) *dPsiLiq_dPsi0(1:nSoil)
- dq_dHydStateBelow(0:nSoil-1) = dq_dHydStateBelow(0:nSoil-1)*dPsiLiq_dPsi0(1:nSoil)
- dq_dHydStateLayerSurfVec(1:nSoil) = dq_dHydStateLayerSurfVec(1:nSoil)*dPsiLiq_dPsi0(1:nSoil)
-
- endif ! if calculating the liquid flux through soil
-
- ! *****
- ! * CALCULATE THE GROUNDWATER FLOW...
- ! ************************************
-
- ! check if computing soil hydrology
- if(nSoilOnlyHyd>0)then
-
- ! set baseflow fluxes to zero if the baseflow routine is not used
- if(local_ixGroundwater/=qbaseTopmodel)then
- ! (diagnostic variables in the data structures)
- scalarExfiltration = 0._dp ! exfiltration from the soil profile (m s-1)
- mLayerColumnOutflow(:) = 0._dp ! column outflow from each soil layer (m3 s-1)
- ! (variables needed for the numerical solution)
- mLayerBaseflow(:) = 0._dp ! baseflow from each soil layer (m s-1)
-
- ! topmodel-ish shallow groundwater
- else ! local_ixGroundwater==qbaseTopmodel
-
- ! check the derivative matrix is sized appropriately
- if(size(dBaseflow_dMatric,1)/=nSoil .or. size(dBaseflow_dMatric,2)/=nSoil)then
- message=trim(message)//'expect dBaseflow_dMatric to be nSoil x nSoil'
- err=20; return
- endif
-
- ! compute the baseflow flux
- call groundwatr(&
- ! input: model control
- nSnow, & ! intent(in): number of snow layers
- nSoil, & ! intent(in): number of soil layers
- nLayers, & ! intent(in): total number of layers
- firstFluxCall, & ! intent(in): logical flag to compute index of the lowest saturated layer
- ! input: state and diagnostic variables
- mLayerdTheta_dPsi, & ! intent(in): derivative in the soil water characteristic w.r.t. matric head in each layer (m-1)
- mLayerMatricHeadLiqTrial, & ! intent(in): liquid water matric potential (m)
- mLayerVolFracLiqTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of liquid water (-)
- mLayerVolFracIceTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of ice (-)
- ! input: data structures
- attr_data, & ! intent(in): model attributes
- mpar_data, & ! intent(in): model parameters
- prog_data, & ! intent(in): model prognostic variables for a local HRU
- diag_data, & ! intent(in): model diagnostic variables for a local HRU
- flux_data, & ! intent(inout): model fluxes for a local HRU
- ! output
- ixSaturation, & ! intent(inout) index of lowest saturated layer (NOTE: only computed on the first iteration)
- mLayerBaseflow, & ! intent(out): baseflow from each soil layer (m s-1)
- dBaseflow_dMatric, & ! intent(out): derivative in baseflow w.r.t. matric head (s-1)
- err,cmessage) ! intent(out): error control
- if(err/=0)then
- message=trim(message)//trim(cmessage)
- print*, message
- return
- endif
- endif ! computing baseflow flux
-
- ! compute total baseflow from the soil zone (needed for mass balance checks)
- scalarSoilBaseflow = sum(mLayerBaseflow)
-
- ! compute total runoff
- scalarTotalRunoff = scalarSurfaceRunoff + scalarSoilDrainage + scalarSoilBaseflow
-
- endif ! if computing soil hydrology
+ ! define forcing for the soil domain for the case of no snow layers
+ ! NOTE: in case where nSnowOnlyHyd==0 AND snow layers exist, then scalarRainPlusMelt is taken from the previous flux evaluation
+ if(nSnow==0)then !no snow layers
+ scalarRainPlusMelt = (scalarThroughfallRain + scalarCanopyLiqDrainage)/iden_water & ! liquid flux from the canopy (m s-1)
+ + drainageMeltPond/iden_water ! melt of the snow without a layer (m s-1)
- ! *****
- ! (7) CALCULATE FLUXES FOR THE DEEP AQUIFER...
- ! ********************************************
-
- ! check if computing aquifer fluxes
- if(ixAqWat/=integerMissing)then
-
- ! identify modeling decision
- if(local_ixGroundwater==bigBucket)then
- ! compute fluxes for the big bucket
- call bigAquifer(&
- ! input: state variables and fluxes
- scalarAquiferStorageTrial, & ! intent(in): trial value of aquifer storage (m)
- scalarCanopyTranspiration, & ! intent(in): canopy transpiration (kg m-2 s-1)
- scalarSoilDrainage, & ! intent(in): soil drainage (m s-1)
- ! input: pre-computed derivatives
- dCanopyTrans_dCanWat, & ! intent(in): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
- dCanopyTrans_dTCanair, & ! intent(in): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTCanopy, & ! intent(in): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
- dCanopyTrans_dTGround, & ! intent(in): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
- ! input: diagnostic variables and parameters
- mpar_data, & ! intent(in): model parameter structure
- diag_data, & ! intent(in): diagnostic variable structure
+ if(ixVegHyd/=integerMissing)then
+ ! save canopy derivatives
+ above_soilLiqFluxDeriv = scalarCanopyLiqDeriv/iden_water ! derivative in (throughfall + drainage) w.r.t. canopy liquid water
+ above_soildLiq_dTk = dCanLiq_dTcanopy ! derivative of canopy liquid storage w.r.t. temperature
+ above_soilFracLiq = scalarFracLiqVeg ! fraction of liquid water in canopy (-)
+ else
+ above_soilLiqFluxDeriv = 0._rkind
+ above_soildLiq_dTk = 0._rkind
+ above_soilFracLiq = 0._rkind
+ endif
+ else ! snow layers, take from previous flux calculation
+ above_soilLiqFluxDeriv = iLayerLiqFluxSnowDeriv(nSnow) ! derivative in vertical liquid water flux at bottom snow layer interface
+ above_soildLiq_dTk = mLayerdTheta_dTk(nSnow) ! derivative in volumetric liquid water content in bottom snow layer w.r.t. temperature
+ above_soilFracLiq = mLayerFracLiqSnow(nSnow) ! fraction of liquid water in bottom snow layer (-)
+ endif ! snow layers or not
+
+ endif ! if calculating the liquid flux through snow
+
+ ! *****
+ ! * CALCULATE THE LIQUID FLUX THROUGH SOIL...
+ ! ********************************************
+
+ ! check the need to calculate the liquid flux through soil
+ if(nSoilOnlyHyd>0)then
+
+ ! calculate the liquid flux through soil
+ call soilLiqFlx(&
+ ! input: model control
+ nSoil, & ! intent(in): number of soil layers
+ firstSplitOper, & ! intent(in): flag indicating first flux call in a splitting operation
+ (scalarSolution .and. .not.firstFluxCall), & ! intent(in): flag to indicate the scalar solution
+ .true., & ! intent(in): flag indicating if derivatives are desired
+ ! input: trial state variables
+ mLayerTempTrial(nSnow+1:nLayers), & ! intent(in): trial temperature at the current iteration (K)
+ mLayerMatricHeadLiqTrial(1:nSoil), & ! intent(in): liquid water matric potential (m)
+ mLayerVolFracLiqTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of liquid water (-)
+ mLayerVolFracIceTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of ice (-)
+ ! input: pre-computed deriavatives
+ mLayerdTheta_dTk(nSnow+1:nLayers), & ! intent(in): derivative in volumetric liquid water content w.r.t. temperature (K-1)
+ dPsiLiq_dTemp(1:nSoil), & ! intent(in): derivative in liquid water matric potential w.r.t. temperature (m K-1)
+ dCanopyTrans_dCanWat, & ! intent(in): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
+ dCanopyTrans_dTCanair, & ! intent(in): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTCanopy, & ! intent(in): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTGround, & ! intent(in): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
+ above_soilLiqFluxDeriv, & ! intent(in): derivative in layer above soil (canopy or snow) liquid flux w.r.t. liquid water
+ above_soildLiq_dTk, & ! intent(in): derivative of layer above soil (canopy or snow) liquid flux w.r.t. temperature
+ above_soilFracLiq, & ! intent(in): fraction of liquid water layer above soil (canopy or snow) (-)
+ ! input: fluxes
+ scalarCanopyTranspiration, & ! intent(in): canopy transpiration (kg m-2 s-1)
+ scalarGroundEvaporation, & ! intent(in): ground evaporation (kg m-2 s-1)
+ scalarRainPlusMelt, & ! intent(in): rain plus melt (m s-1)
+ ! input-output: data structures
+ mpar_data, & ! intent(in): model parameters
+ indx_data, & ! intent(in): model indices
+ prog_data, & ! intent(inout): model prognostic variables for a local HRU
+ diag_data, & ! intent(inout): model diagnostic variables for a local HRU
+ flux_data, & ! intent(inout): model fluxes for a local HRU
+ ! output: diagnostic variables for surface runoff
+ scalarMaxInfilRate, & ! intent(inout): maximum infiltration rate (m s-1)
+ scalarInfilArea, & ! intent(inout): fraction of unfrozen area where water can infiltrate (-)
+ scalarFrozenArea, & ! intent(inout): fraction of area that is considered impermeable due to soil ice (-)
+ scalarSurfaceRunoff, & ! intent(inout): surface runoff (m s-1)
+ ! output: diagnostic variables for model layers
+ mLayerdTheta_dPsi, & ! intent(inout): derivative in the soil water characteristic w.r.t. psi (m-1)
+ mLayerdPsi_dTheta, & ! intent(inout): derivative in the soil water characteristic w.r.t. theta (m)
+ dHydCond_dMatric, & ! intent(inout): derivative in hydraulic conductivity w.r.t matric head (s-1)
! output: fluxes
- scalarAquiferTranspire, & ! intent(out): transpiration loss from the aquifer (m s-1)
- scalarAquiferRecharge, & ! intent(out): recharge to the aquifer (m s-1)
- scalarAquiferBaseflow, & ! intent(out): total baseflow from the aquifer (m s-1)
- dBaseflow_dAquifer, & ! intent(out): change in baseflow flux w.r.t. aquifer storage (s-1)
+ scalarInfiltration, & ! intent(inout): surface infiltration rate (m s-1) -- controls on infiltration only computed for iter==1
+ iLayerLiqFluxSoil, & ! intent(inout): liquid fluxes at layer interfaces (m s-1)
+ mLayerTranspire, & ! intent(inout): transpiration loss from each soil layer (m s-1)
+ mLayerHydCond, & ! intent(inout): hydraulic conductivity in each layer (m s-1)
+ ! output: derivatives in fluxes w.r.t. state variables -- matric head or volumetric lquid water -- in the layer above and layer below (m s-1 or s-1)
+ dq_dHydStateAbove, & ! intent(inout): derivatives in the flux w.r.t. matric head in the layer above (s-1)
+ dq_dHydStateBelow, & ! intent(inout): derivatives in the flux w.r.t. matric head in the layer below (s-1)
+ dq_dHydStateLayerSurfVec, & ! intent(inout): derivative in surface infiltration w.r.t. hydrology state in above soil snow or canopy and every soil layer (m s-1 or s-1)
+ ! output: derivatives in fluxes w.r.t. energy state variables -- now just temperature -- in the layer above and layer below (m s-1 K-1)
+ dq_dNrgStateAbove, & ! intent(inout): derivatives in the flux w.r.t. temperature in the layer above (m s-1 K-1)
+ dq_dNrgStateBelow, & ! intent(inout): derivatives in the flux w.r.t. temperature in the layer below (m s-1 K-1)
+ dq_dNrgStateLayerSurfVec, & ! intent(inout): derivative in surface infiltration w.r.t. energy state in above soil snow or canopy and every soil layer (m s-1 K-1)
! output: derivatives in transpiration w.r.t. canopy state variables
- dAquiferTrans_dTCanair, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy air temperature
- dAquiferTrans_dTCanopy, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy temperature
- dAquiferTrans_dTGround, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. ground temperature
- dAquiferTrans_dCanWat, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy total water
+ mLayerdTrans_dTCanair, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy air temperature
+ mLayerdTrans_dTCanopy, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy temperature
+ mLayerdTrans_dTGround, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. ground temperature
+ mLayerdTrans_dCanWat, & ! intent(inout): derivatives in the soil layer transpiration flux w.r.t. canopy total water
! output: error control
- err,cmessage) ! intent(out): error control
- if(err/=0)then
- message=trim(message)//trim(cmessage)
- print*, message
- return
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
+
+ ! calculate net liquid water fluxes for each soil layer (s-1)
+ do iLayer=1,nSoil
+ mLayerLiqFluxSoil(iLayer) = -(iLayerLiqFluxSoil(iLayer) - iLayerLiqFluxSoil(iLayer-1))/mLayerDepth(iLayer+nSnow)
+ end do
+
+ ! calculate the soil control on infiltration
+ if(nSnow==0) then
+ ! * case of infiltration into soil
+ if(scalarMaxInfilRate > scalarRainPlusMelt)then ! infiltration is not rate-limited
+ scalarSoilControl = (1._rkind - scalarFrozenArea)*scalarInfilArea
+ else
+ scalarSoilControl = 0._rkind ! (scalarRainPlusMelt exceeds maximum infiltration rate
+ endif
+ else
+ ! * case of infiltration into snow
+ scalarSoilControl = 1._rkind
endif
- ! compute total runoff (overwrite previously calculated value before considering aquifer)
- scalarTotalRunoff = scalarSurfaceRunoff + scalarAquiferBaseflow
-
- ! if no aquifer, then fluxes are zero
- else
- scalarAquiferTranspire = 0._dp ! transpiration loss from the aquifer (m s-1)
- scalarAquiferRecharge = 0._dp ! recharge to the aquifer (m s-1)
- scalarAquiferBaseflow = 0._dp ! total baseflow from the aquifer (m s-1)
- dBaseflow_dAquifer = 0._dp ! change in baseflow flux w.r.t. aquifer storage (s-1)
- end if ! no aquifer
-
- endif ! if computing aquifer fluxes
+ ! compute drainage from the soil zone (needed for mass balance checks and in aquifer recharge)
+ scalarSoilDrainage = iLayerLiqFluxSoil(nSoil)
+
+ ! expand derivatives to the total water matric potential
+ ! NOTE: arrays are offset because computing derivatives in interface fluxes, at the top and bottom of the layer respectively
+ if(globalPrintFlag) print*, 'dPsiLiq_dPsi0(1:nSoil) = ', dPsiLiq_dPsi0(1:nSoil)
+ dq_dHydStateAbove(1:nSoil) = dq_dHydStateAbove(1:nSoil) *dPsiLiq_dPsi0(1:nSoil)
+ dq_dHydStateBelow(0:nSoil-1) = dq_dHydStateBelow(0:nSoil-1)*dPsiLiq_dPsi0(1:nSoil)
+ dq_dHydStateLayerSurfVec(1:nSoil) = dq_dHydStateLayerSurfVec(1:nSoil)*dPsiLiq_dPsi0(1:nSoil)
+
+ endif ! if calculating the liquid flux through soil
+
+ ! *****
+ ! * CALCULATE THE GROUNDWATER FLOW...
+ ! ************************************
+
+ ! check if computing soil hydrology
+ if(nSoilOnlyHyd>0)then
+
+ ! set baseflow fluxes to zero if the topmodel baseflow routine is not used
+ if(local_ixGroundwater/=qbaseTopmodel)then
+ ! (diagnostic variables in the data structures)
+ scalarExfiltration = 0._rkind ! exfiltration from the soil profile (m s-1)
+ mLayerColumnOutflow(:) = 0._rkind ! column outflow from each soil layer (m3 s-1)
+ ! (variables needed for the numerical solution)
+ mLayerBaseflow(:) = 0._rkind ! baseflow from each soil layer (m s-1)
+
+ ! topmodel-ish shallow groundwater
+ else ! local_ixGroundwater==qbaseTopmodel
+
+ ! check the derivative matrix is sized appropriately
+ if(size(dBaseflow_dMatric,1)/=nSoil .or. size(dBaseflow_dMatric,2)/=nSoil)then
+ message=trim(message)//'expect dBaseflow_dMatric to be nSoil x nSoil'
+ err=20; return
+ endif
+
+ ! compute the baseflow flux
+ call groundwatr(&
+ ! input: model control
+ nSnow, & ! intent(in): number of snow layers
+ nSoil, & ! intent(in): number of soil layers
+ nLayers, & ! intent(in): total number of layers
+ firstFluxCall, & ! intent(in): logical flag to compute index of the lowest saturated layer
+ ! input: state and diagnostic variables
+ mLayerdTheta_dPsi, & ! intent(in): derivative in the soil water characteristic w.r.t. matric head in each layer (m-1)
+ mLayerMatricHeadLiqTrial, & ! intent(in): liquid water matric potential (m)
+ mLayerVolFracLiqTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of liquid water (-)
+ mLayerVolFracIceTrial(nSnow+1:nLayers), & ! intent(in): volumetric fraction of ice (-)
+ ! input: data structures
+ attr_data, & ! intent(in): model attributes
+ mpar_data, & ! intent(in): model parameters
+ prog_data, & ! intent(in): model prognostic variables for a local HRU
+ diag_data, & ! intent(in): model diagnostic variables for a local HRU
+ flux_data, & ! intent(inout): model fluxes for a local HRU
+ ! output
+ ixSaturation, & ! intent(inout) index of lowest saturated layer (NOTE: only computed on the first iteration)
+ mLayerBaseflow, & ! intent(out): baseflow from each soil layer (m s-1)
+ dBaseflow_dMatric, & ! intent(out): derivative in baseflow w.r.t. matric head (s-1)
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
+ endif ! computing baseflow flux
+
+ ! compute total baseflow from the soil zone (needed for mass balance checks)
+ scalarSoilBaseflow = sum(mLayerBaseflow)
+
+ ! compute total runoff
+ ! (Note: scalarSoilBaseflow is zero if topmodel is not used)
+ ! (Note: scalarSoilBaseflow may need to re-envisioned in topmodel formulation if parts of it flow into neighboring soil rather than exfiltrate)
+ scalarTotalRunoff = scalarSurfaceRunoff + scalarSoilDrainage + scalarSoilBaseflow
+
+ endif ! if computing soil hydrology
+
+
+ ! *****
+ ! (7) CALCULATE FLUXES FOR THE DEEP AQUIFER...
+ ! ********************************************
+
+ ! check if computing aquifer fluxes
+ if(ixAqWat/=integerMissing)then
+
+ ! identify modeling decision
+ if(local_ixGroundwater==bigBucket)then
+
+ ! compute fluxes for the big bucket
+ call bigAquifer(&
+ ! input: state variables and fluxes
+ scalarAquiferStorageTrial, & ! intent(in): trial value of aquifer storage (m)
+ scalarCanopyTranspiration, & ! intent(in): canopy transpiration (kg m-2 s-1)
+ scalarSoilDrainage, & ! intent(in): soil drainage (m s-1)
+ ! input: pre-computed derivatives
+ dCanopyTrans_dCanWat, & ! intent(in): derivative in canopy transpiration w.r.t. canopy total water content (s-1)
+ dCanopyTrans_dTCanair, & ! intent(in): derivative in canopy transpiration w.r.t. canopy air temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTCanopy, & ! intent(in): derivative in canopy transpiration w.r.t. canopy temperature (kg m-2 s-1 K-1)
+ dCanopyTrans_dTGround, & ! intent(in): derivative in canopy transpiration w.r.t. ground temperature (kg m-2 s-1 K-1)
+ ! input: diagnostic variables and parameters
+ mpar_data, & ! intent(in): model parameter structure
+ diag_data, & ! intent(in): diagnostic variable structure
+ ! output: fluxes
+ scalarAquiferTranspire, & ! intent(out): transpiration loss from the aquifer (m s-1)
+ scalarAquiferRecharge, & ! intent(out): recharge to the aquifer (m s-1)
+ scalarAquiferBaseflow, & ! intent(out): total baseflow from the aquifer (m s-1)
+ dBaseflow_dAquifer, & ! intent(out): change in baseflow flux w.r.t. aquifer storage (s-1)
+ ! output: derivatives in transpiration w.r.t. canopy state variables
+ dAquiferTrans_dTCanair, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy air temperature
+ dAquiferTrans_dTCanopy, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy temperature
+ dAquiferTrans_dTGround, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. ground temperature
+ dAquiferTrans_dCanWat, & ! intent(inout): derivatives in the aquifer transpiration flux w.r.t. canopy total water
+ ! output: error control
+ err,cmessage) ! intent(out): error control
+ if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
+
+ ! compute total runoff (overwrite previously calculated value before considering aquifer).
+ ! (Note: SoilDrainage goes into aquifer, not runoff)
+ scalarTotalRunoff = scalarSurfaceRunoff + scalarAquiferBaseflow
+
+ ! if no aquifer, then fluxes are zero
+ else
+ scalarAquiferTranspire = 0._rkind ! transpiration loss from the aquifer (m s-1)
+ scalarAquiferRecharge = 0._rkind ! recharge to the aquifer (m s-1)
+ scalarAquiferBaseflow = 0._rkind ! total baseflow from the aquifer (m s-1)
+ dBaseflow_dAquifer = 0._rkind ! change in baseflow flux w.r.t. aquifer storage (s-1)
+ end if ! no aquifer
+
+ endif ! if computing aquifer fluxes
+
+ ! *****
+ ! (X) WRAP UP...
+ ! *************
+
+ ! define model flux vector for the vegetation sub-domain
+ if(ixCasNrg/=integerMissing) fluxVec(ixCasNrg) = scalarCanairNetNrgFlux/canopyDepth
+ if(ixVegNrg/=integerMissing) fluxVec(ixVegNrg) = scalarCanopyNetNrgFlux/canopyDepth
+ if(ixVegHyd/=integerMissing) fluxVec(ixVegHyd) = scalarCanopyNetLiqFlux ! NOTE: solid fluxes are handled separately
+
+ ! populate the flux vector for energy
+ if(nSnowSoilNrg>0)then
+ do concurrent (iLayer=1:nLayers,ixSnowSoilNrg(iLayer)/=integerMissing) ! (loop through non-missing energy state variables in the snow+soil domain)
+ fluxVec( ixSnowSoilNrg(iLayer) ) = mLayerNrgFlux(iLayer)
+ end do ! looping through non-missing energy state variables in the snow+soil domain
+ endif
- ! *****
- ! (X) WRAP UP...
- ! *************
-
- ! define model flux vector for the vegetation sub-domain
- if(ixCasNrg/=integerMissing) fluxVec(ixCasNrg) = scalarCanairNetNrgFlux/canopyDepth
- if(ixVegNrg/=integerMissing) fluxVec(ixVegNrg) = scalarCanopyNetNrgFlux/canopyDepth
- if(ixVegHyd/=integerMissing) fluxVec(ixVegHyd) = scalarCanopyNetLiqFlux ! NOTE: solid fluxes are handled separately
-
- ! populate the flux vector for energy
- if(nSnowSoilNrg>0)then
- do concurrent (iLayer=1:nLayers,ixSnowSoilNrg(iLayer)/=integerMissing) ! (loop through non-missing energy state variables in the snow+soil domain)
- fluxVec( ixSnowSoilNrg(iLayer) ) = mLayerNrgFlux(iLayer)
- end do ! looping through non-missing energy state variables in the snow+soil domain
- endif
-
- ! populate the flux vector for hydrology
- ! NOTE: ixVolFracWat and ixVolFracLiq can also include states in the soil domain, hence enable primary variable switching
- if(nSnowSoilHyd>0)then ! check if any hydrology states exist
- do iLayer=1,nLayers
- if(ixSnowSoilHyd(iLayer)/=integerMissing)then ! check if a given hydrology state exists
- select case( layerType(iLayer) )
- case(iname_snow); fluxVec( ixSnowSoilHyd(iLayer) ) = mLayerLiqFluxSnow(iLayer)
- case(iname_soil); fluxVec( ixSnowSoilHyd(iLayer) ) = mLayerLiqFluxSoil(iLayer-nSnow)
- case default; err=20; message=trim(message)//'expect layerType to be either iname_snow or iname_soil'; return
- end select
- endif ! if a given hydrology state exists
- end do ! looping through non-missing energy state variables in the snow+soil domain
- endif ! if any hydrology states exist
-
- ! compute the flux vector for the aquifer
- if(ixAqWat/=integerMissing) fluxVec(ixAqWat) = scalarAquiferTranspire + scalarAquiferRecharge - scalarAquiferBaseflow
-
- ! set the first flux call to false
- firstFluxCall=.false.
+ ! populate the flux vector for hydrology
+ ! NOTE: ixVolFracWat and ixVolFracLiq can also include states in the soil domain, hence enable primary variable switching
+ if(nSnowSoilHyd>0)then ! check if any hydrology states exist
+ do iLayer=1,nLayers
+ if(ixSnowSoilHyd(iLayer)/=integerMissing)then ! check if a given hydrology state exists
+ select case( layerType(iLayer) )
+ case(iname_snow); fluxVec( ixSnowSoilHyd(iLayer) ) = mLayerLiqFluxSnow(iLayer)
+ case(iname_soil); fluxVec( ixSnowSoilHyd(iLayer) ) = mLayerLiqFluxSoil(iLayer-nSnow)
+ case default; err=20; message=trim(message)//'expect layerType to be either iname_snow or iname_soil'; return
+ end select
+ endif ! if a given hydrology state exists
+ end do ! looping through non-missing energy state variables in the snow+soil domain
+ endif ! if any hydrology states exist
+
+ ! compute the flux vector for the aquifer
+ if(ixAqWat/=integerMissing) fluxVec(ixAqWat) = scalarAquiferTranspire + scalarAquiferRecharge - scalarAquiferBaseflow
+
+ ! set the first flux call to false
+ firstFluxCall=.false.
! end association to variables in the data structures
end associate
@@ -961,32 +942,32 @@ end subroutine computFlux
! public subroutine soilCmpres: compute soil compressibility (-) and its derivative w.r.t matric head (m-1)
! **********************************************************************************************************
subroutine soilCmpres(&
- ! input:
- ixRichards, & ! intent(in): choice of option for Richards' equation
- ixBeg,ixEnd, & ! intent(in): start and end indices defining desired layers
- mLayerMatricHead, & ! intent(in): matric head at the start of the time step (m)
- mLayerMatricHeadTrial, & ! intent(in): trial value of matric head (m)
- mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each soil layer (-)
- mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice content in each soil layer (-)
- specificStorage, & ! intent(in): specific storage coefficient (m-1)
- theta_sat, & ! intent(in): soil porosity (-)
- ! output:
- compress, & ! intent(out): compressibility of the soil matrix (-)
- dCompress_dPsi, & ! intent(out): derivative in compressibility w.r.t. matric head (m-1)
- err,message) ! intent(out): error code and error message
+ ! input:
+ ixRichards, & ! intent(in): choice of option for Richards' equation
+ ixBeg,ixEnd, & ! intent(in): start and end indices defining desired layers
+ mLayerMatricHead, & ! intent(in): matric head at the start of the time step (m)
+ mLayerMatricHeadTrial, & ! intent(in): trial value of matric head (m)
+ mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each soil layer (-)
+ mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice content in each soil layer (-)
+ specificStorage, & ! intent(in): specific storage coefficient (m-1)
+ theta_sat, & ! intent(in): soil porosity (-)
+ ! output:
+ compress, & ! intent(out): compressibility of the soil matrix (-)
+ dCompress_dPsi, & ! intent(out): derivative in compressibility w.r.t. matric head (m-1)
+ err,message) ! intent(out): error code and error message
implicit none
! input:
integer(i4b),intent(in) :: ixRichards ! choice of option for Richards' equation
integer(i4b),intent(in) :: ixBeg,ixEnd ! start and end indices defining desired layers
- real(dp),intent(in) :: mLayerMatricHead(:) ! matric head at the start of the time step (m)
- real(dp),intent(in) :: mLayerMatricHeadTrial(:) ! trial value for matric head (m)
- real(dp),intent(in) :: mLayerVolFracLiqTrial(:) ! trial value for volumetric fraction of liquid water (-)
- real(dp),intent(in) :: mLayerVolFracIceTrial(:) ! trial value for volumetric fraction of ice (-)
- real(dp),intent(in) :: specificStorage ! specific storage coefficient (m-1)
- real(dp),intent(in) :: theta_sat(:) ! soil porosity (-)
+ real(rkind),intent(in) :: mLayerMatricHead(:) ! matric head at the start of the time step (m)
+ real(rkind),intent(in) :: mLayerMatricHeadTrial(:) ! trial value for matric head (m)
+ real(rkind),intent(in) :: mLayerVolFracLiqTrial(:) ! trial value for volumetric fraction of liquid water (-)
+ real(rkind),intent(in) :: mLayerVolFracIceTrial(:) ! trial value for volumetric fraction of ice (-)
+ real(rkind),intent(in) :: specificStorage ! specific storage coefficient (m-1)
+ real(rkind),intent(in) :: theta_sat(:) ! soil porosity (-)
! output:
- real(dp),intent(inout) :: compress(:) ! soil compressibility (-)
- real(dp),intent(inout) :: dCompress_dPsi(:) ! derivative in soil compressibility w.r.t. matric head (m-1)
+ real(rkind),intent(inout) :: compress(:) ! soil compressibility (-)
+ real(rkind),intent(inout) :: dCompress_dPsi(:) ! derivative in soil compressibility w.r.t. matric head (m-1)
integer(i4b),intent(out) :: err ! error code
character(*),intent(out) :: message ! error message
! local variables
@@ -997,34 +978,36 @@ subroutine soilCmpres(&
! (only compute for the mixed form of Richards' equation)
if(ixRichards==mixdform)then
do iLayer=1,size(mLayerMatricHead)
- if(iLayer>=ixBeg .and. iLayer<=ixEnd)then
+ if(iLayer>=ixBeg .and. iLayer<=ixEnd)then
! compute the derivative for the compressibility term (m-1)
dCompress_dPsi(iLayer) = specificStorage*(mLayerVolFracLiqTrial(iLayer) + mLayerVolFracIceTrial(iLayer))/theta_sat(iLayer)
! compute the compressibility term (-)
compress(iLayer) = (mLayerMatricHeadTrial(iLayer) - mLayerMatricHead(iLayer))*dCompress_dPsi(iLayer)
- endif
+ endif
end do
else
- compress(:) = 0._dp
- dCompress_dPsi(:) = 0._dp
+ compress(:) = 0._rkind
+ dCompress_dPsi(:) = 0._rkind
end if
end subroutine soilCmpres
+
+
! **********************************************************************************************************
! public subroutine soilCmpres: compute soil compressibility (-) and its derivative w.r.t matric head (m-1)
! **********************************************************************************************************
subroutine soilCmpresSundials(&
- ! input:
- ixRichards, & ! intent(in): choice of option for Richards' equation
- ixBeg,ixEnd, & ! intent(in): start and end indices defining desired layers
- mLayerMatricHeadPrime, & ! intent(in): matric head at the start of the time step (m)
- mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each soil layer (-)
- mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice content in each soil layer (-)
- specificStorage, & ! intent(in): specific storage coefficient (m-1)
- theta_sat, & ! intent(in): soil porosity (-)
- ! output:
- compress, & ! intent(out): compressibility of the soil matrix (-)
- dCompress_dPsi, & ! intent(out): derivative in compressibility w.r.t. matric head (m-1)
- err,message) ! intent(out): error code and error message
+ ! input:
+ ixRichards, & ! intent(in): choice of option for Richards' equation
+ ixBeg,ixEnd, & ! intent(in): start and end indices defining desired layers
+ mLayerMatricHeadPrime, & ! intent(in): matric head at the start of the time step (m)
+ mLayerVolFracLiqTrial, & ! intent(in): trial value for the volumetric liquid water content in each soil layer (-)
+ mLayerVolFracIceTrial, & ! intent(in): trial value for the volumetric ice content in each soil layer (-)
+ specificStorage, & ! intent(in): specific storage coefficient (m-1)
+ theta_sat, & ! intent(in): soil porosity (-)
+ ! output:
+ compress, & ! intent(out): compressibility of the soil matrix (-)
+ dCompress_dPsi, & ! intent(out): derivative in compressibility w.r.t. matric head (m-1)
+ err,message) ! intent(out): error code and error message
implicit none
! input:
integer(i4b),intent(in) :: ixRichards ! choice of option for Richards' equation
@@ -1048,15 +1031,15 @@ subroutine soilCmpresSundials(&
if(ixRichards==mixdform)then
do iLayer=1,size(mLayerMatricHeadPrime)
if(iLayer>=ixBeg .and. iLayer<=ixEnd)then
- ! compute the derivative for the compressibility term (m-1)
- dCompress_dPsi(iLayer) = specificStorage*(mLayerVolFracLiqTrial(iLayer) + mLayerVolFracIceTrial(iLayer))/theta_sat(iLayer)
- ! compute the compressibility term (-)
- compress(iLayer) = mLayerMatricHeadPrime(iLayer) * dCompress_dPsi(iLayer)
+ ! compute the derivative for the compressibility term (m-1)
+ dCompress_dPsi(iLayer) = specificStorage*(mLayerVolFracLiqTrial(iLayer) + mLayerVolFracIceTrial(iLayer))/theta_sat(iLayer)
+ ! compute the compressibility term (-)
+ compress(iLayer) = mLayerMatricHeadPrime(iLayer) * dCompress_dPsi(iLayer)
endif
end do
else
- compress(:) = 0._dp
- dCompress_dPsi(:) = 0._dp
+ compress(:) = 0._rkind
+ dCompress_dPsi(:) = 0._rkind
end if
end subroutine soilCmpresSundials