ldndc::PhysiologyPNET Class Reference

Vegetation model PNET. More...

#include <models/physiology/pnet/pnet.h>

Inherits ldndc::MBE_LegacyModel.

Public Member Functions
lerr_t	solve ()
lerr_t	sleep ()

Private Member Functions
void	update_nitrogen_concentration (MoBiLE_Plant *)
void	YearlyReset (MoBiLE_Plant *)
void	Phenology (MoBiLE_Plant *, int, double)
void	Photosyn (MoBiLE_Plant *)
void	AllocateMo (MoBiLE_Plant *)
void	AllocateYr (MoBiLE_Plant *)
void	N_uptake (MoBiLE_Plant *, double &)
void	Update (MoBiLE_Plant *, double, double)
void	TransformPhysiology (MoBiLE_Plant *)
size_t	get_na_max_pnet ()

Static Private Attributes
static const double	TMIN_NUP = 6.0
static const double	TMAX_NUP = 25.0
static const double	KM_NO3 = 1.0
static const double	KM_NH4 = 1.0
static const double	SRATEMAX = 0.01
static const double	C = 0.02
static const double	NUPTAIR = 0.0
static const double	GROWMAX = 0.95

Detailed Description

Vegetation model PNET.

Author

Ruediger Grote

Member Function Documentation

AllocateMo()

void ldndc::PhysiologyPNET::AllocateMo ( MoBiLE_Plant *  )

private

AllocateMo

AllocateYr()

void ldndc::PhysiologyPNET::AllocateYr ( MoBiLE_Plant *  p )

private

AllocateYr

NOTE: the shift in carbon between compartments disturbs the nitrogen balance at this day!!!

{
    double fratio( (*p)->QWODFOLMIN() / (1.0 + (*p)->QWODFOLMIN()));
    double TotalC( WoodC + BudC + FolMass * cbm::CCDM);
    WoodC  = TotalC * fratio;
    BudC   = TotalC - WoodC - FolMass * cbm::CCDM;
 
    double mSapSum( 0.0);
    double mBudSum( 0.0);
//    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//    {
        mSapSum += p->mSap;
        mBudSum += p->mBud;
//    }
 
    // check for consistency of bud biomass
    transfDWbs = 0.0;
    transfDWsb = 0.0;
    double  scale_ay( cbm::KG_IN_G / cbm::CCDM);  // gC in kgDW
//    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//    {
        // shift from buds to sapwood tissue
        double freeBudDW( 0.0);
        double requBudDW(0.0);
        if (mBudSum > 0.0)
        {
            double fm_vt( p->mBud / mBudSum);
            freeBudDW = std::max(0.0, p->mBud - BudC * scale_ay * fm_vt);
            requBudDW = std::max(0.0, BudC * scale_ay * fm_vt - p->mBud);
        }
        else
        {
            requBudDW = BudC * scale_ay / double( m_veg->size());
        }
 
        transfDWbs += freeBudDW;
 
        // shift from sapwood to buds
        double freeSapDW( 0.0);
        if (mSapSum > 0.0)
        {
            freeSapDW = std::max(0.0, p->mSap - WoodC * scale_ay * p->mSap / mSapSum);
        }
        transfDWsb += std::min(freeSapDW, requBudDW);
//    }
 
    // new bud and wood compartments
    BudC  += (transfDWsb - transfDWbs) / scale_ay;
    WoodC += (transfDWbs - transfDWsb) / scale_ay;
    if (BudC < 0.0)  BudC  = 0.0;
    if (WoodC < 0.0) WoodC = 0.0;
 
    // determination of maximum foliage growth next year
    double FolMassMaxOld( FolMassMax);
    FolMassMax = FolMass + BudC / cbm::CCDM; 
 
    // limitation of maximum foliage growth 
    double vegRelGrow( 0.0);
//    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//    {
        vegRelGrow += GROWMAX * p->f_area;
//    }
 
    // limitation of maximum foliage growth by maximum biomass (rg 03.02.11); including height and density impact (rg 19.10.11)
    double FolMassOpt( 0.0);
//    for ( TreeIterator w = m_veg->groupbegin< species::wood >();
//            w != m_veg->groupend< species::wood >(); ++w)
//    {
        FolMassOpt += m_pf.get_optimum_foliage_biomass_trees( p) * cbm::G_IN_KG;
//    }
 
    // 15.05.11 rg maximum density effect considered
    if ( FolMassMax > std::min(FolMassOpt, FolMassMaxOld * (1.0 + vegRelGrow * LightEffMin) * 1.0/p->f_area))
    {
        double excessDW( FolMassMax - std::min(FolMassOpt, FolMassMaxOld * (1.0 + vegRelGrow * LightEffMin)/p->f_area));
        BudC       -= (excessDW * cbm::CCDM);
        FolMassMax  = std::min(FolMassOpt, FolMassMaxOld * (1.0 + vegRelGrow * LightEffMin)/p->f_area);
        transfDWbs  += excessDW;
    }
 
    // determination of minimum foliage in the next year
    if (FolReten > 0.0)
    {
        FolMassMin = FolMassMax * (1.0 - 1.0 / FolReten);
    }
    else
    {
        FolMassMin = 0.0;
    }
 
    // stand death
    if (FolMassMax < 0.0)
    {
//        for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//        {
            p->mCor += p->mSap;
            p->mSap  = 0.0;
            p->mBud  = 0.0;
            p->mFrt  = 0.0;
//        }
    }
}

get_na_max_pnet()

size_t ldndc::PhysiologyPNET::get_na_max_pnet ( )

private

get_na_max_pnet

{
    size_t nb_ageclasses_total( 0.0);
    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
    {
        nb_ageclasses_total += (*vt)->nb_ageclasses();
    }
    return (int)ceil(sqrt((double)nb_ageclasses_total));
}

N_uptake()

void ldndc::PhysiologyPNET::N_uptake ( MoBiLE_Plant *  p, double &  _n_uptake  )

private

N_uptake

