diff --git a/.gitignore b/.gitignore index da1517d60..4397b172c 100644 --- a/.gitignore +++ b/.gitignore @@ -42,3 +42,4 @@ newdoc/_build *~ +/Debug/ diff --git a/.travis.yml b/.travis.yml index 47f9b9620..e1511dddb 100644 --- a/.travis.yml +++ b/.travis.yml @@ -25,7 +25,7 @@ before_install: install: - cd plugins/optim/; tar xvf soplex-1.7.1.tar && cd soplex-1.7.1 && make install && cd ../../.. - - ./waf configure --with-optim --with-soplex=$PWD/plugins/optim/soplex-1.7.1 --with-affine + - ./waf configure --with-optim --with-soplex=$PWD/plugins/optim/soplex-1.7.1 --with-affine --with-optim - sudo ./waf install - export PKG_CONFIG_PATH=/usr/local/share/pkgconfig/ - cd tests/ && make utest && cd .. @@ -35,5 +35,5 @@ install: script: # - __build__/examples/optimizer04 benchs/benchs-optim/coconutbenchmark-library1/ex8_5_2.bch acidhc4 compo smearsumrel 1.e-8 1.e-8 100 1 - cd tests/ && ./utest && cd .. -# - cd plugins/optim/tests/ && ./nonreg nonreg-travis-soplex-filib.res && cd ../../.. + - cd plugins/optim/tests/ && ./nonreg nonreg-travis-soplex-filib.res && cd ../../.. diff --git a/examples/test_ampl.cpp b/examples/test_ampl.cpp index c4c2c327c..773ddacbe 100644 --- a/examples/test_ampl.cpp +++ b/examples/test_ampl.cpp @@ -28,13 +28,13 @@ int main(int argc, char** argv) { ExtendedSystem ext_sys(*sys,1.e-6); cout<<"Extended System"<< endl; cout<(); +} + +void CellCostMaxPFub::set_optim_data(Cell& c, System& sys) { + c.get().compute_pf(*sys.goal,c.box); +} + +double CellCostMaxPFub::cost(const Cell& c) const { + const OptimData *data = &(c.get()); + if (data) { + return -data->pf.ub(); + } else { + ibex_error("CellCostMaxPFub::cost : invalid cost"); return POS_INFINITY; + } +} + } // end namespace ibex diff --git a/plugins/optim/src/strategy/ibex_CellCostFunc.h b/plugins/optim/src/strategy/ibex_CellCostFunc.h index 7d350de78..2ac26b264 100644 --- a/plugins/optim/src/strategy/ibex_CellCostFunc.h +++ b/plugins/optim/src/strategy/ibex_CellCostFunc.h @@ -293,5 +293,26 @@ class CellCostPFub: public CellCostFunc { }; +// ----------------------------------------------------------------------------------------------------------------------------------- + + +class CellCostMaxPFub: public CellCostFunc { +public: + CellCostMaxPFub(); + + /** + * \brief Add OptimData + */ + virtual void add_backtrackable(Cell& root); + + /** + * \brief Set "pf" in OptimData in the cell + */ + virtual void set_optim_data(Cell& c, System& sys); + + /** The "cost" of a element. */ + virtual double cost(const Cell& c) const; + +}; } // end namespace ibex #endif // __IBEX_CELL_COST_FUNC__ diff --git a/plugins/optim/src/strategy/ibex_CellDoubleHeap.h b/plugins/optim/src/strategy/ibex_CellDoubleHeap.h index 09c3e99ee..32fdc5107 100644 --- a/plugins/optim/src/strategy/ibex_CellDoubleHeap.h +++ b/plugins/optim/src/strategy/ibex_CellDoubleHeap.h @@ -103,10 +103,15 @@ inline Cell* CellDoubleHeap::pop() { return DoubleHeap::pop inline Cell* CellDoubleHeap::top() const { return DoubleHeap::top(); } inline std::ostream& CellDoubleHeap::print(std::ostream& os) const - { os << "==============================================================================\n"; - os << " first heap " << " size " << heap1->size() << " top " << heap1->top()->box << std::endl; - os << " second heap " << " size " << heap2->size() << " top " << heap2->top()->box ; - return os << std::endl; + { + os << "==============================================================================\n"; + if (empty()){ + return os << " EMPTY heap" << std::endl; + } else { + os << " first heap " << " size " << heap1->size() << " top " << heap1->top()->box << std::endl; + os << " second heap " << " size " << heap2->size() << " top " << heap2->top()->box ; + return os << std::endl; + } } //{return DoubleHeap::print( os);} diff --git a/plugins/optim/src/strategy/ibex_Multipliers.h b/plugins/optim/src/strategy/ibex_Multipliers.h index 0aabe4b9d..e05bf350b 100644 --- a/plugins/optim/src/strategy/ibex_Multipliers.h +++ b/plugins/optim/src/strategy/ibex_Multipliers.h @@ -54,9 +54,18 @@ class Multipliers : public Backtrackable { std::pair down(); IntervalVector lambda; + protected: + /** + * \brief Constructor by copy. + */ + explicit Multipliers(const Multipliers& e); + + /** + * \brief Create a copy + */ + Backtrackable* copy() const { return new Multipliers(*this);}; - Multipliers(const Multipliers&); }; diff --git a/plugins/optim/src/strategy/ibex_Optim.cpp b/plugins/optim/src/strategy/ibex_Optim.cpp new file mode 100644 index 000000000..c55109336 --- /dev/null +++ b/plugins/optim/src/strategy/ibex_Optim.cpp @@ -0,0 +1,193 @@ +// I B E X +// File : ibex_Optim.cpp +// Author : Gilles Chabert, Bertrand Neveu +// Copyright : Ecole des Mines de Nantes (France) +// License : See the LICENSE file +// Created : May 14, 2012 +// Last Update : December 24, 2012 +//============================================================================ + +#include "ibex_Optim.h" +#include +#include "ibex_Timer.h" + +namespace ibex { + +const double Optim::default_prec = 1e-07; +const double Optim::default_goal_rel_prec = 1e-07; +const double Optim::default_goal_abs_prec = 1e-07; +const int Optim::default_sample_size = 10; +const double Optim::default_equ_eps = 1e-08; +const double Optim::default_loup_tolerance = 0.1; + + +Optim::Optim (int n, CellDoubleHeap* buffer, double prec, + double goal_rel_prec, double goal_abs_prec, int sample_size) : + n(n), buffer(buffer), + prec(prec), goal_rel_prec(goal_rel_prec), goal_abs_prec(goal_abs_prec), + sample_size(sample_size), trace(false), + timeout(1e08), time(0), + loup(POS_INFINITY),uplo(NEG_INFINITY), + loup_point(n), nb_cells(0), + loup_changed(false), initial_loup(POS_INFINITY), + uplo_of_epsboxes(POS_INFINITY) { + + if (trace) std::cout.precision(12); +} + +Optim::~Optim() { +} + + +// compute the value ymax (decreasing the loup with the precision) +// the heap and the current box are contracted with y <= ymax +double Optim::compute_ymax() { + double ymax= loup - goal_rel_prec*fabs(loup); + if (loup - goal_abs_prec < ymax) + ymax = loup - goal_abs_prec; + return ymax;} + + + +void Optim::update_uplo() { + double new_uplo=POS_INFINITY; + + if (! buffer->empty()) { + new_uplo= buffer->minimum(); + if (new_uplo > loup) { + std::cout << " loup = " << loup << " new_uplo=" << new_uplo << std::endl; + ibex_error("Optim: new_uplo>loup (please report bug)"); + } + if (new_uplo < uplo) { + std::cout << "uplo= " << uplo << " new_uplo=" << new_uplo << std::endl; + ibex_error("Optim: new_uplo uplo) { + //std::cout << "uplo update=" << uplo << std::endl; + uplo = new_uplo; + } + } + else uplo= uplo_of_epsboxes; + } + else if (buffer->empty() && loup != POS_INFINITY) { + // empty buffer : new uplo is set to ymax (loup - precision) if a loup has been found + new_uplo=compute_ymax(); // not new_uplo=loup, because constraint y <= ymax was enforced + // std::cout << " new uplo buffer empty " << new_uplo << " uplo " << uplo << std::endl; + + double m = (new_uplo< uplo_of_epsboxes)? new_uplo : uplo_of_epsboxes; + if (uplo < m) uplo = m; // warning: hides the field "m" of the class + // note: we always have uplo <= uplo_of_epsboxes but we may have uplo > new_uplo, because + // ymax is strictly lower than the loup. + } + +} + + + +void Optim::update_uplo_of_epsboxes(double ymin) { + + // the current box cannot be bisected. ymin is a lower bound of the objective on this box + // uplo of epsboxes can only go down, but not under uplo : it is an upperbound for uplo, + //that indicates a lowerbound for the objective in all the small boxes + // found by the precision criterion + assert (uplo_of_epsboxes >= uplo); + assert(ymin >= uplo); + if (uplo_of_epsboxes > ymin) { + uplo_of_epsboxes = ymin; + if (trace) { + std::cout << "uplo_of_epsboxes:" << uplo_of_epsboxes << " uplo " << uplo << std::endl; + } + } +} + + + +void Optim::report() { + + if (timeout >0 && time >=timeout ) { + std::cout << "time limit " << timeout << "s. reached " << std::endl; + } + // No solution found and optimization stopped with empty buffer before the required precision is reached => means infeasible problem + if (buffer->empty() && uplo_of_epsboxes == POS_INFINITY && (loup==POS_INFINITY || (loup==initial_loup && goal_abs_prec==0 && goal_rel_prec==0))) { + std::cout << " infeasible problem " << std::endl; + std::cout << " cpu time used " << time << "s." << std::endl; + std::cout << " number of cells " << nb_cells << std::endl; + } + + else { + std::cout << " best bound in: [" << uplo << "," << loup << "]" << std::endl; + + double rel_prec; + + if (loup==POS_INFINITY) + rel_prec= POS_INFINITY; + else + rel_prec=(loup-uplo)/(fabs (loup))-1.e-15; + + double abs_prec=loup-uplo-1.e-15; + + std::cout << " Relative precision obtained on objective function: " << rel_prec << " " << + (rel_prec <= goal_rel_prec? " [passed]" : " [failed]") << " " << goal_rel_prec << std::endl; + + std::cout << " Absolute precision obtained on objective function: " << abs_prec << " " << + (abs_prec <= goal_abs_prec? " [passed]" : " [failed]") << " " << goal_abs_prec << std::endl; + if (uplo_of_epsboxes != NEG_INFINITY && uplo_of_epsboxes != POS_INFINITY) + std::cout << " precision on variable domains obtained " << prec << " " << " uplo_of_epsboxes " << uplo_of_epsboxes << std::endl; + else if (uplo_of_epsboxes == NEG_INFINITY) + std::cout << " small boxes with negative infinity objective : objective not bound " << std::endl; + if (loup==initial_loup) + std::cout << " no feasible point found " << std::endl; + else + std::cout << " best feasible point " << loup_point << std::endl; + + + std::cout << " cpu time used " << time << "s." << std::endl; + std::cout << " number of cells " << nb_cells << std::endl; + } + /* // statistics on upper bounding + if (trace) { + std::cout << " nbrand " << nb_rand << " nb_inhc4 " << nb_inhc4 << " nb simplex " << nb_simplex << std::endl; + std::cout << " diam_rand " << diam_rand << " diam_inhc4 " << diam_inhc4 << " diam_simplex " << diam_simplex << std::endl; + } + */ +} +/* minimal report for benchmarking */ +void Optim::time_cells_report() { + if (timeout >0 && time >=timeout ) { + std::cout << "timeout" << timeout << " " << uplo << " " << loup << " ";} + else + std::cout << time << " " ; + std::cout << nb_cells << std::endl; +} + + +void Optim::report_perf() { + + double rel_prec; + if (loup==POS_INFINITY) + rel_prec= POS_INFINITY; + else + rel_prec=(loup-uplo)/(fabs(loup))-1.e-15; + + double abs_prec=loup-uplo-1.e-15; + + std::cout << ( ((rel_prec <= goal_rel_prec)|| + (abs_prec <= goal_abs_prec)|| + ((buffer->empty() && uplo_of_epsboxes == POS_INFINITY && loup==POS_INFINITY))|| + (uplo<-1.e300) + )? " T & " : " F & " ); + + std::cout << uplo << " & " << loup << " & "; + std::cout << time << " "<< std::endl ; +} + +void Optim::time_limit_check () { + Timer::stop(); + time += Timer::VIRTUAL_TIMELAPSE(); + if (timeout >0 && time >=timeout ) throw TimeOutException(); + Timer::start(); +} + +} // end namespace ibex diff --git a/plugins/optim/src/strategy/ibex_Optim.h b/plugins/optim/src/strategy/ibex_Optim.h new file mode 100644 index 000000000..c00b5ddb0 --- /dev/null +++ b/plugins/optim/src/strategy/ibex_Optim.h @@ -0,0 +1,222 @@ +//============================================================================ +// I B E X +// File : ibex_Optimizer.h +// Author : Gilles Chabert, Bertrand Neveu +// Copyright : Ecole des Mines de Nantes (France) +// License : See the LICENSE file +// Created : May 14, 2012 +// Last Update : May 14, 2012 +//============================================================================ + +#ifndef __IBEX_OPTIM__ +#define __IBEX_OPTIM__ + +#include "ibex_Interval.h" +#include "ibex_CellDoubleHeap.h" + +namespace ibex { + +/** + * \ingroup strategy + * + * \brief Optimizer. + * + */ +class Optim { +public: + + Optim( int n, CellDoubleHeap* buffer, double prec=default_prec, + double goal_rel_prec=default_goal_rel_prec, double goal_abs_prec=default_goal_abs_prec, + int sample_size=default_sample_size); + + /** + * \brief Delete *this. + */ + virtual ~Optim(); + + /** + * \brief Return status of the optimizer + */ + typedef enum {SUCCESS, INFEASIBLE, NO_FEASIBLE_FOUND, UNBOUNDED_OBJ, TIME_OUT} Status; + + /** + * \brief Run the optimization. + * + * \param init_box The initial box + * \param obj_init_bound (optional) can be set when an initial upper bound of the objective minimum is known a priori. + * This bound can be obtained, e.g., by a local solver. This is equivalent to (but more practical + * than) adding a constraint f(x)<=obj_init_bound. + * + * \return SUCCESS If the global minimum (with respect to the precision required) has been found. + * In particular, at least one feasible point has been found, less than obj_init_bound, and in the time limit. + * + * INFEASIBLE if no feasible point exist less than obj_init_bound. In particular, the function returns INFEASIBLE + * if the initial bound "obj_init_bound" is LESS than the true minimum (this case is only possible if + * goal_abs_prec and goal_rel_prec are 0). In the latter case, there may exist feasible points. + * + * NO_FEASIBLE_FOUND if no feasible point could be found less than obj_init_bound. Contrary to INFEASIBLE, + * infeasibility is not proven here. Warning: this return value is sensitive to the goal_abs_prec and + * goal_rel_prec parameters. The upperbounding makes the optimizer only looking for points less than + * min { (1-goal_rel_prec)*obj_init_bound, obj_init_bound - goal_abs_prec }. + * + * UNBOUNDED_OBJ the objective function seems unbounded (tends to -oo). + * + * TIMEOUT time is out. + */ + virtual Status optimize(const IntervalVector& init_box, double obj_init_bound=POS_INFINITY)=0; + + /** + * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). + * + * Information provided: + *
  • interval of the cost [uplo,loup] + *
  • the best feasible point found + *
  • total running time + *
  • total number of cells created during the exploration + *
