bsearch_newton< real_t > Struct Template Reference

#include <system_descartes1d.hpp>

List of all members.

Public Member Functions

template<class In > bsearch_newton (const In &bzrep, unsigned sz)
~bsearch_newton ()
void reach (real_t *lbzrep, real_t &a, real_t &b, const real_t &eps)

Public Attributes

real_t * m_data
real_t * m_mons
unsigned m_sz

Detailed Description

template<class real_t>
struct mmx::solvers::bsearch_newton< real_t >

Definition at line 88 of file system_descartes1d.hpp.

Constructor & Destructor Documentation

bsearch_newton	(	const In &	bzrep,
		unsigned	sz
	)			`[inline]`

Definition at line 95 of file system_descartes1d.hpp.

References mmx::sparse::copy(), bsearch_newton< real_t >::m_data, and bsearch_newton< real_t >::m_mons.

00095                                                      : m_sz(sz) 
00096     {
00097       m_data = new real_t[ 2*sz ];
00098       std::copy( bzrep, bzrep + sz, m_data );
00099       m_mons = m_data + sz;
00100       std::copy( m_data, m_data + sz, m_mons );
00101       bernstein::bzenv<real_t>::_default_->toMonoms( m_mons, sz );
00102     };

~bsearch_newton ( ) [inline]

Definition at line 103 of file system_descartes1d.hpp.

References bsearch_newton< real_t >::m_data.

00103 { delete[] m_data; };

Member Function Documentation

void reach	(	real_t *	lbzrep,
		real_t &	a,
		real_t &	b,
		const real_t &	eps
	)			`[inline]`

Definition at line 106 of file system_descartes1d.hpp.

References mmx::upoldse_::dhorner(), mmx::upoldse_::horner(), bsearch_newton< real_t >::m_mons, and bsearch_newton< real_t >::m_sz.

00107     {
00108       real_t m;
00109       if ( lbzrep[m_sz-1] > lbzrep[0] ) 
00110         do
00111           {
00112             real_t p,dp,x;
00113             /* évaluation du polynôme en m (p) et de sa derivée (dp) */
00114             upoldse_::dhorner( p, dp, m_mons, m_sz, m = (a+b)/2.0 );
00115             /* étape de bissection                                   */
00116             if ( p < 0  ) a = m; else b = m; 
00117             if ( dp > eps ) /* si la valeur de la derivée n'est pas trop faible */
00118               {
00119                 /* on fait une itération de la méthode de newton */
00120                 x = m - p/dp;
00121                 /* si le résultat est dans l'intervalle de la bissection */
00122                 if ( x > a && x < b )
00123                   {
00124                     /* évaluation du polynôme en ce point */
00125                     real_t xp = upoldse_::horner( m_mons, m_sz, x ); 
00126                     /* correction d'une des deux bornes */
00127                     if ( xp < 0 ) a = x; else b = x;
00128                     /* recherche de la seconde          */
00129                     if ( xp < 0 )
00130                       {
00131                         real_t step = eps;
00132                         while( xp < 0 && x < b )
00133                           {
00134                             x += step;
00135                             xp = upoldse_::horner( m_mons, m_sz, x );
00136                             step *= 2;
00137                           };
00138                         if ( x < b ) b = x;
00139                       }
00140                     else
00141                       {
00142                         real_t step = eps;
00143                         while( xp > 0 && x > a )
00144                           {
00145                             x -= step;
00146                             xp = upoldse_::horner( m_mons, m_sz, x );
00147                             step *= 2;
00148                           };
00149                         if ( x > a ) a = x;
00150                       };
00151                   };
00152               };
00153           }
00154         while( b-a > eps );
00155       else
00156         do
00157           {
00158             real_t p,dp,x;
00159             upoldse_::dhorner( p, dp, m_mons, m_sz, m = (a+b)/2.0 );
00160             if ( p < 0  ) b = m; else a = m; 
00161             if ( dp > eps )
00162               {
00163                 x = m - p/dp;        
00164                 if ( x > a && x < b )
00165                   {
00166                     real_t xp = upoldse_::horner( m_mons, m_sz, x );
00167                     if ( xp < 0 ) b = x; else a = x;
00168                     if ( xp < 0 )
00169                       {
00170                         real_t step = eps/2;
00171                         while( xp < 0 && x > a )
00172                           {
00173                             x -= step;
00174                             xp = upoldse_::horner( m_mons, m_sz, x );
00175                             step *= 2;
00176                           };
00177                         if ( x > a  ) a = x;
00178                       }
00179                     else
00180                       {
00181                         real_t step = eps/2;
00182                         while( xp > 0 && x < b )
00183                           {
00184                             x += step;
00185                             xp = upoldse_::horner( m_mons, m_sz, x );
00186                             step *= 2;
00187                           };
00188                         if ( x < b ) b = x;
00189                       };
00190                   };
00191               };
00192           }
00193         while( b-a > eps ); 
00194     };