{
    bool physio(false);
 
    // nitrogen demand
    double nopt( 0.0);
    double nact( 0.0);
 
    double nFolOpt_vt( m_pf.optimum_nitrogen_content_foliage( p));
    double nBudOpt_vt( m_pf.optimum_nitrogen_content_buds( p));
        
    nopt += (nFolOpt_vt + p->mFrt * (*p)->NC_FINEROOTS_MAX()
             + p->mSap * (*p)->NCSAPOPT() + nBudOpt_vt);              // rg 11.06.10 (end of changes)
    
    nact += (p->mFol * p->ncFol + p->mFrt * p->ncFrt
             + p->mSap * p->ncSap + p->mBud * p->ncBud);
 
    double n_demand( cbm::bound_min(0.0, nopt - nact));
 
    if (physio == false)
    {
        // standard DNDC procedure
 
        // total nitrogen availability
        double  n_avail( 0.0);
        for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
        {
            if (sc_.depth_sl[sl] <= p->rooting_depth)
            {
                n_avail += (*p)->EXPL_NO3() * no3_sl[sl] +
                           (*p)->EXPL_NH4() * nh4_sl[sl];
                      // + (*p)->EXPL_DON() * don_sl[sl];
            }
        }
 
        // nitrogen uptake
        double  n_fact( 0.0);
        double const  transpiration = wc_.accumulated_wateruptake_sl.sum() - accumulated_transpiration_old;
        accumulated_transpiration_old = wc_.accumulated_wateruptake_sl.sum();
        if (( n_avail > 0.0) && ( transpiration > 0.0))
        {
            n_fact = std::max(0.0, std::min(1.0, n_demand / n_avail));
        }
 
        for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
        {
            n_avail = no3_sl[sl] + nh4_sl[sl]; // don_sl[sl];
 
            double  layer_ava_no3( 0.0);
            double  layer_ava_nh4( 0.0);
            if ( cbm::flt_less_equal( sc_.depth_sl[sl], p->rooting_depth) &&
                 cbm::flt_greater_zero( n_avail) &&
                 cbm::flt_greater( wc_.wc_sl[sl], wc_.wc_wp_sl[sl]))
            {
                layer_ava_no3 = (*p)->EXPL_NO3() * no3_sl[sl] * n_fact;
                layer_ava_nh4 = (*p)->EXPL_NH4() * nh4_sl[sl] * n_fact;
            }
 
//            for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//            {
            double const mfrt_veg_sl( ph_.mfrt_sl(m_veg, sl));
            if ( cbm::flt_greater_zero( mfrt_veg_sl))
            {
                double  fm_vt( p->fFrt_sl[sl] * p->mFrt / mfrt_veg_sl);
                _n_uptake += layer_ava_no3 * fm_vt;
                _n_uptake += layer_ava_nh4 * fm_vt;
 
                ph_.accumulated_no3_uptake_sl[sl] += layer_ava_no3 * fm_vt;
                ph_.accumulated_nh4_uptake_sl[sl] += layer_ava_nh4 * fm_vt;
            }
//            }
        }
    }
    else
    {
        // physiological calculation (Arjan)
 
        double no3_mg,nh4_mg,usNH4,usNO3 /*don_mg*/;
        double fact_t,fact_no3,fact_nh4 /*fact_don*/;
        double maxUptake_NO3,maxUptake_NH4,uptakeN /*maxUptake_DON*/;
 
//        for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//        {
            for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
            {
                double no3_uptake_sl( 0.0);
                double nh4_uptake_sl( 0.0);
                
                if ( cbm::flt_greater_zero( n_demand))
                {
                    // temperature dependency of root uptake
                    fact_t = m_pf.nitrogen_uptake_temperature_dependency(
                                                                     TMIN_NUP,
                                                                     TMAX_NUP,
                                                                     mc_.nd_temp_sl[sl]);
 
                    // nitrogen availability factor for nh4 and no3
                    no3_mg   = no3_sl[sl] * cbm::MG_IN_KG / (wc_.wc_sl[sl]* cbm::MM_IN_M * sc_.h_sl[sl]);
                    nh4_mg   = nh4_sl[sl] * cbm::MG_IN_KG / (wc_.wc_sl[sl]* cbm::MM_IN_M * sc_.h_sl[sl]);
                  //don_mg   = don_sl[sl] * cbm::MG_IN_KG / (wc_.wc_sl[sl]* cbm::MM_IN_M * sc_.h_sl[sl]);
                    fact_no3 = no3_mg / (KM_NO3 + no3_mg);
                    fact_nh4 = nh4_mg / (KM_NH4 + nh4_mg);
                  //fact_don = don_mg / (KM_DON + don_mg);
 
                    // maximum specific uptake for nitrogen (rg 01.06.11 - no change in case CFCROPT = 0.0)
                    usNH4 = (*p)->US_NH4() * (1.0-(*p)->CFCROPT()) + (*p)->US_NH4MYC() * (*p)->CFCROPT();
                    usNO3 = (*p)->US_NO3() * (1.0-(*p)->CFCROPT()) + (*p)->US_NO3MYC() * (*p)->CFCROPT();
                  //usDON = (*p)->US_DON() * (1.0-(*p)->CFCROPT()) + (*p)->US_DONMYC() * (*p)->CFCROPT();
 
 
                    // maximum uptake rate
                    maxUptake_NO3 = usNO3 * p->fFrt_sl[sl] * p->mFrt * fact_no3 * fact_t;
                    maxUptake_NH4 = usNH4 * p->fFrt_sl[sl] * p->mFrt * fact_nh4 * fact_t;
                  //maxUptake_DON = usDON * p->fFrt_sl[sl] * p->mFrt * fact_don * fact_t;
                    uptakeN = std::min(n_demand, maxUptake_NO3 + maxUptake_NH4 /* + maxUptake_DON */ );
 
                    // soil layer specific uptake
                    if ( uptakeN > 0.0)
                    {
                        no3_uptake_sl = std::min( no3_sl[sl], uptakeN * maxUptake_NO3 / (maxUptake_NO3 + maxUptake_NH4/* + maxUptake_DON */));
                        nh4_uptake_sl = std::min( nh4_sl[sl], uptakeN * maxUptake_NH4 / (maxUptake_NO3 + maxUptake_NH4/* + maxUptake_DON */));
                      //ph_.nd_don_uptake_sl[sl] = std::min( don_sl[sl], uptakeN * maxUptake_DON / (maxUptake_NO3 + maxUptake_NH4/* + maxUptake_DON */));
 
                    }
 
                    // updating N demand
                    n_demand -= no3_uptake_sl + nh4_uptake_sl;
                }
 
                // species specific uptake
                _n_uptake += no3_uptake_sl;
                _n_uptake += nh4_uptake_sl;
 
                ph_.accumulated_no3_uptake_sl[sl] += no3_uptake_sl;
                ph_.accumulated_nh4_uptake_sl[sl] += nh4_uptake_sl;
 
                nh4_sl[sl] -= nh4_uptake_sl;
                no3_sl[sl] -= no3_uptake_sl;
            }
//        }
    }
}

Phenology()

void ldndc::PhysiologyPNET::Phenology ( MoBiLE_Plant *  p, int  _GrowthPhase, double  _DWater  )

private

Phenology

