Commit 0c06be4a8fc631f9375ff0b767f0c279dce6d2d9

Authored by Olivier Mullier
1 parent d2eceb6d

premier version qui marche

1 #include "integrale.h" 1 #include "integrale.h"
2 2
3 -integrale::integrale(ibex::Interval inf, ibex::Interval sup, ibex::Function* f, ibex::Variable* var) : _inf_bound(inf), _sup_bound(sup), _f(f) 3 +integrale::integrale(const ibex::Interval & inf, const ibex::Interval & sup, ibex::Function* f, ibex::Variable* var) : _inf_bound(inf), _sup_bound(sup), _f(f)
4 { 4 {
5 _zeros_inf = nullptr; 5 _zeros_inf = nullptr;
6 _zeros_sup = nullptr; 6 _zeros_sup = nullptr;
  7 + _result = nullptr;
7 } 8 }
8 9
9 10
10 - std::list<zero> * integrale::construct_zeros_intern(const ibex::Interval & inf_bound) 11 + std::list<zero*> * integrale::construct_zeros_intern(const ibex::Interval & inf_bound)
11 { 12 {
12 #ifdef __DEBUG__ 13 #ifdef __DEBUG__
13 __SHOW_FCT_NAME__ 14 __SHOW_FCT_NAME__
14 #endif 15 #endif
15 //ATTENTION pas de preuve d'unicité 16 //ATTENTION pas de preuve d'unicité
16 - std::list<zero> * res = new std::list<zero>(); 17 + std::list<zero*> * res = new std::list<zero*>();
17 18
18 ibex::CtcNewton newton(*_f); 19 ibex::CtcNewton newton(*_f);
19 20
@@ -47,9 +48,11 @@ integrale::integrale(ibex::Interval inf, ibex::Interval sup, ibex::Function* f, @@ -47,9 +48,11 @@ integrale::integrale(ibex::Interval inf, ibex::Interval sup, ibex::Function* f,
47 } 48 }
48 } 49 }
49 if( (_f->eval_vector(temp)[0].ub() < 0.) ){ 50 if( (_f->eval_vector(temp)[0].ub() < 0.) ){
50 - res->push_back(zero(sols[i][0], -1)); 51 + zero *current = new zero(sols[i][0], -1);
  52 + res->push_back(current);
51 } else { 53 } else {
52 - res->push_back(zero(sols[i][0], 1)); 54 + zero *current = new zero(sols[i][0], 1);
  55 + res->push_back(current);
53 } 56 }
54 } 57 }
55 return res; 58 return res;
@@ -68,16 +71,17 @@ unsigned int integrale::construct_zeros_sup() @@ -68,16 +71,17 @@ unsigned int integrale::construct_zeros_sup()
68 return _zeros_sup->size(); 71 return _zeros_sup->size();
69 } 72 }
70 73
71 -ibex::Interval integrale::compute_integral(const double& inf_bound, const double& sup_bound) 74 +ibex::Interval integrale::compute_integral(const double inf, const double sup)
72 { 75 {
73 ibex::Interval res(0.0); 76 ibex::Interval res(0.0);
74 - double current = inf_bound; 77 + double current = inf;
75 ibex::IntervalVector current_vec(1); 78 ibex::IntervalVector current_vec(1);
76 - while (current < sup_bound){ 79 + while (current < sup){
77 current_vec[0] = ibex::Interval(current, current + __EPS__); 80 current_vec[0] = ibex::Interval(current, current + __EPS__);
78 res += _f->eval_vector(current_vec)[0] *__EPS__; 81 res += _f->eval_vector(current_vec)[0] *__EPS__;
79 current += __EPS__; 82 current += __EPS__;
80 } 83 }
  84 +
81 return res; 85 return res;
82 } 86 }
83 87
@@ -86,10 +90,11 @@ unsigned int integrale::construct_zeros() @@ -86,10 +90,11 @@ unsigned int integrale::construct_zeros()
86 return construct_zeros_inf() + construct_zeros_sup(); 90 return construct_zeros_inf() + construct_zeros_sup();
87 } 91 }
88 92
89 -integrale_sol integrale::compute_integral(const zero& inf, const zero& sup) 93 +integrale_sol integrale::compute_integral(zero *inf, zero *sup)
90 { 94 {
91 - ibex::Interval res_int = compute_integral(inf.get_x().lb(), sup.get_x().ub());  
92 - return integrale_sol(inf, sup, res_int); 95 + ibex::Interval res_inter = compute_integral(inf->get_x().lb(), sup->get_x().ub());
  96 + integrale_sol res(inf, sup, res_inter);
  97 + return res;
93 } 98 }
94 99
95 100
@@ -100,27 +105,31 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_inf(const ib @@ -100,27 +105,31 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_inf(const ib
100 construct_zeros_inf(); 105 construct_zeros_inf();
101 } 106 }
102 107
103 - std::list<zero>::iterator it = _zeros_inf->begin();  
104 - zero bound_inf(ibex::Interval(range.lb()), 0); 108 + std::list<zero*>::iterator it = _zeros_inf->begin();
  109 + zero *bound_inf = new zero(ibex::Interval(range.lb()), 0);
