WebSVN - public iLand - Diff - Rev 798 and 802


Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/production3pg.cpp':
/********************************************************************************************
**    iLand - an individual based forest landscape and disturbance model
**    http://iland.boku.ac.at
**    Copyright (C) 2009-  Werner Rammer, Rupert Seidl
**
**    This program is free software: you can redistribute it and/or modify
**    it under the terms of the GNU General Public License as published by
**    the Free Software Foundation, either version 3 of the License, or
**    (at your option) any later version.
**
**    This program is distributed in the hope that it will be useful,
**    but WITHOUT ANY WARRANTY; without even the implied warranty of
**    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
**    GNU General Public License for more details.
**
**    You should have received a copy of the GNU General Public License
**    along with this program.  If not, see <http://www.gnu.org/licenses/>.
********************************************************************************************/

#include "global.h"
#include "production3pg.h"

#include "resourceunit.h"
#include "species.h"
#include "speciesresponse.h"
#include "model.h"

Production3PG::Production3PG()
{
    mResponse=0;
    mEnvYear = 0.;
}

/**
  This is based on the utilizable photosynthetic active radiation.
  @sa http://iland.boku.ac.at/primary+production
  The resulting radiation is MJ/m2       */
inline double Production3PG::calculateUtilizablePAR(const int month) const
{
    // calculate the available radiation. This is done at SpeciesResponse-Level (SpeciesResponse::calculate())
    // see Equation (3)
    // multiplicative approach: responses are averaged one by one and multiplied on a monthly basis
//    double response = mResponse->absorbedRadiation()[month] *
//                      mResponse->vpdResponse()[month] *
//                      mResponse->soilWaterResponse()[month] *
//                      mResponse->tempResponse()[month];
    // minimum approach: for each day the minimum aof vpd, temp, soilwater is calculated, then averaged for each month
    //double response = mResponse->absorbedRadiation()[month] *
    //                  mResponse->minimumResponses()[month];
    double response = mResponse->utilizableRadiation()[month];

    return response;
}
/** calculate the alphac (=photosynthetic efficiency) for the given month.
   this is based on a global efficiency, and modified per species.
   epsilon is in gC/MJ Radiation
  */
inline double Production3PG::calculateEpsilon(const int month) const
{
    double epsilon = Model::settings().epsilon; // maximum radiation use efficiency
    epsilon *= mResponse->nitrogenResponse() *
               mResponse->co2Response()[month];
    return epsilon;
}

inline double Production3PG::abovegroundFraction() const
{
    double utilized_frac = 1.;
    if (Model::settings().usePARFractionBelowGroundAllocation) {
        // the Landsberg & Waring formulation takes into account the fraction of utilizeable to total radiation (but more complicated)
        // we originally used only nitrogen and added the U_utilized/U_radiation
 
        utilized_frac = mResponse->totalUtilizeableRadiation() / mResponse->yearlyRadiation();
    }
    double harsh =  1 - 0.8/(1 + 2.5 * mResponse->nitrogenResponse() * utilized_frac);
    return harsh;
}

void Production3PG::clear()
{
    for (int i=0;i<12;i++) {
        mGPP[i] = 0.; mUPAR[i]=0.;
    }
    mEnvYear = 0.;
}

/** calculate the stand-level NPP
  @ingroup core
  Standlevel (i.e ResourceUnit-level) production (NPP) following the 3PG approach from Landsberg and Waring.
  @sa http://iland.boku.ac.at/primary+production */
