solver_mv_monomial< FT, POL > Class Template Reference

#include <solver_mv_monomial.hpp>

List of all members.

Public Member Functions

Detailed Description

template<class FT, class POL>
class mmx::realroot::solver_mv_monomial< FT, POL >

Definition at line 66 of file solver_mv_monomial.hpp.

Constructor & Destructor Documentation

solver_mv_monomial ( double  e = 1e-7  )  [inline]

Constructor.

Definition at line 78 of file solver_mv_monomial.hpp.

00079     {
00080       J=NULL;
00081       eps=e;
00082     };

Member Function Documentation

Seq< std::vector<FT> > approximate ( Seq< POL > &  p,
Seq< Interval< C > > &  dom 
) [inline]

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().

00324     {
00325       BOX * sys= new BOX(p, dom.size() );
00326       
00327        /*Global data*/
00328        // free(J); //crash
00329       for (unsigned i=0;i<S0.size();++i) 
00330         free(S0[i]);
00331       S0=Seq<POL *>();
00332       for (unsigned i=0;i<sys->nbpol();++i) 
00333         S0 << new POL( sys->system(i) );            
00334 
00335 
00336 
00337       
00338 
00339       J= jacobian(S0);
00340       /*end Global data*/
00341       
00342       sys->restrict(dom);
00343       unsigned v; //, dim(dom.size() );
00344       Seq<FT> d;
00345       BOX * l, * b;
00346       Seq<C> t;
00347       //FT ev;
00348       Seq<BOX> s(solve_system(*sys) );
00349       
00350       Seq< std::vector<FT> > a;
00351       
00352       for (unsigned i=0;i<s.size(); ++i )
00353       {
00354         b= new BOX(s[i]);
00355         
00356         while ( b->template width<FT>(v) > eps )
00357         {
00358           t= b->middle();
00359           if ( b->is_root(t) ) 
00360           {
00361             d= b->template point<FT>(t);
00362             a << d.rep();
00363           }
00364           
00365           l= new BOX( *b ) ;
00366           l->shift_box( t[v] , v );
00367           
00368           if ( C_INCLUDE ) 
00369           {   
00370             free(b);
00371             b=l; 
00372             continue;
00373           }
00374           else 
00375           {   
00376             free(l);
00377             b->contract_box(t[i],v);
00378             b->reverse_and_shift_box(v);
00379             b->reverse_box(v);
00380           }
00381         }
00382         d= b->template domain<FT>();
00383         a << d.rep();
00384         free(b);
00385       }
00386       return a;
00387     }

