Skip to content

GitLab

Commit 9ef154a8 authored by Ronan Trepos's avatar Ronan Trepos
Browse files

all Pea modules: update 

Integration of work of M. Chanis on Pea.
parent e02e8470
Expand all Hide whitespace changes
Inline Side-by-side
AZODYN/data/Breteniere_1991_2018.csv 0 → 100644
View file @ 9ef154a8
This diff is collapsed.
AZODYN/exp/AzodynPea.vpz
View file @ 9ef154a8
This diff is collapsed.
AZODYN/src/NitrogenFixationPea.cpp
View file @ 9ef154a8
......@@ -64,13 +64,13 @@ public:
Profnod=pp.Profnod_opt;
} else if(ps.da_C0<=pp.damax_C0){
Profnod=((ps.da_C0-ps.da_C0_opt)*(pp.Profnod_opt-pp.Profnod_min)
/(ps.da_C0_opt - pp.damax_C0)+pp.Profnod_opt);
/(ps.da_C0_opt - pp.damax_C0))+pp.Profnod_opt;
} else{
Profnod=pp.Profnod_min;
}
// Quantite d'eau dans la couche avec nodosites(mm)
if (j< pp.date_S){
if (pp.begin_date + j < pp.date_S){
RU_nod = 0.0 ;
} else {
RU_nod = RU_C0()*Profnod()/ps.ep_C0;
......@@ -85,7 +85,7 @@ public:
}
//calcul du nombre de jours ou la couche avec nodosites est en stress hydrique par rapport a la fixation
if(j< pp.date_FLO){
if(pp.begin_date + j < pp.date_FLO){
NJstressH_nod = 0.0 ;
} else if (NJstressH_nod(-1) >= pp.NJstress_nod_seuil) {
NJstressH_nod=NJstressH_nod(-1);
......@@ -100,7 +100,7 @@ public:
//l'inhibition de la fixation due a la presence d'azote mineral dans le sol)
if (Stade() < VEGETATIF){
pNac_fix = 0.0;
} else if (j < pp.date_FLO) {
} else if (pp.begin_date + j < pp.date_FLO) {
pNac_fix = std::max(0.0,std::min(
pp.coefNsol_fix1* QNaccess_C0() + pp.a_pNfix,
pp.pNfix_max));
......@@ -115,10 +115,11 @@ public:
if (STapLEVEE() < pp.STfixmax_seuil){
dQNfixmax_tot = 0.0 ;
} else if (Stade() < DRG_FSLA){
dQNfixmax_tot = pp.coef_tot_a * (MS()-MS(-1))*10.0*
(pp.coef_fixmax1*STapLEVEE() + pp.coef_fixmax3);
dQNfixmax_tot = std::max(0.0,pp.coef_tot_a * (MS()-MS(-1))*10.0*
(pp.coef_fixmax1*STapLEVEE() + pp.coef_fixmax3));
} else if (Stade() < FSLA_MP){
dQNfixmax_tot = (MS()-MS(-1))*10.0*pp.coef_fixmax2 *pp.coef_tot_a;
dQNfixmax_tot =std::max(0.0,(MS()-MS(-1))*10.0*pp.coef_fixmax2
*pp.coef_tot_a);
} else {
dQNfixmax_tot=0.0;
}
......@@ -129,19 +130,12 @@ public:
//quantite totale d'azote fixe le jour j (kgN/ha)
if (j < pp.date_FLO){
if (stressH_nod() == 1) {
dQNfix_tot = 0.0;
} else {
dQNfix_tot = std::min(dQNfix_besoin(),dQNfixmax_tot());
}
} else if (NJstressH_nod(-1) >= pp.NJstress_nod_seuil){
dQNfix_tot = 0.0;
} else {
if (pp.begin_date + j < pp.date_FLO or NJstressH_nod(-1) <= pp.NJstress_nod_seuil ){
dQNfix_tot = std::min(dQNfix_besoin(),dQNfixmax_tot());
} else {
dQNfix_tot = 0.0;
}
//quantite d'azote dans les parties aeriennes issu de la fixation jusqu'au jour j (kgN/ha)
if (Stade() < VEGETATIF){
QNfix = 0.0;}
......@@ -151,7 +145,7 @@ public:
//Nitrogen derived from fixation: proportion d'azote fixe sur l'azote
//total accumule dans les parties aeriennes jusqu'au jour j (kgN/ha)
//total accumule dans les parties aeriennes jusqu'au jour j (%)
if (Stade() < VEGETATIF){
Ndfa = 0.0;
} else {
......
AZODYN/src/NitrogenPlantPea.cpp
View file @ 9ef154a8
......@@ -4,8 +4,7 @@
#include <vle/DiscreteTime.hpp>
#include <Phenology.hpp>
#include <ParametersPlant.hpp>
#include <ParametersPlantPea.hpp>
#include<iostream>
#include<cstdlib>
#include<time.h>
......@@ -29,6 +28,7 @@ public:
Stade.init(this, "Stade", events);
MS.init(this, "MS", events);
MS_DRG.init(this,"MS_DRG", events);
dQNfix_tot.init(this, "dQNfix_tot", events);
Tmoy.init(this, "Tmoy", events);
QN_offre.init(this, "QN_offre", events);
......@@ -46,8 +46,9 @@ public:
dQNabs_C0.init(this, "dQNabs_C0", events);
dQNabs_C1.init(this, "dQNabs_C1", events);
tNa.init(this, "tNa", events);
tNa_DF.init(this, "tNa_DF", events);
tNveg.init(this, "tNveg", events);
tNveg_MP.init(this, "tNveg_MP", events);
}
virtual ~NitrogenPlantPea()
......@@ -55,7 +56,7 @@ public:
virtual void compute(const vd::Time& /*j*/)
{
//Quantite d azote totale par la culture le jour j (kg N/ha)
//Quantite d azote totale accumulee par la culture le jour j (kg N/ha)
if(Stade() < VEGETATIF){
QN_tot=0.0;
......@@ -75,19 +76,19 @@ public:
}
// besoin de la culture le jour j (kg N/ha)
if (Stade() < VEGETATIF){
if (Stade() < VEGETATIF or Stade () ==MP){
dQN_besoin = 0.0 ;
} else if(Stade() <DRG_FSLA){
if((MS()/100)<pp.MS_Nmax_seuil) {
dQN_besoin= (pp.coef_tot_a*MS()*10.0*pp.tNmax_Nmax/100)-
(pp.coef_tot_a*MS(-1)*10.0*pp.tNmax_Nmax/100) ;
dQN_besoin= std::max(0.0,(pp.coef_tot_a*MS()*10.0*pp.tNmax_Nmax/100)-
(pp.coef_tot_a*MS(-1)*10.0*pp.tNmax_Nmax/100)) ;
} else {
dQN_besoin= pp.coef_tot_a*MS()*10.0*pp.coef_Nmax/100*std::pow((MS()/100),pp.pente_Nmax) -
pp.coef_tot_a*MS(-1)*10.0*pp.coef_Nmax/100*std::pow((MS(-1)/100),pp.pente_Nmax);
dQN_besoin= std::max(0.0,pp.coef_tot_a*MS()*10.0*pp.coef_Nmax/100*std::pow((MS()/100),pp.pente_Nmax) -
pp.coef_tot_a*MS(-1)*10.0*pp.coef_Nmax/100*std::pow((MS(-1)/100),pp.pente_Nmax));
}
} else {
dQN_besoin= pp.coef_tot_a*MS()*10.0*(QN_tot()/pp.coef_tot_a)/(MS()*10.0)-
pp.coef_tot_a*MS(-1)*10.0*(QN_tot()/pp.coef_tot_a)/(MS(-1)*10.0);
dQN_besoin= std::max(0.0,pp.coef_tot_a*MS()*10.0*pp.coef_Nmax/100*std::pow(MS_DRG()*pp.pMSapDRG/100,pp.pente_Nmax) -
pp.coef_tot_a*MS(-1)*10.0*pp.coef_Nmax/100*std::pow(MS_DRG()*pp.pMSapDRG/100,pp.pente_Nmax));
}
......@@ -109,7 +110,7 @@ public:
}
//quantites d'azote absorbees dans C0
//quantites d'azote absorbees dans C0 (kg/ha)
if (QN_offre() == 0.0) {
dQNabs_C0 = 0;
......@@ -117,7 +118,7 @@ public:
dQNabs_C0 = QNaccess_C0()/(QNaccess_C0() + QNaccess_C1()) * dQNabs();
}
//quantites d'azote absorbees dans C1
//quantites d'azote absorbees dans C1 (kg/ha)
if (QN_offre() == 0.0) {
dQNabs_C1 = 0;
......@@ -125,7 +126,7 @@ public:
dQNabs_C1 = QNaccess_C1()/(QNaccess_C0() + QNaccess_C1()) * dQNabs();
}
//teneur en azote aerien
//teneur en azote aerien (%)
if(Stade() < VEGETATIF){
tNa = 0;
......@@ -134,32 +135,60 @@ public:
} else {
tNa = 100*(QN_tot()/pp.coef_tot_a)/(MS()*10);
}
//Teneur en azote des parties aériennes à début floraison
if(Stade() < DF_DRG) {
tNa_DF = 0 ;
} else if (Stade(-1) == VEGETATIF and Stade() == DF_DRG) {
tNa_DF = tNa() ;
} else {
tNa_DF = tNa_DF(-1) ;
}
//teneur en azote de la partie vegetative, partie aerienne(sans les graines) (%)
//teneur en azote de la partie vegétative (sans les graines)
if(Stade() < VEGETATIF){
tNveg = 0;
} else if (Stade() < DRG_FSLA) {
} else if (Stade() <= DRG_FSLA) {
tNveg = tNa();
} else {
tNveg = 100*QNveg()/((MS()-MSG())*10);
} else {
tNveg = std::min(100*QNveg() /((MS()-MSG())*10 ),tNveg(-1));
}
//teneur en azote des parties aériennes végétatives à maturité physiologique (%)
if(Stade() < MP){
tNveg_MP = 0;
} else if (Stade(-1) == FSLA_MP and Stade() == MP) {
tNveg_MP = tNveg();
} else {
tNveg_MP = tNveg_MP(-1);
}
// Quantite d azote critique le jour j (kg N/ha)
if(Stade() < VEGETATIF){
QNc=0.0;
} else if (MS()/100<pp.MS_Nc_seuil){
QNc= MS()*10.0* (pp.tNc_Nc/100);
if (Stade() < VEGETATIF ){
dQN_besoin = 0.0 ;
} else if(Stade() <DRG_FSLA){
if((MS()/100)<pp.MS_Nmax_seuil) {
QNc= MS()*10.0* (pp.tNc_Nc/100);
} else {
QNc= MS()*10.0*pp.coef_Nc/100*std::pow(MS()/100,pp.pente_Nc);
}
} else {
QNc= MS()*10.0*pp.coef_Nc/100*std::pow((MS()/100),pp.pente_Nc);
QNc= MS()*10.0*pp.coef_Nc/100*std::pow((MS_DRG()*pp.pMSapDRG/100),pp.pente_Nc);
}
}
private:
ParametersPlant pp;
ParametersPlantPea pp ;
/*Sync*/ Var Stade;
/*Sync*/ Var MS;
......@@ -171,6 +200,7 @@ private:
/*Sync*/ Var densRac;
/*Sync*/ Var QNveg;
/*Sync*/ Var MSG;
/*Sync*/ Var MS_DRG ;
Var dQN_besoin;
Var dQN_demande;
......@@ -181,6 +211,8 @@ private:
Var dQNabs_C1;
Var tNa;
Var tNveg;
Var tNa_DF ;
Var tNveg_MP ;
};
} // namespace
......
AZODYN/src/ParametersPlantPea.hpp
View file @ 9ef154a8
......@@ -142,6 +142,9 @@ public:
//proteines
double coef_prot_gr;
// Poids de mille grain moyen
double PMGmoy_var;
void initialiser( const vle::devs::InitEventList& events ){
ParametersPlant::initialiser(events);
// Paramètres liés à la phénologie
......@@ -181,8 +184,9 @@ public:
coef_vN = Utils::extractDouble(events, "coef_vN");
pNfsla = Utils::extractDouble(events, "pNfsla");
coef_vcrgr1 = Utils::extractDouble(events, "coef_vcrgr1");
coef_vcrgr2 = Utils::extractDouble(events, "coef_vcrgr2");
coef_prot_gr = Utils::extractDouble(events, "coef_prot_gr");
coef_vcrgr2 = Utils::extractDouble(events, "coef_vcrgr2");
coef_prot_gr = Utils::extractDouble(events, "coef_prot_gr");
PMGmoy_var = Utils::extractDouble(events, "PMGmoy_var");
}
};
......
AZODYN/src/PhenologyPea.cpp
View file @ 9ef154a8
......@@ -28,8 +28,6 @@ public:
Tmin.init(this, "Tmin", events);
Tmax.init(this, "Tmax", events);
tNveg.init(this, "tNveg", events);
P1G.init(this, "P1G", events);
IndicMP.init(this, "IndicMP", events);
Tmoy.init(this, "Tmoy", events);
......@@ -48,9 +46,11 @@ public:
STapFSLA.init(this, "STapFSLA", events);
Nj_LEVEE.init(this, "Nj_LEVEE", events);
Nj_DRG.init(this, "NjDRG", events);
Nj_DRG.init(this, "Nj_DRG", events);
Nj_FSLA.init(this, "Nj_FSLA", events);
Nj_MP.init(this, "Nj_MP", events);
Nj_S_DF.init(this, "Nj_S_DF", events);
Nj_S_MP.init(this, "Nj_S_MP", events);
}
......@@ -60,14 +60,14 @@ public:
virtual void compute(const vd::Time& j)
{
// Calcul de la tempeature moyenne au jour j
// Calcul de la temperature moyenne au jour j
Tmoy= (Tmin()+Tmax())/2;
// Tempetaure moyenne en prenant en compte un Tempeature de base
// pour le developement de la plante
// Degrés efficaces pour le developement de la plante,
//en prenant en compte un Tempeature de base
Tmoy_base = Tmoy() - pp.Tbase;
// Somme de temperature apres semis
// Temps thermique depuis le semis
if (pp.begin_date + j < pp.date_S) {
STapS = 0.0;
} else {
......@@ -85,8 +85,7 @@ public:
STapDRG(-1) >= (pp.DRG_ET*(pp.NETpot-1))) {
Stade = FSLA_MP;
} else if ((Stade(-1) == FSLA_MP) and (
(tNveg(-1)<(pp.tNstruct/100)) or
(P1G(-1)>=pp.PMGmax_var/1000) or
(IndicMP(-1)==1) or
(STapFSLA(-1) >= pp.STmax_FSLA_MP))){
Stade=MP;
......@@ -101,7 +100,7 @@ public:
date_LEVEE = date_LEVEE(-1);
}
//calcul date RG
//calcul date DRG
if (Stade(-1) == DF_DRG and Stade() == DRG_FSLA) {
date_DRG = pp.begin_date + j;
} else {
......@@ -121,28 +120,28 @@ public:
} else {
date_MP = date_MP(-1);
}
// Somme de temperature apres levee
// Temps thermique depuis la levee
if (Stade() < VEGETATIF) {
STapLEVEE = 0.0;
} else {
STapLEVEE = STapLEVEE(-1) + Tmoy_base();
}
// Somme de temperature apres debut floraison
// Temps thermique depuis debut floraison
if (Stade() < DF_DRG) {
STapDF = 0.0;
} else {
STapDF = STapDF(-1) + Tmoy_base();
}
// Somme de temperature apres DRG
// Temps thermique depuis DRG
if (Stade() < DRG_FSLA) {
STapDRG = 0.0;
} else {
STapDRG = STapDRG(-1) + Tmoy_base();
}
// Somme de temperature apres FSLA
// Temps thermique depuis FSLA
if (Stade() < FSLA_MP) {
STapFSLA = 0.0;
} else {
......@@ -181,6 +180,16 @@ public:
if(Nj_MP()==7){
date_MP_prev=pp.begin_date + j;
}
//Nombre de jours entre le semis et la date de début floraison
if (Stade() >= DF_DRG){
Nj_S_DF= pp.date_FLO - pp.date_S;
}
//Nombre de jours entre le semis et la date de de maturité physiologique
if (Stade() >= MP){
Nj_S_MP= date_MP() - pp.date_S;
}
}
......@@ -194,14 +203,10 @@ private:
ParametersPlantPea pp;
/*Sync*/ Var Tmin;
/*Sync*/ Var Tmax;
/*Nosync*/ Var tNveg;
/*Nosync*/ Var P1G;
/*Nosync*/ Var IndicMP;
Var Tmoy;
Var Tmoy_base;
Var Stade;
......@@ -219,6 +224,8 @@ private:
Var Nj_DRG;
Var Nj_FSLA;
Var Nj_MP;
Var Nj_S_DF;
Var Nj_S_MP;
};
} // namespace AZODYN
......
AZODYN/src/PlantGrowthPea.cpp
View file @ 9ef154a8
......@@ -51,7 +51,14 @@ public:
MS_DF.init(this, "MS_DF", events);
MS_DRG.init(this, "MS_DRG", events);
MS_FSLA.init(this, "MS_FSLA", events);
MS_MP.init(this, "MS_MP", events);
MS_noGEL.init(this, "MS_noGEL", events);
NjstressN_Eb.init (this, "NjstressN_Eb",events) ;
stressN_cumul.init (this, "stressN_cumul",events) ;
NjstressH.init(this, "NjstressH", events);
stressH_cumul.init(this, "stressH_cumul", events);
NjstressT.init(this, "NjstressT", events);
stressT_cumul.init(this, "stressT_cumul", events);
}
virtual ~PlantGrowthPea()
......@@ -59,7 +66,7 @@ public:
virtual void compute(const vd::Time& j)
{
// calcul de la surface foliaire (LAI) le jour j
// calcul de l'indice foliaire (LAI) le jour j
if((Stade() < VEGETATIF) or (IndicMP(-1)==1)){
LAI=0.0;
......@@ -80,51 +87,50 @@ public:
Ei = pp.Eimax*(1-std::exp(-pp.k*LAI()));
}
//Efficience de conversion potentielle du rayonnement le jour j (sans stress N et H)
// Efficience de conversion potentielle du rayonnement le jour j
// (sans stress N et H) en g/MJ
if ((Stade() < VEGETATIF) or (IndicMP(-1)==1)){
Ebpot = 0.0;
Ebpot = 0;
} else if (pp.begin_date + j< pp.date_FLO) {
Ebpot = pp.Ebv * stressT(-1);
Ebpot = pp.Ebv;
} else if ( Stade() >= DRG_FSLA and tNveg(-1)<pp.tNveg_seuil){
Ebpot = std::max(0.0,std::min(stressT(-1) * pp.Ebr,
stressT(-1)*pp.Ebr*(tNveg(-1) - pp.tNstruct)
Ebpot = std::max(0.0,std::min(pp.Ebr,
pp.Ebr*(tNveg(-1) - pp.tNstruct)
/(pp.tNveg_seuil - pp.tNstruct)));
} else {
Ebpot = (pp.Ebv +((pp.Ebr - pp.Ebv)
*(STapDF()/pp.ST_DF_DRG)))*stressT(-1);
*(STapDF()/pp.ST_DF_DRG)));
}
// ------ fonctionnement de la plante avec prise en compte des stress N et H --------------------------------
// efficience de conversion du rayonnement
// efficience de conversion du rayonnement (g/MJ)
if (Stade() < VEGETATIF) {
if (Stade() < VEGETATIF ) {
Eb=0.0;
} else {
Eb=Ebpot()*stressH(-1)*stressN_Eb(-1);
Eb=Ebpot()*stressT(-1)*stressH(-1)*stressN_Eb(-1);
}
// matiere seche accumulee depuis l initilisation jusqu au jour j
if ((Stade() < VEGETATIF) or (IndicMP(-1)==1)){
if (Stade() < VEGETATIF){
MS= 0.0;
} else if(Stade(-1)==SOL_NU and Stade()==VEGETATIF) {
if (pp.MS_init == -1) {//MS_init not available
MS= (pp.densite_semis * pp.PMGmax_var/1000) / pp.coef_MS_init;
MS= (pp.densite_semis * pp.PMGmoy_var/1000) / pp.coef_MS_init;
} else {
MS = pp.MS_init;
}
} else {
MS=MS(-1)+ RG()/100.0*pp.Ec*Ei()*Eb()*(stressGEL(-1));
} else if (Ei()==0.0) {
MS=MS(-1)*stressGEL(-1);
} else{
MS=(MS(-1)+ RG()/100.0*pp.Ec*Ei()*Eb())*(stressGEL(-1));
}
// matiere seche accumulee depuis l initilisation jusqu'au jour j sans prise en compte du gel
MS_noGEL= MS()/stressGEL(-1);
// matiere seche accumulee depuis l initilisation jusqu a DF
// matiere seche accumulee depuis l initilisation jusqu a DF (g/m²)
if (pp.begin_date + j< pp.date_FLO) {
MS_DF= 0.0;
......@@ -134,7 +140,7 @@ public:
MS_DF=MS_DF(-1);
}
// matiere seche accumulee depuis l initilisation jusqu a DRG
// matiere seche accumulee depuis l initilisation jusqu a DRG (g/m²)
if (Stade() < DRG_FSLA) {
MS_DRG= 0.0;
......@@ -144,7 +150,7 @@ public:
MS_DRG=MS_DRG(-1);
}
// matiere seche accumulee depuis l initilisation jusqu a FSLA
// matiere seche accumulee depuis l initilisation jusqu a FSLA (g/m²)
if (Stade() < FSLA_MP) {
MS_FSLA= 0.0;
......@@ -154,6 +160,72 @@ public:
MS_FSLA=MS_FSLA(-1);
}
// matiere seche à maturité physiologique (g/m²)
if (Stade() < MP) {
MS_MP= 0.0;
} else if((Stade(-1)== FSLA_MP) and (Stade() == MP)) {
MS_MP= MS();
} else {
MS_MP=MS_MP(-1);
}
//Variables intégratives de stress sur la totalité du cycle de culture
//Nombre de jour de stressN_Eb : Cumul des jours où stressN_Eb<1
if (Stade () < VEGETATIF ) {
NjstressN_Eb = 0;
} else if (Stade () >= MP or stressN_Eb () ==1 ) {
NjstressN_Eb = NjstressN_Eb(-1);
} else if (stressN_Eb ()< 1) {
NjstressN_Eb = NjstressN_Eb(-1) + 1 ;
}
//Intensité cumulée de stressN_Eb sur tout le cycle: cumul de 1-stressN_Eb
if (Stade () < VEGETATIF ) {
stressN_cumul = 0;
} else if (Stade () >= MP ) {
stressN_cumul = stressN_cumul(-1) ;
} else {
stressN_cumul = stressN_cumul(-1) + (1-stressN_Eb()) ;
}
//Nombre de jour de stressH : Cumul des jours où stressH<1
if (Stade () < VEGETATIF ) {
NjstressH = 0;
} else if (Stade () >= MP or stressH () ==1 ) {
NjstressH = NjstressH(-1);
} else if (stressH ()< 1) {
NjstressH = NjstressH(-1) + 1 ;
}
//Intensité cumulée de stressH sur tout le cycle: cumul de 1-stressH
if (Stade () < VEGETATIF ) {
stressH_cumul = 0;
}else if (Stade () >= MP ) {
stressH_cumul = stressH_cumul(-1) ;
}else {
stressH_cumul = stressH_cumul(-1) + (1-stressH()) ;
}
//Nombre de jour de stressT : Cumul des jours où stressT<0.8
if (Stade () < VEGETATIF ) {
NjstressT = 0;
}else if (Stade () >= MP or stressT () >= 0.8 ) {
NjstressT = NjstressT(-1);
}else if (stressT ()< 0.8) {
NjstressT = NjstressT(-1) + 1 ;
}
//Intensité cumulée de stressT sur tout le cycle, jours ou stressT<0.8: cumul de 1-stressH
if (Stade () < VEGETATIF ) {
stressT_cumul = 0;
}else if (Stade () >= MP or stressT () >= 0.8 ) {
stressT_cumul = stressT_cumul(-1) ;
}else if (stressT ()< 0.8){
stressT_cumul = stressT_cumul(-1) + (1-stressT()) ;
}
}
......@@ -184,7 +256,14 @@ private:
Var MS_DF;
Var MS_DRG;
Var MS_FSLA;
Var MS_MP;
Var MS_noGEL;
Var NjstressN_Eb ;
Var stressN_cumul;
Var NjstressH;
Var stressH_cumul;
Var NjstressT;