double Production3PG::calculate()
{
    Q_ASSERT(mResponse!=0);
    // Radiation: sum over all days of each month with foliage
    double year_raw_gpp = 0.;
    clear();
    double utilizable_rad, epsilon;
    // conversion from gC to kg Biomass: C/Biomass=0.5
    const double gC_to_kg_biomass = 1. / (biomassCFraction * 1000.);
    for (int i=0;i<12;i++) {
 
        utilizable_rad = calculateUtilizablePAR(i); // utilizable radiation of the month ... (MJ/m2)
        epsilon = calculateEpsilon(i); // ... photosynthetic efficiency ... (gC/MJ)
        mUPAR[i] = utilizable_rad ;
        mGPP[i] =utilizable_rad * epsilon * gC_to_kg_biomass; // ... results in GPP of the month kg Biomass/m2 (converted from gC/m2)
        year_raw_gpp += mGPP[i]; // kg Biomass/m2
    }

    // calculate f_env,yr: see http://iland.boku.ac.at/sapling+growth+and+competition
    double f_sum = 0.;
    for (int i=0;i<12;i++)
        f_sum += mGPP[i] / gC_to_kg_biomass; // == uAPar * epsilon_eff

    //  the factor f_ref: parameter that scales response values to the range 0..1 (1 for best growth conditions) (species parameter)
    const double perf_factor = mResponse->species()->saplingGrowthParameters().referenceRatio;
    // f_env,yr=(uapar*epsilon_eff) / (APAR * epsilon_0 * fref)
    mEnvYear = f_sum / (Model::settings().epsilon * mResponse->yearlyRadiation() * perf_factor);
    if (mEnvYear > 1.) {
        qDebug() << "ERROR: fEnvYear > 1 for " << mResponse->species()->id() << mEnvYear << "f_sum, epsilon, yearlyRad, perf_factor" <<  f_sum << Model::settings().epsilon <<  mResponse->yearlyRadiation() << perf_factor;
        mEnvYear = 1.;
    }

    // calculate fraction for belowground biomass
    mRootFraction = 1. - abovegroundFraction();

    // global value set?
    double dbg = GlobalSettings::instance()->settings().paramValue("gpp_per_year",0);
    if (dbg) {
        year_raw_gpp = dbg;
        mRootFraction = 0.4;
    }

    // year GPP/rad: kg Biomass/m2
    mGPPperArea = year_raw_gpp;
    return mGPPperArea; // yearly GPP in kg Biomass/m2
}
 

Subversion Repositories public iLand

(root)/src/core/production3pg.cpp @ 1222 - Rev 798 → 802

