Frost.cpp 6.47 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
// @@tagdynamic@@
// @@tagdepends: vle.discrete-time @@endtagdepends

#include <vle/DiscreteTime.hpp>

#include <ParametersPlant.hpp>
#include<iostream>
#include<cstdlib>
#include<time.h>

namespace vd = vle::devs;
namespace vv = vle::value;

namespace AZODYN {

using namespace vle::utils;
using namespace vle::discrete_time;

class Frost : public DiscreteTimeDyn
{
public:
    Frost(const vd::DynamicsInit& atom,
          const vd::InitEventList& events)
: DiscreteTimeDyn(atom, events)
    {
        pp.initialiser(events);

        Tmoy.init(this, "Tmoy", events);
        Tmin.init(this, "Tmin", events);
        Tmax.init(this, "Tmax", events);
        STapS.init(this, "STapS", events);
        STapLEVEE.init(this, "STapLEVEE", events);

        ST_S_i.init(this, "ST_S_i", events);
        ST_S_f.init(this, "ST_S_f", events);
        Rmax.init(this, "Rmax", events);
        Rpot.init(this, "Rpot", events);
        dR.init(this, "dR", events);
        R.init(this, "R", events);
        stressGEL.init(this, "stressGEL", events);
        dDes.init(this, "dDes", events);
Ronan Trepos's avatar
Ronan Trepos committed
42
43
44
45
        Des_cumul.init(this, "Des_cumul", events);
        dindiceGEL.init(this, "dindiceGEL", events);
        indiceGEL_cumul.init(this, "indiceGEL_cumul", events);
        indiceGELmax.init(this, "indiceGELmax", events);
46
        NjGel.init(this, "NjGel", events);
Ronan Trepos's avatar
Ronan Trepos committed
47
48
49
50
51
        stressGEL_cumul.init(this, "stressGEL_cumul", events);
        dEnd.init(this, "dEnd", events);
        End_cumul.init(this, "End_cumul", events);
        pEnd.init(this, "pEnd", events);
        pEndmax.init(this, "pEndmax", events);
52
53
54
55
56
57
58
59
60
61
62
    }

    virtual ~Frost()
    {}