+ */ + void report(); + + /** + * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). + * + * Information provided: + *
  • interval of the cost [uplo,loup] in case of termination due to timelimit + *
  • total running time + *
  • total number of cells created during the exploration + *
+ */ + + void time_cells_report(); + + /** + * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). + * + * Information provided: + *
  • interval of the cost [uplo,loup] + *
  • total running time + *
+ */ + void report_perf(); + + + /** Number of variables. */ + const int n; + + /** Cell buffers. + Two buffers are used for node selection. the first one corresponds to minimize the minimum of the objective estimate, + the second one to minimize another criterion (by default the maximum of the objective estimate). + The second one is chosen at each node with a probability critpr/100 (default value critpr=50) + */ + CellDoubleHeap *buffer; + + /** + * \brief Index of the goal variable y in the extended box. + * + */ + + /** Precision (bisection control) */ + const double prec; + + /** Relative precision on the objective */ + const double goal_rel_prec; + + /** Absolute precision on the objective */ + const double goal_abs_prec; + + /** Number of samples used to update the loup */ + const int sample_size; + + /** Trace activation flag. + * The value can be fixed by the user. By default: 0 nothing is printed + 1 for printing each better found feasible point + 2 for printing each handled node */ + int trace; + + /** + * \brief Time limit. + * + * Maximum CPU time used by the strategy. + * This parameter allows to bound time consumption. + * The value can be fixed by the user. + */ + double timeout; + + /* Remember running time of the last exploration */ + double time; + + void time_limit_check(); + + /** Default bisection precision: 1e-07 */ + static const double default_prec; + + /** Default goal relative precision */ + static const double default_goal_rel_prec; + + /** Default goal absolute precision */ + static const double default_goal_abs_prec; + + /** Default sample size */ + static const int default_sample_size; + + /** Default epsilon applied to equations */ + static const double default_equ_eps; + + /** Default tolerance increase ratio for the pseudo-loup. */ + static const double default_loup_tolerance; + /** + * \brief The "loup" (lowest upper bound of the criterion) + * + * In rigor mode, represents the real-loup (not the pseudo-loup). + */ + double loup; + + /** The "uplo" (uppermost lower bound of the criterion) */ + double uplo; + + /** The point satisfying the constraints corresponding to the loup */ + Vector loup_point; + + /** Number of cells put into the heap (which passed through the contractors) */ + int nb_cells; + +protected: + + /** + * \brief Update the uplo of non bisectable boxes + */ + void update_uplo_of_epsboxes(double ymin); + + /** + * \brief Update the uplo + */ + void update_uplo(); + + /** + * \brief Computes and returns the value ymax (the loup decreased with the precision) + * the heap and the current box are actually contracted with y <= ymax + * + */ + double compute_ymax (); + + + bool loup_changed; + + /** + * \brief The bound on the objective given by the user, +oo otherwise. + * + * Used to see if at least a loup-point has been found. + * + */ + double initial_loup; + + + /** Lower bound of the small boxes taken by the precision */ + double uplo_of_epsboxes; + + +}; + + +} // end namespace ibex +#endif // __IBEX_OPTIM__ diff --git a/plugins/optim/src/strategy/ibex_OptimData.h b/plugins/optim/src/strategy/ibex_OptimData.h index 38923c4b0..7c995af25 100644 --- a/plugins/optim/src/strategy/ibex_OptimData.h +++ b/plugins/optim/src/strategy/ibex_OptimData.h @@ -68,7 +68,16 @@ class OptimData : public Backtrackable { protected: - OptimData(const OptimData&); + /** + * \brief Constructor by copy. + */ + explicit OptimData(const OptimData& e); + + /** + * \brief Create a copy + */ + Backtrackable* copy() const { return new OptimData(*this);}; + }; } // end namespace ibex diff --git a/plugins/optim/src/strategy/ibex_Optimizer.cpp b/plugins/optim/src/strategy/ibex_Optimizer.cpp index 795de24cc..83efd04e1 100644 --- a/plugins/optim/src/strategy/ibex_Optimizer.cpp +++ b/plugins/optim/src/strategy/ibex_Optimizer.cpp @@ -32,13 +32,6 @@ using namespace std; namespace ibex { -const double Optimizer::default_prec = 1e-07; -const double Optimizer::default_goal_rel_prec = 1e-07; -const double Optimizer::default_goal_abs_prec = 1e-07; -const int Optimizer::default_sample_size = 10; -const double Optimizer::default_equ_eps = 1e-08; -const double Optimizer::default_loup_tolerance = 0.1; - void Optimizer::write_ext_box(const IntervalVector& box, IntervalVector& ext_box) { int i2=0; for (int i=0; i0? new CtcUnion(sys) : NULL; // ============================================================= - if (trace) cout.precision(12); - // objshaver= new CtcOptimShaving (*new CtcHC4 (ext_sys.ctrs,0.1,true),20,1,1.e-11); - int niter=100; - // int niter=1000; - if (niter < 3*n) niter=3*n; //==================================== #ifdef _IBEX_WITH_NOLP_ mylp = NULL; #else + + int niter=100; + // int niter=1000; + if (niter < 3*n) niter=3*n; //lr = new LinearRelaxCombo(sys, LinearRelaxCombo::XNEWTON); //mylp = new LinearSolver(sys.nb_var,sys.nb_ctr,niter); lr = new LinearRelaxCombo(ext_sys, LinearRelaxCombo::XNEWTON); @@ -120,25 +110,18 @@ Optimizer::~Optimizer() { delete is_inside; } - buffer.flush(); if (equs) delete equs; if (df) delete df; delete mylp; delete lr; - delete &buffer.cost1(); - delete &buffer.cost2(); + buffer->flush(); + delete &(buffer->cost1()); + delete &(buffer->cost2()); + delete buffer; // delete &(objshaver->ctc); // delete objshaver; } -// compute the value ymax (decreasing the loup with the precision) -// the heap and the current box are contracted with y <= ymax -double Optimizer::compute_ymax() { - double ymax= loup - goal_rel_prec*fabs(loup); - if (loup - goal_abs_prec < ymax) - ymax = loup - goal_abs_prec; - return ymax;} - // launch Hansen test bool Optimizer::update_real_loup() { @@ -230,42 +213,6 @@ double minimum (double a, double b) { else return b; } -void Optimizer::update_uplo() { - double new_uplo=POS_INFINITY; - - if (! buffer.empty()) { - new_uplo= buffer.minimum(); - if (new_uplo > loup) { - cout << " loup = " << loup << " new_uplo=" << new_uplo << endl; - ibex_error("optimizer: new_uplo>loup (please report bug)"); - } - if (new_uplo < uplo) { - cout << "uplo= " << uplo << " new_uplo=" << new_uplo << endl; - ibex_error("optimizer: new_uplo uplo) { - //cout << "uplo update=" << uplo << endl; - uplo = new_uplo; - } - } - else uplo= uplo_of_epsboxes; - } - else if (buffer.empty() && loup != POS_INFINITY) { - // empty buffer : new uplo is set to ymax (loup - precision) if a loup has been found - new_uplo=compute_ymax(); // not new_uplo=loup, because constraint y <= ymax was enforced - // cout << " new uplo buffer empty " << new_uplo << " uplo " << uplo << endl; - - double m = minimum(new_uplo, uplo_of_epsboxes); - if (uplo < m) uplo = m; // warning: hides the field "m" of the class - // note: we always have uplo <= uplo_of_epsboxes but we may have uplo > new_uplo, because - // ymax is strictly lower than the loup. - } - -} - - /* contract the box of the cell c , try to find a new loup ; push the cell in the 2 heaps or if the contraction makes the box empty, delete the cell. @@ -281,10 +228,10 @@ void Optimizer::handle_cell(Cell& c, const IntervalVector& init_box ){ // objshaver->contract(c.box); // we know cost1() does not require OptimData - buffer.cost2().set_optim_data(c,sys); + buffer->cost2().set_optim_data(c,sys); // the cell is put into the 2 heaps - buffer.push(&c); + buffer->push(&c); nb_cells++; } @@ -405,7 +352,7 @@ void Optimizer::firstorder_contract(IntervalVector& box, const IntervalVector& i Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj_init_bound) { loup=obj_init_bound; pseudo_loup=obj_init_bound; - buffer.contract(loup); + //buffer.contract(loup); // ??JN il y a un flush just apres?? uplo=NEG_INFINITY; uplo_of_epsboxes=POS_INFINITY; @@ -416,7 +363,7 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj nb_rand=0; diam_rand=0; - buffer.flush(); + buffer->flush(); Cell* root=new Cell(IntervalVector(n+1)); @@ -426,7 +373,7 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj bsc.add_backtrackable(*root); // add data "pu" and "pf" (if required) - buffer.cost2().add_backtrackable(*root); + buffer->cost2().add_backtrackable(*root); // add data required by optimizer + Fritz John contractor root->add(); @@ -444,21 +391,21 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj update_uplo(); try { - while (!buffer.empty()) { + while (!buffer->empty()) { // if (trace >= 2) cout << " buffer " << buffer << endl; - if (trace >= 2) buffer.print(cout); + if (trace >= 2) buffer->print(cout); // cout << "buffer size " << buffer.size() << " " << buffer2.size() << endl; loup_changed=false; - Cell *c = buffer.top(); + Cell *c = buffer->top(); try { - pair boxes=bsc.bisect(*c); + //pair boxes=bsc.bisect(*c); + //pair new_cells=c->bisect(boxes.first,boxes.second); + pair new_cells=bsc.bisect_cell(*c); - pair new_cells=c->bisect(boxes.first,boxes.second); - - buffer.pop(); + buffer->pop(); delete c; // deletes the cell. handle_cell(*new_cells.first, init_box); @@ -477,7 +424,7 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj double ymax=compute_ymax(); - buffer.contract(ymax); + buffer->contract(ymax); //cout << " now buffer is contracted and min=" << buffer.minimum() << endl; @@ -493,7 +440,7 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj } catch (NoBisectableVariableException& ) { update_uplo_of_epsboxes((c->box)[ext_sys.goal_var()].lb()); - buffer.pop(); + buffer->pop(); delete c; // deletes the cell. //if (trace>=1) cout << "epsilon-box found: uplo cannot exceed " << uplo_of_epsboxes << endl; update_uplo(); // the heap has changed -> recalculate the uplo @@ -518,108 +465,4 @@ Optimizer::Status Optimizer::optimize(const IntervalVector& init_box, double obj return SUCCESS; } -void Optimizer::update_uplo_of_epsboxes(double ymin) { - - // the current box cannot be bisected. ymin is a lower bound of the objective on this box - // uplo of epsboxes can only go down, but not under uplo : it is an upperbound for uplo, - //that indicates a lowerbound for the objective in all the small boxes - // found by the precision criterion - assert (uplo_of_epsboxes >= uplo); - assert(ymin >= uplo); - if (uplo_of_epsboxes > ymin) { - uplo_of_epsboxes = ymin; - if (trace) { - cout << "uplo_of_epsboxes:" << setprecision(12) << uplo_of_epsboxes << " uplo " << uplo << endl; - } - } -} - - - -void Optimizer::report() { - - if (timeout >0 && time >=timeout ) { - cout << "time limit " << timeout << "s. reached " << endl; - } - // No solution found and optimization stopped with empty buffer before the required precision is reached => means infeasible problem - if (buffer.empty() && uplo_of_epsboxes == POS_INFINITY && (loup==POS_INFINITY || (loup==initial_loup && goal_abs_prec==0 && goal_rel_prec==0))) { - cout << " infeasible problem " << endl; - cout << " cpu time used " << time << "s." << endl; - cout << " number of cells " << nb_cells << endl; - } - - else { - cout << " best bound in: [" << uplo << "," << loup << "]" << endl; - - double rel_prec; - - if (loup==POS_INFINITY) - rel_prec= POS_INFINITY; - else - rel_prec=(loup-uplo)/(fabs (loup))-1.e-15; - - double abs_prec=loup-uplo-1.e-15; - - cout << " Relative precision obtained on objective function: " << rel_prec << " " << - (rel_prec <= goal_rel_prec? " [passed]" : " [failed]") << " " << goal_rel_prec << endl; - - cout << " Absolute precision obtained on objective function: " << abs_prec << " " << - (abs_prec <= goal_abs_prec? " [passed]" : " [failed]") << " " << goal_abs_prec << endl; - if (uplo_of_epsboxes != NEG_INFINITY && uplo_of_epsboxes != POS_INFINITY) - cout << " precision on variable domains obtained " << prec << " " << " uplo_of_epsboxes " << uplo_of_epsboxes << endl; - else if (uplo_of_epsboxes == NEG_INFINITY) - cout << " small boxes with negative infinity objective : objective not bound " << endl; - if (loup==initial_loup) - cout << " no feasible point found " << endl; - else - cout << " best feasible point " << loup_point << endl; - - - cout << " cpu time used " << time << "s." << endl; - cout << " number of cells " << nb_cells << endl; - } - /* // statistics on upper bounding - if (trace) { - cout << " nbrand " << nb_rand << " nb_inhc4 " << nb_inhc4 << " nb simplex " << nb_simplex << endl; - cout << " diam_rand " << diam_rand << " diam_inhc4 " << diam_inhc4 << " diam_simplex " << diam_simplex << endl; - } - */ -} -/* minimal report for benchmarking */ -void Optimizer::time_cells_report() { - if (timeout >0 && time >=timeout ) { - cout << "timeout" << timeout << " " << uplo << " " << loup << " ";} - else - cout << time << " " ; - cout << nb_cells << endl; -} - - -void Optimizer::report_perf() { - - double rel_prec; - if (loup==POS_INFINITY) - rel_prec= POS_INFINITY; - else - rel_prec=(loup-uplo)/(fabs(loup))-1.e-15; - - double abs_prec=loup-uplo-1.e-15; - - cout << ( ((rel_prec <= goal_rel_prec)|| - (abs_prec <= goal_abs_prec)|| - ((buffer.empty() && uplo_of_epsboxes == POS_INFINITY && loup==POS_INFINITY))|| - (uplo<-1.e300) - )? " T & " : " F & " ); - - cout << uplo << " & " << loup << " & "; - cout << time << " "<< endl ; -} - -void Optimizer::time_limit_check () { - Timer::stop(); - time += Timer::VIRTUAL_TIMELAPSE(); - if (timeout >0 && time >=timeout ) throw TimeOutException(); - Timer::start(); -} - } // end namespace ibex diff --git a/plugins/optim/src/strategy/ibex_Optimizer.h b/plugins/optim/src/strategy/ibex_Optimizer.h index 346dd3601..c30406fac 100644 --- a/plugins/optim/src/strategy/ibex_Optimizer.h +++ b/plugins/optim/src/strategy/ibex_Optimizer.h @@ -11,6 +11,7 @@ #ifndef __IBEX_OPTIMIZER_H__ #define __IBEX_OPTIMIZER_H__ +#include "ibex_Optim.h" #include "ibex_Bsc.h" #include "ibex_CtcHC4.h" #include "ibex_Ctc3BCid.h" @@ -25,6 +26,7 @@ #include "ibex_PdcHansenFeasibility.h" #include "ibex_Random.h" #include "ibex_LinearRelaxCombo.h" +#include "ibex_IntervalVector.h" namespace ibex { @@ -40,7 +42,7 @@ namespace ibex { * \remark In all the comments of this class, "loup" means "lowest upper bound" of the criterion f * and "uplo" means "uppermost lower bound" of the criterion. */ -class Optimizer { +class Optimizer : public Optim { public: /** * \brief Create an optimizer. @@ -78,11 +80,6 @@ class Optimizer { */ virtual ~Optimizer(); - /** - * \brief Return status of the optimizer - */ - typedef enum {SUCCESS, INFEASIBLE, NO_FEASIBLE_FOUND, UNBOUNDED_OBJ, TIME_OUT} Status; - /** * \brief Run the optimization. * @@ -109,40 +106,6 @@ class Optimizer { */ Status optimize(const IntervalVector& init_box, double obj_init_bound=POS_INFINITY); - /** - * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). - * - * Information provided: - *
  • interval of the cost [uplo,loup] - *
  • the best feasible point found - *
  • total running time - *
  • total number of cells created during the exploration - *
