My Project
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Member Functions | Private Attributes
rootContainer Class Reference

complex root finder for univariate polynomials based on laguers algorithm More...

#include <mpr_numeric.h>

Public Types

enum  rootType {
  none , cspecial , cspecialmu , det ,
  onepoly
}
 

Public Member Functions

 rootContainer ()
 
 ~rootContainer ()
 
void fillContainer (number *_coeffs, number *_ievpoint, const int _var, const int _tdg, const rootType _rt, const int _anz)
 
bool solver (const int polishmode=PM_NONE)
 
poly getPoly ()
 
gmp_complexoperator[] (const int i)
 
gmp_complexevPointCoord (const int i)
 
gmp_complexgetRoot (const int i)
 
bool swapRoots (const int from, const int to)
 
int getAnzElems ()
 
int getLDim ()
 
int getAnzRoots ()
 

Private Member Functions

 rootContainer (const rootContainer &v)
 
bool laguer_driver (gmp_complex **a, gmp_complex **roots, bool polish=true)
 Given the degree tdg and the tdg+1 complex coefficients ad[0..tdg] (generated from the number coefficients coeffs[0..tdg]) of the polynomial this routine successively calls "laguer" and finds all m complex roots in roots[0..tdg]. More...
 
bool isfloat (gmp_complex **a)
 
void divlin (gmp_complex **a, gmp_complex x, int j)
 
void divquad (gmp_complex **a, gmp_complex x, int j)
 
void solvequad (gmp_complex **a, gmp_complex **r, int &k, int &j)
 
void sortroots (gmp_complex **roots, int r, int c, bool isf)
 
void sortre (gmp_complex **r, int l, int u, int inc)
 
void laguer (gmp_complex **a, int m, gmp_complex *x, int *its, bool type)
 Given the degree m and the m+1 complex coefficients a[0..m] of the polynomial, and given the complex value x, this routine improves x by Laguerre's method until it converges, within the achievable roundoff limit, to a root of the given polynomial. More...
 
void computefx (gmp_complex **a, gmp_complex x, int m, gmp_complex &f0, gmp_complex &f1, gmp_complex &f2, gmp_float &ex, gmp_float &ef)
 
void computegx (gmp_complex **a, gmp_complex x, int m, gmp_complex &f0, gmp_complex &f1, gmp_complex &f2, gmp_float &ex, gmp_float &ef)
 
void checkimag (gmp_complex *x, gmp_float &e)
 

Private Attributes

int var
 
int tdg
 
number * coeffs
 
number * ievpoint
 
rootType rt
 
gmp_complex ** theroots
 
int anz
 
bool found_roots
 

Detailed Description

complex root finder for univariate polynomials based on laguers algorithm

Definition at line 65 of file mpr_numeric.h.

Member Enumeration Documentation

◆ rootType

Enumerator
none 
cspecial 
cspecialmu 
det 
onepoly 

Definition at line 68 of file mpr_numeric.h.

Constructor & Destructor Documentation

◆ rootContainer() [1/2]

rootContainer::rootContainer ( )

Definition at line 264 of file mpr_numeric.cc.

265{
266 rt=none;
267
268 coeffs= NULL;
269 ievpoint= NULL;
270 theroots= NULL;
271
272 found_roots= false;
273}
rootType rt
Definition: mpr_numeric.h:137
gmp_complex ** theroots
Definition: mpr_numeric.h:139
number * ievpoint
Definition: mpr_numeric.h:136
The main handler for Singular numbers which are suitable for Singular polynomials.
#define NULL
Definition: omList.c:12

◆ ~rootContainer()

rootContainer::~rootContainer ( )

Definition at line 277 of file mpr_numeric.cc.