    virtual void compute(const vd::Time& j)
    {

        //---------------------------------//
        
        //Calcul de la somme de température entre le semis et le stade initial ou on l'augmente le Rmax
63
        ST_S_i = STapS() - STapLEVEE() + pp.NbFi * pp.phyllo;
64
65
      
        //Calcul de la somme de température entre le semis et le stade final ou on l'augmente plus le Rmax    
66
        ST_S_f = STapS() - STapLEVEE()  + pp.NbFf * pp.phyllo;
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90

        //Calcul de Rmax
        if(STapS()<= ST_S_i()) {
            Rmax = pp.R_Coleo;
        } else if(STapS()<=ST_S_f()){
            Rmax =((pp.Rmax_var-pp.R_Coleo)*STapS()+
                    (pp.R_Coleo*ST_S_f()-pp.Rmax_var*ST_S_i()))/
                            (ST_S_f()-ST_S_i());
        } else{
            Rmax = pp.Rmax_var;
        }

        //Calcul de Rpot
        if(pp.begin_date + j<=pp.date_S){
            Rpot = pp.Rmin;
        } else if((Tmoy() <  pp.Tmax_R) and (Tmoy() < pp.Tmin_R)){
            Rpot = Rmax();
        } else if ((Tmoy() <  pp.Tmax_R) and (Tmoy() >= pp.Tmin_R)){
            Rpot =(((pp.Rmin - Rmax()) * Tmoy()) +
                    ((Rmax() * pp.Tmax_R) - (pp.Rmin * pp.Tmin_R)))/(pp.Tmax_R-pp.Tmin_R);
        } else {
            Rpot = pp.Rmin;
        }

Ronan Trepos's avatar
Ronan Trepos committed
91
92
        //Calcul de de la variation journaliere de la résistance (degres-jour)
        if(pp.begin_date + j<=pp.date_S){
93
            dR = 0;
Ronan Trepos's avatar
Ronan Trepos committed
94
        } else if(Rpot() < R(-1)){
95
96
            dR =(Rpot()-pp.Rmin)/pp.Njend;
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
97
            dR =(pp.Rmin-Rmax()) * Tmoy()/(20.0 * pp.Njdes);
98
99
100
101
102
103
104
105
106
107
108
109
110
        }

        //Calcul de R la resistance
        if(pp.begin_date + j< pp.date_S){
            R=pp.Rmin;
        } else if(R(-1)==Rpot()){
            R=R(-1);
        } else if(R(-1)>Rpot()){
            R=std::max((R(-1)+dR()),Rpot());
        } else{
            R=std::min((R(-1)+dR()),Rpot());
        }

Ronan Trepos's avatar
Ronan Trepos committed
111
112
113
114
115
116
117
118
        //Calcul du stress ou ecart journalier entre la
        //resistance et la temperature minimale dindiceGEL
        if(pp.begin_date + j < pp.date_S ) {
            dindiceGEL = 0;
        } else {
            dindiceGEL = std::max(0.0,  R() - Tmin());
        }

119
        //Calcul de la destruction de la biomasse stressGEL (de 0: pas de stress à 1 : stress TOTAL)
Ronan Trepos's avatar
Ronan Trepos committed
120
        if((STapLEVEE() == 0) or (pp.begin_date + j  >= pp.date_fin_GEL)) {
121
122
123
124
125
126
            stressGEL = 1;
        } else  if (Tmin()>= R()){
            stressGEL=1;
        } else if((Tmin())<= R()- pp.GELmax_seuil){
            stressGEL = 0.0;
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
127
            stressGEL= 1.0-(1/ pp.GELmax_seuil)*dindiceGEL();
128
129
        }

Ronan Trepos's avatar
Ronan Trepos committed
130
131
132
133
134
135
136
137
138
        //Intensité cumulée de stressGEL entre la levée et la date_fin_GEL: cumul de 1-stressGEL
        if (STapLEVEE() == 0 ) {
            stressGEL_cumul = 0;
        } else {
          stressGEL_cumul = stressGEL_cumul(-1) + (1-stressGEL()) ;		
        }
		
        //Calcul du desendurcissement journalier dDes
		if (R() > R(-1)){
139
140
141
142
143
            dDes = R() - R(-1);
        } else {
            dDes = 0;
        }

Ronan Trepos's avatar
Ronan Trepos committed
144
145
146
147
148
149
150
        //Calcul du desendurcissement cumule

        Des_cumul = dDes() + Des_cumul(-1);        

       //Calcul de l'endurcissement journalier dEnd (°C)
        if (R() < R(-1)){
            dEnd = R() - R(-1);
151
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
152
            dEnd = 0;
153
154
        }

Ronan Trepos's avatar
Ronan Trepos committed
155
156
157
158
159
160
161
        //Calcul de l'endurcissement cumulé End_cumul (°C)

        End_cumul = dEnd() + End_cumul(-1);    

        //Pourcentage d'endurcissement pEnd
        if(pp.begin_date + j< pp.date_S){
             pEnd= 0;
162
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
163
             pEnd= (R()-pp.Rmin) / (Rmax() -pp.Rmin)*100;
164
        }
Ronan Trepos's avatar
Ronan Trepos committed
165
166
167
168
169
170
171
172
173
        //Pourcentage d'endurcissement maximal, pEndmax 
        if(pp.begin_date + j< pp.date_S){
             pEndmax= 0;
        } else if (pEnd()>pEndmax(-1)){
             pEndmax= pEnd();
        } else {
             pEndmax=pEndmax(-1);
        }
        
174
        //Calcul du stress cumulé Gel
Ronan Trepos's avatar
Ronan Trepos committed
175
176
        if(dindiceGEL()==0) {
            indiceGEL_cumul = indiceGEL_cumul(-1);
177
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
178
            indiceGEL_cumul = dindiceGEL() + indiceGEL_cumul(-1);
179
180
        }

Ronan Trepos's avatar
Ronan Trepos committed
181
182
183
        //Calcul de l'écart maimum journalier indiceGELmax
        if(dindiceGEL()> indiceGELmax(-1)){
            indiceGELmax = dindiceGEL();
184
        } else {
Ronan Trepos's avatar
Ronan Trepos committed
185
            indiceGELmax = indiceGELmax(-1);
186
187
188
        }

        //Calcul du nombre de jours de gel NjGel
Ronan Trepos's avatar
Ronan Trepos committed
189
        if(dindiceGEL()>0) {
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
            NjGel = NjGel(-1) + 1;
        } else {
            NjGel=NjGel(-1);
        }
    }

private:

    ParametersPlant pp;

    //double Njdes;

    /*Sync*/ Var Tmoy;
    /*Sync*/ Var Tmin;
    /*Sync*/ Var Tmax;
    /*Sync*/ Var STapS;
    /*Sync*/ Var STapLEVEE;

    Var ST_S_i;
    Var ST_S_f;
    Var Rmax;
    Var Rpot;
    Var dR;
    Var R;
    Var stressGEL;
    Var dDes;
Ronan Trepos's avatar
Ronan Trepos committed
216
217
218
219
    Var Des_cumul;
    Var dindiceGEL;
    Var indiceGEL_cumul;
    Var indiceGELmax;
220
    Var NjGel;
Ronan Trepos's avatar
Ronan Trepos committed
221
222
223
224
225
    Var stressGEL_cumul;
    Var dEnd;
    Var End_cumul;
    Var pEnd ; 
    Var pEndmax;
226
227
228
229
230
231
};

} // namespace AZODYN

DECLARE_DYNAMICS(AZODYN::Frost)