105 integrale_sol integ0 = compute_integral(bound_inf, *(it)); 110 integrale_sol integ0 = compute_integral(bound_inf, *(it));
106 res_tmp.push_back(integ0); 111 res_tmp.push_back(integ0);
107 112
108 while (it != _zeros_inf->end()){ 113 while (it != _zeros_inf->end()){
109 if (*it != _zeros_inf->back()){ 114 if (*it != _zeros_inf->back()){
110 - std::list<zero>::iterator it_first = it; 115 + std::list<zero*>::iterator it_first = it;
111 it++; 116 it++;
112 - std::list<zero>::iterator it_second = it; 117 + std::list<zero*>::iterator it_second = it;
113 integrale_sol integ = compute_integral(*it_first, *it_second); 118 integrale_sol integ = compute_integral(*it_first, *it_second);
114 - for(auto elem : res_tmp){  
115 - elem+=integ; 119 + std::list<integrale_sol>::iterator it_tmp = res_tmp.begin();
  120 + while (it_tmp != res_tmp.end()){
  121 + *it_tmp = *it_tmp + integ;
  122 + it_tmp++;
116 } 123 }
117 res_tmp.push_back(integ); 124 res_tmp.push_back(integ);
118 125
119 } else { 126 } else {
120 - zero bound(ibex::Interval(range.ub()), 0); 127 + zero *bound = new zero(ibex::Interval(range.ub()), 0);
121 integrale_sol integ = compute_integral(*it, bound); 128 integrale_sol integ = compute_integral(*it, bound);
122 - for(auto elem : res_tmp){  
123 - elem+=integ; 129 + std::list<integrale_sol>::iterator it_tmp = res_tmp.begin();
  130 + while (it_tmp != res_tmp.end()){
  131 + (*it_tmp) = (*it_tmp) + integ;
  132 + it_tmp++;
124 } 133 }
125 it++; 134 it++;
126 res_tmp.push_back(integ); 135 res_tmp.push_back(integ);
@@ -151,39 +160,39 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib @@ -151,39 +160,39 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib
151 construct_zeros_sup(); 160 construct_zeros_sup();
152 } 161 }
153 162
154 - std::list<zero>::iterator it = _zeros_sup->begin();  
155 - zero bound_sup(ibex::Interval(range.ub()), 0);  
156 - std::cout << "inf bound: " << *it << endl;  
157 - std::cout << "sup bound: " << bound_sup << endl; 163 + std::list<zero*>::iterator it = _zeros_sup->begin();
  164 + zero * bound_sup = new zero(ibex::Interval(range.ub()), 0);
  165 +// std::cout << "inf bound: " << *it << endl;
  166 +// std::cout << "sup bound: " << bound_sup << endl;
158 integrale_sol integ0 = compute_integral((*it), bound_sup); 167 integrale_sol integ0 = compute_integral((*it), bound_sup);
159 res_temp.push_back(integ0); 168 res_temp.push_back(integ0);
160 169
161 while (it != _zeros_sup->end()){ 170 while (it != _zeros_sup->end()){
162 if (*it != _zeros_sup->back()){ 171 if (*it != _zeros_sup->back()){
163 - std::list<zero>::iterator it_second = it; 172 + std::list<zero*>::iterator it_second = it;
164 it++; 173 it++;
165 - std::list<zero>::iterator it_first = it; 174 + std::list<zero*>::iterator it_first = it;
166 // std::cout << "inf bound: " << *it_first << endl; 175 // std::cout << "inf bound: " << *it_first << endl;
167 // std::cout << "sup bound: " << *it_second << endl; 176 // std::cout << "sup bound: " << *it_second << endl;
168 integrale_sol integ = compute_integral(*it_first, *it_second); 177 integrale_sol integ = compute_integral(*it_first, *it_second);
169 - for(auto elem : res_temp){  
170 -// cout << "elem : " << elem << endl;  
171 -// cout << "integ : " << integ << endl;  
172 - elem+=integ;  
173 -// cout << "new elem : " << elem << endl; 178 + std::list<integrale_sol>::iterator it_tmp = res_temp.begin();
  179 + while (it_tmp != res_temp.end()){
  180 + *it_tmp = *it_tmp + integ;
  181 + it_tmp++;
174 } 182 }
175 res_temp.push_back(integ); 183 res_temp.push_back(integ);
176 184
177 } else { 185 } else {
178 - zero bound(ibex::Interval(range.lb()), 0); 186 +// cout << "on passe ici" << endl;
  187 + zero * bound = new zero(ibex::Interval(range.lb()), 0);
179 // std::cout << "inf bound: " << bound << endl; 188 // std::cout << "inf bound: " << bound << endl;
180 // std::cout << "sup bound: " << *it << endl; 189 // std::cout << "sup bound: " << *it << endl;
181 integrale_sol integ = compute_integral(bound, *it); 190 integrale_sol integ = compute_integral(bound, *it);
182 - for(auto elem : res_temp){  
183 -// cout << "elem : " << elem << endl;  
184 -// cout << "integ : " << integ << endl;  
185 - elem+=integ;  
186 -// cout << "new elem : " << elem << endl; 191 +// cout << integ << endl;
  192 + std::list<integrale_sol>::iterator it_tmp = res_temp.begin();
  193 + while (it_tmp != res_temp.end()){
  194 + *it_tmp = *it_tmp + integ;
  195 + it_tmp++;
187 } 196 }
188 it++; 197 it++;
189 res_temp.push_back(integ); 198 res_temp.push_back(integ);
@@ -196,7 +205,7 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib @@ -196,7 +205,7 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib
196 integrale_sol min = *(res_temp.begin()); 205 integrale_sol min = *(res_temp.begin());
197 integrale_sol max = *(res_temp.begin()); 206 integrale_sol max = *(res_temp.begin());
198 for (auto elem : res_temp){ 207 for (auto elem : res_temp){
199 - cout << "\tcandidate : " << elem << endl; 208 +// cout << "\tcandidate : " << elem << endl;
200 if (elem.get_res_integrale().lb() < min.get_res_integrale().lb()){ 209 if (elem.get_res_integrale().lb() < min.get_res_integrale().lb()){
201 min = elem; 210 min = elem;
202 } 211 }
@@ -207,6 +216,38 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib @@ -207,6 +216,38 @@ std::pair&lt;integrale_sol, integrale_sol&gt; integrale::compute_integral_sup(const ib
207 return pair<integrale_sol, integrale_sol>(min, max); 216 return pair<integrale_sol, integrale_sol>(min, max);
208 } 217 }
209 218
  219 +std::pair<integrale_sol, integrale_sol> integrale::compute_integral()
  220 +{
  221 + std::pair<integrale_sol, integrale_sol> res = compute_integral_inf(_inf_bound);
  222 +// std::cout << "in inf : " << endl;
  223 +// std::cout << "\tmin : " << res.first << endl;
  224 +// std::cout << "\tmax : " << res.second << endl;
  225 + zero * inf = res.first.get_sup_bound();
  226 +
  227 + std::pair<integrale_sol, integrale_sol> res_sup = compute_integral_sup(_sup_bound);
  228 +// std::cout << "in sup : " << endl;
  229 +// std::cout << "\tmin : " << res_sup.first << endl;
  230 +// std::cout << "\tmax : " << res_sup.second << endl;
  231 +
  232 + zero * sup = res_sup.first.get_inf_bound();
  233 +
  234 + ibex::Interval mid = compute_integral(_inf_bound.ub(), _sup_bound.lb());
  235 + cout << *inf << endl;
  236 + cout << *sup << endl;
  237 + integrale_sol integ_mid(inf, sup, mid);
  238 +// cout << "integrale mid : " << integ_mid << endl;
  239 + res.first+= integ_mid;
  240 + res.second+= integ_mid;
  241 +// std::cout << "mid : " << mid << endl;
  242 + res.first+= res_sup.first;
  243 + res.second+= res_sup.second;
  244 + if (_result == nullptr){
  245 + _result = new std::pair<integrale_sol, integrale_sol>(res);
  246 + } else {
  247 + *_result = res;
  248 + }
  249 + return res;
  250 +}
210 251
211 252
212 253
@@ -24,24 +24,35 @@ private: @@ -24,24 +24,35 @@ private:
24 ibex::Function * _f; 24 ibex::Function * _f;
25 ibex::Variable * _x_var; 25 ibex::Variable * _x_var;
26 26
27 - std::list<zero> * _zeros_inf;  
28 - std::list<zero> * _zeros_sup; 27 + std::list<zero*> * _zeros_inf;
  28 + std::list<zero*> * _zeros_sup;
29 29
  30 + std::pair<integrale_sol, integrale_sol> * _result;
30 31
  32 + std::list<zero*> * construct_zeros_intern(const ibex::Interval & inf_bound);
31 33
32 - std::list<zero> * construct_zeros_intern(const ibex::Interval & inf_bound); 34 + integrale_sol compute_integral(zero* inf, zero* sup);
  35 +
  36 + std::list<ibex::Interval> compute_integral(const double & inf_bound, const ibex::Interval & sup_bound);
  37 +
  38 + ibex::Interval compute_integral(const double inf_bound, const double sup_bound);
  39 +
  40 + std::pair<integrale_sol, integrale_sol> compute_integral_inf(const ibex::Interval & range);
  41 +
  42 + std::pair<integrale_sol, integrale_sol> compute_integral_sup(const ibex::Interval & range);
33 43
34 public: 44 public:
35 //CONSTRUCTORS 45 //CONSTRUCTORS
36 - integrale(ibex::Interval inf, ibex::Interval sup, ibex::Function * f, ibex::Variable * var); 46 + integrale(const ibex::Interval& inf, const ibex::Interval& sup, ibex::Function * f, ibex::Variable * var);
37 47
38 //GETTERS 48 //GETTERS
39 inline ibex::Interval get_inf_bound() const {return _inf_bound;}; 49 inline ibex::Interval get_inf_bound() const {return _inf_bound;};
40 inline ibex::Interval get_sup_bound() const {return _sup_bound;}; 50 inline ibex::Interval get_sup_bound() const {return _sup_bound;};
41 inline ibex::Function * get_f() const {return _f;}; 51 inline ibex::Function * get_f() const {return _f;};
42 inline ibex::Variable * get_x_var() const {return _x_var;}; 52 inline ibex::Variable * get_x_var() const {return _x_var;};
43 - inline std::list<zero> * get_zeros_inf(){return _zeros_inf;};  
44 - inline std::list<zero> * get_zeros_sup(){return _zeros_sup;}; 53 + inline std::list<zero*> * get_zeros_inf(){return _zeros_inf;};
  54 + inline std::list<zero*> * get_zeros_sup(){return _zeros_sup;};
  55 + inline std::pair<integrale_sol, integrale_sol> * get_result(){return _result;};
45 56
46 unsigned int construct_zeros_inf(); 57 unsigned int construct_zeros_inf();
47 58
@@ -49,15 +60,7 @@ public: @@ -49,15 +60,7 @@ public:
49 60
50 unsigned int construct_zeros(); 61 unsigned int construct_zeros();
51 62
52 - integrale_sol compute_integral(const zero & inf, const zero & sup);  
53 -  
54 - std::list<ibex::Interval> compute_integral(const double & inf_bound, const ibex::Interval & sup_bound);  
55 -  
56 - ibex::Interval compute_integral(const double & inf_bound, const double & sup_bound);  
57 -  
58 - std::pair<integrale_sol, integrale_sol> compute_integral_inf(const ibex::Interval & range);  
59 -  
60 - std::pair<integrale_sol, integrale_sol> compute_integral_sup(const ibex::Interval & range); 63 + std::pair<integrale_sol, integrale_sol> compute_integral();
61 }; 64 };
62 65
63 66
src/integrale_sol.cpp
@@ -18,10 +18,10 @@ @@ -18,10 +18,10 @@
18 18
19 #include "integrale_sol.h" 19 #include "integrale_sol.h"
20 20
21 -integrale_sol::integrale_sol(const zero& inf, const zero& sup, const ibex::Interval& res) : _res_integrale(res) 21 +integrale_sol::integrale_sol(zero* inf, zero* sup, const ibex::Interval& res) : _res_integrale(res)
22 { 22 {
23 - _inf = new zero(inf);  
24 - _sup = new zero(sup); 23 + _inf = (inf);
  24 + _sup = (sup);
25 _min_candidate = false; 25 _min_candidate = false;
26 _max_candidate = false; 26 _max_candidate = false;
27 if (res.ub() < 0.){ 27 if (res.ub() < 0.){
@@ -39,12 +39,12 @@ integrale_sol &amp; integrale_sol::operator+=(const integrale_sol &amp;rhs){ @@ -39,12 +39,12 @@ integrale_sol &amp; integrale_sol::operator+=(const integrale_sol &amp;rhs){
39 39
40 integrale_sol operator+(const integrale_sol& lhs, const integrale_sol& rhs) 40 integrale_sol operator+(const integrale_sol& lhs, const integrale_sol& rhs)
41 { 41 {
42 - if (lhs.get_sup_bound()->get_x().lb() == rhs.get_sup_bound()->get_x().lb()){  
43 - cout << "ici" << endl;  
44 - return integrale_sol(*lhs.get_inf_bound(), *rhs.get_sup_bound(), (lhs.get_res_integrale() + rhs.get_res_integrale())); 42 + if (lhs.get_sup_bound() == rhs.get_inf_bound()){
  43 +// cout << "ici" << endl;
  44 + return integrale_sol(&*lhs.get_inf_bound(), &*rhs.get_sup_bound(), (lhs.get_res_integrale() + rhs.get_res_integrale()));
45 } else { 45 } else {
46 - cout << "la" << endl;  
47 - return integrale_sol(*rhs.get_inf_bound(), *lhs.get_sup_bound(), (lhs.get_res_integrale() + rhs.get_res_integrale())); 46 +// cout << "la" << endl;
  47 + return integrale_sol(&*rhs.get_inf_bound(), &*lhs.get_sup_bound(), (lhs.get_res_integrale() + rhs.get_res_integrale()));
48 } 48 }
49 49
50 } 50 }
src/integrale_sol.h
@@ -34,7 +34,7 @@ private: @@ -34,7 +34,7 @@ private:
34 bool _max_candidate; 34 bool _max_candidate;
35 35
36 public: 36 public:
37 - integrale_sol(const zero & inf, const zero & sup, const ibex::Interval & res); 37 + integrale_sol(zero *inf, zero *sup, const ibex::Interval & res);
38 38
39 inline zero * get_inf_bound() const {return _inf;}; 39 inline zero * get_inf_bound() const {return _inf;};
40 inline zero * get_sup_bound() const {return _sup;}; 40 inline zero * get_sup_bound() const {return _sup;};
@@ -15,21 +15,14 @@ void example(){ @@ -15,21 +15,14 @@ void example(){
15 integrale example(inf, sup, &f, &x); 15 integrale example(inf, sup, &f, &x);
16 16
17 cout << example << endl; 17 cout << example << endl;
18 -  
19 example.construct_zeros(); 18 example.construct_zeros();
20 -// pair<integrale_sol, integrale_sol> res = example.compute_integral_inf(inf);  
21 -// cout << "minimum: " << res.first << endl;  
22 -// cout << "maximum: " << res.second << endl;  
23 -// for (auto elem : *(example.get_zeros_inf())){  
24 -// cout << elem << endl;  
25 -// }  
26 - cout << "======================================================" << endl;  
27 - pair<integrale_sol, integrale_sol> res2 = example.compute_integral_sup(sup);  
28 - cout << "minimum: " << res2.first << endl;  
29 - cout << "maximum: " << res2.second << endl;  
30 - for (auto elem : *(example.get_zeros_sup())){ 19 + for (auto elem : *example.get_zeros_sup()){
31 cout << elem << endl; 20 cout << elem << endl;
32 } 21 }
  22 + example.compute_integral();
  23 + cout << "minimum: " << example.get_result()->first << endl;
  24 + cout << "maximum: " << example.get_result()->second << endl;
  25 +
33 } 26 }
34 27
35 28