278{
279 int i;
280 // free coeffs, ievpoint
281 if ( ievpoint != NULL )
282 {
283 for ( i=0; i < anz+2; i++ ) nDelete( ievpoint + i );
284 omFreeSize( (void *)ievpoint, (anz+2) * sizeof( number ) );
285 }
286
287 for ( i=0; i <= tdg; i++ )
288 if (coeffs[i]!=NULL) nDelete( coeffs + i );
289 omFreeSize( (void *)coeffs, (tdg+1) * sizeof( number ) );
290
291 // theroots loeschen
292 for ( i=0; i < tdg; i++ ) delete theroots[i];
293 omFreeSize( (void *) theroots, (tdg)*sizeof(gmp_complex*) );
294
295 //mprPROTnl("~rootContainer()");
296}
int i
Definition: cfEzgcd.cc:132
gmp_complex numbers based on
Definition: mpr_complex.h:179
#define nDelete(n)
Definition: numbers.h:16
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260

◆ rootContainer() [2/2]

rootContainer::rootContainer ( const rootContainer v)
private

Member Function Documentation

◆ checkimag()

void rootContainer::checkimag ( gmp_complex x,
gmp_float e 
)
private

Definition at line 617 of file mpr_numeric.cc.

618{
619 if(abs(x->imag())<abs(x->real())*e)
620 {
621 x->imag(0.0);
622 }
623}
Rational abs(const Rational &a)
Definition: GMPrat.cc:436
Variable x
Definition: cfModGcd.cc:4082

◆ computefx()

void rootContainer::computefx ( gmp_complex **  a,
gmp_complex  x,
int  m,
gmp_complex f0,
gmp_complex f1,
gmp_complex f2,
gmp_float ex,
gmp_float ef 
)
private

Definition at line 801 of file mpr_numeric.cc.

804{
805 int k;
806
807 f0= *a[m];
808 ef= abs(f0);
809 f1= gmp_complex(0.0);
810 f2= f1;
811 ex= abs(x);
812
813 for ( k= m-1; k >= 0; k-- )
814 {
815 f2 = ( x * f2 ) + f1;
816 f1 = ( x * f1 ) + f0;
817 f0 = ( x * f0 ) + *a[k];
818 ef = abs( f0 ) + ( ex * ef );
819 }
820}
int m
Definition: cfEzgcd.cc:128
int k
Definition: cfEzgcd.cc:99

◆ computegx()

void rootContainer::computegx ( gmp_complex **  a,
gmp_complex  x,
int  m,
gmp_complex f0,
gmp_complex f1,
gmp_complex f2,
gmp_float ex,
gmp_float ef 
)
private

Definition at line 822 of file mpr_numeric.cc.

825{
826 int k;
827
828 f0= *a[0];
829 ef= abs(f0);
830 f1= gmp_complex(0.0);
831 f2= f1;
832 ex= abs(x);
833
834 for ( k= 1; k <= m; k++ )
835 {
836 f2 = ( x * f2 ) + f1;
837 f1 = ( x * f1 ) + f0;
838 f0 = ( x * f0 ) + *a[k];
839 ef = abs( f0 ) + ( ex * ef );
840 }
841}

◆ divlin()

void rootContainer::divlin ( gmp_complex **  a,
gmp_complex  x,
int  j 
)
private

Definition at line 638 of file mpr_numeric.cc.

639{
640 int i;
641 gmp_float o(1.0);
642
643 if (abs(x)<o)
644 {
645 for (i= j-1; i > 0; i-- )
646 *a[i] += (*a[i+1]*x);
647 for (i= 0; i < j; i++ )
648 *a[i] = *a[i+1];
649 }
650 else
651 {
652 gmp_complex y(o/x);
653 for (i= 1; i < j; i++)
654 *a[i] += (*a[i-1]*y);
655 }
656}
const CanonicalForm int const CFList const Variable & y
Definition: facAbsFact.cc:53
int j
Definition: facHensel.cc:110

◆ divquad()

void rootContainer::divquad ( gmp_complex **  a,
gmp_complex  x,
int  j 
)
private

Definition at line 658 of file mpr_numeric.cc.

