|
realroot_doc 0.1.1
|
|
realroot_doc 0.1.1
|
#include <solver_mv_monomial.hpp>
Definition at line 66 of file solver_mv_monomial.hpp.
| solver_mv_monomial | ( | double | e = 1e-7 | ) | [inline] |
Definition at line 295 of file solver_mv_monomial.hpp.
References BOX, mmx::realroot::jacobian(), POL, Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().
Referenced by solver< C, MCFapproximate >::solve(), and solver< C, M >::solve().
{
homography_mv<C> h(d);
BOX sys= BOX(p,h);
/*Global data*/
free(J);
for (unsigned i=0;i<S0.size();++i)
free(S0[i]);
S0=Seq<POL *>();
for (unsigned i=0;i<sys->nbpol();++i)
S0 << new POL( sys->system(i) );
J= jacobian(S0);
/*end Global data*/
Seq<FT> t;
Seq< BOX> s(solve_system(sys) );
Seq< std::vector<FT> > a;
// ... inf ?
// free(J);
// for (unsigned i=0;i<S0.size();++i)
// free(S0[i]);
return a;
}
Definition at line 323 of file solver_mv_monomial.hpp.
References BOX, C_INCLUDE, mmx::realroot::jacobian(), POL, Seq< C, R >::rep(), Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().
{
BOX * sys= new BOX(p, dom.size() );
/*Global data*/
// free(J); //crash
for (unsigned i=0;i<S0.size();++i)
free(S0[i]);
S0=Seq<POL *>();
for (unsigned i=0;i<sys->nbpol();++i)
S0 << new POL( sys->system(i) );
J= jacobian(S0);
/*end Global data*/
sys->restrict(dom);
unsigned v; //, dim(dom.size() );
Seq<FT> d;
BOX * l, * b;
Seq<C> t;
//FT ev;
Seq<BOX> s(solve_system(*sys) );
Seq< std::vector<FT> > a;
for (unsigned i=0;i<s.size(); ++i )
{
b= new BOX(s[i]);
while ( b->template width<FT>(v) > eps )
{
t= b->middle();
if ( b->is_root(t) )
{
d= b->template point<FT>(t);
a << d.rep();
}
l= new BOX( *b ) ;
l->shift_box( t[v] , v );
if ( C_INCLUDE )
{
free(b);
b=l;
continue;
}
else
{
free(l);
b->contract_box(t[i],v);
b->reverse_and_shift_box(v);
b->reverse_box(v);
}
}
d= b->template domain<FT>();
a << d.rep();
free(b);
}
return a;
}
| bool check_root | ( | const Seq< double > | t, |
| const double | eps | ||
| ) | [inline] |
Definition at line 390 of file solver_mv_monomial.hpp.
References mmx::abs(), mmx::assign(), DPOL, mmx::eval(), and Seq< C, R >::size().
Referenced by solver_mv_monomial< FT, POL >::solve_system().
{
DPOL p(0);
double ev;
for (unsigned j=0; j!=S0.size(); ++j)
{
let::assign(p, *S0[j]);
tensor::eval( ev , p.rep() ,
t , t.size() );
//std::cout<<"check: "<< ev<<std::endl;
if ( abs(ev) > 1e-10 )
return false;
}
std::cout<<"Root on split: "<< t <<std::endl;
return true;
};
Definition at line 266 of file solver_mv_monomial.hpp.
References BOX, mmx::realroot::jacobian(), POL, Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().
{
BOX * sys= new BOX(p, dom.size() );
/*Global data*/
free(J);
for (unsigned i=0;i<S0.size();++i)
free(S0[i]);
S0=Seq<POL *>();
for (unsigned i=0;i<sys->nbpol();++i)
S0 << new POL( sys->system(i) );
J= jacobian(S0);
/*end Global data*/
sys->restrict(dom);
Seq< BOX> s(solve_system( *sys) );
Seq< std::vector<Interval<FT> > > a;
for (unsigned i=0; i<s.size(); ++i)
a << s[i].template domain<FT>().rep();
// free(J);
// for (unsigned i=0;i<S0.size();++i)
// free(S0[i]);
return a;
}
Definition at line 238 of file solver_mv_monomial.hpp.
References BOX, mmx::realroot::jacobian(), POL, Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().
Referenced by solver< C, MCFisolate >::solve().
{
homography_mv<C> h(d);
BOX sys= BOX(p,h);
/*Global data*/
free(J);
for (unsigned i=0;i<S0.size();++i)
free(S0[i]);
S0=Seq<POL *>();
for (unsigned i=0;i<sys.nbpol();++i)
S0 << new POL( sys.system(i) );
J= jacobian(S0);
/*end Global data*/
Seq< BOX> s(solve_system(sys) );
Seq< std::vector<Interval<FT> > > a;
for (unsigned i=0; i<s.size(); ++i)
a << s[i].template domain<FT>().rep();
// free(J);
// for (unsigned i=0;i<S0.size();++i)
// free(S0[i]);
return a;
}
Solve routine.
Definition at line 85 of file solver_mv_monomial.hpp.
References ALLBOXES, BOX, C_EXCLUDE, C_INCLUDE, solver_mv_monomial< FT, POL >::check_root(), Seq< C, R >::erase(), mmx::realroot::exclude2(), N_ITER, Seq< C, R >::size(), mmx::sparse::subs(), and VERBOSE.
Referenced by solver_mv_monomial< FT, POL >::approximate(), and solver_mv_monomial< FT, POL >::isolate().
{
unsigned c=0,cand=0, i, it=0, subs=0, ver=0;
BOX * b = NULL;
BOX * r = NULL;
Seq<BOX> sol;
Seq<BOX> psol;
bool red, out;
STACK boxes;
FT v(0), bx;
bx= sys.template volume <FT>();
if (VERBOSE) sys.print();
boxes.push( new BOX(sys) );
while ( !boxes.empty() )
{
it++;
b = boxes.top() ;
boxes.pop();
/*reduce the domain */
out= false;
red= true;
while ( !out && red )
{
if ( C_EXCLUDE )
{
out=true;
if (VERBOSE) {
//std::cout<<"- Excluded:"<<std::endl;
//b->print();
}
if (ALLBOXES) //FOR AXEL
{
v+= b->template volume<FT>();
r = new BOX( *b ) ;
sol << (*r);
free(r);
}
c++;
}
else
{
red= b->reduce_domain();
if (VERBOSE && red) {
//std::cout<<"- Reduced:"<<std::endl;
//b->print();
}
//red= false;
}
}
if ( out )
{
free(b);
continue;
}
/*check for root */
if ( C_INCLUDE )
{
if (VERBOSE) {
std::cout<<"- Solution found:"<<std::endl;
b->print();
}
ver++ ;
sol << (*b);
free(b);
continue;
}
/*Subdivide*/
if ( it > N_ITER )
{
cand++;
std::cout<<"MAX iters"<<" ("<<N_ITER<<") "
<<"reached!" << std::endl;
b->print();
sol << (*b);
free(b);
}
else
{
if ( b->template width<double>(i) > eps )
{
subs++;
if (check_root( b->subdiv_point(i),eps) )
{
psol <<BOX( *b, i);
//free(b);
//continue;
//sol[sol.size()-1].print();
ver++;
}
b->subdivide (i,boxes, b->middle(i) );
//b->subdivide (i,boxes);
//b->subdivide (i,boxes, 1 );
//b->subdivide (boxes);
free(b);
}
else
{
if ( C_EXCLUDE ){
if (ALLBOXES) sol << (*b); //FOR AXEL
}else
{
cand++;
sol << (*b);
if (VERBOSE) {
std::cout<<"- EPS reached:"<<std::endl;
//b->print();
}
}
free(b);
}
}
}/*while*/
unsigned j=0;
//if (0)
if ( !ALLBOXES && S0.size()==2)
while (j<sol.size())
{
if ( exclude2( &sol[j],J) )
{
sol.erase(j);
cand--;
}
else
{ //std::cout <<"td="<<sol[j].td<<std::endl;
++j;
}
}
if (VERBOSE) {
std::cout<< "Iterations= "<< it <<std::endl;
std::cout<< "Excluded = "<< c <<std::endl;
std::cout<< "Verified = "<< ver <<std::endl;
std::cout<< "Subdivs = "<< subs <<std::endl;
if (ALLBOXES)
std::cout<< "Reduced volume= "<<
as<double>( 100*(bx-v)/bx )<< "\%" <<std::endl;
std::cout<< "#stack= "<< cand <<std::endl;
}
sol<< psol;
return sol;
}
