1 | #ifndef MPR_NUMERIC_H |
---|
2 | #define MPR_NUMERIC_H |
---|
3 | /**************************************** |
---|
4 | * Computer Algebra System SINGULAR * |
---|
5 | ****************************************/ |
---|
6 | |
---|
7 | /* $Id: mpr_numeric.h,v 1.3 1999-06-29 09:03:46 wenk Exp $ */ |
---|
8 | |
---|
9 | /* |
---|
10 | * ABSTRACT - multipolynomial resultants - numeric stuff |
---|
11 | * ( root finder, vandermonde system solver, simplex ) |
---|
12 | */ |
---|
13 | |
---|
14 | //-> include & define stuff |
---|
15 | #include "numbers.h" |
---|
16 | #include "mpr_global.h" |
---|
17 | #include "mpr_complex.h" |
---|
18 | |
---|
19 | // define polish mode when finding roots |
---|
20 | #define PM_NONE 0 |
---|
21 | #define PM_POLISH 1 |
---|
22 | #define PM_CORRUPT 2 |
---|
23 | //<- |
---|
24 | |
---|
25 | //-> vandermonde system solver |
---|
26 | /** |
---|
27 | * vandermonde system solver for interpolating polynomials from their values |
---|
28 | */ |
---|
29 | class vandermonde |
---|
30 | { |
---|
31 | public: |
---|
32 | vandermonde( const long _cn, const long _n, |
---|
33 | const long _maxdeg, number *_p, const bool _homog = true ); |
---|
34 | ~vandermonde(); |
---|
35 | |
---|
36 | /** Solves the Vandermode linear system |
---|
37 | * \sum_{i=1}^{n} x_i^k-1 w_i = q_k, k=1,..,n. |
---|
38 | * Any computations are done using type number to get high pecision results. |
---|
39 | * @param q n-tuple of results (right hand of equations) |
---|
40 | * @return w n-tuple of coefficients of resulting polynomial, lowest deg first |
---|
41 | */ |
---|
42 | number * interpolateDense( const number * q ); |
---|
43 | |
---|
44 | poly numvec2poly(const number * q ); |
---|
45 | |
---|
46 | private: |
---|
47 | void init(); |
---|
48 | |
---|
49 | private: |
---|
50 | long n; // number of variables |
---|
51 | long cn; // real number of coefficients of poly to interpolate |
---|
52 | long maxdeg; // degree of the polynomial to interpolate |
---|
53 | long l; // max number of coefficients in poly of deg maxdeg = (maxdeg+1)^n |
---|
54 | |
---|
55 | number *p; // evaluation point |
---|
56 | number *x; // coefficients, determinend by init() from *p |
---|
57 | |
---|
58 | bool homog; |
---|
59 | }; |
---|
60 | //<- |
---|
61 | |
---|
62 | //-> rootContainer |
---|
63 | /** |
---|
64 | * complex root finder for univariate polynomials based on laguers algorithm |
---|
65 | */ |
---|
66 | class rootContainer |
---|
67 | { |
---|
68 | public: |
---|
69 | enum rootType { none, cspecial, cspecialmu, det, onepoly }; |
---|
70 | |
---|
71 | rootContainer(); |
---|
72 | ~rootContainer(); |
---|
73 | |
---|
74 | void fillContainer( number *_coeffs, number *_ievpoint, |
---|
75 | const int _var, const int _tdg, |
---|
76 | const rootType _rt, const int _anz ); |
---|
77 | |
---|
78 | bool solver( const int polishmode= PM_NONE ); |
---|
79 | |
---|
80 | poly getPoly(); |
---|
81 | |
---|
82 | //gmp_complex & operator[] ( const int i ); |
---|
83 | inline gmp_complex & operator[] ( const int i ) { |
---|
84 | return *theroots[i]; |
---|
85 | } |
---|
86 | gmp_complex & evPointCoord( const int i ); |
---|
87 | |
---|
88 | inline gmp_complex * getRoot( const int i ) { |
---|
89 | return theroots[i]; |
---|
90 | } |
---|
91 | |
---|
92 | bool swapRoots( const int from, const int to ); |
---|
93 | |
---|
94 | int getAnzElems() { return anz; } |
---|
95 | int getLDim() { return anz; } |
---|
96 | int getAnzRoots() { return tdg; } |
---|
97 | |
---|
98 | private: |
---|
99 | rootContainer( const rootContainer & v ); |
---|
100 | |
---|
101 | /** Given the degree tdg and the tdg+1 complex coefficients ad[0..tdg] |
---|
102 | * (generated from the number coefficients coeffs[0..tdg]) of the polynomial |
---|
103 | * this routine successively calls "laguer" and finds all m complex roots in |
---|
104 | * roots[0..tdg]. The bool var "polish" should be input as "true" if polishing |
---|
105 | * (also by "laguer") is desired, "false" if the roots will be subsequently |
---|
106 | * polished by other means. |
---|
107 | */ |
---|
108 | bool laguer_driver( gmp_complex ** a, gmp_complex ** roots, bool polish = true ); |
---|
109 | |
---|
110 | /** Given the degree m and the m+1 complex coefficients a[0..m] of the |
---|
111 | * polynomial, and given the complex value x, this routine improves x by |
---|
112 | * Laguerre's method until it converges, within the achievable roundoff limit, |
---|
113 | * to a root of the given polynomial. The number of iterations taken is |
---|
114 | * returned at its. |
---|
115 | */ |
---|
116 | void laguer(gmp_complex ** a, int m, gmp_complex * x, int * its); |
---|
117 | |
---|
118 | int var; |
---|
119 | int tdg; |
---|
120 | |
---|
121 | number * coeffs; |
---|
122 | number * ievpoint; |
---|
123 | rootType rt; |
---|
124 | |
---|
125 | gmp_complex ** theroots; |
---|
126 | |
---|
127 | int anz; |
---|
128 | bool found_roots; |
---|
129 | }; |
---|
130 | //<- |
---|
131 | |
---|
132 | //-> class rootArranger |
---|
133 | class rootArranger |
---|
134 | { |
---|
135 | public: |
---|
136 | rootArranger( rootContainer ** _roots, |
---|
137 | rootContainer ** _mu, |
---|
138 | const int _howclean = PM_CORRUPT ); |
---|
139 | ~rootArranger() {} |
---|
140 | |
---|
141 | void solve_all(); |
---|
142 | void arrange(); |
---|
143 | |
---|
144 | lists listOfRoots( const unsigned int oprec ); |
---|
145 | |
---|
146 | const bool success() { return found_roots; } |
---|
147 | |
---|
148 | private: |
---|
149 | rootArranger( const rootArranger & ); |
---|
150 | |
---|
151 | rootContainer ** roots; |
---|
152 | rootContainer ** mu; |
---|
153 | |
---|
154 | int howclean; |
---|
155 | int rc,mc; |
---|
156 | bool found_roots; |
---|
157 | }; |
---|
158 | //<- |
---|
159 | |
---|
160 | //-> simplex computation |
---|
161 | // (used by sparse matrix construction) |
---|
162 | #define SIMPLEX_EPS 1.0e-12 |
---|
163 | |
---|
164 | void simplx( mprfloat **a, int m, int n, int m1, int m2, int m3, int *icase, int izrov[], int iposv[] ); |
---|
165 | //<- |
---|
166 | |
---|
167 | #endif MPR_NUMERIC_H |
---|
168 | |
---|
169 | // local Variables: *** |
---|
170 | // folded-file: t *** |
---|
171 | // compile-command-1: "make installg" *** |
---|
172 | // compile-command-2: "make install" *** |
---|
173 | // End: *** |
---|