659{
660 int i;
661 gmp_float o(1.0),p(x.real()+x.real()),
662 q((x.real()*x.real())+(x.imag()*x.imag()));
663
664 if (abs(x)<o)
665 {
666 *a[j-1] += (*a[j]*p);
667 for (i= j-2; i > 1; i-- )
668 *a[i] += ((*a[i+1]*p)-(*a[i+2]*q));
669 for (i= 0; i < j-1; i++ )
670 *a[i] = *a[i+2];
671 }
672 else
673 {
674 p = p/q;
675 q = o/q;
676 *a[1] += (*a[0]*p);
677 for (i= 2; i < j-1; i++)
678 *a[i] += ((*a[i-1]*p)-(*a[i-2]*q));
679 }
680}
int p
Definition: cfModGcd.cc:4078

◆ evPointCoord()

gmp_complex & rootContainer::evPointCoord ( const int  i)

Definition at line 388 of file mpr_numeric.cc.

389{
390 if (! ((i >= 0) && (i < anz+2) ) )
391 WarnS("rootContainer::evPointCoord: index out of range");
392 if (ievpoint == NULL)
393 WarnS("rootContainer::evPointCoord: ievpoint == NULL");
394
395 if ( (rt == cspecialmu) && found_roots ) // FIX ME
396 {
397 if ( ievpoint[i] != NULL )
398 {
399 gmp_complex *tmp= new gmp_complex();
400 *tmp= numberToComplex(ievpoint[i], currRing->cf);
401 return *tmp;
402 }
403 else
404 {
405 Warn("rootContainer::evPointCoord: NULL index %d",i);
406 }
407 }
408
409 // warning
410 Warn("rootContainer::evPointCoord: Wrong index %d, found_roots %s",i,found_roots?"true":"false");
411 gmp_complex *tmp= new gmp_complex();
412 return *tmp;
413}
#define Warn
Definition: emacs.cc:77
#define WarnS
Definition: emacs.cc:78
gmp_complex numberToComplex(number num, const coeffs r)
Definition: mpr_complex.h:312
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:13

◆ fillContainer()

void rootContainer::fillContainer ( number *  _coeffs,
number *  _ievpoint,
const int  _var,
const int  _tdg,
const rootType  _rt,
const int  _anz 
)

Definition at line 300 of file mpr_numeric.cc.

303{
304 int i;
305 number nn= nInit(0);
306 var=_var;
307 tdg=_tdg;
308 coeffs=_coeffs;
309 rt=_rt;
310 anz=_anz;
311
312 for ( i=0; i <= tdg; i++ )
313 {
314 if ( nEqual(coeffs[i],nn) )
315 {
316 nDelete( &coeffs[i] );
317 coeffs[i]=NULL;
318 }
319 }
320 nDelete( &nn );
321
322 if ( rt == cspecialmu && _ievpoint ) // copy ievpoint
323 {
324 ievpoint= (number *)omAlloc( (anz+2) * sizeof( number ) );
325 for (i=0; i < anz+2; i++) ievpoint[i]= nCopy( _ievpoint[i] );
326 }
327
328 theroots= NULL;
329 found_roots= false;
330}
#define nEqual(n1, n2)
Definition: numbers.h:20
#define nCopy(n)
Definition: numbers.h:15
#define nInit(i)
Definition: numbers.h:24
#define omAlloc(size)
Definition: omAllocDecl.h:210

◆ getAnzElems()

int rootContainer::getAnzElems ( )
inline

Definition at line 95 of file mpr_numeric.h.

95{ return anz; }

◆ getAnzRoots()

int rootContainer::getAnzRoots ( )
inline

Definition at line 97 of file mpr_numeric.h.

97{ return tdg; }

◆ getLDim()

int rootContainer::getLDim ( )
inline

Definition at line 96 of file mpr_numeric.h.

96{ return anz; }

◆ getPoly()

poly rootContainer::getPoly ( )

Definition at line 334 of file mpr_numeric.cc.