{
    if ( _GrowthPhase == 1)
    {
        if ( _DWater >= ((*p)->H2OREF_A()*0.5))
        {
            p->growing_degree_days += cbm::bound_min( 0.0, mc_.nd_airtemperature);
        }
 
        if (( p->growing_degree_days >= (*p)->GDDFOLSTART()) && ( lclock_ref().yearday() >= year_start_) && ( _DWater >= (*p)->H2OREF_A()*0.5))
        {
            p->dvsFlushOld = p->dvsFlush;
            p->dvsFlush = cbm::bound(0.0, (p->growing_degree_days - (*p)->GDDFOLSTART()) / ((*p)->GDDFOLEND() - (*p)->GDDFOLSTART()), 1.0);
            double const delta_dvsFlush( p->dvsFlush - p->dvsFlushOld);
 
            if ( p->dvsFlushOld < 0.99)
            {
                mature = false;
            }
 
            double gFolC( std::min(BudCStart * delta_dvsFlush, BudC));
            double const OldFolMass( FolMass);
            FolMass     += (gFolC / cbm::CCDM);
            FolProdCMo   = ((FolMass - OldFolMass) * cbm::CCDM);
            FolGRespMo   = (FolProdCMo * (*p)->FYIELD());
        }
        else
        {
            FolProdCMo = 0.0;
            FolGRespMo = 0.0;
        }
    }
    else if ( _GrowthPhase == 2)
    {
        FolLitSum = 0.0;
 
        if (( PosCBalMass < FolMass) && mature
            && (( lclock_ref().yearday() > (unsigned int)((*p)->SENESCSTART()+year_start_-1)) || ( lclock_ref().yearday() < year_start_))
            && ( lclock_ref().yearday() != year_end_))
        {
            // rg 28.03.12 stand density already considered in annual update of foliage values
            double const FolMassNew( std::max( PosCBalMass, FolMassMin));
            if (FolMassNew < FolMass)
            {
                FolLitSum = FolMass - FolMassNew;
            }
 
            if ( FolMassNew > 0.000001)
            {
                FolMass = std::min( FolMass, FolMassNew);
            }
            else
            {
                FolMass = 0.0;
            }
        }
 
        // This function accounts for foliage not lost during the year because carbon exchange was never negative
        if ( lclock_ref().yearday() == year_end_)
        {
            FolLitSum = std::max( 0.0, FolMass - FolMassMin);
            FolMass   = std::min( FolMass, FolMassMin);
        }
    }
    else
    {
        KLOGFATAL( "unknown growth phase");
    }
}

Photosyn()

void ldndc::PhysiologyPNET::Photosyn ( MoBiLE_Plant *  p )

private

Photosyn

NOTE: This module is written for photosynthesis responses to monthly average climate values and should be used with care if applied on an daily time step. In PNET-N-DNDC the available light is reduced by: PAR_f = par / (Lai * 0.1 + 1.0); and the light response is changed to: LightEff = (1.0 - pow(2.666, (-Il / HALFSAT[p]))); If the Lai-impact on LightEff had not been limited to Lai 15, the responses would be equal at Lai app. 22. As it is, the photosynthesis in PNET is always larger (app. twice as high).

CHANGES:

• The assumption of equal foliage biomass distribution is replaced by the explicit foliage distribution calculated in the vegetation structure module.
• Also specific leaf area is taken form the annual calculation and distributed into the calculated layers. Thereby replacing the internally parameterized values.

AMAXA and AMAXB are specific for the PnET module. The parameters define the photosynthesis rate in dependence on nitrogen, assuming a linear relationship. AMAXA is a mere hypothetical photosynthesis rate at nitrogen content = 0. AMAXB is the increase of photosynthesis with increase of nitrogen percentage.

For parametrization, there are three different sources that have defined general values for broadleaved trees. These are Aber et al. 1995 (AMAXA: -46, AMAXB: 71.9, defined from investigations on Quercus rubra). Reich et al. 1995 (AMAXA 0.31*; AMAXB 5.45, defined from oak, maple and birch leaves). Kattge et al. 2009 (AMAXA: 5.73; AMAXB: 29.81 * 10.0 / SLAMIN, defined from a literature review with many different species).

The general practice now is to take the AMAXA from the newest and most elaborated dataset (Kattge et al. = 5.73) and adjust AMAXB according to the specific leaf area value at the top of the canopy. As long as there is no indication of SLAMIN, a standard of 18.1 is used (resulting from the assumption SLAMIN = 16.5).

{
    double Amax(0.0);                     // Foliage nitrogen-related photosynthesis rate (umol m-2 s-1)
    double BaseFolResp(0.0);
    double GrossAmax(0.0);
    double LightEff(0.0);
    double LayerResp(0.0);
    double DTemp(0.0);
    double Il(0.0);                       // Light extinction coefficient (-)
    double LayerGrossPsnRate(0.0);
    double LayerNetPsn(0.0);
    double LayerGrossPsn_i(0.0);
    double fFol_fl(0.0);
    double laiEst(0.0);                   // estimate of LAI based on evenly distributed foliage biomass
    double sla(0.0);                      // Specific foliage area (m2 kg-1)
    double FolMass_fl(0.0);               // total foliage in one layer (gDW m-2)
 
    double ten9 = 1000000000.0;
 
    double const Tnight( (mc_.nd_minimumairtemperature + mc_.nd_airtemperature) * 0.5);
    double const Tday( (mc_.nd_maximumairtemperature + mc_.nd_airtemperature) * 0.5);
 
 
    // relative temperature impact on photosynthesis, integrated over time
    DTemp = (((*p)->PSNTMAX() - Tday) * (Tday - (*p)->PSNTMIN())) / ( pow((((*p)->PSNTMAX() - (*p)->PSNTMIN()) / 2.0), 2));
 
    // considering frost effects empirically
    if (mc_.nd_minimumairtemperature < 6.0 && DTemp > 0.0 && p->growing_degree_days >= (*p)->GDDFOLEND())
    {
        DTemp *= (1.0 - ((6.0 - mc_.nd_minimumairtemperature) / 6.0) * (1.0 / 30.0));
    }
 
    double const DayLength   = ldndc::meteo::daylength( latitude, lclock_ref().yearday()) * cbm::SEC_IN_HR;
    double const NightLength = cbm::SEC_IN_DAY - DayLength;
 
    DTemp = std::max( 0.0, DTemp);
 
    // relative vapour pressure deficit impact on photosynthesis
    Dvpd = 1.0 - (*p)->DVPD1() * ( pow( mc_.nd_watervaporsaturationdeficit,(*p)->DVPD2()));
 
    // Photosynthesis response to foliage N
    Amax = std::max( 0.0, (*p)->AMAXA() + (*p)->AMAXB() * FolNCon);
 
    // respiration in relation to photosynthesis
    BaseFolResp = (*p)->BASEFOLRESPFRAC() * Amax;
 
    // Gross photosynthesis
    GrossAmax = (Amax * (*p)->AMAXFRAC() + BaseFolResp) * Dvpd * DTemp * DayLength * 12.0 / ten9;
    if (GrossAmax < 0.0) GrossAmax = 0.0;
    GrossAmax *= fFrost;
 
    // Light respiration
    DayResp = (BaseFolResp * (pow((*p)->RESPQ10(),((Tday - (*p)->PSNTOPT()) / 10.0)))
              * DayLength * 12.0) / ten9;
 
    // Dark respiration
    NightResp = (BaseFolResp * (pow((*p)->RESPQ10(),((Tnight - (*p)->PSNTOPT()) / 10.0)))
                * NightLength * 12.0) / ten9;
 
    if ( FolMass > 0.0)
    {
        PosCBalMass = FolMass;
 
        // initializations
        laiEst = 0.0;
        double  fiCum( 0.0);
        CanopyNetPsn = 0.0;
        CanopyGrossPsn = 0.0;
        size_t  flMin( 0);
 
        size_t  fl_cnt( m_veg->canopy_layers_used());
        for ( int  fl = (int)(fl_cnt-1);  fl != 0;  --fl)
        {
            LayerGrossPsn[fl] = 0.0;
 
            FolMass_fl = 0.0;   // rg 25.11.10: prevent array mismatch by using FolMass for scaling
            for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
            {
                FolMass_fl += p->fFol_fl[fl] * FolMass / double( m_veg->size());
            }
            if ( FolMass_fl > 0.0)
            {
                flMin = fl;
            }
        }
 
        // number of steps calculated within each canopy layer
        size_t  ifl( cbm::bound_min( (size_t)1, (size_t)( 50.0 / double( fl_cnt - flMin))));
 
 
        // calculations per canopy layer
        for ( int  fl = (int)(fl_cnt-1);  fl >= (int)flMin;  --fl)
        {
            // foliage mass and specific leaf area for each layer
            /* using indicators from outside the module instead equal mass distribution */
            fFol_fl = 0.0; sla = 0.0;
            for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
            {
                fFol_fl += p->fFol_fl[fl];
                sla     += p->sla_fl[fl] / double( m_veg->size());
            }
 
            double  fi(( fFol_fl * FolMass / double( m_veg->size())) / double(ifl));
 
            // photosysnthesis and respiration in each layer fraction
            for ( size_t  i = 0;  i < ifl;  i++)
            {
                // cumlative leaf biomass and area
                fiCum  += fi;
                laiEst += (fi * sla * cbm::KG_IN_G);
 
                // photosynthesis
 
                // photosynthetic active radiation (umol PAR m-2 s-1)
                double const par( mc_.nd_shortwaveradiation_in * cbm::FPAR *
                                  cbm::UMOL_IN_W * cbm::HR_IN_DAY
                                 / ldndc::meteo::daylength( latitude, lclock_ref().yearday()));
                Il                 = par * exp(-(*p)->EXT() * laiEst);
                LightEff           = std::max(0.0, (1.0 - exp(-Il * log(2.0) / (*p)->HALFSAT())));
                LayerGrossPsnRate  = GrossAmax * LightEff * DWater;
                LayerGrossPsn_i    = LayerGrossPsnRate * fi;
                LayerGrossPsn[fl] += LayerGrossPsn_i;
                LayerResp          = (DayResp + NightResp) * fi;
                LayerNetPsn        = LayerGrossPsn_i - LayerResp;
                CanopyNetPsn      += LayerNetPsn;
                CanopyGrossPsn    += LayerGrossPsn_i;
 
                if ((LayerNetPsn <= 0.0) && cbm::flt_equal( PosCBalMass, FolMass))
                {
                    PosCBalMass = fiCum;
                }
 
                if (LightEff < LightEffMin)
                {
                    LightEffMin = LightEff;
                }
            }
        }
 
        if ( DTemp > 0.0  &&  p->growing_degree_days > (*p)->GDDFOLEND()  &&  lclock_ref().yearday() < (*p)->SENESCSTART()+year_start_-1)
        {
            PosCBalMassTot += ( PosCBalMass);
            PosCBalMassIx  += 1.0;
        }
    } // End if FolMass > 0.0
    else
    {
        PosCBalMass = 0.0;
        CanopyNetPsn    = 0.0;
        CanopyGrossPsn  = 0.0;
        DayResp         = 0.0;
        NightResp       = 0.0;
    }
 
    CanopyGrossPsnMo = CanopyGrossPsn;
    NetPsnMo = (CanopyGrossPsnMo - (DayResp + NightResp) * FolMass);
}