Rev 798		Rev 802
1	Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/production3pg.cpp':	1	Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/production3pg.cpp':
2	/********************************************************************************************	2	/********************************************************************************************
3	** iLand - an individual based forest landscape and disturbance model	3	** iLand - an individual based forest landscape and disturbance model
4	** http://iland.boku.ac.at	4	** http://iland.boku.ac.at
5	** Copyright (C) 2009- Werner Rammer, Rupert Seidl	5	** Copyright (C) 2009- Werner Rammer, Rupert Seidl
6	**	6	**
7	** This program is free software: you can redistribute it and/or modify	7	** This program is free software: you can redistribute it and/or modify
8	** it under the terms of the GNU General Public License as published by	8	** it under the terms of the GNU General Public License as published by
9	** the Free Software Foundation, either version 3 of the License, or	9	** the Free Software Foundation, either version 3 of the License, or
10	** (at your option) any later version.	10	** (at your option) any later version.
11	**	11	**
12	** This program is distributed in the hope that it will be useful,	12	** This program is distributed in the hope that it will be useful,
13	** but WITHOUT ANY WARRANTY; without even the implied warranty of	13	** but WITHOUT ANY WARRANTY; without even the implied warranty of
14	** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the	14	** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15	** GNU General Public License for more details.	15	** GNU General Public License for more details.
16	**	16	**
17	** You should have received a copy of the GNU General Public License	17	** You should have received a copy of the GNU General Public License
18	** along with this program. If not, see <http://www.gnu.org/licenses/>.	18	** along with this program. If not, see <http://www.gnu.org/licenses/>.
19	********************************************************************************************/	19	********************************************************************************************/
20		20
21	#include "global.h"	21	#include "global.h"
22	#include "production3pg.h"	22	#include "production3pg.h"
23		23
24	#include "resourceunit.h"	24	#include "resourceunit.h"
25	#include "species.h"	25	#include "species.h"
26	#include "speciesresponse.h"	26	#include "speciesresponse.h"
27	#include "model.h"	27	#include "model.h"
28		28
29	Production3PG::Production3PG()	29	Production3PG::Production3PG()
30	{	30	{
31	mResponse=0;	31	mResponse=0;
32	mEnvYear = 0.;	32	mEnvYear = 0.;
33	}	33	}
34		34
35	/**	35	/**
36	This is based on the utilizable photosynthetic active radiation.	36	This is based on the utilizable photosynthetic active radiation.
37	@sa http://iland.boku.ac.at/primary+production	37	@sa http://iland.boku.ac.at/primary+production
38	The resulting radiation is MJ/m2 */	38	The resulting radiation is MJ/m2 */
39	inline double Production3PG::calculateUtilizablePAR(const int month) const	39	inline double Production3PG::calculateUtilizablePAR(const int month) const
40	{	40	{
41	// calculate the available radiation. This is done at SpeciesResponse-Level (SpeciesResponse::calculate())	41	// calculate the available radiation. This is done at SpeciesResponse-Level (SpeciesResponse::calculate())
42	// see Equation (3)	42	// see Equation (3)
43	// multiplicative approach: responses are averaged one by one and multiplied on a monthly basis	43	// multiplicative approach: responses are averaged one by one and multiplied on a monthly basis
44	// double response = mResponse->absorbedRadiation()[month] *	44	// double response = mResponse->absorbedRadiation()[month] *
45	// mResponse->vpdResponse()[month] *	45	// mResponse->vpdResponse()[month] *
46	// mResponse->soilWaterResponse()[month] *	46	// mResponse->soilWaterResponse()[month] *
47	// mResponse->tempResponse()[month];	47	// mResponse->tempResponse()[month];
48	// minimum approach: for each day the minimum aof vpd, temp, soilwater is calculated, then averaged for each month	48	// minimum approach: for each day the minimum aof vpd, temp, soilwater is calculated, then averaged for each month
49	//double response = mResponse->absorbedRadiation()[month] *	49	//double response = mResponse->absorbedRadiation()[month] *
50	// mResponse->minimumResponses()[month];	50	// mResponse->minimumResponses()[month];
51	double response = mResponse->utilizableRadiation()[month];	51	double response = mResponse->utilizableRadiation()[month];
52		52
53	return response;	53	return response;
54	}	54	}
55	/** calculate the alphac (=photosynthetic efficiency) for the given month.	55	/** calculate the alphac (=photosynthetic efficiency) for the given month.
56	this is based on a global efficiency, and modified per species.	56	this is based on a global efficiency, and modified per species.
57	epsilon is in gC/MJ Radiation	57	epsilon is in gC/MJ Radiation
58	*/	58	*/
59	inline double Production3PG::calculateEpsilon(const int month) const	59	inline double Production3PG::calculateEpsilon(const int month) const
60	{	60	{
61	double epsilon = Model::settings().epsilon; // maximum radiation use efficiency	61	double epsilon = Model::settings().epsilon; // maximum radiation use efficiency
62	epsilon = mResponse->nitrogenResponse()	62	epsilon = mResponse->nitrogenResponse()
63	mResponse->co2Response()[month];	63	mResponse->co2Response()[month];
64	return epsilon;	64	return epsilon;
65	}	65	}
66		66
67	inline double Production3PG::abovegroundFraction() const	67	inline double Production3PG::abovegroundFraction() const
68	{	68	{
69	double utilized_frac = 1.;	69	double utilized_frac = 1.;
70	if (Model::settings().usePARFractionBelowGroundAllocation) {	70	if (Model::settings().usePARFractionBelowGroundAllocation) {
71	// the Landsberg & Waring formulation takes into account the fraction of utilizeable to total radiation (but more complicated)	71	// the Landsberg & Waring formulation takes into account the fraction of utilizeable to total radiation (but more complicated)
72	// we originally used only nitrogen and added the U_utilized/U_radiation	72	// we originally used only nitrogen and added the U_utilized/U_radiation
73	utilized_frac = mResponse->totalUtilizedRadiation() / mResponse->yearlyRadiation();	-
-		73	utilized_frac = mResponse->totalUtilizeableRadiation() / mResponse->yearlyRadiation();
74	}	74	}
75	double harsh = 1 - 0.8/(1 + 2.5 * mResponse->nitrogenResponse() * utilized_frac);	75	double harsh = 1 - 0.8/(1 + 2.5 * mResponse->nitrogenResponse() * utilized_frac);
76	return harsh;	76	return harsh;
77	}	77	}
78		78
79	void Production3PG::clear()	79	void Production3PG::clear()
80	{	80	{
81	for (int i=0;i<12;i++) {	81	for (int i=0;i<12;i++) {
82	mGPP[i] = 0.; mUPAR[i]=0.;	82	mGPP[i] = 0.; mUPAR[i]=0.;
83	}	83	}
84	mEnvYear = 0.;	84	mEnvYear = 0.;
85	}	85	}
86		86
87	/** calculate the stand-level NPP	87	/** calculate the stand-level NPP
88	@ingroup core	88	@ingroup core
89	Standlevel (i.e ResourceUnit-level) production (NPP) following the 3PG approach from Landsberg and Waring.	89	Standlevel (i.e ResourceUnit-level) production (NPP) following the 3PG approach from Landsberg and Waring.
90	@sa http://iland.boku.ac.at/primary+production */	90	@sa http://iland.boku.ac.at/primary+production */
91	double Production3PG::calculate()	91	double Production3PG::calculate()
92	{	92	{
93	Q_ASSERT(mResponse!=0);	93	Q_ASSERT(mResponse!=0);
94	// Radiation: sum over all days of each month with foliage	94	// Radiation: sum over all days of each month with foliage
95	double year_raw_gpp = 0.;	95	double year_raw_gpp = 0.;
96	clear();	96	clear();
97	double utilizable_rad, epsilon;	97	double utilizable_rad, epsilon;
98	// conversion from gC to kg Biomass: C/Biomass=0.5	98	// conversion from gC to kg Biomass: C/Biomass=0.5
99	const double gC_to_kg_biomass = 1. / (biomassCFraction * 1000.);	99	const double gC_to_kg_biomass = 1. / (biomassCFraction * 1000.);
100	for (int i=0;i<12;i++) {	100	for (int i=0;i<12;i++) {
101	utilizable_rad = calculateUtilizablePAR(i); // utilizable radiation of the month times ... (MJ/m2)	-
-		101	utilizable_rad = calculateUtilizablePAR(i); // utilizable radiation of the month ... (MJ/m2)
102	epsilon = calculateEpsilon(i); // ... photosynthetic efficiency ... (gC/MJ)	102	epsilon = calculateEpsilon(i); // ... photosynthetic efficiency ... (gC/MJ)
103	mUPAR[i] = utilizable_rad ;	103	mUPAR[i] = utilizable_rad ;
104	mGPP[i] =utilizable_rad * epsilon * gC_to_kg_biomass; // ... results in GPP of the month kg Biomass/m2 (converted from gC/m2)	104	mGPP[i] =utilizable_rad * epsilon * gC_to_kg_biomass; // ... results in GPP of the month kg Biomass/m2 (converted from gC/m2)
105	year_raw_gpp += mGPP[i]; // kg Biomass/m2	105	year_raw_gpp += mGPP[i]; // kg Biomass/m2
106	}	106	}
107		107
108	// calculate f_env,yr: see http://iland.boku.ac.at/sapling+growth+and+competition	108	// calculate f_env,yr: see http://iland.boku.ac.at/sapling+growth+and+competition
109	double f_sum = 0.;	109	double f_sum = 0.;
110	for (int i=0;i<12;i++)	110	for (int i=0;i<12;i++)
111	f_sum += mGPP[i] / gC_to_kg_biomass; // == uAPar * epsilon_eff	111	f_sum += mGPP[i] / gC_to_kg_biomass; // == uAPar * epsilon_eff
112		112
113	// the factor f_ref: parameter that scales response values to the range 0..1 (1 for best growth conditions) (species parameter)	113	// the factor f_ref: parameter that scales response values to the range 0..1 (1 for best growth conditions) (species parameter)
114	const double perf_factor = mResponse->species()->saplingGrowthParameters().referenceRatio;	114	const double perf_factor = mResponse->species()->saplingGrowthParameters().referenceRatio;
115	// f_env,yr=(uaparepsilon_eff) / (APAR epsilon_0 * fref)	115	// f_env,yr=(uaparepsilon_eff) / (APAR epsilon_0 * fref)
116	mEnvYear = f_sum / (Model::settings().epsilon * mResponse->yearlyRadiation() * perf_factor);	116	mEnvYear = f_sum / (Model::settings().epsilon * mResponse->yearlyRadiation() * perf_factor);
117	if (mEnvYear > 1.) {	117	if (mEnvYear > 1.) {
118	qDebug() << "ERROR: fEnvYear > 1 for " << mResponse->species()->id() << mEnvYear << "f_sum, epsilon, yearlyRad, perf_factor" << f_sum << Model::settings().epsilon << mResponse->yearlyRadiation() << perf_factor;	118	qDebug() << "ERROR: fEnvYear > 1 for " << mResponse->species()->id() << mEnvYear << "f_sum, epsilon, yearlyRad, perf_factor" << f_sum << Model::settings().epsilon << mResponse->yearlyRadiation() << perf_factor;
119	mEnvYear = 1.;	119	mEnvYear = 1.;
120	}	120	}
121		121
122	// calculate fraction for belowground biomass	122	// calculate fraction for belowground biomass
123	mRootFraction = 1. - abovegroundFraction();	123	mRootFraction = 1. - abovegroundFraction();
124		124
125	// global value set?	125	// global value set?
126	double dbg = GlobalSettings::instance()->settings().paramValue("gpp_per_year",0);	126	double dbg = GlobalSettings::instance()->settings().paramValue("gpp_per_year",0);
127	if (dbg) {	127	if (dbg) {
128	year_raw_gpp = dbg;	128	year_raw_gpp = dbg;
129	mRootFraction = 0.4;	129	mRootFraction = 0.4;
130	}	130	}
131		131
132	// year GPP/rad: kg Biomass/m2	132	// year GPP/rad: kg Biomass/m2
133	mGPPperArea = year_raw_gpp;	133	mGPPperArea = year_raw_gpp;
134	return mGPPperArea; // yearly GPP in kg Biomass/m2	134	return mGPPperArea; // yearly GPP in kg Biomass/m2
135	}	135	}
136		136