335{
336 int i;
337
338 poly result= NULL;
339 poly ppos;
340
341 if ( (rt == cspecial) || ( rt == cspecialmu ) )
342 {
343 for ( i= tdg; i >= 0; i-- )
344 {
345 if ( coeffs[i] )
346 {
347 poly p= pOne();
348 //pSetExp( p, var+1, i);
349 pSetExp( p, 1, i);
350 pSetCoeff( p, nCopy( coeffs[i] ) );
351 pSetm( p );
352 if (result)
353 {
354 ppos->next=p;
355 ppos=ppos->next;
356 }
357 else
358 {
359 result=p;
360 ppos=p;
361 }
362
363 }
364 }
365 if (result!=NULL) pSetm( result );
366 }
367
368 return result;
369}
return result
Definition: facAbsBiFact.cc:75
#define pSetm(p)
Definition: polys.h:271
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
Definition: polys.h:31
#define pSetExp(p, i, v)
Definition: polys.h:42
#define pOne()
Definition: polys.h:315

◆ getRoot()

gmp_complex * rootContainer::getRoot ( const int  i)
inline

Definition at line 88 of file mpr_numeric.h.

89 {
90 return theroots[i];
91 }

◆ isfloat()

bool rootContainer::isfloat ( gmp_complex **  a)
private

Definition at line 625 of file mpr_numeric.cc.

626{
627 gmp_float z(0.0);
628 gmp_complex *b;
629 for (int i=tdg; i >= 0; i-- )
630 {
631 b = &(*a[i]);
632 if (!(b->imag()==z))
633 return false;
634 }
635 return true;
636}
CanonicalForm b
Definition: cfModGcd.cc:4103

◆ laguer()

void rootContainer::laguer ( gmp_complex **  a,
int  m,
gmp_complex x,
int *  its,
bool  type 
)
private

Given the degree m and the m+1 complex coefficients a[0..m] of the polynomial, and given the complex value x, this routine improves x by Laguerre's method until it converges, within the achievable roundoff limit, to a root of the given polynomial.

The number of iterations taken is returned at its.

Definition at line 550 of file mpr_numeric.cc.

551{
552 int iter,j;
553 gmp_float zero(0.0),one(1.0),deg(m);
554 gmp_float abx_g, err_g, fabs;
555 gmp_complex dx, x1, b, d, f, g, h, sq, gp, gm, g2;
556 gmp_float frac_g[MR+1] = { 0.0, 0.5, 0.25, 0.75, 0.125, 0.375, 0.625, 0.875, 1.0 };
557
558 gmp_float epss(0.1);
559 mpf_pow_ui(*epss._mpfp(),*epss.mpfp(),gmp_output_digits);
560
561 for ( iter= 1; iter <= MAXIT; iter++ )
562 {
564 *its=iter;
565 if (type)
566 computefx(a,*x,m,b,d,f,abx_g,err_g);
567 else
568 computegx(a,*x,m,b,d,f,abx_g,err_g);
569 err_g *= epss; // EPSS;
570
571 fabs = abs(b);
572 if (fabs <= err_g)
573 {
574 if ((fabs==zero) || (abs(d)==zero)) return;
575 *x -= (b/d); // a last newton-step
576 goto ende;
577 }
578
579 g= d / b;
580 g2 = g * g;
581 h= g2 - (((f+f) / b ));
582 sq= sqrt(( ( h * deg ) - g2 ) * (deg - one));
583 gp= g + sq;
584 gm= g - sq;
585 if (abs(gp)<abs(gm))
586 {
587 dx = deg/gm;
588 }
589 else
590 {
591 if((gp.real()==zero)&&(gp.imag()==zero))
592 {
593 dx.real(cos((mprfloat)iter));
594 dx.imag(sin((mprfloat)iter));
595 dx = dx*(one+abx_g);
596 }
597 else
598 {
599 dx = deg/gp;
600 }
601 }
602 x1= *x - dx;
603
604 if (*x == x1) goto ende;
605
606 j = iter%MMOD;
607 if (j==0) j=MT;
608 if ( j % MT ) *x= x1;
609 else *x -= ( dx * frac_g[ j / MT ] );
610 }
611
612 *its= MAXIT+1;
613ende:
614 checkimag(x,epss);
615}
g
Definition: cfModGcd.cc:4090
CanonicalForm gp
Definition: cfModGcd.cc:4102
FILE * f
Definition: checklibs.c:9
gmp_float imag() const
Definition: mpr_complex.h:235
gmp_float real() const
Definition: mpr_complex.h:234
void computegx(gmp_complex **a, gmp_complex x, int m, gmp_complex &f0, gmp_complex &f1, gmp_complex &f2, gmp_float &ex, gmp_float &ef)
Definition: mpr_numeric.cc:822
void computefx(gmp_complex **a, gmp_complex x, int m, gmp_complex &f0, gmp_complex &f1, gmp_complex &f2, gmp_float &ex, gmp_float &ef)
Definition: mpr_numeric.cc:801
void checkimag(gmp_complex *x, gmp_float &e)
Definition: mpr_numeric.cc:617
CFFListIterator iter
Definition: facAbsBiFact.cc:53
STATIC_VAR Poly * h
Definition: janet.cc:971
EXTERN_VAR size_t gmp_output_digits
Definition: mpr_base.h:115
gmp_float sin(const gmp_float &a)
Definition: mpr_complex.cc:333
gmp_float sqrt(const gmp_float &a)
Definition: mpr_complex.cc:327
gmp_float cos(const gmp_float &a)
Definition: mpr_complex.cc:338
#define mprSTICKYPROT(msg)
Definition: mpr_global.h:54
#define ST_ROOTS_LG
Definition: mpr_global.h:82
double mprfloat
Definition: mpr_global.h:17
#define MT
Definition: mpr_numeric.cc:258
#define MMOD
Definition: mpr_numeric.cc:259
#define MR
Definition: mpr_numeric.cc:257
#define MAXIT
Definition: mpr_numeric.cc:260