sleep()

lerr_t ldndc::PhysiologyPNET::sleep ( )

inline

Note

it is sadly wrong to assume that you can put me to sleep during times when multiple species exist: in general, i cannot handle the resulting slot layout. if you do, i also have no way of knowing .. so beware!

{ return  LDNDC_ERR_OK; }

solve()

lerr_t ldndc::PhysiologyPNET::solve

(

)

Note

species loops are in general incorrect: however, we only allow one single species with this module. if the species is not alive we return immediately.

@fixme The following two variables need to re-initialized after land use change. However, currently only arrays can be re-initialized in the treedyn module. Therefore, pnet can not be run for re-planting within the year after harvest (possible work-around: harvest at 31.12. of the year).

{
    if ( m_veg->size() == 0)
    {
        return  LDNDC_ERR_OK;
    }
    else if ( m_veg->size() > 1)
    {
        KLOGERROR( "can't handle more than one species :(  bye bye.");
        return  LDNDC_ERR_FAIL;
    }
 
    step_init();
 
    for ( TreeIterator w = m_veg->groupbegin< species::wood >();
            w != m_veg->groupend< species::wood >(); ++w)
    {
        MoBiLE_Plant *p = *w;
 
        // variables init
        FolLitSum = 0.0;        // litterfall (gDW m-2)
 
        tra = 0.0;              // transpiration during a single time step (mm)
        uptNAir = 0.0;          // uptake of nitrogen by canopy processes (kg N m-2)
 
        Dvpd = 0.0;             // Gradient on vapour pressure deficit (?)
        DayResp = 0.0;          // Respiration during daytime (?)
        NightResp = 0.0;        // Respiration during nighttime (?)
        CanopyNetPsn = 0.0;     // Canopy net photosynthesis (?)
        CanopyGrossPsn = 0.0;   // Canopy gross photosynthesis (?)
        FolProdCMo = 0.0;       // foliage carbon production over the integration timestep (day or month) (?)
        FolGRespMo = 0.0;       // foliage growth respiration over the integration timestep (day or month) (?)
        NetPsnMo = 0.0;         // net photosynthesis over the integration timestep (day or month) (?)
        CanopyGrossPsnMo = 0.0; // canopy gross photosynthesis over the integration timestep (day or month) (?)
        RootCAddMo = 0.0;       // fine root growth over the integration timestep (day or month)
 
        // - allocation and mass variables
        qwodfolact = 0.0;       // relation between wood and foliage biomass at the start of the year
        transfDWbs = 0.0;       // annual mass transfer from reserves (buds) to (sap)wood (gC m-2)
        transfDWsb = 0.0;       // annual mass transfer from (sap)wood to reserves (buds) (gC m-2)
        WoodC = 0.0;            // living wood carbon (gC m-2)
        BudC = 0.0;             // reserve carbon for future foliage (gC m-2)
        RootC = 0.0;            // fine root carbon (gC m-2)
        FolMass = 0.0;          // foliage biomass (gDW m-2)
        FolMassMax = 0.0;       // maximum summer foliar biomass (gDW m-2) (not really a parameter)
        FolMassMin = 0.0;       // minimum summer foliar biomass (gDW m-2) (not really a parameter)
        FolReten = 0.0;         // foliage retention time (yr) (not really a parameter)
 
        WoodMRespMo = 0.0;      // Daily maintanance respiration of the wood
        WoodGRespMo = 0.0;      // Daily growth respiration of the wood
        RootMRespMo = 0.0;      // Daily maintanance root respiration
        RootGRespMo = 0.0;      // Daily root growth respiration
        WoodProdCMo = 0.0;      // Daily wood production
        RootProdCMo = 0.0;      // Daily carbon loss form the root system other than respiration;
 
        for ( size_t  fl = 0;  fl < m_setup.canopylayers();  ++fl)
        {
            LayerGrossPsn[fl] = 0.0;
        }
 
        // other variables
        fFrost = 1.0;
        mature = false;
 
        // - at the beginning of each year
        YearlyReset( p);
 
        static double const scale( cbm::G_IN_KG * cbm::CCDM);
 
        // - at the beginning of each year
        if (( lclock_ref().cycles() == 0u) || ( lclock_ref().yearday() == year_start_))
        {
            LightEffMin    = 1.0;
            PosCBalMass    = 0.0;
            PosCBalMassTot = 0.0;
            PosCBalMassIx  = 0.0;
 
            // PnET specific re-initialization of sapwood
            double const mSapIni( p->mSap);
            double const mCorIni( p->mCor);
            p->mSap  = std::min( mSapIni + mCorIni, w->QWODFOLMIN() * ( p->mFol + p->mBud));
            p->mCor += (mSapIni - p->mSap);
 
            if ((p->mFol + p->mBud) > 0.0)
            {
                qwodfolact = p->mSap / (p->mFol + p->mBud);
            }
            else
            {
                qwodfolact = 0.0;
            }
 
            if (p->mCor > 0.0)
            {
                p->ncCor = (mCorIni * p->ncCor + (mSapIni - p->mSap) * p->ncSap) / p->mCor;
            }
            else
            {
                p->ncCor = 0.0;
            }
 
            // carbon pools in PnET initialized from biomass in different compartments
            WoodCStart = (qwodfolact * p->mBud * scale); // virtual carbon pool for this year wood growth
            BudCStart  = (p->mBudStart * scale);         // carbon for this year foliage growth
        }
 
        // transformation of units
        RootC      = (p->mFrt * scale);
        WoodC      = (p->mSap * scale);
        BudC       = (p->mBud * scale);
        FolMass    = (p->mFol * cbm::G_IN_KG);
        FolMassMax = (p->mFolMax * cbm::G_IN_KG);
        FolMassMin = (p->mFolMin * cbm::G_IN_KG);
 
        if (( FolMassMax - FolMassMin) > 0.0)
        {
            FolReten = FolMassMax / (FolMassMax - FolMassMin);
        }
 
        balance_check(p, 1);
 
        // variables from the start of the time step
        mFolOld_vt[p->slot] = p->mFol;
        mFrtOld_vt[p->slot] = p->mFrt;
        mBudOld_vt[p->slot] = p->mBud;
        mSapOld_vt[p->slot] = p->mSap;
        mCorOld_vt[p->slot] = p->mCor;
 
        // consideration of nitrogen supply for photosynthesis
        update_nitrogen_concentration( p);
 
        transfDWbs  = 0.0;
        transfDWsb  = 0.0;
        mature = true;
 
        // average relative water availability (drought) during an integration step (mm)
        DroughtStress  droughtstress;
        droughtstress.fh2o_ref = w->H2OREF_A();
        (*w)->f_h2o = m_state->get_fh2o( *w, m_sipar->CP_WP(), (*w)->f_area);
        DWater = droughtstress.linear_threshold( (*w)->f_h2o);
 
        Phenology( p, 1, DWater);
 
        Photosyn( p);
 
        AllocateMo( p);
 
        Phenology( p, 2, DWater);
 
        if ( lclock_ref().yearday() == year_end_)
        {
            AllocateYr( p);
            PosCBalMass = FolMassMin * (*w)->f_area;
        }
 
        TransformPhysiology( p);
 
        double uptNSoil( 0.0);
        N_uptake( p, uptNSoil);
 
        Update( p, uptNSoil, uptNAir);
 
        wc_.accumulated_potentialtranspiration += potential_wood_transpiration(
                                                                               sl_,
                                                                               mc_.nd_watervaporsaturationdeficit,
                                                                               cbm::sum( p->d_carbonuptake_fl, p->nb_foliagelayers()),
                                                                               p->f_area,
                                                                               w->WUECMAX(),
                                                                               w->WUECMIN(),
                                                                               sc_.h_sl,
                                                                               wc_.wc_sl,
                                                                               wc_.wc_wp_sl,
                                                                               wc_.wc_fc_sl);
 
        balance_check( p, 1);
    }
 
    lerr_t rc_vs = CalcVegStructure();
    if ( rc_vs)
    {
        KLOGERROR("Simulation of vegetation structure not successful");
        return  LDNDC_ERR_FAIL;
    }
 
    step_out();
 
    return  LDNDC_ERR_OK;
}

References ldndc::potential_wood_transpiration().

Here is the call graph for this function:

TransformPhysiology()

void ldndc::PhysiologyPNET::TransformPhysiology ( MoBiLE_Plant *  p )

private

TransformPhysiology

{
    double mSum(0.0),maintResp(0.0),fResp(0.0);
    double dcTot(0.0),dcSapC(0.0),dcBudC(0.0),dcFrtC(0.0),dBudY(0.0),dSapY(0.0);
    double fm_vt(0.0),fwood(0.0);
 
    // distribution of net carbon gain into wood and foliage reserves for next year
    /* root growth is already calculated as a fraction of gross gain, so it has to
       be substracted first */
    dcTot = (CanopyGrossPsnMo
             - (CanopyGrossPsnMo - NetPsnMo) - WoodMRespMo - RootMRespMo
             - FolGRespMo - WoodGRespMo - RootGRespMo);
 
    // Reduction of respiration if total growth is negative
    /* This has been newly introduced to account for the large biomass increases
       particular in Eucalyptus plantation under severe drought */
//    double massLoss = 0.0;
    if (dcTot < 0.0)
    {
        maintResp = WoodMRespMo + RootMRespMo + (CanopyGrossPsnMo - NetPsnMo);
        if (maintResp > 0.0)
        {
            fResp = (maintResp + dcTot) / maintResp;
        }
 
        if (fResp < 0.0)
        {
            fResp     = 0.0;
//            massLoss -= dcTot;
        }
 
        WoodMRespMo *= fResp;
        RootMRespMo *= fResp;
 
        NetPsnMo = CanopyGrossPsnMo - fResp * (CanopyGrossPsnMo - NetPsnMo);
        dcTot = (CanopyGrossPsnMo
                 - (CanopyGrossPsnMo - NetPsnMo) - WoodMRespMo - RootMRespMo
                 - FolGRespMo - WoodGRespMo - RootGRespMo);
 
        // if we have negative growth carbon is taken from sap wood
        // here we make sure that negative carbon growth is not exceeding available carbon
        if (dcTot < 0.0 && std::abs(dcTot) > p->mSap)
        {
            FolGRespMo    = 0.0;
            WoodMRespMo    = 0.0;
            WoodGRespMo    = 0.0;
            RootMRespMo    = 0.0;
            RootGRespMo    = 0.0;
 
            dcTot = (CanopyGrossPsnMo
                     - (CanopyGrossPsnMo - NetPsnMo) - WoodMRespMo - RootMRespMo
                     - FolGRespMo - WoodGRespMo - RootGRespMo);
        }
    }
 
    // wood and reserve change
    /* If growth is negative, it is taken only from free available tissue. This is
       assumed to be the whole labile reserve pool and 20 percent of woody tissue. */
    if (dcTot > 0.0 && p->dvsFlushOld >= 1.0)    // positive gain, after flushing
    {
        fwood = (*p)->QWODFOLMIN() / (1.0 + (*p)->QWODFOLMIN());
    }
    else if (BudC + WoodC * 0.2 > 0.0 && p->dvsFlushOld >= 1.0)    // negative gain, after flushing
    {
        fwood = WoodC * 0.2 / (BudC + WoodC * 0.2);
    }
    else    // before flushing
    {
        fwood = 1.0;
    }
 
    dcFrtC = std::max(0.0, std::min(dcTot, RootCAddMo));  // root growth has priority, but is limited by availability
    dcSapC = (dcTot-dcFrtC) * fwood;
    dcBudC = (dcTot-dcFrtC) * (1.0 - fwood);
 
    mSum = 0.0;
//    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//    {
        mSum += (p->mFol + p->mFrt + p->mSap + p->mBud);
        p->mFolMax = 0.0;
        p->mFolMin = 0.0;
//    }
 
//    for ( PlantIterator vt = m_veg->begin(); vt != m_veg->end(); ++vt)
//    {
        // biomass fraction (and unit conversion)
        if (mSum > 0.0)
        {
            fm_vt = ((p->mFol + p->mFrt + p->mSap + p->mBud) / mSum);
 
            if (fm_vt > 1.0)
            {
                fm_vt = 1.0;
            }
 
            fm_vt *= cbm::KG_IN_G;
        }
        else
        {
            fm_vt = 0.0;
        }
 
        // growth
        p->d_dcFol = FolProdCMo * fm_vt; // shift from buds to foliage
        mBudLoss_vt[p->slot]     = p->d_dcFol / cbm::CCDM;
        p->d_dcFrt = dcFrtC * fm_vt;     // growth from assimilated carbon
        p->d_dcSap = dcSapC * fm_vt;
        p->d_dcBud = dcBudC * fm_vt;
 
        // respiration
        double const sap_wood_frac( sap_wood_fraction( p->lai_pot, p->qsfa, p->height_at_canopy_start,
                                                       p->height_max, p->rooting_depth));
        double const dc_sum(  std::max( 0.0, p->d_dcFol)
                            + std::max( 0.0, p->d_dcFrt)
                            + std::max( 0.0, p->d_dcSap)
                            + std::max( 0.0, p->d_dcBud));
        double const dc_below( std::max( 0.0, p->d_dcFrt) + std::max( 0.0, p->d_dcSap * sap_wood_frac));
 
        p->d_rFol    = (CanopyGrossPsnMo - NetPsnMo) * fm_vt;
        p->d_rFrt    = RootMRespMo * fm_vt;
        p->d_rSap    = WoodMRespMo * fm_vt;
        p->d_rSapBelow = (p->d_rSap * sap_wood_frac);
        p->d_rBud    = 0.0;
        p->d_rRes    = p->d_rFol + p->d_rFrt + p->d_rSap + p->d_rBud;
        p->d_rTra    = 0.0;
        p->d_rGro    = (WoodGRespMo + RootGRespMo + FolGRespMo) * fm_vt;
        p->d_rGroBelow = ((dc_sum > 0.0) ? p->d_rGro * dc_below / dc_sum : 0.0);
 
        // senescence
        p->d_sFol = FolLitSum * fm_vt;
        for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
        {
            p->d_sFrt_sl[sl] += RootProdCMo / cbm::CCDM * fm_vt * p->fFrt_sl[sl];
        }
 
        mSapLoss_vt[p->slot] = 0.0;
        if ((( p->mSap * 0.99999)) > ((p->mFol + p->mBud) * (*p)->QWODFOLMIN()))
        {
            mSapLoss_vt[p->slot] = std::min( SRATEMAX * p->mSap, p->mSap - (p->mFol + p->mBud) * (*p)->QWODFOLMIN());
            if (mSapLoss_vt[p->slot] < 0.0)
            {
                mSapLoss_vt[p->slot] = 0.0;
            }
        }
        p->d_sWoodAbove = mSapLoss_vt[p->slot] * p->f_branch;
 
        // senescence due to annual shift between biomasses
        dBudY = transfDWsb * fm_vt;
        dSapY = transfDWbs * fm_vt;
 
        // rg 19.11.10 security
        if (dBudY < 0.0)
        {
            dBudY = std::max(-1 * (p->mCor + mSapLoss_vt[p->slot] * (1.0-p->f_branch)), transfDWsb * fm_vt);
        }
        if (dSapY < 0.0)
        {
            dSapY = std::max(-1 * (p->mSap + (p->d_dcSap / cbm::CCDM - mSapLoss_vt[p->slot])), transfDWbs * fm_vt);
        }
 
        // biomass update
        p->mFol  = FolMass * fm_vt;
        p->mFrt += ((p->d_dcFrt / cbm::CCDM - cbm::sum( p->d_sFrt_sl, sl_.soil_layer_cnt())));
        p->mSap += ((p->d_dcSap / cbm::CCDM - mSapLoss_vt[p->slot]) - dBudY);
        p->mBud += ((p->d_dcBud / cbm::CCDM - mBudLoss_vt[p->slot]) + dBudY - dSapY);
        p->mCor += (mSapLoss_vt[p->slot] * (1.0-p->f_branch) + dSapY);
 
        // foliage biomass in age classes     // rg 08.06.10 (new)
        p->mFol_na[0] += p->d_dcFol;
        if ( get_na_max_pnet() > 0)
        {
            p->mFol_na[get_na_max_pnet()-1] -= p->d_sFol;
            for ( size_t  na = 0;  na < get_na_max_pnet();  ++na)
            {
                if ( p->mFol_na[na] < 0.0)
                {
                    p->mFol_na[na] = 0.0;
                }
            }
        }
 
        // control values for next year
        p->mBudStart = p->mBud;
        p->mFolMax  += (FolMassMax * fm_vt);
        p->mFolMin  += (FolMassMin * fm_vt);
 
        /* Note that if one of the following equations applies, the carbon balance is not closed any more */
        if ( p->mFrt < 0.0)
        {
            if (( p->mFrt > 1e-7) || (-1 * p->mFrt > 1e-7))
            {
                KLOGWARN(" mFrt-C-leak in PNET  [day=",lclock_ref().now().to_string(),",mFrt=",p->mFrt,"]");
            }
            p->mFrt = 0.0;
        }
 
        if ( p->mSap < 0.0)
        {
            if (p->mSap > 1e-7 || -1 * p->mSap > 1e-7)
            {
                KLOGWARN(" mSap-C-leak in PNET  [day=",lclock_ref().now().to_string(),",mSap=",p->mSap,"]");
            }
            p->mSap = 0.0;
        }
 
        if ( p->mBud < 0.0)
        {
            if (p->mBud > 1e-7 || -1 * p->mBud > 1e-7)
            {
                KLOGWARN(" mBud-C-leak in PNET  [day=",lclock_ref().now().to_string(),",mBud=",p->mBud,"]");
            }
            p->mBud   = 0.0;
        }
 
        if (( p->growing_degree_days >= (*p)->GDDFOLSTART()) && ( p->dEmerg == 0))
        {
            p->dEmerg = lclock_ref().yearday();
        }
 
 
        if ( mature
            && ( p->f_area > 0.0)
            && (lclock_ref().yearday() > (unsigned int)((*p)->SENESCSTART()+year_start_-1)
                || lclock_ref().yearday() < year_start_)
            && lclock_ref().yearday() != lclock_ref().days_in_year())
        {
            // rg 08.06.10 (no senescence at the last day of the year)
            p->dvsMort = 1.0 - (FolMass / p->f_area - FolMassMin) / ( std::max(FolMass / p->f_area, FolMassMax) - FolMassMin);
        }
        else
        {
            p->dvsMort = 0.0;
        }
 
 
        if (p->dvsMort < 0.0)
        {
            p->dvsMort = 0.0;
        }
        else if (p->dvsMort > 1.0)
        {
            p->dvsMort = 1.0;
        }
 
        size_t  fl_cnt( p->nb_foliagelayers());
        m_pf.update_foliage_layer_biomass_and_lai( p, fl_cnt, 0.001);
 
        // fraction of photosynthesis is used with smaller-timestep models
        for (  size_t  fl = 0;  fl < fl_cnt;  ++fl)
        {
            p->d_carbonuptake_fl[fl] = LayerGrossPsn[fl] * cbm::KG_IN_G;
        }
//    }
}

