Merge branch 'leafSenescence' of git@gitlab.com:rootmodels/OpenSimRoot.git into FixesForFarquhar (7da02bf3) · Commits · rootmodels / OpenSimRoot

OpenSimRoot/InputFiles/templates/plantTemplate.IncludeLeafSenescence.xml

+38 −5

Original line number	Diff line number	Diff line
		<SimulationModelIncludeFile xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:noNamespaceSchemaLocation="../../scripts/XML/SimulaXMLSchema.xsd">

		<!--

		This adds leaf senescence to the model. You need to specify next to the leafAreaExpantionRate a leafAreaSenescenceRate which
		typically is the expansion rate shifted in time, where the time shift the lifeexpectency of a leaf is. Of course you are
		welcome to make the senescenc rate dependent on other factors, using for example multipliers.

		This will not reduce the leafArea, or the leafDryWeight. That is because carbon allocation to leaves must be based on the increase
		in leaf area, that is the production of new leaves, without the leave senescence rate.


		This is made transparent by putting in the 'green' version of dry weight and leaf area.

		The leaf area index is corrected, and with it the photo synthesis and the transpiration.

		note that this requires historical information on leave area, so the garbage collection must be turned off. The code should do this automatically.


		-->


		<SimulaDirective path="/shootTemplate">
		<SimulaVariable
		name="senescedLeafArea"
		unit="cm2"
		function="leafSenescence" />
		<!--maxTimeStep="0.1">
		0. </SimulaVariable-->
		<SimulaVariable
		name="senescedLeafDryWeight"
		unit="g"
		function="leafSenescence" />
		<!--maxTimeStep="0.1">
		0. </SimulaVariable-->
		<SimulaDerivative
		name="greenLeafArea"
		unit="cm2"
		@@ -25,6 +42,22 @@
		<SimulaDirective path="/plants">
		<SimulaVariable name="senescedLeafArea" function="sumOverAllPlantShoots" unit="cm2" />
		</SimulaDirective>
		</SimulationModelIncludeFile>


		<SimulaDirective path="/plantTemplate">
		<SimulaDerivative name="greenShootDryWeight" function="useFormula">
		<SimulaConstant name="formula" type="string"> shootDryWeight - senescedLeafDryWeight </SimulaConstant>
		<SimulaBase name="variablePaths">
		<SimulaConstant name="senescedLeafDryWeight" type="string"> plantPosition/shoot/senescedLeafDryWeight </SimulaConstant>
		</SimulaBase>
		</SimulaDerivative>
		<SimulaDerivative name="greenPlantDryWeight" function="useFormula">
		<SimulaConstant name="formula" type="string"> plantDryWeight - senescedLeafDryWeight </SimulaConstant>
		<SimulaBase name="variablePaths">
		<SimulaConstant name="senescedLeafDryWeight" type="string"> plantPosition/shoot/senescedLeafDryWeight </SimulaConstant>
		</SimulaBase>
		</SimulaDerivative>
		</SimulaDirective>

		</SimulationModelIncludeFile>

OpenSimRoot/src/modules/PlantNutrients/OptimalNutrientContent.cpp

+40 −4

Original line number	Diff line number	Diff line
		@@ -170,13 +170,14 @@ DerivativeBase * newInstantiationOptimalNutrientContent(SimulaDynamic* const pSD
		}

		ShootOptimalNutrientContent::ShootOptimalNutrientContent(SimulaDynamic* pSD) :
		DerivativeBase(pSD), conc(nullptr), dw(nullptr), refTime(0) { //leafs or stem
		DerivativeBase(pSD), conc(nullptr), dw(nullptr), dws(nullptr), refTime(0) { //leafs or stem
		std::string s(pSD->getName());
		std::size_t n(s.size() - 22);
		std::string sp(s.substr(0, n));
		std::string som(s.substr(0, n + 7));
		//dryweight
		dw = pSD->getParent(2)->getChild(sp + "DryWeight", "g");
		if(sp=="leaf") dws = pSD->getParent(2)->existingChild("senescedLeafDryWeight", "g");
		//concentration
		conc = pSD->existingSibling(som + "NutrientConcentration", pSD->getUnit() / "g");
		if(!conc){
		@@ -192,6 +193,11 @@ ShootOptimalNutrientContent::ShootOptimalNutrientContent(SimulaDynamic* pSD) :
		}
		void ShootOptimalNutrientContent::calculate(const Time &t, double &result) {
		dw->get(t, result);
		if(dws){
		double senescence;
		dws->get(t, senescence);
		result-=senescence;
		}
		double c;
		conc->get(t-refTime, c);
		result *= c;
		@@ -204,7 +210,7 @@ DerivativeBase * newInstantiationShootOptimalNutrientContent(SimulaDynamic* cons
		return new ShootOptimalNutrientContent(pSD);
		}

		ActualNutrientContent::ActualNutrientContent(SimulaDynamic* pSD):TotalBaseLabeled(pSD), organSize(nullptr) {
		ActualNutrientContent::ActualNutrientContent(SimulaDynamic* pSD):TotalBaseLabeled(pSD), organSize(nullptr), leafSenescence(nullptr) {
		std::string s(pSD->getName());
		std::size_t n(s.size() - 21);
		std::string type = s.substr(s.size() - 7, s.size());
		@@ -228,6 +234,17 @@ ActualNutrientContent::ActualNutrientContent(SimulaDynamic* pSD):TotalBaseLabele
		} else{
		msg::error("ActualNutrientContent: Unknown unit");
		}
		if(organ=="leaf"){
		SimulaBase *p1 = pSD->getParent()->existingSibling("senescedLeafDryWeight", "g");
		SimulaBase *p2 = pSD->getParent()->existingSibling("senescedLeafArea", "cm2");
		if (p1 && pSD->getUnit() == "uMol/cm2"){
		leafSenescence=p1;
		} else if (p2 && pSD->getUnit() == "uMol/g"){
		leafSenescence=p2;
		} else{
		if(p1 \|\| p2) msg::error("ActualNutrientContent: Unit not compatible with leaf senescence, use area or mass based");
		}
		}
		} else{
		msg::error("ActualNutrientContent:: Unexpected tag name");
		}
		@@ -269,6 +286,11 @@ void ActualNutrientContent::calculate(const Time &t, double &result) {
		if (organSize){
		double oSize;
		organSize->get(t, oSize);
		if(leafSenescence){
		double sen;
		leafSenescence->get(t,sen);
		oSize-=sen; //todo, assuming senesced leaves contain no nutrients.
		}
		if (oSize > 0) {
		result = result/oSize;
		} else{
		@@ -285,7 +307,7 @@ DerivativeBase * newInstantiationActualNutrientContent(SimulaDynamic* const pSD)
		return new ActualNutrientContent(pSD);
		}

		LeafNutrientContent::LeafNutrientContent(SimulaDynamic* pSD):DerivativeBase(pSD), leafSize(nullptr), leafArea(nullptr), leafPartArea(nullptr) {
		LeafNutrientContent::LeafNutrientContent(SimulaDynamic* pSD):DerivativeBase(pSD), leafSize(nullptr), leafArea(nullptr), leafPartArea(nullptr), leafSenescence(nullptr) {
		std::string s(pSD->getName());
		std::string name = s.substr(0, 6); // name = sunlit or shaded
		if (name == "sunlit" \|\| name == "shaded"){
		@@ -309,6 +331,15 @@ LeafNutrientContent::LeafNutrientContent(SimulaDynamic* pSD):DerivativeBase(pSD)
		} else{
		msg::error("leafNutrientContent: Unknown unit");
		}
		SimulaBase *p1 = pSD->getParent()->existingSibling("senescedLeafDryWeight", "g");
		SimulaBase *p2 = pSD->getParent()->existingSibling("senescedLeafArea", "cm2");
		if (p1 && pSD->getUnit() == "uMol/cm2"){
		leafSenescence=p1;
		} else if (p2 && pSD->getUnit() == "uMol/g"){
		leafSenescence=p2;
		} else{
		if(p1 \|\| p2) msg::error("ActualNutrientContent: Unit not compatible with leaf senescence, use area or mass based");
		}
		} else{
		msg::error("leafNutrientContent:: Unexpected tag name");
		}
		@@ -358,6 +389,11 @@ void LeafNutrientContent::calculate(const Time &t, double &result) {
		if (leafSize){
		double lSize;
		leafSize->get(t, lSize);
		if(leafSenescence){
		double sen;
		leafSenescence->get(t,sen);
		lSize-=sen; //todo, assuming senesced leaves contain no nutrients.
		}
		if (lSize > 0){
		result = result/lSize;
		} else{

OpenSimRoot/src/modules/PlantNutrients/OptimalNutrientContent.hpp

+3 −3

Original line number	Diff line number	Diff line
		@@ -49,7 +49,7 @@ public:
		std::string getName() const;
		protected:
		void calculate(const Time &t, double &var);
		SimulaBase conc, dw;
		SimulaBase conc, dw, *dws;
		Time refTime;
		};
		class ActualNutrientContent: public TotalBaseLabeled {
		@@ -58,7 +58,7 @@ public:
		std::string getName() const;
		protected:
		void calculate(const Time &t, double &var);
		SimulaBase totalUptake, plantOptimalContent, plantMinimalContent, organOptimalContent, organMinimalContent, organSize;
		SimulaBase totalUptake, plantOptimalContent, plantMinimalContent, organOptimalContent, organMinimalContent, organSize, *leafSenescence;
		};
		class LeafNutrientContent: public DerivativeBase {
		public:
		@@ -66,7 +66,7 @@ public:
		std::string getName() const;
		protected:
		void calculate(const Time &t, double &var);
		SimulaBase totalUptake, plantOptimalContent, plantMinimalContent, leafOptimalContent, leafMinimalContent, leafSize, leafArea, leafPartArea;
		SimulaBase totalUptake, plantOptimalContent, plantMinimalContent, leafOptimalContent, leafMinimalContent, leafSize, leafArea, leafPartArea, *leafSenescence;
		};

		#endif

OpenSimRoot/src/modules/ShootModule/LeafSenescence.cpp

+19 −4

Original line number	Diff line number	Diff line
		@@ -29,13 +29,25 @@ protected:

		void calculate(const Time &t, double& r){
		senescenceRate->get(t-plantingTime,r);

		// make sure r is not so large that leafArea-senescedLeafArea becomes negative
		double rla, la, sl;
		leafArea->get(t,la);//current leaf area
		leafArea->getRate(t,rla);//current leaf area
		senescedLeafArea->get(t,sl);//current amount of dead leaves

		if( (la+rla-sl)<r) r=0.999*(la+rla-sl); //0.999 for estimation errors in rla causing slight overshoot of senescence and thereby a small negative value.
		if(r<0) r=0;
		if(SLA){
		Time st(lt); //time at which the leaf was created, to get appropriate sla.
		double sl,sla;
		senescedLeafArea->get(t,sl);//current amount of senesced leafs
		//std::cout<<std::endl<<sl<<" "<<st<<" "<<lt<<" "<<t;
		if(la>sl){
		leafArea->getTime(sl,st,lt,t);
		lt=st;
		}else{
		st=t;
		}
		double sla;
		SLA->get(st-plantingTime,sla);
		if(SLAisLAR_){
		r=sla;//convert from cm2 to g; cm2g/cm2=g
		@@ -55,6 +67,8 @@ public:
		SimulaBase *parameters(GETSHOOTPARAMETERS(plantType));

		senescenceRate=parameters->getChild("leafAreaSenescenceRate","cm2/day");
		leafArea = pSD->getSibling("leafArea","cm2");
		senescedLeafArea = pSD->getSibling("senescedLeafArea","cm2");

		//check if unit given in input file agrees with this function
		if(pSD->getUnit()=="g") {
		@@ -67,6 +81,7 @@ public:
		if(SLA->getUnit()!="g/cm2") msg::error("LeafArea: Unkown unit for "+SLA->getPath()+". Use g/cm2 or cm2/g");
		}
		leafArea = pSD->getSibling("leafArea","cm2");
		leafArea->garbageCollectionOff(); //need to turn it of as we will need historic data.
		senescedLeafArea = pSD->getSibling("senescedLeafArea","cm2");
		}else if (pSD->getUnit()=="cm2"){
		//do nothing