◆ laguer_driver()

bool rootContainer::laguer_driver ( gmp_complex **  a,
gmp_complex **  roots,
bool  polish = true 
)
private

Given the degree tdg and the tdg+1 complex coefficients ad[0..tdg] (generated from the number coefficients coeffs[0..tdg]) of the polynomial this routine successively calls "laguer" and finds all m complex roots in roots[0..tdg].

The bool var "polish" should be input as "true" if polishing (also by "laguer") is desired, "false" if the roots will be subsequently polished by other means.

Definition at line 467 of file mpr_numeric.cc.

468{
469 int i,j,k,its;
470 gmp_float zero(0.0);
471 gmp_complex x(0.0),o(1.0);
472 bool ret= true, isf=isfloat(a), type=true;
473
474 gmp_complex ** ad= (gmp_complex**)omAlloc( (tdg+1)*sizeof(gmp_complex*) );
475 for ( i=0; i <= tdg; i++ ) ad[i]= new gmp_complex( *a[i] );
476
477 k = 0;
478 i = tdg;
479 j = i-1;
480 while (i>2)
481 {
482 // run laguer alg
483 x = zero;
484 laguer(ad, i, &x, &its, type);
485 if ( its > MAXIT )
486 {
487 type = !type;
488 x = zero;
489 laguer(ad, i, &x, &its, type);
490 }
491
493 if ( its > MAXIT )
494 { // error
495 WarnS("Laguerre solver: Too many iterations!");
496 ret= false;
497 goto theend;
498 }
499 if ( polish )
500 {
501 laguer( a, tdg, &x, &its , type);
502 if ( its > MAXIT )
503 { // error
504 WarnS("Laguerre solver: Too many iterations in polish!");
505 ret= false;
506 goto theend;
507 }
508 }
509 if ((!type)&&(!((x.real()==zero)&&(x.imag()==zero)))) x = o/x;
510 if (x.imag() == zero)
511 {
512 *roots[k] = x;
513 k++;
514 divlin(ad,x,i);
515 i--;
516 }
517 else
518 {
519 if(isf)
520 {
521 *roots[j] = x;
522 *roots[j-1]= gmp_complex(x.real(),-x.imag());
523 j -= 2;
524 divquad(ad,x,i);
525 i -= 2;
526 }
527 else
528 {
529 *roots[j] = x;
530 j--;
531 divlin(ad,x,i);
532 i--;
533 }
534 }
535 type = !type;
536 }
537 solvequad(ad,roots,k,j);
538 sortroots(roots,k,j,isf);
539
540theend:
541 mprSTICKYPROT("\n");
542 for ( i=0; i <= tdg; i++ ) delete ad[i];
543 omFreeSize( (void *) ad, (tdg+1)*sizeof( gmp_complex* ));
544
545 return ret;
546}
void laguer(gmp_complex **a, int m, gmp_complex *x, int *its, bool type)
Given the degree m and the m+1 complex coefficients a[0..m] of the polynomial, and given the complex ...
Definition: mpr_numeric.cc:550
void divlin(gmp_complex **a, gmp_complex x, int j)
Definition: mpr_numeric.cc:638
void sortroots(gmp_complex **roots, int r, int c, bool isf)
Definition: mpr_numeric.cc:735
void divquad(gmp_complex **a, gmp_complex x, int j)
Definition: mpr_numeric.cc:658
void solvequad(gmp_complex **a, gmp_complex **r, int &k, int &j)
Definition: mpr_numeric.cc:682
bool isfloat(gmp_complex **a)
Definition: mpr_numeric.cc:625
#define ST_ROOTS_LGSTEP
Definition: mpr_global.h:80