References ldndc::sap_wood_fraction().

Here is the call graph for this function:

Update()

void ldndc::PhysiologyPNET::Update ( MoBiLE_Plant *  _vt, double  _upt1_vt, double  _upt2_vt  )

private

Update

{
    // component demand
    double nFolOpt_vt( m_pf.optimum_nitrogen_content_foliage( _vt));
    double nBudOpt_vt( m_pf.optimum_nitrogen_content_buds( _vt));
 
    // optimum nitrogen contents
    double ncFolOpt( (_vt->mFol > 0.0) ? nFolOpt_vt/_vt->mFol : 0.0);
    double ncBudOpt( (_vt->mBud > 0.0) ? nBudOpt_vt/_vt->mBud : 0.0);
 
 
    // component demand of compartments
    double demFol( std::max(0.0, (ncFolOpt - _vt->ncFol) * mFolOld_vt[_vt]
                            + ncFolOpt * (_vt->mFol - mFolOld_vt[_vt])));
    double demFrt( std::max(0.0, ((*_vt)->NC_FINEROOTS_MAX() - _vt->ncFrt) * mFrtOld_vt[_vt]
                            + (*_vt)->NC_FINEROOTS_MAX() * (_vt->mFrt - mFrtOld_vt[_vt])));
    double demSap( std::max(0.0, ((*_vt)->NCSAPOPT() - _vt->ncSap) * mSapOld_vt[_vt]
                            + (*_vt)->NCSAPOPT() * (_vt->mSap - mSapOld_vt[_vt])));
    double demBud( std::max(0.0, (ncBudOpt - _vt->ncBud) * mBudOld_vt[_vt]
                            + ncBudOpt * (_vt->mBud - mBudOld_vt[_vt])));
    double demTot = (demFol + demFrt + demSap + demBud);
 
 
    // accounting for biomass shifts at the end of a year
    double transfBud( (transfDWbs > 0.0) ?
                     std::max(0.0,
                              (mBudOld_vt[_vt] - mBudLoss_vt[_vt->slot] + _vt->d_dcBud/cbm::CCDM
                               - _vt->mBud) * _vt->ncBud)
                     : 0.0);
 
    double transfSap( (transfDWsb > 0.0) ?
                     std::max(0.0,
                              (mSapOld_vt[_vt] - mSapLoss_vt[_vt->slot] + _vt->d_dcSap/cbm::CCDM
                               - _vt->mSap) * _vt->ncSap)
                     : 0.0);
 
 
    // element losses // rg 08.12.10 simplified
    // the minimum function is needed to guarantee mass conservation.
    // there is no more nitrogen loss allowed than present at the beginning of the time step
    double xlossFol( std::min(_vt->d_sFol, mFolOld_vt[_vt]) * _vt->ncFol);
    double xlossFrt( std::min(cbm::sum(_vt->d_sFrt_sl, sl_.soil_layer_cnt()), mFrtOld_vt[_vt]) * _vt->ncFrt);
    double xlossBud( std::min(mBudLoss_vt[_vt->slot], mBudOld_vt[_vt]) * _vt->ncBud);
    double xlossSap( std::min(mSapLoss_vt[_vt->slot], mSapOld_vt[_vt]) * _vt->ncSap);
 
 
    // nitrogen uptake
    double xdist_up( _upt1_vt + _upt2_vt);
 
 
    // nitrogen retention
    // note: elements are principally neither retranslocated from sapwood
    // when turned into corewood nor from fine roots (Gordon and Jackson 2000)
    _vt->d_n_retention = ( xlossFol * (*_vt)->FRET_N());
 
 
    // total nitrogen availability for distribution
    double xdistTot( xdist_up + _vt->d_n_retention);
 
 
    // correct nitrogen retention in case exceeding demand
    if (demTot < xdistTot)
    {
        // reduced nitrogen retention
        double const xret_reduced( std::max(demTot - xdist_up, 0.0));
        _vt->d_n_retention = xret_reduced;
        xdistTot = (xdist_up + xret_reduced);
    }
 
 
    // components from uptake and retranslocation to be distributed between plant tissues
    // by internal translocation from sapwood to corewood
    double xdistCor( xlossSap * (1.0 - _vt->f_branch));
 
 
    // by internal translocation from buds to foliage
    double xlitBud( 0.0);
    double xfolBud( xlossBud);
    if (xlossBud < demFol)
    {
        demFol -= xlossBud;
        demTot -= xlossBud;
    }
    else
    {
        xlitBud = xlossBud;
        xfolBud = 0.0;
    }
 
 
    // nitrogen distribution to components
    double xdistFol( 0.0);
    double xdistFrt( 0.0);
    double xdistSap( 0.0);
    double xdistBud( 0.0);
    double xrest( 0.0);
 
    if (demTot > xdistTot)    // any available element is taken up
    {
        xdistFol = xdistTot * demFol / demTot;
        xdistFrt = xdistTot * demFrt / demTot;
        xdistSap = xdistTot * demSap / demTot;
        xdistBud = xdistTot * demBud / demTot;
        xrest    = 0.0;
    }
    else if (xdistTot > 0.0)    // every element demand is saturated
    {
        xdistFol = demFol;
        xdistFrt = demFrt;
        xdistSap = demSap;
        xdistBud = demBud;
        xrest    = xdistTot - demTot;
    }
 
 
    // not needed elements are distributed to litter
    if ( xrest > 0.0)
    {
        double const m_tot( mFolOld_vt[_vt]+mBudOld_vt[_vt]+mFrtOld_vt[_vt]+mSapOld_vt[_vt]);
        if (cbm::flt_equal_zero( m_tot))
        {
            KLOGWARN("No more living biomass for species ", _vt->name());
        }
        else
        {
            _vt->d_nLitFol += (xrest * mFolOld_vt[_vt] / m_tot);
            _vt->d_nLitWoodAbove += (xrest * mSapOld_vt[_vt] / m_tot);
            _vt->d_nLitBud += (xrest * mBudOld_vt[_vt] / m_tot);
            for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
            {
                _vt->d_nLitFrt_sl[sl] += (xrest * mFrtOld_vt[_vt] / m_tot * _vt->fFrt_sl[sl]);
            }
        }
    }
 
    _vt->d_nLitBud += xlitBud;
    _vt->d_nLitFol += (xlossFol - _vt->d_n_retention);
    _vt->d_nLitWoodAbove += (xlossSap * _vt->f_branch);
    for ( size_t  sl = 0;  sl < sl_.soil_layer_cnt();  ++sl)
    {
        _vt->d_nLitFrt_sl[sl] += xlossFrt * _vt->fFrt_sl[sl];
    }
 
    // new concentrations
    if (_vt->mFol > 1.0e-9)
    {
        _vt->ncFol = std::max(0., mFolOld_vt[_vt] * _vt->ncFol + xdistFol + xfolBud - xlossFol) / _vt->mFol;
    }
    else
    {
        _vt->ncFol = 0.0;
    }
 
    if (_vt->mBud > 1.0e-9)
    {
        _vt->ncBud = std::max(0., mBudOld_vt[_vt] * _vt->ncBud + xdistBud - xlossBud - transfBud + transfSap) / _vt->mBud;
    }
    else
    {
        _vt->ncBud = 0.0;
    }
 
    if (_vt->mFrt > 1.0e-9)
    {
        _vt->ncFrt = std::max(0., mFrtOld_vt[_vt] * _vt->ncFrt + xdistFrt - xlossFrt) / _vt->mFrt;
    }
    else
    {
        _vt->ncFrt = 0.0;
    }
 
    if (_vt->mCor > 1.0e-9)
    {
        _vt->ncCor = std::max(0., mCorOld_vt[_vt] * _vt->ncCor + xdistCor + transfBud) / _vt->mCor;
    }
    else
    {
        _vt->ncCor = 0.0;
    }
 
    if (_vt->mSap > 1.0e-9)
    {
        _vt->ncSap = std::max(0., mSapOld_vt[_vt] * _vt->ncSap + xdistSap - xlossSap - transfSap) / _vt->mSap;
        /*
         // additional nitrogen allocation from core wood to reach minimum nitrogen concentration
         if (_vt->ncSap < cMin && _vt->mSap > 0.0 && _vt->mCor > 0.0){                   // rg 11.06.10 (avoid divide by zero problem)
         transfer = Min(_vt->mCor * _vt->ncCor, _vt->mSap * (cMin - _vt->ncSap));
         _vt->ncSap  = (_vt->mSap * _vt->ncSap + transfer) / _vt->mSap;
         _vt->ncCor  = (_vt->mCor * _vt->ncCor - transfer) / _vt->mCor;
         }
         */
    }
    else
    {
        _vt->ncSap = 0.0;
    }
}