- */ - void report(); - - /** - * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). - * - * Information provided: - *
  • interval of the cost [uplo,loup] in case of termination due to timelimit - *
  • total running time - *
  • total number of cells created during the exploration - *
- */ - - void time_cells_report(); - - /** - * \brief Displays on standard output a report of the last call to #optimize(const IntervalVector&). - * - * Information provided: - *
  • interval of the cost [uplo,loup] - *
  • total running time - *
- */ - void report_perf(); - /** * \brief The original system * @@ -158,8 +121,6 @@ class Optimizer { */ NormalizedSystem sys; - /** Number of variables. */ - const int n; /** Number of constraints. */ const int m; @@ -192,30 +153,6 @@ class Optimizer { /** Bisector. */ Bsc& bsc; - /** Cell buffers. - Two buffers are used for node selection. the first one corresponds to minimize the minimum of the objective estimate, - the second one to minimize another criterion (by default the maximum of the objective estimate). - The second one is chosen at each node with a probability critpr/100 (default value critpr=50) - */ - CellDoubleHeap buffer; - - /** - * \brief Index of the goal variable y in the extended box. - * - */ - - /** Precision (bisection control) */ - const double prec; - - /** Relative precision on the objective */ - const double goal_rel_prec; - - /** Absolute precision on the objective */ - const double goal_abs_prec; - - /** Number of samples used to update the loup */ - const int sample_size; - /** Flag for applying monotonicity analysis. * The value can be fixed by the user. By default: true. */ bool mono_analysis_flag; @@ -225,56 +162,6 @@ class Optimizer { * The value can be fixed by the user. By default: true. */ bool in_HC4_flag; - /** Trace activation flag. - * The value can be fixed by the user. By default: 0 nothing is printed - 1 for printing each better found feasible point - 2 for printing each handled node */ - int trace; - - /** - * \brief Time limit. - * - * Maximum CPU time used by the strategy. - * This parameter allows to bound time consumption. - * The value can be fixed by the user. - */ - double timeout; - - /* Remember running time of the last exploration */ - double time; - - void time_limit_check(); - - /** Default bisection precision: 1e-07 */ - static const double default_prec; - - /** Default goal relative precision */ - static const double default_goal_rel_prec; - - - - - - - - /** Default goal absolute precision */ - static const double default_goal_abs_prec; - - /** Default sample size */ - static const int default_sample_size; - - /** Default epsilon applied to equations */ - static const double default_equ_eps; - - /** Default tolerance increase ratio for the pseudo-loup. */ - static const double default_loup_tolerance; - /** - * \brief The "loup" (lowest upper bound of the criterion) - * - * In rigor mode, represents the real-loup (not the pseudo-loup). - */ - double loup; - /** * \brief The pseudo-loup. * @@ -282,18 +169,9 @@ class Optimizer { */ double pseudo_loup; - /** The "uplo" (uppermost lower bound of the criterion) */ - double uplo; - - /** The point satisfying the constraints corresponding to the loup */ - Vector loup_point; - /** Rigor mode: the box satisfying the constraints corresponding to the loup */ IntervalVector loup_box; - /** Number of cells put into the heap (which passed through the contractors) */ - int nb_cells; - protected: /** * \brief Return an upper bound of f(x). @@ -352,16 +230,6 @@ class Optimizer { */ bool update_entailed_ctr(const IntervalVector& box); - /** - * \brief Update the uplo of non bisectable boxes - */ - void update_uplo_of_epsboxes(double ymin); - - /** - * \brief Update the uplo - */ - void update_uplo(); - /** * \brief Main procedure for updating the loup. @@ -511,24 +379,7 @@ class Optimizer { */ Function* df; - /** - * \brief Computes and returns the value ymax (the loup decreased with the precision) - * the heap and the current box are actually contracted with y <= ymax - * - */ - double compute_ymax (); - - bool loup_changed; - - /** - * \brief The bound on the objective given by the user, +oo otherwise. - * - * Used to see if at least a loup-point has been found. - * - */ - double initial_loup; - - Ctc3BCid* objshaver; + //Ctc3BCid* objshaver; private: @@ -542,9 +393,6 @@ class Optimizer { /** Inner contractor (for the negation of g) */ CtcUnion* is_inside; - /** Lower bound of the small boxes taken by the precision */ - double uplo_of_epsboxes; - /** Currently entailed constraints */ EntailedCtr* entailed; diff --git a/plugins/optim/tests/TestCell2.cpp b/plugins/optim/tests/TestCell2.cpp new file mode 100644 index 000000000..44f216484 --- /dev/null +++ b/plugins/optim/tests/TestCell2.cpp @@ -0,0 +1,265 @@ +//============================================================================ +// I B E X +// File : TestCell2.cpp +// Author : Jordan Ninin +// Copyright : ENSTA Bretagne (France) +// License : See the LICENSE file +// Created : Oct 10, 2016 +// Last Update : Oct 10, 2016 +//============================================================================ + +#include "TestCell2.h" +#include "ibex_Cell.h" +#include "ibex_EntailedCtr.h" +#include "ibex_Bsc.h" +#include "ibex_LargestFirst.h" +#include "ibex_SystemFactory.h" +#include "ibex_NormalizedSystem.h" +#include "ibex_Cell.h" +#include "ibex_Multipliers.h" +#include "ibex_OptimData.h" + +//using namespace std; + +namespace ibex { + + +void TestCell2::test01() { + IntervalVector box (2, Interval(-1,1)); + Cell * c = new Cell(box); + + Cell * copy =new Cell(*c); + check(copy->box,c->box); + + LargestFirst bsc; + std::pair new_cells = bsc.bisect_cell(*c); + + + check(new_cells.first->box|new_cells.second->box,c->box); + delete c; + + check(new_cells.first->box|new_cells.second->box,copy->box); + delete new_cells.first; + delete new_cells.second; + delete copy; + +} + +void TestCell2::test02() { + IntervalVector box (2, Interval(-1,1)); + Cell * root = new Cell(box); + + Variable x,y; + Function func(x,y,pow(cos(y)+cos(2*y+x),2)); + Function c1f(x,y,y-x*(x+6.28) ); + NumConstraint c1(c1f,LEQ); + Function c2f(x,y,y-x*(x-6.28)); + NumConstraint c2(c2f,LEQ); + + SystemFactory fac; + fac.add_var(x); + fac.add_var(y); + fac.add_ctr(c1); + fac.add_ctr(c2); + NormalizedSystem sys(fac); + + root->add(); + EntailedCtr * entailed = &root->get(); + entailed->init_root(sys,sys); + entailed->original(1) = true; + entailed->normalized(0) = true; + + root->add(); + Multipliers * multi = &root->get(); + multi->init_root(1,1,1); + + + Cell * copy =new Cell(*root); + check(copy->box,root->box); + CPPUNIT_ASSERT(entailed->original(1)); + CPPUNIT_ASSERT(entailed->normalized(0)); + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + check(copy->get().lambda,multi->lambda); + + LargestFirst bsc; + std::pair new_cells = bsc.bisect_cell(*root); + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + check(new_cells.first->get().lambda,multi->lambda); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + check(new_cells.second->get().lambda,multi->lambda); + + + check(new_cells.first->box|new_cells.second->box,root->box); + delete root; + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + check(new_cells.first->get().lambda,copy->get().lambda); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + check(new_cells.second->get().lambda,copy->get().lambda); + + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + + delete copy; + delete new_cells.first; + delete new_cells.second; + +} + + +void TestCell2::test03() { + IntervalVector box (2, Interval(-1,1)); + Cell * root = new Cell(box); + + Variable x,y; + Function func(x,y,pow(cos(y)+cos(2*y+x),2)); + Function c1f(x,y,y-x*(x+6.28) ); + NumConstraint c1(c1f,LEQ); + Function c2f(x,y,y-x*(x-6.28)); + NumConstraint c2(c2f,LEQ); + + SystemFactory fac; + fac.add_var(x); + fac.add_var(y); + fac.add_ctr(c1); + fac.add_ctr(c2); + NormalizedSystem sys(fac); + + root->add(); + EntailedCtr * entailed = &root->get(); + entailed->init_root(sys,sys); + entailed->original(1) = true; + entailed->normalized(0) = true; + + root->add(); + Multipliers * multi = &root->get(); + multi->init_root(1,1,1); + + + root->add(); + OptimData * data = &root->get(); + data->pf= Interval(5,6); + data->pu=10; + + + Cell * copy =new Cell(*root); + + check(copy->box,root->box); + CPPUNIT_ASSERT(entailed->original(1)); + CPPUNIT_ASSERT(entailed->normalized(0)); + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + check(copy->get().lambda,multi->lambda); + check(copy->get().pf,data->pf); + + + LargestFirst bsc; + std::pair new_cells = bsc.bisect_cell(*root); + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + check(new_cells.first->get().lambda,multi->lambda); + check(new_cells.first->get().pf,data->pf); + CPPUNIT_ASSERT(new_cells.first->get().pu==10); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + check(new_cells.second->get().lambda,multi->lambda); + check(new_cells.second->get().pf,data->pf); + CPPUNIT_ASSERT(new_cells.second->get().pu==10); + + + data->pf= Interval(1,10); + data->pu=144; + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + check(new_cells.first->get().lambda,copy->get().lambda); + check(new_cells.first->get().pf,copy->get().pf); + CPPUNIT_ASSERT(new_cells.first->get().pu==10); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + check(new_cells.second->get().lambda,copy->get().lambda); + check(new_cells.second->get().pf,copy->get().pf); + CPPUNIT_ASSERT(new_cells.second->get().pu==10); + + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + CPPUNIT_ASSERT(copy->get().pf==Interval(5,6)); + CPPUNIT_ASSERT(copy->get().pu==10); + + delete root; + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + check(new_cells.first->get().lambda,copy->get().lambda); + check(new_cells.first->get().pf,copy->get().pf); + CPPUNIT_ASSERT(new_cells.first->get().pu==10); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + check(new_cells.second->get().lambda,copy->get().lambda); + check(new_cells.second->get().pf,copy->get().pf); + CPPUNIT_ASSERT(new_cells.second->get().pu==10); + + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + CPPUNIT_ASSERT(copy->get().pf==Interval(5,6)); + CPPUNIT_ASSERT(copy->get().pu==10); + + delete copy; + delete new_cells.first; + delete new_cells.second; + +} + + + +} // end namespace + + + diff --git a/plugins/optim/tests/TestCell2.h b/plugins/optim/tests/TestCell2.h new file mode 100644 index 000000000..fbad64bd9 --- /dev/null +++ b/plugins/optim/tests/TestCell2.h @@ -0,0 +1,42 @@ +/* ============================================================================ + * I B E X - CtcCell2 Tests + * ============================================================================ + * Copyright : ENSTA Bretagne (FRANCE) + * License : This program can be distributed under the terms of the GNU LGPL. + * See the file COPYING.LESSER. + * + * Author(s) : Jordan Ninin + * Created : Oct 10, 2016 + * ---------------------------------------------------------------------------- */ + +#ifndef __TEST_CELL2_H__ +#define __TEST_CELL2_H__ + +#include +#include + +#include "utils.h" + +namespace ibex { + +class TestCell2 : public CppUnit::TestFixture { + +public: + + CPPUNIT_TEST_SUITE(TestCell2); + CPPUNIT_TEST(test01); + CPPUNIT_TEST(test02); + CPPUNIT_TEST(test03); + CPPUNIT_TEST_SUITE_END(); + + void test01(); + void test02(); + void test03(); + +}; + +CPPUNIT_TEST_SUITE_REGISTRATION(TestCell2); + +} // namespace ibex + +#endif // __TEST_CELL_H__ diff --git a/plugins/optim/tests/TestCellHeap.cpp b/plugins/optim/tests/TestCellHeap.cpp index 3e67f99aa..6fca40beb 100644 --- a/plugins/optim/tests/TestCellHeap.cpp +++ b/plugins/optim/tests/TestCellHeap.cpp @@ -54,7 +54,7 @@ void TestCellHeap::test02() { h2.push(cell); } - h2.pop(); h2.pop(); + delete h2.pop(); delete h2.pop(); CPPUNIT_ASSERT(h2.size()==nb-2); } @@ -220,6 +220,7 @@ void TestCellHeap::test_D03() { delete h1.pop(); delete h1.pop(); CPPUNIT_ASSERT(h1.size()==((unsigned int)(nb-2))); + h1.flush(); } @@ -292,6 +293,8 @@ void TestCellHeap::test_D05() { h1.flush(); CPPUNIT_ASSERT(h1.size()==0); + h2.flush(); + CPPUNIT_ASSERT(h2.size()==0); } diff --git a/plugins/optim/tests/TestDoubleHeap.cpp b/plugins/optim/tests/TestDoubleHeap.cpp deleted file mode 100644 index 7de6c9eff..000000000 --- a/plugins/optim/tests/TestDoubleHeap.cpp +++ /dev/null @@ -1,168 +0,0 @@ -/* ============================================================================ - * I B E X - CellHeap Tests - * ============================================================================ - * Copyright : Ecole des Mines de Nantes (FRANCE) - * License : This program can be distributed under the terms of the GNU LGPL. - * See the file COPYING.LESSER. - * - * Author(s) : Jordan Ninin - * Created : Mar 2, 2014 - * ---------------------------------------------------------------------------- */ - -#include "TestDoubleHeap.h" -#include "ibex_DoubleHeap.h" - -using namespace std; - -namespace ibex { - - - -class TestCostFunc1 : public CostFunc { -public: - TestCostFunc1() { } ; - double cost(const Interval& c) const { return c.diam();}; -}; - -class TestCostFunc2 : public CostFunc { -public: - TestCostFunc2() { } ; - double cost(const Interval& c) const { return c.lb();}; -}; - -class TestCostFunc3 : public CostFunc { -public: - TestCostFunc3(): loup(10) { } ; - double cost(const Interval& c) const { return c.ub()*loup;}; - void set_loup(double lb) { loup = lb; } - -private: - double loup; -}; - - -void TestDoubleHeap::test01() { - - int nb= 10; - TestCostFunc1 costf1; - TestCostFunc2 costf2; - - DoubleHeap h(costf1,false,costf2,false,50); - - for (int i=1; i<=nb ;i++) { - if ((i%2)==1) h.push(new Interval(i,2*i)); - else h.push(new Interval(i,20*i)); - } - - CPPUNIT_ASSERT(h.minimum1()==1); - CPPUNIT_ASSERT(h.minimum2()==1); - - h.pop1(); - h.pop1(); - CPPUNIT_ASSERT(h.minimum1()==5); - CPPUNIT_ASSERT(h.minimum2()==2); - - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==5); - CPPUNIT_ASSERT(h.minimum2()==4); - - h.contract(20); - CPPUNIT_ASSERT(h.minimum1()==5); - CPPUNIT_ASSERT(h.minimum2()==5); - - h.pop2(); - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==9); - CPPUNIT_ASSERT(h.minimum2()==9); - - CPPUNIT_ASSERT(h.size()==1); - - h.push(new Interval(10,11)); - h.push(new Interval(12,14)); - h.push(new Interval(12,12)); - h.push(new Interval(12,13.1)); - - CPPUNIT_ASSERT(h.minimum1()==0); - CPPUNIT_ASSERT(h.minimum2()==9); - CPPUNIT_ASSERT(h.size()==5); - - h.pop1(); - CPPUNIT_ASSERT(h.minimum1()==1); - CPPUNIT_ASSERT(h.minimum2()==9); - CPPUNIT_ASSERT(h.size()==4); - - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==1); - CPPUNIT_ASSERT(h.minimum2()==10); - CPPUNIT_ASSERT(h.size()==3); - - h.flush(); - CPPUNIT_ASSERT(h.size()==0); -} - - -void TestDoubleHeap::test02() { - - int nb= 10; - TestCostFunc2 costf2; - TestCostFunc3 costf3; - - DoubleHeap h(costf2,false,costf3,true,50); - - for (int i=1; i<=nb ;i++) { - if ((i%2)==1) h.push(new Interval(i,2*i)); - else h.push(new Interval(i,20*i)); - } - - CPPUNIT_ASSERT(h.minimum1()==1); - CPPUNIT_ASSERT(h.minimum2()==20); - - h.pop2(); - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==2); - CPPUNIT_ASSERT(h.minimum2()==100); - - h.pop1(); - CPPUNIT_ASSERT(h.minimum1()==4); - CPPUNIT_ASSERT(h.minimum2()==100); - - // it is necessary to update the loup. - costf3.set_loup(100); - h.contract(7.5); - CPPUNIT_ASSERT(h.minimum1()==4); - CPPUNIT_ASSERT(h.minimum2()==1000); - CPPUNIT_ASSERT(h.size()==4); - - - h.pop2(); - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==4); - CPPUNIT_ASSERT(h.minimum2()==8000); - CPPUNIT_ASSERT(h.size()==2); - - h.push(new Interval(10,11)); - h.push(new Interval(12,14)); - h.push(new Interval(12,12)); - h.push(new Interval(12,13.1)); - - CPPUNIT_ASSERT(h.minimum1()==4); - CPPUNIT_ASSERT(h.minimum2()==1100); - CPPUNIT_ASSERT(h.size()==6); - - h.pop1(); - CPPUNIT_ASSERT(h.minimum1()==6); - CPPUNIT_ASSERT(h.minimum2()==1100); - CPPUNIT_ASSERT(h.size()==5); - - h.pop2(); - CPPUNIT_ASSERT(h.minimum1()==6); - CPPUNIT_ASSERT(h.minimum2()==1200); - CPPUNIT_ASSERT(h.size()==4); - - h.flush(); - CPPUNIT_ASSERT(h.size()==0); -} - - - -} // end namespace diff --git a/plugins/optim/tests/TestOptimizer.cpp b/plugins/optim/tests/TestOptimizer.cpp index cf9ed59a4..f9552ab48 100644 --- a/plugins/optim/tests/TestOptimizer.cpp +++ b/plugins/optim/tests/TestOptimizer.cpp @@ -19,7 +19,7 @@ using namespace std; namespace ibex { // true minimum is 0. -Optimizer::Status issue50(double init_loup, double prec) { +Optim::Status issue50(double init_loup, double prec) { SystemFactory f; const ExprSymbol& x=ExprSymbol::new_(); f.add_var(x); @@ -34,19 +34,19 @@ Optimizer::Status issue50(double init_loup, double prec) { } void TestOptimizer::issue50_1() { - CPPUNIT_ASSERT(issue50(1e-10, 0.1)==Optimizer::NO_FEASIBLE_FOUND); + CPPUNIT_ASSERT(issue50(1e-10, 0.1)==Optim::NO_FEASIBLE_FOUND); } void TestOptimizer::issue50_2() { - CPPUNIT_ASSERT(issue50(1e-10, 0)==Optimizer::SUCCESS); + CPPUNIT_ASSERT(issue50(1e-10, 0)==Optim::SUCCESS); } void TestOptimizer::issue50_3() { - CPPUNIT_ASSERT(issue50(-1e-10, 0.1)==Optimizer::NO_FEASIBLE_FOUND); + CPPUNIT_ASSERT(issue50(-1e-10, 0.1)==Optim::NO_FEASIBLE_FOUND); } void TestOptimizer::issue50_4() { - CPPUNIT_ASSERT(issue50(-1e-10, 0)==Optimizer::INFEASIBLE); + CPPUNIT_ASSERT(issue50(-1e-10, 0)==Optim::INFEASIBLE); } diff --git a/plugins/optim/tests/TestUnconstrainedLocalSearch.cpp b/plugins/optim/tests/TestUnconstrainedLocalSearch.cpp index 5a7279caa..474440468 100644 --- a/plugins/optim/tests/TestUnconstrainedLocalSearch.cpp +++ b/plugins/optim/tests/TestUnconstrainedLocalSearch.cpp @@ -54,7 +54,7 @@ void TestUnconstrainedLocalSearch::almost_diag() { int max_iter=1000; Vector xk(n); UnconstrainedLocalSearch::ReturnCode ret=o.minimize(x0,xk,eps,max_iter); - cout << " xk=" << xk << endl; + //cout << " xk=" << xk << endl; CPPUNIT_ASSERT(ret==UnconstrainedLocalSearch::SUCCESS); CPPUNIT_ASSERT(almost_eq(IntervalVector(xk),IntervalVector(xsol), 1.e-9)); diff --git a/plugins/optim/tests/nonreg.cpp b/plugins/optim/tests/nonreg.cpp index acc5a571a..ff535521a 100644 --- a/plugins/optim/tests/nonreg.cpp +++ b/plugins/optim/tests/nonreg.cpp @@ -137,7 +137,7 @@ int main (int argc, char** argv) { Optimizer o(p.get_sys(), p.get_ctc(), p.get_bsc(), p.prec, p.goal_rel_prec, p.goal_abs_prec, p.sample_size, p.eq_eps); - Optimizer::Status status=o.optimize(p.get_sys().box); + Optim::Status status=o.optimize(p.get_sys().box); double scaled_time = ratio_perf*time; double time_gain_ratio = (o.time-scaled_time)/scaled_time; @@ -145,11 +145,11 @@ int main (int argc, char** argv) { //cerr << "number of cells=" << o.nb_cells << " time=" << o.time << endl; switch (status) { - case Optimizer::INFEASIBLE : cerr << "FAILED: infeasible"; break; - case Optimizer::NO_FEASIBLE_FOUND : cerr << "FAILED: no feasible point found"; break; - case Optimizer::UNBOUNDED_OBJ : cerr << "FAILED: unbounded objective"; break; - case Optimizer::TIME_OUT : cerr << "FAILED: timeout"; break; - case Optimizer::SUCCESS : { + case Optim::INFEASIBLE : cerr << "FAILED: infeasible"; break; + case Optim::NO_FEASIBLE_FOUND : cerr << "FAILED: no feasible point found"; break; + case Optim::UNBOUNDED_OBJ : cerr << "FAILED: unbounded objective"; break; + case Optim::TIME_OUT : cerr << "FAILED: timeout"; break; + case Optim::SUCCESS : { if (o.loup < lb) { cerr.precision(20); cerr << "FAILED: upper bound (loup=" << o.loup << ") is wrong"; } else if (o.uplo > ub) { cerr.precision(20); cerr << "FAILED: lower bound (uplo=" << o.uplo << ") is wrong"; } else { diff --git a/plugins/wscript b/plugins/wscript index f31353e13..7bad156ed 100644 --- a/plugins/wscript +++ b/plugins/wscript @@ -13,8 +13,9 @@ def get_plugins(path): return map(basename,folders) def configure (conf): - - ################################################################################################## + + + ################################################################################################## # Enable AFFINE plugin if option --with-affine-extended is enable if (conf.options.WITH_AFFINE_EXTEND): conf.options.WITH_AFFINE =True diff --git a/src/arithmetic/ibex_Matrix.cpp b/src/arithmetic/ibex_Matrix.cpp index 6f28d0efd..3bf76cb18 100644 --- a/src/arithmetic/ibex_Matrix.cpp +++ b/src/arithmetic/ibex_Matrix.cpp @@ -134,4 +134,40 @@ std::ostream& operator<<(std::ostream& os, const Matrix& m) { return _displayM(os,m); } +double Matrix::det(){ + if (!is_square()) //TODO: return error message + return 0; + if(_nb_cols ==1 ) + return M[0][0]; + else + return compute_det(); +} + +double Matrix::compute_det(){ + double res(0); + if(_nb_cols ==2) + return M[0][0]*M[1][1]-M[1][0]*M[0][1]; + else { + for(int i=0;i<_nb_rows;i++) { + Matrix sub(_nb_rows-1,_nb_rows-1); + int row_it(0); + for(int j=0;j<_nb_rows;j++) { //create sub matrix + if(j!=i){ + sub.set_row(row_it,(this->row(j).subvector(1,_nb_cols-1))); + row_it++; + } + } + if(i%2 == 0) + res += M[i][0]*sub.compute_det(); + else + res -= M[i][0]*sub.compute_det(); + } + } + return res; +} + +bool Matrix::is_square(){ + return (_nb_cols == _nb_rows); +} // namespace ibex + } // namespace ibex diff --git a/src/arithmetic/ibex_Matrix.h b/src/arithmetic/ibex_Matrix.h index 05a564c8a..868f64464 100644 --- a/src/arithmetic/ibex_Matrix.h +++ b/src/arithmetic/ibex_Matrix.h @@ -155,6 +155,17 @@ class Matrix { */ void put(int row_start_index, int col_start_index, const Vector& v, bool row_vec); + /** + * \brief Compute determinant + */ + double det(); + + + /** + * \brief True if square matrix + */ + bool is_square(); + /** * \brief Set *this to (*this)+m. */ @@ -209,6 +220,11 @@ class Matrix { Matrix(); + /** + * \brief Compute determinant: recursive function + */ + double compute_det(); + int _nb_rows; int _nb_cols; Vector* M; diff --git a/src/bisector/ibex_Bsc.cpp b/src/bisector/ibex_Bsc.cpp index 4ec42cf20..383d7c53c 100644 --- a/src/bisector/ibex_Bsc.cpp +++ b/src/bisector/ibex_Bsc.cpp @@ -12,7 +12,6 @@ #include "ibex_Cell.h" #include "ibex_Exception.h" -using std::pair; namespace ibex { @@ -30,9 +29,6 @@ Bsc::Bsc(const Vector& prec) : _prec(prec) { if (prec[i]<=0) ibex_error("precision must be a nonnegative number"); } -pair Bsc::bisect(Cell& cell) { - return bisect(cell.box); -} void Bsc::add_backtrackable(Cell& root) { root.add(); diff --git a/src/bisector/ibex_Bsc.h b/src/bisector/ibex_Bsc.h index cf0d5dfe0..94e1ab4ab 100644 --- a/src/bisector/ibex_Bsc.h +++ b/src/bisector/ibex_Bsc.h @@ -16,8 +16,6 @@ namespace ibex { -class Cell; - /** * \defgroup bisector Bisectors */ @@ -70,11 +68,6 @@ class Bsc { */ virtual ~Bsc() { } - /** - * \brief Bisect the current box and return the result. - */ - virtual std::pair bisect(const IntervalVector& box)=0; - /** * \brief Bisect the current cell and return the result. * @@ -83,7 +76,12 @@ class Bsc { * Implementation is optional. By default, this function call bisect(cell.box). * See #bisect(const IntervalVector&). */ - virtual std::pair bisect(Cell& cell); + std::pair bisect_cell(const Cell& cell); + + /** + * \brief Bisect the current box and return the result. + */ + virtual std::pair bisect(const IntervalVector& box)=0; /** * Allows to add the backtrackable data required @@ -141,11 +139,23 @@ class BisectedVar : public Backtrackable { } int var; + +protected: + + explicit BisectedVar(const BisectedVar& e) : var(e.var) { } + Backtrackable* copy() const { return new BisectedVar(*this);}; + }; /*============================================ inline implementation ============================================ */ + +inline std::pair Bsc::bisect_cell(const Cell& cell) { + std::pair boxes=this->bisect(cell.box); + return cell.bisect(boxes.first,boxes.second); +} + inline bool Bsc::uniform_prec() const { return _prec.size()==1; } diff --git a/src/bisector/ibex_RoundRobin.cpp b/src/bisector/ibex_RoundRobin.cpp index 8de6e5a8b..4fa6a9839 100644 --- a/src/bisector/ibex_RoundRobin.cpp +++ b/src/bisector/ibex_RoundRobin.cpp @@ -48,12 +48,14 @@ pair RoundRobin::bisect(const IntervalVector& box return bisect(box,i); } -pair RoundRobin::bisect(Cell& cell) { +pair RoundRobin::bisect(Cell& cell) { BisectedVar& v=cell.get(); // the following instruction will update v.var // and the new value of v.var will be copied to child nodes - return bisect(cell.box,v.var); + + pair boxes=this->bisect(cell.box, v.var); + return cell.bisect(boxes.first,boxes.second); } void RoundRobin::add_backtrackable(Cell& root) { diff --git a/src/bisector/ibex_RoundRobin.h b/src/bisector/ibex_RoundRobin.h index 6810bedf7..e458c77f3 100644 --- a/src/bisector/ibex_RoundRobin.h +++ b/src/bisector/ibex_RoundRobin.h @@ -69,7 +69,7 @@ class RoundRobin : public Bsc { * bisected variable and this information is copied to * the subcells data. */ - virtual std::pair bisect(Cell& cell); + virtual std::pair bisect(Cell& cell); /** * \brief Add an instance of #ibex::BisectedVar to the backtrackle data of the root cell. diff --git a/src/set/ibex_Set.cpp b/src/set/ibex_Set.cpp index 00caadcbe..37590ad62 100644 --- a/src/set/ibex_Set.cpp +++ b/src/set/ibex_Set.cpp @@ -232,6 +232,7 @@ class NodeAndDist : public Backtrackable { public: NodeAndDist() : node(NULL), dist(-1) { } + NodeAndDist(SetNode* _node) : node(_node), dist(-1) { } /** @@ -258,6 +259,12 @@ class NodeAndDist : public Backtrackable { SetNode* node; double dist; + +protected: + + explicit NodeAndDist(const NodeAndDist& e) : node(e.node), dist(e.dist) { } // TODO JN: Chabs, can you check it? + Backtrackable* copy() const {return new NodeAndDist(*this);}; + }; /** diff --git a/src/strategy/ibex_Backtrackable.h b/src/strategy/ibex_Backtrackable.h index 9d900287f..8e63eaf1c 100644 --- a/src/strategy/ibex_Backtrackable.h +++ b/src/strategy/ibex_Backtrackable.h @@ -15,6 +15,8 @@ namespace ibex { +class Cell; + /** * \ingroup strategy * @@ -49,6 +51,15 @@ class Backtrackable { * \brief Delete *this. */ virtual ~Backtrackable() { } + +protected: + friend class Cell; + + /** + * \brief Create a copy + */ + virtual Backtrackable* copy() const =0; + }; } // end namespace ibex diff --git a/src/strategy/ibex_Cell.cpp b/src/strategy/ibex_Cell.cpp index e2569f57d..3bd72d560 100644 --- a/src/strategy/ibex_Cell.cpp +++ b/src/strategy/ibex_Cell.cpp @@ -22,14 +22,19 @@ int id_count=0; assert(id_count Cell::bisect(const IntervalVector& left, const IntervalVector& right) { - Cell* cleft = new Cell(left); - Cell* cright = new Cell(right); - for (IBEXMAP(Backtrackable*)::iterator it=data.begin(); it!=data.end(); it++) { - std::pair child_data=it->second->down(); - cleft->data.insert_new(it->first,child_data.first); - cright->data.insert_new(it->first,child_data.second); - } + Cell::Cell(const Cell& e) : box(e.box), id(id_count++) { + for ( IBEXMAP(Backtrackable*)::const_iterator it=e.data.begin(); it!=e.data.end(); it++) { + data.insert_new(it->first, it->second->copy()); + } + } + + + +std::pair Cell::bisect(const IntervalVector& left, const IntervalVector& right) const { + Cell* cleft = new Cell(*this); + Cell* cright = new Cell(*this); + cleft->box = left; + cright->box = right; return std::pair(cleft,cright); } diff --git a/src/strategy/ibex_Cell.h b/src/strategy/ibex_Cell.h index 63b0976ba..b0133ec0d 100644 --- a/src/strategy/ibex_Cell.h +++ b/src/strategy/ibex_Cell.h @@ -45,7 +45,13 @@ class Cell { * * \param box - Box (passed by copy). */ - Cell(const IntervalVector& box); + explicit Cell(const IntervalVector& box); + + /** + * \brief Constructor by copy. + * + */ + explicit Cell (const Cell& e); /** * \brief Bisect this cell. @@ -59,7 +65,7 @@ class Cell { * bisector class can simply bisect a box into two subboxes, the * cell bisection has a default implementation in #ibex::Bsc. */ - std::pair bisect(const IntervalVector& left, const IntervalVector& right); + std::pair bisect(const IntervalVector& left, const IntervalVector& right) const; /** * \brief Delete *this. diff --git a/src/strategy/ibex_EntailedCtr.cpp b/src/strategy/ibex_EntailedCtr.cpp index a331952e5..63eb9477c 100644 --- a/src/strategy/ibex_EntailedCtr.cpp +++ b/src/strategy/ibex_EntailedCtr.cpp @@ -13,7 +13,8 @@ namespace ibex { -EntailedCtr::EntailedCtr() : orig_sys(NULL), norm_sys(NULL) { +EntailedCtr::EntailedCtr() : orig_sys(NULL), norm_sys(NULL), + orig_entailed(NULL), norm_entailed(NULL) { } diff --git a/src/strategy/ibex_EntailedCtr.h b/src/strategy/ibex_EntailedCtr.h index e4a4508cd..e5b599fe8 100644 --- a/src/strategy/ibex_EntailedCtr.h +++ b/src/strategy/ibex_EntailedCtr.h @@ -28,6 +28,8 @@ class EntailedCtr : public Backtrackable { */ EntailedCtr(); + + /** * \brief Mark the constraints as not entailed (by default). * @@ -87,6 +89,11 @@ class EntailedCtr : public Backtrackable { protected: + /** + * \brief Constructor by copy. + */ + explicit EntailedCtr( const EntailedCtr& e); + const System* orig_sys; const NormalizedSystem* norm_sys; @@ -97,9 +104,13 @@ class EntailedCtr : public Backtrackable { bool *orig_entailed; bool *norm_entailed; - EntailedCtr(const EntailedCtr&); - friend std::ostream& operator<<(std::ostream& os, const EntailedCtr&); + + + /** + * \brief Create a copy + */ + Backtrackable* copy() const { return new EntailedCtr(*this); }; }; diff --git a/src/strategy/ibex_Paver.cpp b/src/strategy/ibex_Paver.cpp index 1f9c8c150..4c47ecc03 100644 --- a/src/strategy/ibex_Paver.cpp +++ b/src/strategy/ibex_Paver.cpp @@ -16,7 +16,7 @@ using namespace std; namespace ibex { Paver::Paver(const Array& c, Bsc& b, CellBuffer& buffer) : - capacity(-1), ctc_loop(true), ctc(c), bsc(b), buffer(buffer) { + capacity(-1), timeout(3600), ctc_loop(true), trace(false), ctc(c), bsc(b), buffer(buffer) { assert(ctc.size()>0); } @@ -73,8 +73,7 @@ void Paver::contract(Cell& cell, SubPaving* paving) { void Paver::bisect(Cell& c) { - pair boxes=bsc.bisect(c); - pair new_cells=c.bisect(boxes.first,boxes.second); + pair new_cells=bsc.bisect_cell(c); delete buffer.pop(); buffer.push(new_cells.first); diff --git a/src/strategy/ibex_Solver.cpp b/src/strategy/ibex_Solver.cpp index efe39ec66..0b2975f13 100644 --- a/src/strategy/ibex_Solver.cpp +++ b/src/strategy/ibex_Solver.cpp @@ -77,8 +77,7 @@ bool Solver::next(std::vector& sols) { try { - pair boxes=bsc.bisect(*c); - pair new_cells=c->bisect(boxes.first,boxes.second); + pair new_cells=bsc.bisect_cell(*c); delete buffer.pop(); buffer.push(new_cells.first); diff --git a/src/system/ibex_SystemFactory.cpp b/src/system/ibex_SystemFactory.cpp index 20c96651e..86c116f7b 100644 --- a/src/system/ibex_SystemFactory.cpp +++ b/src/system/ibex_SystemFactory.cpp @@ -108,6 +108,7 @@ void SystemFactory::add_ctr(const ExprCtr& ctr) { ctrs.push_back(new NumConstraint(*new Function(ctr_vars, ctr_expr), ctr.op, true)); } +/* void SystemFactory::add_ctr2(const ExprCtr& ctr) { if (!args) args = new Array(tmp_args); @@ -117,6 +118,7 @@ void SystemFactory::add_ctr2(const ExprCtr& ctr) { ctrs.push_back(new NumConstraint(*new Function(ctr_vars, ctr_expr), ctr.op, true)); } +*/ void SystemFactory::add_ctr(const NumConstraint& ctr) { init_arg_bound(); @@ -206,9 +208,7 @@ System::System(const SystemFactory& fac) : nb_var(0), nb_ctr(0), box(1) { void System::init(const SystemFactory& fac) { - // the field fac.args is initialized upon addition of an objective - // function or a constraint. - if (!fac.args) throw EmptySystemException(); + (int&) nb_var = fac.nb_var; (int&) nb_ctr = fac.ctrs.size(); @@ -218,13 +218,31 @@ void System::init(const SystemFactory& fac) { // f is not initialized here --> see init_f_from_ctrs() below - // =========== init args ============== - args.resize(fac.nb_arg); - varcopy(*fac.args, args); + // the field fac.args is initialized upon addition of an objective + // function or a constraint. + //if (!fac.args) throw EmptySystemException(); + if (!fac.args) { + // =========== init args ============== + args.resize(fac.nb_arg); + Array tmp(fac.tmp_args); + varcopy(tmp, args); + + box.resize(nb_var); + int i=0; + for (typename std::vector::const_iterator it=fac.tmp_bound.begin(); it!=fac.tmp_bound.end(); it++) { + box.put(i,*it); + i+=(*it).size(); + } - box.resize(nb_var); - box=fac.bound_init; + } else { + // =========== init args ============== + args.resize(fac.nb_arg); + varcopy(*fac.args, args); + box.resize(nb_var); + box=fac.bound_init; + + } // =========== init ctrs ============== ctrs.resize(nb_ctr); for (int i=0; i ctrs; -private: - void init_arg_bound(); }; diff --git a/src/tools/ibex_DoubleHeap.h b/src/tools/ibex_DoubleHeap.h index 5c0b60a2f..42caba1aa 100644 --- a/src/tools/ibex_DoubleHeap.h +++ b/src/tools/ibex_DoubleHeap.h @@ -35,6 +35,11 @@ class DoubleHeap { */ DoubleHeap(CostFunc& cost1, bool update_cost1_when_sorting, CostFunc& cost2, bool update_cost2_when_sorting, int critpr=50); + /** copy constructor **/ + explicit DoubleHeap(const DoubleHeap& dhcp); + + DoubleHeap* deepcopy() const; + /** * \brief Flush the buffer. * @@ -42,6 +47,13 @@ class DoubleHeap { */ void flush(); + /** + * \brief Clear the buffer. + * + * All the remaining data will be *removed* without being *deleted* + */ + void clear(); + /** \brief Return the size of the buffer. */ unsigned int size() const; @@ -139,6 +151,7 @@ class DoubleHeap { */ void contract_rec(double new_loup, HeapNode* node, SharedHeap& heap, bool percolate); +private: /** * Erase all the subnodes of node (including itself) in the first heap * and manage the impact on the second heap. @@ -153,6 +166,12 @@ class DoubleHeap { void erase_subnodes(HeapNode* node, bool percolate); std::ostream& print(std::ostream& os) const; + + // for deepcopy only + explicit DoubleHeap(int critproba); + + void delete_subnodes1(HeapNode* node) ; + }; @@ -166,22 +185,81 @@ DoubleHeap::DoubleHeap(CostFunc& cost1, bool update_cost1_when_sorting, Co } +template +DoubleHeap::DoubleHeap(int critproba) : nb_nodes(0), heap1(NULL), heap2(NULL), critpr(critproba), current_heap_id(0) { + +} + +template +DoubleHeap* DoubleHeap::deepcopy() const { + DoubleHeap* copy =new DoubleHeap(this->critpr); + copy->nb_nodes=this->nb_nodes; + copy->current_heap_id= this->current_heap_id; + std::pair *,std::vector*> *> p = this->heap1->deepcopy_sheap(2); + copy->heap1 = p.first; + copy->heap2 = new SharedHeap(this->heap2->costf,this->heap2->update_cost_when_sorting,this->heap2->heap_id); + while(!p.second->empty()) { + copy->heap2->push_elt(p.second->back()); + p.second->pop_back(); + } + delete p.second; + return copy; +} + +template +DoubleHeap::DoubleHeap(const DoubleHeap &dhcp):nb_nodes(dhcp.nb_nodes),heap1(NULL),heap2(NULL),critpr(dhcp.critpr),current_heap_id(dhcp.current_heap_id) { + std::pair *,std::vector*> *> p; + p= dhcp.heap1->copy_sheap(2); + heap1 = p.first; + heap2 = new SharedHeap(dhcp.heap2->costf,dhcp.heap2->update_cost_when_sorting,dhcp.heap2->heap_id); + while(!p.second->empty()) { + heap2->push_elt(p.second->back()); + p.second->pop_back(); + } + delete p.second; +} + template DoubleHeap::~DoubleHeap() { + clear(); if (heap1) delete heap1; if (heap2) delete heap2; + } template void DoubleHeap::flush() { if (nb_nodes>0) { - erase_subnodes(heap1->root,false); + delete_subnodes1(heap1->root); + heap1->nb_nodes=0; + heap1->root=NULL; + + heap2->flush(); + nb_nodes=0; + } +} + +template +void DoubleHeap::delete_subnodes1(HeapNode* node) { + if (node->left) delete_subnodes1(node->left); + if (node->right)delete_subnodes1(node->right); + delete node; +} + + +template +void DoubleHeap::clear() { + if (nb_nodes>0) { + delete_subnodes1(heap1->root); heap1->nb_nodes=0; heap1->root=NULL; + + heap2->clear(); nb_nodes=0; } } + template unsigned int DoubleHeap::size() const { assert(heap1->size()==heap2->size()); @@ -200,7 +278,8 @@ void DoubleHeap::contract(double new_loup1) { heap1->root = copy1->root; heap1->nb_nodes = copy1->size(); nb_nodes = copy1->size(); - copy1->root = NULL; // avoid to delete heap1 with copy1 + copy1->root = NULL; + copy1->nb_nodes=0;// avoid to delete heap1 with copy1 delete copy1; if (heap2->update_cost_when_sorting) heap2->sort(); @@ -241,6 +320,7 @@ void DoubleHeap::erase_subnodes(HeapNode* node, bool percolate) { else heap2->erase_node(node->elt->holder[1]); + node->elt->flush(); delete node->elt; delete node; } diff --git a/src/tools/ibex_Heap.h b/src/tools/ibex_Heap.h index 4059e597d..3252e77eb 100644 --- a/src/tools/ibex_Heap.h +++ b/src/tools/ibex_Heap.h @@ -55,6 +55,9 @@ class Heap { Heap(CostFunc& costf); + Heap(const Heap& original_heap); + + ~Heap() { flush(); } /** * \brief Flush the buffer. * @@ -131,6 +134,17 @@ Heap::Heap(CostFunc& costf) : costf(costf) { } +template +Heap::Heap(const Heap& original_heap) : costf(original_heap.costf) { + for(unsigned i=0;i p; + p.first = elem; + p.second = original_heap.l.at(i).second; + l.push_back(p); + } +} + template bool Heap::operator()(const std::pair& c1, const std::pair& c2) const { return c1.second >= c2.second; diff --git a/src/tools/ibex_SharedHeap.h b/src/tools/ibex_SharedHeap.h index dea41502b..b5cbdd593 100644 --- a/src/tools/ibex_SharedHeap.h +++ b/src/tools/ibex_SharedHeap.h @@ -65,6 +65,20 @@ class SharedHeap { /** \brief Return true if the buffer is empty. */ bool empty() const; + /** + * \brief Flush the buffer. + * + * All the remaining data will be *deleted* + */ + void flush(); + + /** + * \brief Clear the buffer. + * + * All the remaining data will be *removed* without being *deleted* + */ + void clear(); + /** * \brief Return the next box (but does not pop it). * @@ -163,6 +177,9 @@ class SharedHeap { */ HeapNode* get_node(unsigned int i) const; + + std::ostream& print(std::ostream& os) const; + /** * \brief Streams out the heap * @@ -176,10 +193,25 @@ class SharedHeap { */ bool heap_state(); + /** copy constructor + * Copy all the elements and the data + **/ + std::pair *,std::vector*> *> deepcopy_sheap(int nb_crit); + + /** copy constructor + * Copy only the structure without the data + **/ + std::pair *,std::vector*> *> copy_sheap(int nb_crit); + + private: /** Used in the sort function (proceed by recursivity) */ void sort_rec(HeapNode* node, SharedHeap& heap); + + + void flush_subnodes(HeapNode* node); + void clear_subnodes(HeapNode* node); }; @@ -231,6 +263,10 @@ class HeapNode { template friend std::ostream& operator<<(std::ostream& os, const SharedHeap& heap); + /** Copy the heap **/ + void copy_tree(HeapNode * node,std::vector*> * elm_vect,int heap_id,int nb_crit); + void deepcopy_tree(HeapNode * node,std::vector*> * elm_vect,int heap_id,int nb_crit); + }; /** @@ -257,9 +293,18 @@ class HeapElt { /** Create an HeapElt with a data and two criteria */ explicit HeapElt(T* data, double crit_1, double crit_2); + explicit HeapElt(double* crit, HeapNode** holder, T* data); + + /** Copy constructor **/ + HeapElt* copy(int nb_crit); + HeapElt* deepcopy(int nb_crit); + /** Delete the element */ ~HeapElt() ; + /** Delete the data */ + void flush() ; + /** * Compare the criterion of a given heap with the value d. * Return true if the criterion is greater. @@ -296,12 +341,89 @@ SharedHeap::SharedHeap(CostFunc& cost, bool update_cost, int id) : nb_node } +template +std::pair *,std::vector*> *> SharedHeap::deepcopy_sheap(int nb_crit){ + SharedHeap* new_heap = new SharedHeap(costf,update_cost_when_sorting,heap_id); + new_heap->nb_nodes = nb_nodes; + std::vector*> * elm_vect = new std::vector*>() ; + if(root != NULL) { + new_heap->root = new HeapNode(root->elt->deepcopy(nb_crit)); + new_heap->root->elt->holder[new_heap->heap_id] = new_heap->root; + + elm_vect->push_back(new_heap->root->elt); + new_heap->root->deepcopy_tree(root,elm_vect,heap_id,nb_crit); + } + std::pair *,std::vector*> *> rtrn; + rtrn.first = new_heap; + rtrn.second = elm_vect; + return rtrn; +} + + + +template +std::pair *,std::vector*> *> SharedHeap::copy_sheap(int nb_crit){ + SharedHeap* new_heap = new SharedHeap(costf,update_cost_when_sorting,heap_id); + new_heap->nb_nodes = nb_nodes; + std::vector*> * elm_vect = new std::vector*>() ; + if(root != NULL) { + new_heap->root = new HeapNode(root->elt->copy(nb_crit)); + new_heap->root->elt->holder[new_heap->heap_id] = new_heap->root; + + elm_vect->push_back(new_heap->root->elt); + new_heap->root->copy_tree(root,elm_vect,heap_id,nb_crit); + } + std::pair *,std::vector*> *> rtrn; + rtrn.first = new_heap; + rtrn.second = elm_vect; + return rtrn; +} + template SharedHeap::~SharedHeap() { - if (root) delete root; // warning: delete all sub-nodes - root = NULL; + clear(); +} + + +template +void SharedHeap::flush() { + if (nb_nodes>0) { + flush_subnodes(root); + nb_nodes=0; + root= NULL; + } } +template +void SharedHeap::flush_subnodes(HeapNode* node) { + if (node->left) flush_subnodes(node->left); + if (node->right) flush_subnodes(node->right); + node->elt->flush(); + delete node->elt; + delete node; +} + + +template +void SharedHeap::clear() { + if (nb_nodes>0) { + clear_subnodes(root); + nb_nodes=0; + root= NULL; + } +} + +template +void SharedHeap::clear_subnodes(HeapNode* node) { + if (node->left) clear_subnodes(node->left); + if (node->right) clear_subnodes(node->right); + //node->elt->flush(); + delete node->elt; + delete node; +} + + + template inline double SharedHeap::minimum() const { return root->elt->crit[heap_id]; @@ -336,6 +458,7 @@ void SharedHeap::sort() { nb_nodes = heap_tmp->size(); heap_tmp->root = NULL; + heap_tmp->nb_nodes=0; delete heap_tmp; // warning: delete all sub-nodes } @@ -574,6 +697,44 @@ HeapNode::HeapNode(HeapElt* elt, HeapNode* father): elt(elt), right(NUL // if (left) delete left; //} +template +void HeapNode::deepcopy_tree(HeapNode * node,std::vector*> * elm_vect,int heap_id,int nb_crit){ + HeapElt* elem; + if(node->left != NULL){ + elem = node->left->elt->deepcopy(nb_crit);//new HeapElt(*(root->elt),nb_crit); + elm_vect->push_back(elem); + left = new HeapNode(elem,this); + elem->holder[heap_id] = left; + left->deepcopy_tree(node->left,elm_vect,heap_id,nb_crit); + } + if(node->right != NULL){ + elem = node->right->elt->deepcopy(nb_crit);//new HeapElt(*(root->elt),nb_crit); + right = new HeapNode(elem,this); + elm_vect->push_back(elem); + elem->holder[heap_id] = right; + right->deepcopy_tree(node->right,elm_vect,heap_id,nb_crit); + } +} + +template +void HeapNode::copy_tree(HeapNode * node,std::vector*> * elm_vect,int heap_id,int nb_crit){ + HeapElt* elem; + if(node->left != NULL){ + elem = node->left->elt->copy(nb_crit); + elm_vect->push_back(elem); + left = new HeapNode(elem,this); + elem->holder[heap_id] = left; + left->copy_tree(node->left,elm_vect,heap_id,nb_crit); + } + if(node->right != NULL){ + elem = node->right->elt->copy(nb_crit); + right = new HeapNode(elem,this); + elm_vect->push_back(elem); + elem->holder[heap_id] = right; + right->copy_tree(node->right,elm_vect,heap_id,nb_crit); + } +} + template bool HeapNode::is_sup(HeapNode* node, int heap_id) const { return elt->is_sup(node->elt->crit[heap_id],heap_id); @@ -605,6 +766,11 @@ void HeapNode::switch_elt(HeapNode* node, int heap_id) { // } //} +template +HeapElt::HeapElt(double* crit, HeapNode** holder, T* data): data(data), crit(crit), holder(holder) { + +} + template HeapElt::HeapElt(T* data, double crit_1) : data(data), /*nb_heaps(1),*/ crit(new double[1]), holder(new HeapNode*[1]){ crit[0] = crit_1; @@ -619,13 +785,43 @@ HeapElt::HeapElt(T* data, double crit_1, double crit_2) : data(data), /*nb_he holder[1] = NULL; } +template +HeapElt* HeapElt::deepcopy(int nb_crit) { + double * copycrit = new double[nb_crit]; + HeapNode** copyholder = new HeapNode*[nb_crit]; + for(unsigned i=0;icrit)[i]; + copyholder[i] = NULL; + } + T* copydata = new T(*(this->data)); + return new HeapElt(copycrit,copyholder,copydata); +} + +template +HeapElt* HeapElt::copy(int nb_crit) { + double * copycrit = new double[nb_crit]; + HeapNode** copyholder = new HeapNode*[nb_crit]; + for(unsigned i=0;icrit)[i]; + copyholder[i] = NULL; + } + T* copydata = this->data; + return new HeapElt(copycrit,copyholder,copydata); +} + + template HeapElt::~HeapElt() { - if (data) delete data; + //if (data) delete data; delete [] crit; delete [] holder; } +template +void HeapElt::flush() { + if (data) delete data; +} + template bool HeapElt::is_sup(double d, int ind_crit) const { return (crit[ind_crit] > d); @@ -647,15 +843,21 @@ std::ostream& operator<<(std::ostream& os, const HeapNode& node) { template std::ostream& operator<<(std::ostream& os, const SharedHeap& heap) { - if (!heap.root) return os << "(empty heap)"; + return heap.print(os); +} + + +template +std::ostream& SharedHeap::print(std::ostream& os) const{ + if (!root) return os << "(empty heap)"; os << std::endl; std::stack*,int> > s; - s.push(std::pair*,int>(heap.root,0)); + s.push(std::pair*,int>(root,0)); while (!s.empty()) { std::pair*,int> p=s.top(); s.pop(); for (int i=0; ielt->crit[heap.heap_id]) << std::endl; + os << (p.first->elt->crit[heap_id]) << std::endl; if (p.first->right) s.push(std::pair*,int>(p.first->right,p.second+1)); if (p.first->left) s.push(std::pair*,int>(p.first->left,p.second+1)); } diff --git a/tests/Ponts30.h b/tests/Ponts30.h index 21f4e9bda..47615ec5a 100644 --- a/tests/Ponts30.h +++ b/tests/Ponts30.h @@ -19,6 +19,7 @@ namespace ibex { class Ponts30 { public: Ponts30(); + ~Ponts30() { delete f;} Function* f; IntervalVector init_box; diff --git a/tests/TestArith.cpp b/tests/TestArith.cpp index 832357eb0..60ec192ea 100644 --- a/tests/TestArith.cpp +++ b/tests/TestArith.cpp @@ -181,6 +181,16 @@ void TestArith::mulMM01() { CPPUNIT_ASSERT(M==M_expected); } +void TestArith::detMM01() { + double _tab[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + Matrix M1(3,3,_tab); + double _res = 0; + CPPUNIT_ASSERT(M1.det()==_res); + + Matrix M2(1,2); + CPPUNIT_ASSERT(M1.det()==0); //TODO: check error message is returned +} + void TestArith::div01() { check_div(Interval::EMPTY_SET, Interval(0,1), Interval::EMPTY_SET); } void TestArith::div02() { check_div(Interval::ZERO, Interval::ZERO, Interval::EMPTY_SET); } void TestArith::div03() { check_div(Interval(1,2), Interval::ZERO, Interval::EMPTY_SET); } diff --git a/tests/TestArith.h b/tests/TestArith.h index e6e59818e..6adce575b 100644 --- a/tests/TestArith.h +++ b/tests/TestArith.h @@ -433,6 +433,7 @@ class TestArith : public CppUnit::TestFixture { void mul19(); void mulMM01(); + void detMM01(); /* test: * ======= diff --git a/tests/TestCell.cpp b/tests/TestCell.cpp new file mode 100644 index 000000000..6be3cad8b --- /dev/null +++ b/tests/TestCell.cpp @@ -0,0 +1,120 @@ +//============================================================================ +// I B E X +// File : TestCell.cpp +// Author : Jordan Ninin +// Copyright : ENSTA Bretagne (France) +// License : See the LICENSE file +// Created : Oct 10, 2016 +// Last Update : Oct 10, 2016 +//============================================================================ + +#include "TestCell.h" +#include "ibex_Cell.h" +#include "ibex_EntailedCtr.h" +#include "ibex_Bsc.h" +#include "ibex_LargestFirst.h" +#include "ibex_SystemFactory.h" +#include "ibex_NormalizedSystem.h" +#include "ibex_Cell.h" + +//using namespace std; + +namespace ibex { + + +void TestCell::test01() { + IntervalVector box (2, Interval(-1,1)); + Cell * c = new Cell(box); + + Cell * copy =new Cell(*c); + check(copy->box,c->box); + + LargestFirst bsc; + std::pair new_cells = bsc.bisect_cell(*c); + + + check(new_cells.first->box|new_cells.second->box,c->box); + delete c; + + check(new_cells.first->box|new_cells.second->box,copy->box); + +} + +void TestCell::test02() { + IntervalVector box (2, Interval(-1,1)); + Cell * root = new Cell(box); + + Variable x,y; + Function func(x,y,pow(cos(y)+cos(2*y+x),2)); + Function c1f(x,y,y-x*(x+6.28) ); + NumConstraint c1(c1f,LEQ); + Function c2f(x,y,y-x*(x-6.28)); + NumConstraint c2(c2f,LEQ); + + SystemFactory fac; + fac.add_var(x); + fac.add_var(y); + fac.add_ctr(c1); + fac.add_ctr(c2); + NormalizedSystem sys(fac); + + root->add(); + EntailedCtr * entailed = &root->get(); + entailed->init_root(sys,sys); + entailed->original(1) = true; + entailed->normalized(0) = true; + + Cell * copy =new Cell(*root); + check(copy->box,root->box); + CPPUNIT_ASSERT(entailed->original(1)); + CPPUNIT_ASSERT(entailed->normalized(0)); + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + + LargestFirst bsc; + std::pair new_cells = bsc.bisect_cell(*root); + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + + + check(new_cells.first->box|new_cells.second->box,root->box); + delete root; + + + CPPUNIT_ASSERT(new_cells.first->get().original(1)); + CPPUNIT_ASSERT(new_cells.first->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.first->get().original(0)); + CPPUNIT_ASSERT(!new_cells.first->get().normalized(1)); + + CPPUNIT_ASSERT(new_cells.second->get().original(1)); + CPPUNIT_ASSERT(new_cells.second->get().normalized(0)); + CPPUNIT_ASSERT(!new_cells.second->get().original(0)); + CPPUNIT_ASSERT(!new_cells.second->get().normalized(1)); + + CPPUNIT_ASSERT(copy->get().original(1)); + CPPUNIT_ASSERT(copy->get().normalized(0)); + CPPUNIT_ASSERT(!copy->get().original(0)); + CPPUNIT_ASSERT(!copy->get().normalized(1)); + + delete copy; + delete new_cells.first; + delete new_cells.second; + +} + + +} // end namespace + + + diff --git a/tests/TestCell.h b/tests/TestCell.h new file mode 100644 index 000000000..28d750979 --- /dev/null +++ b/tests/TestCell.h @@ -0,0 +1,40 @@ +/* ============================================================================ + * I B E X - CtcCell Tests + * ============================================================================ + * Copyright : ENSTA Bretagne (FRANCE) + * License : This program can be distributed under the terms of the GNU LGPL. + * See the file COPYING.LESSER. + * + * Author(s) : Jordan Ninin + * Created : Oct 10, 2016 + * ---------------------------------------------------------------------------- */ + +#ifndef __TEST_CELL_H__ +#define __TEST_CELL_H__ + +#include +#include + +#include "utils.h" + +namespace ibex { + +class TestCell : public CppUnit::TestFixture { + +public: + + CPPUNIT_TEST_SUITE(TestCell); + CPPUNIT_TEST(test01); + CPPUNIT_TEST(test02); + CPPUNIT_TEST_SUITE_END(); + + void test01(); + void test02(); + +}; + +CPPUNIT_TEST_SUITE_REGISTRATION(TestCell); + +} // namespace ibex + +#endif // __TEST_CELL_H__ diff --git a/tests/TestDoubleHeap.cpp b/tests/TestDoubleHeap.cpp new file mode 100644 index 000000000..cb55b65db --- /dev/null +++ b/tests/TestDoubleHeap.cpp @@ -0,0 +1,219 @@ +/* ============================================================================ + * I B E X - CellHeap Tests + * ============================================================================ + * Copyright : Ecole des Mines de Nantes (FRANCE) + * License : This program can be distributed under the terms of the GNU LGPL. + * See the file COPYING.LESSER. + * + * Author(s) : Jordan Ninin + * Created : Mar 2, 2014 + * ---------------------------------------------------------------------------- */ + +#include "TestDoubleHeap.h" +#include "ibex_DoubleHeap.h" + +using namespace std; + +namespace ibex { + + + +class TestCostFunc1 : public CostFunc { +public: + TestCostFunc1() { } ; + double cost(const Interval& c) const { return c.diam();}; +}; + +class TestCostFunc2 : public CostFunc { +public: + TestCostFunc2() { } ; + double cost(const Interval& c) const { return c.lb();}; +}; + +class TestCostFunc3 : public CostFunc { +public: + TestCostFunc3(): loup(10) { } ; + double cost(const Interval& c) const { return c.ub()*loup;}; + void set_loup(double lb) { loup = lb; } + +private: + double loup; +}; + + +void TestDoubleHeap::test01() { + + int nb= 10; + TestCostFunc1 costf1; + TestCostFunc2 costf2; + + DoubleHeap h(costf1,false,costf2,false,50); + + for (int i=1; i<=nb ;i++) { + if ((i%2)==1) h.push(new Interval(i,2*i)); + else h.push(new Interval(i,20*i)); + } + + CPPUNIT_ASSERT(h.minimum1()==1); + CPPUNIT_ASSERT(h.minimum2()==1); + + delete h.pop1(); + delete h.pop1(); + CPPUNIT_ASSERT(h.minimum1()==5); + CPPUNIT_ASSERT(h.minimum2()==2); + + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==5); + CPPUNIT_ASSERT(h.minimum2()==4); + + h.contract(20); + CPPUNIT_ASSERT(h.minimum1()==5); + CPPUNIT_ASSERT(h.minimum2()==5); + + delete h.pop2(); + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==9); + CPPUNIT_ASSERT(h.minimum2()==9); + + CPPUNIT_ASSERT(h.size()==1); + + h.push(new Interval(10,11)); + h.push(new Interval(12,14)); + h.push(new Interval(12,12)); + h.push(new Interval(12,13.1)); + + CPPUNIT_ASSERT(h.minimum1()==0); + CPPUNIT_ASSERT(h.minimum2()==9); + CPPUNIT_ASSERT(h.size()==5); + + delete h.pop1(); + CPPUNIT_ASSERT(h.minimum1()==1); + CPPUNIT_ASSERT(h.minimum2()==9); + CPPUNIT_ASSERT(h.size()==4); + + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==1); + CPPUNIT_ASSERT(h.minimum2()==10); + CPPUNIT_ASSERT(h.size()==3); + + h.flush(); + CPPUNIT_ASSERT(h.size()==0); +} + + +void TestDoubleHeap::test02() { + + int nb= 10; + TestCostFunc2 costf2; + TestCostFunc3 costf3; + + DoubleHeap h(costf2,false,costf3,true,50); + + for (int i=1; i<=nb ;i++) { + if ((i%2)==1) h.push(new Interval(i,2*i)); + else h.push(new Interval(i,20*i)); + } + + CPPUNIT_ASSERT(h.minimum1()==1); + CPPUNIT_ASSERT(h.minimum2()==20); + + delete h.pop2(); + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==2); + CPPUNIT_ASSERT(h.minimum2()==100); + + delete h.pop1(); + CPPUNIT_ASSERT(h.minimum1()==4); + CPPUNIT_ASSERT(h.minimum2()==100); + + // it is necessary to update the loup. + costf3.set_loup(100); + h.contract(7.5); + CPPUNIT_ASSERT(h.minimum1()==4); + CPPUNIT_ASSERT(h.minimum2()==1000); + CPPUNIT_ASSERT(h.size()==4); + + + delete h.pop2(); + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==4); + CPPUNIT_ASSERT(h.minimum2()==8000); + CPPUNIT_ASSERT(h.size()==2); + + h.push(new Interval(10,11)); + h.push(new Interval(12,14)); + h.push(new Interval(12,12)); + h.push(new Interval(12,13.1)); + + CPPUNIT_ASSERT(h.minimum1()==4); + CPPUNIT_ASSERT(h.minimum2()==1100); + CPPUNIT_ASSERT(h.size()==6); + + delete h.pop1(); + CPPUNIT_ASSERT(h.minimum1()==6); + CPPUNIT_ASSERT(h.minimum2()==1100); + CPPUNIT_ASSERT(h.size()==5); + + delete h.pop2(); + CPPUNIT_ASSERT(h.minimum1()==6); + CPPUNIT_ASSERT(h.minimum2()==1200); + CPPUNIT_ASSERT(h.size()==4); + + h.flush(); + CPPUNIT_ASSERT(h.size()==0); +} + +void TestDoubleHeap::test03() { + int nb= 10; + TestCostFunc2 costf2; + TestCostFunc3 costf3; + + DoubleHeap h(costf2,false,costf3,true,50); + + for (int i=1; i<=nb ;i++) { + if ((i%2)==1) h.push(new Interval(i,2*i)); + else h.push(new Interval(i,20*i)); + } + + DoubleHeap newh(h); + + while (h.size() > 0) { + CPPUNIT_ASSERT(h.top1() == newh.top1()); + CPPUNIT_ASSERT(h.top2() == newh.top2()); + h.pop1(); + delete newh.pop1(); + + } + + +} + +void TestDoubleHeap::test04() { + + int nb= 10; + TestCostFunc2 costf2; + TestCostFunc3 costf3; + + DoubleHeap h(costf2,false,costf3,true,50); + + for (int i=1; i<=nb ;i++) { + if ((i%2)==1) h.push(new Interval(i,2*i)); + else h.push(new Interval(i,20*i)); + } + + DoubleHeap *newh = h.deepcopy(); + + while (h.size() > 0) { + CPPUNIT_ASSERT(*h.top1() == *newh->top1()); + CPPUNIT_ASSERT(*h.top2() == *newh->top2()); + delete h.pop1(); + delete newh->pop1(); + + } + delete newh; + +} + + + +} // end namespace diff --git a/plugins/optim/tests/TestDoubleHeap.h b/tests/TestDoubleHeap.h similarity index 89% rename from plugins/optim/tests/TestDoubleHeap.h rename to tests/TestDoubleHeap.h index cde984be0..bebf5d543 100644 --- a/plugins/optim/tests/TestDoubleHeap.h +++ b/tests/TestDoubleHeap.h @@ -14,32 +14,33 @@ #include #include -#include "utils.h" #include "ibex_DoubleHeap.h" +#include "utils.h" -namespace ibex { +namespace ibex { class TestDoubleHeap : public CppUnit::TestFixture { + public: CPPUNIT_TEST_SUITE(TestDoubleHeap); - - - CPPUNIT_TEST(test01); - CPPUNIT_TEST(test02); + CPPUNIT_TEST(test01); + CPPUNIT_TEST(test02); + CPPUNIT_TEST(test03); + CPPUNIT_TEST(test04); CPPUNIT_TEST_SUITE_END(); - void test01(); void test02(); + void test03(); + void test04(); }; CPPUNIT_TEST_SUITE_REGISTRATION(TestDoubleHeap); +} // namespace ibex - -} // namespace ibex #endif // __TEST_CTC_DOUBLEHEAP_H__ diff --git a/tests/TestExpr.cpp b/tests/TestExpr.cpp index 0adca2dc2..fe3d740ac 100644 --- a/tests/TestExpr.cpp +++ b/tests/TestExpr.cpp @@ -87,6 +87,9 @@ void TestExpr::addxy02() { CPPUNIT_ASSERT(e.dim==Dim::row_vec(3)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::ROW_VECTOR); + delete &x; + delete &y; + delete &e; } void TestExpr::addxy03() { @@ -97,6 +100,9 @@ void TestExpr::addxy03() { CPPUNIT_ASSERT(e.dim==Dim::col_vec(3)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::COL_VECTOR); + delete &x; + delete &y; + delete &e; } void TestExpr::addxy04() { @@ -106,6 +112,9 @@ void TestExpr::addxy04() { CPPUNIT_ASSERT(e.dim==Dim::matrix(2,3)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::MATRIX); + delete &x; + delete &y; + delete &e; } void TestExpr::addxx01() { @@ -133,6 +142,9 @@ void TestExpr::mulxy01() { CPPUNIT_ASSERT(e.dim==Dim::scalar()); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::SCALAR); + delete &x; + delete &y; + delete &e; } @@ -144,6 +156,9 @@ void TestExpr::mulxy02() { CPPUNIT_ASSERT(e.dim==Dim::scalar()); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::SCALAR); + delete &x; + delete &y; + delete &e; } // matrix * matrix @@ -154,6 +169,9 @@ void TestExpr::mulxy03() { CPPUNIT_ASSERT(e.dim==Dim::matrix(2,4)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::MATRIX); + delete &x; + delete &y; + delete &e; } // matrix * vector @@ -164,6 +182,9 @@ void TestExpr::mulxy04() { CPPUNIT_ASSERT(e.dim==Dim::col_vec(2)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::COL_VECTOR); + delete &x; + delete &y; + delete &e; } // scalar * matrix @@ -174,6 +195,9 @@ void TestExpr::mulxy05() { CPPUNIT_ASSERT(e.dim==Dim::matrix(2,3)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::MATRIX); + delete &x; + delete &y; + delete &e; } // scalar * vector @@ -186,6 +210,9 @@ void TestExpr::mulxy06() { CPPUNIT_ASSERT(e.type()==Dim::ROW_VECTOR); CPPUNIT_ASSERT(sameExpr(e,"(x*y)")); + delete &x; + delete &y; + delete &e; } // vector * vector (outer product) @@ -196,6 +223,9 @@ void TestExpr::mulxy07() { CPPUNIT_ASSERT(e.dim==Dim::matrix(2,3)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::MATRIX); + delete &x; + delete &y; + delete &e; } // vector * matrix @@ -206,6 +236,9 @@ void TestExpr::mulxy08() { CPPUNIT_ASSERT(e.dim==Dim::row_vec(4)); CPPUNIT_ASSERT(!e.is_zero()); CPPUNIT_ASSERT(e.type()==Dim::ROW_VECTOR); + delete &x; + delete &y; + delete &e; } @@ -279,6 +312,7 @@ void TestExpr::unaryOp() { CPPUNIT_ASSERT(sameExpr(acosh(x),"acosh(x)")); CPPUNIT_ASSERT(sameExpr(asinh(x),"asinh(x)")); CPPUNIT_ASSERT(sameExpr(atanh(x),"atanh(x)")); + delete &x; } void TestExpr::binaryOp() { @@ -291,6 +325,8 @@ void TestExpr::binaryOp() { CPPUNIT_ASSERT(sameExpr(max(x,y),"max(x,y)")); CPPUNIT_ASSERT(sameExpr(min(x,y),"min(x,y)")); CPPUNIT_ASSERT(sameExpr(atan2(x,y),"atan2(x,y)")); + delete &x; + delete &y; } void TestExpr::cst01() { @@ -306,6 +342,8 @@ void TestExpr::cst01() { const ExprConstant& z=ExprConstant::new_scalar(0); CPPUNIT_ASSERT(z.is_zero()); + delete &c; + delete &z; } void TestExpr::cst02() { @@ -316,6 +354,7 @@ void TestExpr::cst02() { CPPUNIT_ASSERT(!c1.is_zero()); CPPUNIT_ASSERT(c1.type()==Dim::COL_VECTOR); CPPUNIT_ASSERT(c1.get_vector_value()==v); + delete &c1; } void TestExpr::cst03() { @@ -326,6 +365,7 @@ void TestExpr::cst03() { CPPUNIT_ASSERT(!c2.is_zero()); CPPUNIT_ASSERT(c2.type()==Dim::MATRIX); CPPUNIT_ASSERT(c2.get_matrix_value()==M); + delete &c2; } @@ -337,6 +377,7 @@ void TestExpr::cst04() { CPPUNIT_ASSERT(!c1.is_zero()); CPPUNIT_ASSERT(c1.type()==Dim::ROW_VECTOR); CPPUNIT_ASSERT(c1.get_vector_value()==v); + delete &c1; } void TestExpr::cst05() { @@ -352,6 +393,8 @@ void TestExpr::cst05() { const ExprConstant& z2=ExprConstant::new_matrix(M); CPPUNIT_ASSERT(z1.is_zero()); CPPUNIT_ASSERT(z2.is_zero()); + delete &z1; + delete &z2; } void TestExpr::vector01() { @@ -389,6 +432,11 @@ void TestExpr::vector02() { CPPUNIT_ASSERT(v.row_vector()); CPPUNIT_ASSERT(v.type()==Dim::MATRIX); CPPUNIT_ASSERT(sameExpr(v,"(x,y,(x+y),(x+(x+y)))")); + delete &x; + delete &y; + delete &e1; + delete &e2; + delete &v; } void TestExpr::index01() { @@ -443,6 +491,8 @@ void TestExpr::index03() { CPPUNIT_ASSERT(p.first==&x); bool mask[3][4]={{false,false,false,false},{false,false,false,false},{false,true,false,false}}; CPPUNIT_ASSERT(same_mask(3,4,(bool*) mask,p.second)); + delete &x; + delete &e; } void TestExpr::index04() { @@ -456,6 +506,9 @@ void TestExpr::index04() { {false,false, true, true}, {false,false, true, true}}; CPPUNIT_ASSERT(same_mask(3,4,(bool*) mask,p.second)); + delete &x; + delete &tmp; + delete &e; } @@ -535,6 +588,7 @@ void TestExpr::subnodes01() { CPPUNIT_ASSERT(&nodes[3]==&e1); CPPUNIT_ASSERT((&nodes[5]==&x && &nodes[6]==&y) || (&nodes[5]==&y && &nodes[6]==&x)); + } void TestExpr::subnodes02() { @@ -559,6 +613,9 @@ void TestExpr::subnodes02() { CPPUNIT_ASSERT(&nodes2[0]==&e2); CPPUNIT_ASSERT(&nodes2[1]==&x1); CPPUNIT_ASSERT(&nodes2[2]==&x2); + delete &x1; + delete &x2; + delete &y; } @@ -582,6 +639,12 @@ void TestExpr::subnodes03() { CPPUNIT_ASSERT(&nodes[6]==&x3); CPPUNIT_ASSERT(&nodes[7]==&x4); CPPUNIT_ASSERT(&nodes[8]==&x5); + delete &x1; + delete &x2; + delete &x3; + delete &x4; + delete &x5; + delete &y; } void TestExpr::subnodes04() { @@ -605,6 +668,11 @@ void TestExpr::subnodes04() { CPPUNIT_ASSERT(&nodes[4]==&two); CPPUNIT_ASSERT(&nodes[5]==&x1); CPPUNIT_ASSERT(&nodes[6]==&x2); + delete &x1; + delete &x2; + delete &one; + delete &two; + delete &y; } @@ -613,6 +681,8 @@ void TestExpr::bug81() { IntervalVector y(3); const ExprNode& z=x+y; CPPUNIT_ASSERT(z.dim == x.dim); + delete &x; + delete &z; } } // end namespace diff --git a/wscript b/wscript index 38939fc3b..66e836d95 100644 --- a/wscript +++ b/wscript @@ -129,7 +129,7 @@ def configure (conf): # Disable rounding interval if (conf.options.WITH_STANDALONE): conf.env.WITHOUT_ROUNDING =True - + ################################################################################################## # Bison / Flex env.append_unique ("BISONFLAGS", ["--name-prefix=ibex", "--report=all", "--file-prefix=parser"])