◆ operator[]()

gmp_complex & rootContainer::operator[] ( const int  i)
inline

Definition at line 82 of file mpr_numeric.h.

83 {
84 return *theroots[i];
85 }

◆ solvequad()

void rootContainer::solvequad ( gmp_complex **  a,
gmp_complex **  r,
int &  k,
int &  j 
)
private

Definition at line 682 of file mpr_numeric.cc.

683{
684 gmp_float zero(0.0);
685
686 if ((j>k)
687 &&((!(*a[2]).real().isZero())||(!(*a[2]).imag().isZero())))
688 {
689 gmp_complex sq(zero);
690 gmp_complex h1(*a[1]/(*a[2] + *a[2])), h2(*a[0] / *a[2]);
691 gmp_complex disk((h1 * h1) - h2);
692 if (disk.imag().isZero())
693 {
694 if (disk.real()<zero)
695 {
696 sq.real(zero);
697 sq.imag(sqrt(-disk.real()));
698 }
699 else
700 sq = (gmp_complex)sqrt(disk.real());
701 }
702 else
703 sq = sqrt(disk);
704 *r[k+1] = sq - h1;
705 sq += h1;
706 *r[k] = (gmp_complex)0.0-sq;
707 if(sq.imag().isZero())
708 {
709 k = j;
710 j++;
711 }
712 else
713 {
714 j = k;
715 k--;
716 }
717 }
718 else
719 {
720 if (((*a[1]).real().isZero()) && ((*a[1]).imag().isZero()))
721 {
722 WerrorS("precision lost, try again with higher precision");
723 }
724 else
725 {
726 *r[k]= (gmp_complex)0.0-(*a[0] / *a[1]);
727 if(r[k]->imag().isZero())
728 j++;
729 else
730 k--;
731 }
732 }
733}
bool isZero(const CFArray &A)
checks if entries of A are zero
void WerrorS(const char *s)
Definition: feFopen.cc:24

◆ solver()

bool rootContainer::solver ( const int  polishmode = PM_NONE)

Definition at line 437 of file mpr_numeric.cc.