update_nitrogen_concentration()

void ldndc::PhysiologyPNET::update_nitrogen_concentration ( MoBiLE_Plant *  p )

private

update_nitrogen_concentration

{
    double const nFolOpt( m_pf.optimum_nitrogen_content_foliage( p));
 
    FolNCon = 0.0;
    double FolNConOpt( 0.0);
    if (p->mFol > 0.0)
    {
        FolNCon    += (p->ncFol * cbm::PRO_IN_GG * p->mFol / p->mFol);
        FolNConOpt += (nFolOpt / p->mFol * cbm::PRO_IN_GG * p->mFol / p->mFol);
    }
 
    FolNCon = std::min( FolNCon, FolNConOpt);
}

YearlyReset()

void ldndc::PhysiologyPNET::YearlyReset ( MoBiLE_Plant *  p )

private

Resets:

• Growing degree days

{
    if ( lclock_ref().yearday() !=  year_start_)
    {
        return;
    }
 
    p->dEmerg = 0;
    p->growing_degree_days = 0.0;
    p->dvsFlush = 0.0;
    p->dvsFlushOld = 0.0;
    p->dvsWood = 0.0;
 
    /* rg 24.02.11
     * (re-initialization of biomass according foliage longvity based on the expected maximum foliage mFolMax)
     * rg 20.01.12
     * (shift from calculating only once to calculating every year start
     * (account for new initialization after regeneration)
     * rg 28.03.12
     * (considering new nitrogen concentrations after adjustment of biomasses, including an adjustment factor)
     */
    if ( get_na_max_pnet() > 1)
    {
        double mFolNew  = p->mFolMax * (1.0 - 1.0 / double(get_na_max_pnet()));
        double mBudNew  = mFolNew / double(get_na_max_pnet()-1);//assuming no flush of foliage
        double mSapNew  = p->mSap - (mFolNew + mBudNew - p->mFol - p->mBud);
        double dnfol    = (p->mFol - mFolNew) * p->ncFol;
        double ncFolNew = (p->mFol * p->ncFol - dnfol) / mFolNew;
        double ncBudNew = (p->mBud * p->ncBud + dnfol) / mBudNew;
        double ncSapNew = (p->mSap * p->ncSap) / mSapNew;
 
        p->ncFol     = ncFolNew;
        p->ncBud     = ncBudNew;
        p->ncSap     = ncSapNew;
        p->mFol      = mFolNew;
        p->mBud      = mBudNew;
        p->mSap      = mSapNew;
        p->mBudStart = p->mBud;
        p->mFolMin   = p->mFol;
    }
 
    double const  mfol_va( p->mFolMax - p->mFolMin);
    p->mFol_na[0] = 0.0;
    for ( size_t  na = 1;  na < get_na_max_pnet();  ++na)
    {
        p->mFol_na[na] = mfol_va;
    }
}

Member Data Documentation

C

double const ldndc::PhysiologyPNET::C = 0.02

staticprivate

parameter from frost hardiness function

GROWMAX

double const ldndc::PhysiologyPNET::GROWMAX = 0.95

staticprivate

fraction of uncovered area that can be reoccuppied in one year

KM_NH4

double const ldndc::PhysiologyPNET::KM_NH4 = 1.0

staticprivate

??

KM_NO3

double const ldndc::PhysiologyPNET::KM_NO3 = 1.0

staticprivate

??

NUPTAIR

double const ldndc::PhysiologyPNET::NUPTAIR = 0.0

staticprivate

fraction of nitrogen from total deposition that is taken up by the canopy

SRATEMAX

double const ldndc::PhysiologyPNET::SRATEMAX = 0.01

staticprivate

??

TMAX_NUP

double const ldndc::PhysiologyPNET::TMAX_NUP = 25.0

staticprivate

Optimum temperature where maximum N uptake occurs (Gessler et al. 1998)

TMIN_NUP

double const ldndc::PhysiologyPNET::TMIN_NUP = 6.0

staticprivate

Minimum temperature where N uptake is possible (estimated from Gessler et al. 1998, corroborated by Alvarez-Uria and Koerner 2007)