Seq< std::vector<FT> > approximate ( Seq< POL > &  p,
unsigned &  d 
) [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().

00296     {
00297       homography_mv<C> h(d);                            
00298       BOX sys= BOX(p,h);
00299       
00300       /*Global data*/
00301       free(J);
00302       for (unsigned i=0;i<S0.size();++i) 
00303         free(S0[i]);
00304       S0=Seq<POL *>();
00305       for (unsigned i=0;i<sys->nbpol();++i) 
00306         S0 << new POL( sys->system(i) );            
00307       J= jacobian(S0);
00308       /*end Global data*/
00309       
00310       Seq<FT> t;
00311       Seq< BOX> s(solve_system(sys) );
00312       
00313       Seq< std::vector<FT> > a;
00314       
00315       // ... inf ?
00316       
00317 //      free(J);
00318 //      for (unsigned i=0;i<S0.size();++i) 
00319 //        free(S0[i]);
00320       return a;
00321         }

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().

00391     {
00392       DPOL p(0);
00393       double ev;
00394       for (unsigned j=0; j!=S0.size(); ++j)
00395       {
00396         
00397         let::assign(p, *S0[j]);
00398         tensor::eval( ev , p.rep() , 
00399                       t , t.size() );
00400 
00401         //std::cout<<"check: "<< ev<<std::endl;
00402         if ( abs(ev) > 1e-10 )
00403           return false;
00404       } 
00405       std::cout<<"Root on split: "<< t  <<std::endl;
00406       return true;
00407 
00408     };

Seq< std::vector<Interval<FT> > > isolate ( Seq< POL > &  p,
Seq< Interval< C > > &  dom 
) [inline]

Definition at line 266 of file solver_mv_monomial.hpp.

References BOX, mmx::realroot::jacobian(), POL, Seq< C, R >::rep(), Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().

00267         {
00268           BOX * sys= new BOX(p, dom.size() );
00269           
00270           /*Global data*/
00271           free(J);
00272           for (unsigned i=0;i<S0.size();++i) 
00273             free(S0[i]);
00274           S0=Seq<POL *>();
00275           for (unsigned i=0;i<sys->nbpol();++i) 
00276             S0 << new POL( sys->system(i) );        
00277           J= jacobian(S0);
00278           /*end Global data*/
00279 
00280           sys->restrict(dom);
00281           Seq< BOX> s(solve_system( *sys) );
00282       
00283           Seq< std::vector<Interval<FT> > > a;
00284 
00285           for (unsigned i=0; i<s.size(); ++i)
00286             a << s[i].template domain<FT>().rep();
00287           
00288 //          free(J);
00289 //          for (unsigned i=0;i<S0.size();++i) 
00290 //            free(S0[i]);
00291           return a;
00292         }

Seq< std::vector<Interval<FT> > > isolate ( Seq< POL > &  p,
unsigned &  d 
) [inline]

Definition at line 238 of file solver_mv_monomial.hpp.

References BOX, mmx::realroot::jacobian(), POL, Seq< C, R >::rep(), Seq< C, R >::size(), and solver_mv_monomial< FT, POL >::solve_system().

Referenced by solver< C, MCFisolate >::solve().

00239     {
00240       homography_mv<C> h(d);                            
00241       BOX sys= BOX(p,h);
00242       
00243       /*Global data*/
00244       free(J);
00245       for (unsigned i=0;i<S0.size();++i) 
00246         free(S0[i]);
00247       S0=Seq<POL *>();
00248       for (unsigned i=0;i<sys.nbpol();++i) 
00249         S0 << new POL( sys.system(i) );     
00250       J= jacobian(S0);
00251       /*end Global data*/
00252       
00253       Seq< BOX> s(solve_system(sys) );
00254       
00255       Seq< std::vector<Interval<FT> > > a;
00256       
00257       for (unsigned i=0; i<s.size(); ++i)
00258         a << s[i].template domain<FT>().rep();
00259                
00260 //      free(J);
00261 //      for (unsigned i=0;i<S0.size();++i) 
00262 //        free(S0[i]);  
00263       return a;
00264     }

Seq< box_rep<POL> > solve_system ( box_rep< POL > &  sys  )  [inline]

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().

00086     {
00087       unsigned c=0,cand=0, i, it=0, subs=0, ver=0;
00088       BOX * b = NULL;
00089       BOX * r = NULL;
00090       Seq<BOX> sol;
00091       Seq<BOX> psol;
00092       bool red, out;
00093       STACK boxes;
00094 
00095       FT v(0), bx;
00096       bx= sys.template volume <FT>();
00097       if (VERBOSE) sys.print();
00098       
00099       boxes.push( new BOX(sys) );
00100       
00101       while ( !boxes.empty() )
00102       {
00103         it++;
00104         b = boxes.top() ;
00105         boxes.pop();
00106 
00107         /*reduce the domain */
00108         out= false;
00109         red= true;
00110         while ( !out && red )
00111         {
00112           if ( C_EXCLUDE )
00113           { 
00114             out=true;
00115             if (VERBOSE) { 
00116               //std::cout<<"- Excluded:"<<std::endl;
00117               //b->print();
00118             }                       
00119             if (ALLBOXES) //FOR AXEL
00120             {
00121               v+= b->template volume<FT>();
00122               r = new BOX( *b ) ;
00123               sol << (*r); 
00124               free(r);
00125             }
00126             c++;
00127           }
00128           else
00129           {
00130             red= b->reduce_domain(); 
00131             if (VERBOSE && red) { 
00132               //std::cout<<"- Reduced:"<<std::endl;
00133               //b->print();
00134             }            
00135             //red= false;
00136           }
00137         }
00138 
00139         if ( out ) 
00140         {
00141           free(b);
00142           continue;
00143         }               
00144         /*check for root */
00145         if ( C_INCLUDE )
00146         { 
00147           if (VERBOSE) { 
00148             std::cout<<"- Solution found:"<<std::endl;
00149             b->print();
00150           }
00151           
00152           ver++ ;
00153           sol << (*b);
00154           free(b);
00155           continue; 
00156         }
00157         
00158         /*Subdivide*/
00159         if ( it > N_ITER )
00160         {
00161           cand++;
00162           std::cout<<"MAX iters"<<" ("<<N_ITER<<") "
00163                    <<"reached!" << std::endl;
00164           b->print();
00165           sol << (*b);
00166           free(b);
00167         }
00168         else
00169         {
00170           if ( b->template width<double>(i) > eps )
00171           {     
00172             subs++;
00173           
00174             if (check_root( b->subdiv_point(i),eps) )
00175             {
00176               psol <<BOX( *b, i);
00177               //free(b);
00178               //continue; 
00179               //sol[sol.size()-1].print();
00180               ver++;
00181             }
00182 
00183             b->subdivide (i,boxes, b->middle(i) );
00184             //b->subdivide (i,boxes);
00185             //b->subdivide (i,boxes,   1 );
00186             //b->subdivide (boxes);
00187             free(b); 
00188           }
00189           else 
00190           {
00191             if ( C_EXCLUDE  ){
00192               if (ALLBOXES) sol << (*b); //FOR AXEL
00193             }else
00194             {   
00195               cand++;
00196               sol << (*b);
00197               if (VERBOSE) { 
00198                 std::cout<<"- EPS reached:"<<std::endl;
00199                 //b->print();
00200               }
00201             }
00202             free(b);
00203           }
00204         }               
00205       }/*while*/
00206       
00207       unsigned j=0;
00208       //if (0)
00209       if ( !ALLBOXES && S0.size()==2)
00210       while (j<sol.size())
00211       {
00212         if ( exclude2( &sol[j],J) )
00213         {
00214           sol.erase(j);
00215           cand--;
00216         }
00217         else 
00218         { //std::cout  <<"td="<<sol[j].td<<std::endl;
00219           ++j; 
00220         }
00221       }
00222 
00223       if (VERBOSE) {
00224       std::cout<< "Iterations= "<< it           <<std::endl;    
00225       std::cout<< "Excluded  = "<< c            <<std::endl;    
00226       std::cout<< "Verified  = "<< ver          <<std::endl;    
00227       std::cout<< "Subdivs   = "<< subs         <<std::endl;
00228       if (ALLBOXES)
00229         std::cout<< "Reduced volume= "<< 
00230           as<double>(  100*(bx-v)/bx )<< "\%"     <<std::endl;  
00231       std::cout<< "#stack=     "<< cand         <<std::endl;    
00232       }
00233       sol<< psol;
00234 
00235       return sol;
00236     }

The documentation for this class was generated from the following file:
Generated on 6 Dec 2012 for realroot by  doxygen 1.6.1