438{
439 int i;
440
441 // there are maximal tdg roots, so *roots ranges form 0 to tdg-1.
443 for ( i=0; i < tdg; i++ ) theroots[i]= new gmp_complex();
444
445 // copy the coefficients of type number to type gmp_complex
446 gmp_complex **ad= (gmp_complex**)omAlloc( (tdg+1)*sizeof(gmp_complex*) );
447 for ( i=0; i <= tdg; i++ )
448 {
449 ad[i]= new gmp_complex();
450 if ( coeffs[i] ) *ad[i] = numberToComplex( coeffs[i], currRing->cf );
451 }
452
453 // now solve
454 found_roots= laguer_driver( ad, theroots, polishmode != 0 );
455 if (!found_roots)
456 WarnS("rootContainer::solver: No roots found!");
457
458 // free memory
459 for ( i=0; i <= tdg; i++ ) delete ad[i];
460 omFreeSize( (void *) ad, (tdg+1)*sizeof(gmp_complex*) );
461
462 return found_roots;
463}
bool laguer_driver(gmp_complex **a, gmp_complex **roots, bool polish=true)
Given the degree tdg and the tdg+1 complex coefficients ad[0..tdg] (generated from the number coeffic...
Definition: mpr_numeric.cc:467

◆ sortre()

void rootContainer::sortre ( gmp_complex **  r,
int  l,
int  u,
int  inc 
)
private

Definition at line 754 of file mpr_numeric.cc.

755{
756 int pos,i;
757 gmp_complex *x,*y;
758
759 pos = l;
760 x = r[pos];
761 for (i=l+inc; i<=u; i+=inc)
762 {
763 if (r[i]->real()<x->real())
764 {
765 pos = i;
766 x = r[pos];
767 }
768 }
769 if (pos>l)
770 {
771 if (inc==1)
772 {
773 for (i=pos; i>l; i--)
774 r[i] = r[i-1];
775 r[l] = x;
776 }
777 else
778 {
779 y = r[pos+1];
780 for (i=pos+1; i+1>l; i--)
781 r[i] = r[i-2];
782 if (x->imag()>y->imag())
783 {
784 r[l] = x;
785 r[l+1] = y;
786 }
787 else
788 {
789 r[l] = y;
790 r[l+1] = x;
791 }
792 }
793 }
794 else if ((inc==2)&&(x->imag()<r[l+1]->imag()))
795 {
796 r[l] = r[l+1];
797 r[l+1] = x;
798 }
799}
int l
Definition: cfEzgcd.cc:100

◆ sortroots()

void rootContainer::sortroots ( gmp_complex **  roots,
int  r,
int  c,
bool  isf 
)
private

Definition at line 735 of file mpr_numeric.cc.

736{
737 int j;
738
739 for (j=0; j<r; j++) // the real roots
740 sortre(ro, j, r, 1);
741 if (c>=tdg) return;
742 if (isf)
743 {
744 for (j=c; j+2<tdg; j+=2) // the complex roots for a real poly
745 sortre(ro, j, tdg-1, 2);
746 }
747 else
748 {
749 for (j=c; j+1<tdg; j++) // the complex roots for a general poly
750 sortre(ro, j, tdg-1, 1);
751 }
752}
void sortre(gmp_complex **r, int l, int u, int inc)
Definition: mpr_numeric.cc:754

◆ swapRoots()

bool rootContainer::swapRoots ( const int  from,
const int  to 
)

Definition at line 417 of file mpr_numeric.cc.

418{
419 if ( found_roots && ( from >= 0) && ( from < tdg ) && ( to >= 0) && ( to < tdg ) )
420 {
421 if ( to != from )
422 {
423 gmp_complex tmp( *theroots[from] );
424 *theroots[from]= *theroots[to];
425 *theroots[to]= tmp;
426 }
427 return true;
428 }
429
430 // warning
431 Warn(" rootContainer::changeRoots: Wrong index %d, %d",from,to);
432 return false;
433}

Field Documentation

◆ anz

int rootContainer::anz
private

Definition at line 141 of file mpr_numeric.h.

◆ coeffs

number* rootContainer::coeffs
private

Definition at line 135 of file mpr_numeric.h.

◆ found_roots

bool rootContainer::found_roots
private

Definition at line 142 of file mpr_numeric.h.

◆ ievpoint

number* rootContainer::ievpoint
private

Definition at line 136 of file mpr_numeric.h.

◆ rt

rootType rootContainer::rt
private

Definition at line 137 of file mpr_numeric.h.

◆ tdg

int rootContainer::tdg
private

Definition at line 133 of file mpr_numeric.h.

◆ theroots

gmp_complex** rootContainer::theroots
private

Definition at line 139 of file mpr_numeric.h.

◆ var

int rootContainer::var
private

Definition at line 132 of file mpr_numeric.h.


The documentation for this class was generated from the following files: