root / pobysoPythonSage / src / sageSLZ / sagePolynomialOperations.sage @ 266
Historique | Voir | Annoter | Télécharger (51,7 ko)
1 | 166 | storres | load("/home/storres/recherche/arithmetique/pobysoPythonSage/src/sageSLZ/sageMatrixOperations.sage") |
---|---|---|---|
2 | 166 | storres | #load(str('/home/storres/recherche/arithmetique/pobysoPythonSage/src/sageSLZ/sageMatrixOperations.sage')) |
3 | 246 | storres | sys.stderr.write("sagePolynomialOperations loading...\n") |
4 | 202 | storres | |
5 | 106 | storres | def spo_add_polynomial_coeffs_to_matrix_row(poly, |
6 | 83 | storres | knownMonomials, |
7 | 83 | storres | protoMatrixRows, |
8 | 83 | storres | columnsWidth=0): |
9 | 80 | storres | """ |
10 | 106 | storres | For a given polynomial , |
11 | 80 | storres | add the coefficients of the protoMatrix (a list of proto matrix rows). |
12 | 80 | storres | Coefficients are added to the protoMatrix row in the order imposed by the |
13 | 80 | storres | monomials discovery list (the knownMonomials list) built as construction |
14 | 80 | storres | goes on. |
15 | 83 | storres | As a bonus, data can be printed out for a visual check. |
16 | 106 | storres | poly : the polynomial; in argument; |
17 | 106 | storres | knownMonomials : the list of the already known monomials; will determine |
18 | 106 | storres | the order of the coefficients appending to a row; in-out |
19 | 106 | storres | argument (new monomials may be discovered and then |
20 | 106 | storres | appended the the knowMonomials list); |
21 | 80 | storres | protoMatrixRows: a list of lists, each one holding the coefficients of the |
22 | 106 | storres | monomials of a polynomial; in-out argument: a new row is |
23 | 106 | storres | added at each call; |
24 | 80 | storres | columnWith : the width, in characters, of the displayed column ; if 0, |
25 | 106 | storres | do not display anything; in argument. |
26 | 80 | storres | """ |
27 | 106 | storres | pMonomials = poly.monomials() |
28 | 106 | storres | pCoefficients = poly.coefficients() |
29 | 80 | storres | # We have started with the smaller degrees in the first variable. |
30 | 80 | storres | pMonomials.reverse() |
31 | 80 | storres | pCoefficients.reverse() |
32 | 80 | storres | # New empty proto matrix row. |
33 | 80 | storres | protoMatrixRowCoefficients = [] |
34 | 80 | storres | # We work according to the order of the already known monomials |
35 | 80 | storres | # No known monomials yet: add the pMonomials to knownMonomials |
36 | 80 | storres | # and add the coefficients to the proto matrix row. |
37 | 80 | storres | if len(knownMonomials) == 0: |
38 | 80 | storres | for pmIdx in xrange(0, len(pMonomials)): |
39 | 80 | storres | knownMonomials.append(pMonomials[pmIdx]) |
40 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[pmIdx]) |
41 | 80 | storres | if columnsWidth != 0: |
42 | 80 | storres | monomialAsString = str(pCoefficients[pmIdx]) + " " + \ |
43 | 80 | storres | str(pMonomials[pmIdx]) |
44 | 80 | storres | print monomialAsString, " " * \ |
45 | 80 | storres | (columnsWidth - len(monomialAsString)), |
46 | 80 | storres | # There are some known monomials. We search for them in pMonomials and |
47 | 80 | storres | # add their coefficients to the proto matrix row. |
48 | 80 | storres | else: |
49 | 80 | storres | for knownMonomialIndex in xrange(0,len(knownMonomials)): |
50 | 80 | storres | # We lazily use an exception here since pMonomials.index() function |
51 | 80 | storres | # may fail throwing the ValueError exception. |
52 | 80 | storres | try: |
53 | 80 | storres | indexInPmonomials = \ |
54 | 80 | storres | pMonomials.index(knownMonomials[knownMonomialIndex]) |
55 | 80 | storres | if columnsWidth != 0: |
56 | 80 | storres | monomialAsString = str(pCoefficients[indexInPmonomials]) + \ |
57 | 80 | storres | " " + str(knownMonomials[knownMonomialIndex]) |
58 | 80 | storres | print monomialAsString, " " * \ |
59 | 80 | storres | (columnsWidth - len(monomialAsString)), |
60 | 155 | storres | # Add the coefficient to the proto matrix row and delete the |
61 | 80 | storres | # known monomial from the current pMonomial list |
62 | 168 | storres | # (and the corresponding coefficient as well). |
63 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[indexInPmonomials]) |
64 | 80 | storres | del pMonomials[indexInPmonomials] |
65 | 80 | storres | del pCoefficients[indexInPmonomials] |
66 | 80 | storres | # The knownMonomials element is not in pMonomials |
67 | 80 | storres | except ValueError: |
68 | 80 | storres | protoMatrixRowCoefficients.append(0) |
69 | 80 | storres | if columnsWidth != 0: |
70 | 80 | storres | monomialAsString = "0" + " "+ \ |
71 | 80 | storres | str(knownMonomials[knownMonomialIndex]) |
72 | 80 | storres | print monomialAsString, " " * \ |
73 | 80 | storres | (columnsWidth - len(monomialAsString)), |
74 | 80 | storres | # End for knownMonomialKey loop. |
75 | 80 | storres | # We now append the remaining monomials of pMonomials to knownMonomials |
76 | 80 | storres | # and the corresponding coefficients to proto matrix row. |
77 | 80 | storres | for pmIdx in xrange(0, len(pMonomials)): |
78 | 80 | storres | knownMonomials.append(pMonomials[pmIdx]) |
79 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[pmIdx]) |
80 | 80 | storres | if columnsWidth != 0: |
81 | 80 | storres | monomialAsString = str(pCoefficients[pmIdx]) + " " \ |
82 | 80 | storres | + str(pMonomials[pmIdx]) |
83 | 80 | storres | print monomialAsString, " " * \ |
84 | 80 | storres | (columnsWidth - len(monomialAsString)), |
85 | 80 | storres | # End for pmIdx loop. |
86 | 80 | storres | # Add the new list row elements to the proto matrix. |
87 | 80 | storres | protoMatrixRows.append(protoMatrixRowCoefficients) |
88 | 80 | storres | if columnsWidth != 0: |
89 | 80 | storres | |
90 | 83 | storres | # End spo_add_polynomial_coeffs_to_matrix_row |
91 | 80 | storres | |
92 | 186 | storres | def spo_float_poly_of_float_to_rat_poly_of_rat(polyOfFloat): |
93 | 186 | storres | """ |
94 | 186 | storres | Create a polynomial over the rationals from a polynomial over |
95 | 186 | storres | a RealField. |
96 | 186 | storres | Important warning: default Sage behavior is to convert coefficients |
97 | 186 | storres | using continued fractions instead of making a simple conversion |
98 | 186 | storres | with a powers of two at denominators (and possible simplification). |
99 | 186 | storres | Hence conversion is not exact but with a relative error around 10^(-521). |
100 | 186 | storres | """ |
101 | 186 | storres | ratPolynomialRing = QQ[str(polyOfFloat.variables()[0])] |
102 | 186 | storres | return(ratPolynomialRing(polyOfFloat)) |
103 | 186 | storres | # End spo_float_poly_of_float_to_rat_poly_of_rat. |
104 | 186 | storres | |
105 | 186 | storres | def spo_float_poly_of_float_to_rat_poly_of_rat_pow_two(polyOfFloat): |
106 | 186 | storres | """ |
107 | 186 | storres | Create a polynomial over the rationals from a polynomial over |
108 | 186 | storres | a RealField where all denominators are |
109 | 186 | storres | powers of two. |
110 | 186 | storres | Allows for exact conversions (and lcm computation of the coefficients |
111 | 186 | storres | denominator). |
112 | 186 | storres | """ |
113 | 186 | storres | polyVariable = polyOfFloat.variables()[0] |
114 | 186 | storres | RPR = QQ[str(polyVariable)] |
115 | 186 | storres | polyCoeffs = polyOfFloat.coefficients() |
116 | 186 | storres | #print polyCoeffs |
117 | 186 | storres | polyExponents = polyOfFloat.exponents() |
118 | 186 | storres | #print polyExponents |
119 | 186 | storres | polyDenomPtwoCoeffs = [] |
120 | 186 | storres | for coeff in polyCoeffs: |
121 | 186 | storres | polyDenomPtwoCoeffs.append(sno_float_to_rat_pow_of_two_denom(coeff)) |
122 | 186 | storres | #print "Converted coefficient:", sno_float_to_rat_pow_of_two_denom(coeff), |
123 | 186 | storres | #print type(sno_float_to_rat_pow_of_two_denom(coeff)) |
124 | 186 | storres | ratPoly = RPR(0) |
125 | 186 | storres | #print type(ratPoly) |
126 | 186 | storres | ## !!! CAUTION !!! Do not use the RPR(coeff * polyVariagle^exponent) |
127 | 186 | storres | # construction. |
128 | 186 | storres | # The coefficient becomes plainly wrong when exponent == 0. |
129 | 186 | storres | # No clue as to why. |
130 | 186 | storres | for coeff, exponent in zip(polyDenomPtwoCoeffs, polyExponents): |
131 | 186 | storres | ratPoly += coeff * RPR(polyVariable^exponent) |
132 | 186 | storres | return ratPoly |
133 | 186 | storres | # End slz_float_poly_of_float_to_rat_poly_of_rat. |
134 | 186 | storres | |
135 | 186 | storres | |
136 | 109 | storres | def spo_get_coefficient_for_monomial(monomialsList, coefficientsList, monomial): |
137 | 109 | storres | """ |
138 | 109 | storres | Get, for a polynomial, the coefficient for a given monomial. |
139 | 109 | storres | The polynomial is given as two lists (monomials and coefficients as |
140 | 109 | storres | return by the respective methods ; indexes of the two lists must match). |
141 | 109 | storres | If the monomial is not found, 0 is returned. |
142 | 109 | storres | """ |
143 | 109 | storres | monomialIndex = 0 |
144 | 109 | storres | for mono in monomialsList: |
145 | 109 | storres | if mono == monomial: |
146 | 109 | storres | return coefficientsList[monomialIndex] |
147 | 109 | storres | monomialIndex += 1 |
148 | 109 | storres | return 0 |
149 | 109 | storres | # End spo_get_coefficient_for_monomial. |
150 | 109 | storres | |
151 | 109 | storres | |
152 | 111 | storres | def spo_expression_as_string(powI, boundI, powT, boundT, powP, powN): |
153 | 80 | storres | """ |
154 | 80 | storres | Computes a string version of the i^k + t^l + p^m + N^n expression for |
155 | 80 | storres | output. |
156 | 80 | storres | """ |
157 | 80 | storres | expressionAsString ="" |
158 | 80 | storres | if powI != 0: |
159 | 111 | storres | expressionAsString += str(iBound^powI) + " i^" + str(powI) |
160 | 80 | storres | if powT != 0: |
161 | 80 | storres | if len(expressionAsString) != 0: |
162 | 80 | storres | expressionAsString += " * " |
163 | 111 | storres | expressionAsString += str(tBound^powT) + " t^" + str(powT) |
164 | 80 | storres | if powP != 0: |
165 | 80 | storres | if len(expressionAsString) != 0: |
166 | 80 | storres | expressionAsString += " * " |
167 | 80 | storres | expressionAsString += "p^" + str(powP) |
168 | 105 | storres | if (powN) != 0 : |
169 | 80 | storres | if len(expressionAsString) != 0: |
170 | 80 | storres | expressionAsString += " * " |
171 | 105 | storres | expressionAsString += "N^" + str(powN) |
172 | 80 | storres | return(expressionAsString) |
173 | 80 | storres | # End spo_expression_as_string. |
174 | 80 | storres | |
175 | 87 | storres | def spo_norm(poly, p=2): |
176 | 81 | storres | """ |
177 | 81 | storres | Behaves more or less (no infinity defined) as the norm for the |
178 | 81 | storres | univariate polynomials. |
179 | 107 | storres | Quoting Sage documentation: |
180 | 107 | storres | "Definition: For integer p, the p-norm of a polynomial is the pth root of |
181 | 81 | storres | the sum of the pth powers of the absolute values of the coefficients of |
182 | 107 | storres | the polynomial." |
183 | 87 | storres | |
184 | 81 | storres | """ |
185 | 87 | storres | # TODO: check the arguments (for p see below).. |
186 | 81 | storres | norm = 0 |
187 | 87 | storres | # For infinity norm. |
188 | 87 | storres | if p == Infinity: |
189 | 87 | storres | for coefficient in poly.coefficients(): |
190 | 87 | storres | coefficientAbs = coefficient.abs() |
191 | 87 | storres | if coefficientAbs > norm: |
192 | 87 | storres | norm = coefficientAbs |
193 | 87 | storres | return norm |
194 | 87 | storres | # TODO: check here the value of p |
195 | 107 | storres | # p must be a positive integer >= 1. |
196 | 107 | storres | if p < 1 or (not p in ZZ): |
197 | 94 | storres | return None |
198 | 87 | storres | # For 1 norm. |
199 | 87 | storres | if p == 1: |
200 | 87 | storres | for coefficient in poly.coefficients(): |
201 | 87 | storres | norm += coefficient.abs() |
202 | 87 | storres | return norm |
203 | 87 | storres | # For other norms |
204 | 81 | storres | for coefficient in poly.coefficients(): |
205 | 103 | storres | norm += coefficient.abs()^p |
206 | 87 | storres | return pow(norm, 1/p) |
207 | 81 | storres | # end spo_norm |
208 | 81 | storres | |
209 | 100 | storres | def spo_polynomial_to_proto_matrix(p, alpha, N, columnsWidth=0): |
210 | 74 | storres | """ |
211 | 83 | storres | From a (bivariate) polynomial and some other parameters build a proto |
212 | 87 | storres | matrix (an array of "rows") to be converted into a "true" matrix and |
213 | 83 | storres | eventually by reduced by fpLLL. |
214 | 102 | storres | The matrix is such as those found in Boneh-Durphee and Stehlé. |
215 | 74 | storres | |
216 | 83 | storres | Parameters |
217 | 83 | storres | ---------- |
218 | 87 | storres | p: the (bivariate) polynomial; |
219 | 87 | storres | pRing: the ring over which p is defined; |
220 | 74 | storres | alpha: |
221 | 74 | storres | N: |
222 | 83 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
223 | 83 | storres | printed in colums of columnsWitdth width. |
224 | 74 | storres | """ |
225 | 100 | storres | pRing = p.parent() |
226 | 77 | storres | knownMonomials = [] |
227 | 77 | storres | protoMatrixRows = [] |
228 | 92 | storres | polynomialsList = [] |
229 | 74 | storres | pVariables = p.variables() |
230 | 123 | storres | #print "In spo...", p, p.variables() |
231 | 74 | storres | iVariable = pVariables[0] |
232 | 76 | storres | tVariable = pVariables[1] |
233 | 87 | storres | polynomialAtPower = pRing(1) |
234 | 87 | storres | currentPolynomial = pRing(1) |
235 | 74 | storres | pIdegree = p.degree(pVariables[0]) |
236 | 74 | storres | pTdegree = p.degree(pVariables[1]) |
237 | 87 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
238 | 105 | storres | nAtAlpha = N^alpha |
239 | 105 | storres | nAtPower = nAtAlpha |
240 | 92 | storres | polExpStr = "" |
241 | 74 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
242 | 74 | storres | for pPower in xrange(0, alpha + 1): |
243 | 76 | storres | # pPower == 0 is a special case. We introduce all the monomials but one |
244 | 78 | storres | # in i and those in t necessary to be able to introduce |
245 | 76 | storres | # p. We arbitrary choose to introduce the highest degree monomial in i |
246 | 76 | storres | # with p. We also introduce all the mixed i^k * t^l monomials with |
247 | 77 | storres | # k < p.degree(i) and l <= p.degree(t). |
248 | 78 | storres | # Mixed terms introduction is necessary here before we start "i shifts" |
249 | 78 | storres | # in the next iteration. |
250 | 74 | storres | if pPower == 0: |
251 | 78 | storres | # Notice that i^pIdegree is excluded as the bound of the xrange is |
252 | 78 | storres | # pIdegree |
253 | 74 | storres | for iPower in xrange(0, pIdegree): |
254 | 74 | storres | for tPower in xrange(0, pTdegree + 1): |
255 | 77 | storres | if columnsWidth != 0: |
256 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
257 | 76 | storres | tPower, |
258 | 76 | storres | pPower, |
259 | 105 | storres | alpha-pPower) |
260 | 92 | storres | print "->", polExpStr |
261 | 74 | storres | currentExpression = iVariable^iPower * \ |
262 | 91 | storres | tVariable^tPower * nAtAlpha |
263 | 78 | storres | # polynomialAtPower == 1 here. Next line should be commented |
264 | 78 | storres | # out but it does not work! Some conversion problem? |
265 | 91 | storres | currentPolynomial = pRing(currentExpression) |
266 | 106 | storres | polynomialsList.append(currentPolynomial) |
267 | 74 | storres | pMonomials = currentPolynomial.monomials() |
268 | 74 | storres | pCoefficients = currentPolynomial.coefficients() |
269 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
270 | 83 | storres | pCoefficients, |
271 | 83 | storres | knownMonomials, |
272 | 83 | storres | protoMatrixRows, |
273 | 83 | storres | columnsWidth) |
274 | 78 | storres | # End tPower. |
275 | 78 | storres | # End for iPower. |
276 | 77 | storres | else: # pPower > 0: (p^1..p^alpha) |
277 | 78 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
278 | 77 | storres | # polynomial. |
279 | 77 | storres | # This step could technically be fused as the first iteration |
280 | 77 | storres | # of the next loop (with iPower starting at 0). |
281 | 77 | storres | # We set it apart for clarity. |
282 | 77 | storres | if columnsWidth != 0: |
283 | 105 | storres | polExpStr = spo_expression_as_string(0, 0, pPower, alpha-pPower) |
284 | 92 | storres | print "->", polExpStr |
285 | 77 | storres | currentPolynomial = polynomialAtPower * nAtPower |
286 | 106 | storres | polynomialsList.append(currentPolynomial) |
287 | 77 | storres | pMonomials = currentPolynomial.monomials() |
288 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
289 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
290 | 83 | storres | pCoefficients, |
291 | 83 | storres | knownMonomials, |
292 | 83 | storres | protoMatrixRows, |
293 | 83 | storres | columnsWidth) |
294 | 77 | storres | |
295 | 77 | storres | # The i^iPower * p^pPower polynomials: they add i^k monomials to |
296 | 77 | storres | # p^pPower up to k < pIdegree * pPower. This only introduces i^k |
297 | 77 | storres | # monomials since mixed terms (that were introduced at a previous |
298 | 77 | storres | # stage) are only shifted to already existing |
299 | 77 | storres | # ones. p^pPower is "shifted" to higher degrees in i as far as |
300 | 77 | storres | # possible, one step short of the degree in i of p^(pPower+1) . |
301 | 77 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
302 | 77 | storres | for iPower in xrange(1, pIdegree): |
303 | 87 | storres | if columnsWidth != 0: |
304 | 92 | storres | polExpStr = spo_expression_as_string(iPower, \ |
305 | 87 | storres | 0, \ |
306 | 87 | storres | pPower, \ |
307 | 87 | storres | alpha) |
308 | 92 | storres | print "->", polExpStr |
309 | 77 | storres | currentExpression = i^iPower * nAtPower |
310 | 87 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
311 | 106 | storres | polynomialsList.append(currentPolynomial) |
312 | 77 | storres | pMonomials = currentPolynomial.monomials() |
313 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
314 | 87 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, \ |
315 | 87 | storres | pCoefficients, \ |
316 | 87 | storres | knownMonomials, \ |
317 | 87 | storres | protoMatrixRows, \ |
318 | 83 | storres | columnsWidth) |
319 | 77 | storres | # End for iPower |
320 | 77 | storres | # We want now to introduce a t * p^pPower polynomial. But before |
321 | 77 | storres | # that we must introduce some mixed monomials. |
322 | 77 | storres | # This loop is no triggered before pPower == 2. |
323 | 78 | storres | # It introduces a first set of high i degree mixed monomials. |
324 | 77 | storres | for iPower in xrange(1, pPower): |
325 | 77 | storres | tPower = pPower - iPower + 1 |
326 | 77 | storres | if columnsWidth != 0: |
327 | 92 | storres | polExpStr = spo_expression_as_string(iPower * pIdegree, |
328 | 77 | storres | tPower, |
329 | 77 | storres | 0, |
330 | 77 | storres | alpha) |
331 | 92 | storres | print "->", polExpStr |
332 | 91 | storres | currentExpression = i^(iPower * pIdegree) * t^tPower * nAtAlpha |
333 | 87 | storres | currentPolynomial = pRing(currentExpression) |
334 | 106 | storres | polynomialsList.append(currentPolynomial) |
335 | 77 | storres | pMonomials = currentPolynomial.monomials() |
336 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
337 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
338 | 83 | storres | pCoefficients, |
339 | 83 | storres | knownMonomials, |
340 | 83 | storres | protoMatrixRows, |
341 | 83 | storres | columnsWidth) |
342 | 77 | storres | # End for iPower |
343 | 78 | storres | # |
344 | 78 | storres | # This is the mixed monomials main loop. It introduces: |
345 | 77 | storres | # - the missing mixed monomials needed before the |
346 | 78 | storres | # t^l * p^pPower * N^(alpha-pPower) polynomial; |
347 | 78 | storres | # - the t^l * p^pPower * N^(alpha-pPower) itself; |
348 | 78 | storres | # - for each of i^k * t^l * p^pPower * N^(alpha-pPower) polynomials: |
349 | 78 | storres | # - the the missing mixed monomials needed polynomials, |
350 | 78 | storres | # - the i^k * t^l * p^pPower * N^(alpha-pPower) itself. |
351 | 78 | storres | # The t^l * p^pPower * N^(alpha-pPower) is introduced when |
352 | 78 | storres | # |
353 | 77 | storres | for iShift in xrange(0, pIdegree): |
354 | 77 | storres | # When pTdegree == 1, the following loop only introduces |
355 | 77 | storres | # a single new monomial. |
356 | 77 | storres | #print "++++++++++" |
357 | 77 | storres | for outerTpower in xrange(1, pTdegree + 1): |
358 | 77 | storres | # First one high i degree mixed monomial. |
359 | 77 | storres | iPower = iShift + pPower * pIdegree |
360 | 77 | storres | if columnsWidth != 0: |
361 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
362 | 77 | storres | outerTpower, |
363 | 77 | storres | 0, |
364 | 77 | storres | alpha) |
365 | 92 | storres | print "->", polExpStr |
366 | 91 | storres | currentExpression = i^iPower * t^outerTpower * nAtAlpha |
367 | 87 | storres | currentPolynomial = pRing(currentExpression) |
368 | 106 | storres | polynomialsList.append(currentPolynomial) |
369 | 77 | storres | pMonomials = currentPolynomial.monomials() |
370 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
371 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
372 | 83 | storres | pCoefficients, |
373 | 83 | storres | knownMonomials, |
374 | 83 | storres | protoMatrixRows, |
375 | 83 | storres | columnsWidth) |
376 | 77 | storres | #print "+++++" |
377 | 78 | storres | # At iShift == 0, the following innerTpower loop adds |
378 | 78 | storres | # duplicate monomials, since no extra i^l * t^k is needed |
379 | 78 | storres | # before introducing the |
380 | 77 | storres | # i^iShift * t^outerPpower * p^pPower * N^(alpha-pPower) |
381 | 77 | storres | # polynomial. |
382 | 77 | storres | # It introduces smaller i degree monomials than the |
383 | 77 | storres | # one(s) added previously (no pPower multiplication). |
384 | 77 | storres | # Here the exponent of t decreases as that of i increases. |
385 | 78 | storres | # This conditional is not entered before pPower == 1. |
386 | 78 | storres | # The innerTpower loop does not produce anything before |
387 | 78 | storres | # pPower == 2. We keep it anyway for other configuration of |
388 | 78 | storres | # p. |
389 | 77 | storres | if iShift > 0: |
390 | 77 | storres | iPower = pIdegree + iShift |
391 | 77 | storres | for innerTpower in xrange(pPower, 1, -1): |
392 | 77 | storres | if columnsWidth != 0: |
393 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
394 | 77 | storres | innerTpower, |
395 | 77 | storres | 0, |
396 | 77 | storres | alpha) |
397 | 77 | storres | currentExpression = \ |
398 | 91 | storres | i^(iPower) * t^(innerTpower) * nAtAlpha |
399 | 87 | storres | currentPolynomial = pRing(currentExpression) |
400 | 106 | storres | polynomialsList.append(currentPolynomial) |
401 | 77 | storres | pMonomials = currentPolynomial.monomials() |
402 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
403 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
404 | 77 | storres | pCoefficients, |
405 | 77 | storres | knownMonomials, |
406 | 77 | storres | protoMatrixRows, |
407 | 77 | storres | columnsWidth) |
408 | 77 | storres | iPower += pIdegree |
409 | 77 | storres | # End for innerTpower |
410 | 77 | storres | # End of if iShift > 0 |
411 | 78 | storres | # When iShift == 0, just after each of the |
412 | 78 | storres | # p^pPower * N^(alpha-pPower) polynomials has |
413 | 78 | storres | # been introduced (followed by a string of |
414 | 78 | storres | # i^k * p^pPower * N^(alpha-pPower) polynomials) a |
415 | 78 | storres | # t^l * p^pPower * N^(alpha-pPower) is introduced here. |
416 | 78 | storres | # |
417 | 77 | storres | # Eventually, the following section introduces the |
418 | 105 | storres | # i^iShift * t^outerTpower * p^iPower * N^(alpha-pPower) |
419 | 77 | storres | # polynomials. |
420 | 77 | storres | if columnsWidth != 0: |
421 | 92 | storres | polExpStr = spo_expression_as_string(iShift, |
422 | 77 | storres | outerTpower, |
423 | 77 | storres | pPower, |
424 | 105 | storres | alpha-pPower) |
425 | 92 | storres | print "->", polExpStr |
426 | 77 | storres | currentExpression = i^iShift * t^outerTpower * nAtPower |
427 | 105 | storres | currentPolynomial = pRing(currentExpression) * \ |
428 | 105 | storres | polynomialAtPower |
429 | 106 | storres | polynomialsList.append(currentPolynomial) |
430 | 77 | storres | pMonomials = currentPolynomial.monomials() |
431 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
432 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
433 | 83 | storres | pCoefficients, |
434 | 83 | storres | knownMonomials, |
435 | 83 | storres | protoMatrixRows, |
436 | 83 | storres | columnsWidth) |
437 | 77 | storres | # End for outerTpower |
438 | 77 | storres | #print "++++++++++" |
439 | 77 | storres | # End for iShift |
440 | 77 | storres | polynomialAtPower *= p |
441 | 77 | storres | nAtPower /= N |
442 | 77 | storres | # End for pPower loop |
443 | 92 | storres | return ((protoMatrixRows, knownMonomials, polynomialsList)) |
444 | 83 | storres | # End spo_polynomial_to_proto_matrix |
445 | 81 | storres | |
446 | 111 | storres | def spo_polynomial_to_polynomials_list_2(p, alpha, N, iBound, tBound, |
447 | 111 | storres | columnsWidth=0): |
448 | 105 | storres | """ |
449 | 112 | storres | Badly out of sync code: check with versions 3 or 4. |
450 | 112 | storres | |
451 | 106 | storres | From p, alpha, N build a list of polynomials... |
452 | 106 | storres | TODO: clean up the comments below! |
453 | 106 | storres | |
454 | 105 | storres | From a (bivariate) polynomial and some other parameters build a proto |
455 | 105 | storres | matrix (an array of "rows") to be converted into a "true" matrix and |
456 | 105 | storres | eventually by reduced by fpLLL. |
457 | 105 | storres | The matrix is based on a list of polynomials that are built in a way |
458 | 105 | storres | that one and only monomial is added at each new polynomial. Among the many |
459 | 105 | storres | possible ways to build this list we pick one strongly dependent on the |
460 | 105 | storres | structure of the polynomial and of the problem. |
461 | 105 | storres | We consider here the polynomials of the form: |
462 | 105 | storres | a_k*i^k + a_(k-1)*i^(k-1) + ... + a_1*i + a_0 - t |
463 | 105 | storres | The values of i and t are bounded and we eventually look for (i_0,t_0) |
464 | 105 | storres | pairs such that: |
465 | 105 | storres | a_k*i_0^k + a_(k-1)*i_0^(k-1) + ... + a_1*i_0 + a_0 = t_0 |
466 | 105 | storres | Hence, departing from the procedure in described in Boneh-Durfee, we will |
467 | 105 | storres | not use "t-shifts" but only "i-shifts". |
468 | 105 | storres | |
469 | 105 | storres | Parameters |
470 | 105 | storres | ---------- |
471 | 105 | storres | p: the (bivariate) polynomial; |
472 | 105 | storres | pRing: the ring over which p is defined; |
473 | 105 | storres | alpha: |
474 | 105 | storres | N: |
475 | 105 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
476 | 105 | storres | printed in colums of columnsWitdth width. |
477 | 105 | storres | """ |
478 | 105 | storres | pRing = p.parent() |
479 | 105 | storres | polynomialsList = [] |
480 | 105 | storres | pVariables = p.variables() |
481 | 105 | storres | iVariable = pVariables[0] |
482 | 105 | storres | tVariable = pVariables[1] |
483 | 105 | storres | polynomialAtPower = pRing(1) |
484 | 105 | storres | currentPolynomial = pRing(1) |
485 | 105 | storres | pIdegree = p.degree(iVariable) |
486 | 105 | storres | pTdegree = p.degree(tVariable) |
487 | 105 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
488 | 105 | storres | nAtAlpha = N^alpha |
489 | 105 | storres | nAtPower = nAtAlpha |
490 | 105 | storres | polExpStr = "" |
491 | 105 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
492 | 105 | storres | for pPower in xrange(0, alpha + 1): |
493 | 105 | storres | # pPower == 0 is a special case. We introduce all the monomials in i |
494 | 105 | storres | # up to i^pIdegree. |
495 | 105 | storres | if pPower == 0: |
496 | 105 | storres | # Notice who iPower runs up to i^pIdegree. |
497 | 105 | storres | for iPower in xrange(0, pIdegree + 1): |
498 | 105 | storres | # No t power is taken into account as we limit our selves to |
499 | 105 | storres | # degree 1 in t and make no "t-shifts". |
500 | 105 | storres | if columnsWidth != 0: |
501 | 111 | storres | polExpStr = spo_expression_as_string(iPower, |
502 | 111 | storres | iBound, |
503 | 105 | storres | 0, |
504 | 111 | storres | tBound, |
505 | 105 | storres | 0, |
506 | 105 | storres | alpha) |
507 | 105 | storres | print "->", polExpStr |
508 | 111 | storres | currentExpression = iVariable^iPower * nAtAlpha * iBound^iPower |
509 | 105 | storres | # polynomialAtPower == 1 here. Next line should be commented |
510 | 105 | storres | # out but it does not work! Some conversion problem? |
511 | 105 | storres | currentPolynomial = pRing(currentExpression) |
512 | 105 | storres | polynomialsList.append(currentPolynomial) |
513 | 105 | storres | # End for iPower. |
514 | 105 | storres | else: # pPower > 0: (p^1..p^alpha) |
515 | 105 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
516 | 105 | storres | # polynomial. This is also where the t^pPower monomials shows up for |
517 | 105 | storres | # the first time. |
518 | 105 | storres | if columnsWidth != 0: |
519 | 111 | storres | polExpStr = spo_expression_as_string(0, iBound, 0, tBound, \ |
520 | 111 | storres | pPower, alpha-pPower) |
521 | 105 | storres | print "->", polExpStr |
522 | 105 | storres | currentPolynomial = polynomialAtPower * nAtPower |
523 | 105 | storres | polynomialsList.append(currentPolynomial) |
524 | 106 | storres | # Exit when pPower == alpha |
525 | 106 | storres | if pPower == alpha: |
526 | 110 | storres | return polynomialsList |
527 | 105 | storres | # This is where the "i-shifts" take place. Mixed terms, i^k * t^l |
528 | 105 | storres | # (that were introduced at a previous |
529 | 105 | storres | # stage or are introduced now) are only shifted to already existing |
530 | 105 | storres | # ones with the notable exception of i^iPower * t^pPower, which |
531 | 105 | storres | # must be manually introduced. |
532 | 105 | storres | # p^pPower is "shifted" to higher degrees in i as far as |
533 | 105 | storres | # possible, up to of the degree in i of p^(pPower+1). |
534 | 105 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
535 | 105 | storres | for iPower in xrange(1, pIdegree + 1): |
536 | 105 | storres | # The i^iPower * t^pPower monomial. Notice the alpha exponent |
537 | 105 | storres | # for N. |
538 | 105 | storres | internalIpower = iPower |
539 | 105 | storres | for tPower in xrange(pPower,0,-1): |
540 | 105 | storres | if columnsWidth != 0: |
541 | 111 | storres | polExpStr = spo_expression_as_string(internalIpower, |
542 | 111 | storres | iBound, |
543 | 111 | storres | tPower, |
544 | 111 | storres | tBound, |
545 | 111 | storres | 0, |
546 | 105 | storres | alpha) |
547 | 105 | storres | print "->", polExpStr |
548 | 111 | storres | currentExpression = i^internalIpower * t^tPower * \ |
549 | 111 | storres | nAtAlpha * iBound^internalIpower * \ |
550 | 111 | storres | tBound^tPower |
551 | 111 | storres | |
552 | 105 | storres | currentPolynomial = pRing(currentExpression) |
553 | 105 | storres | polynomialsList.append(currentPolynomial) |
554 | 105 | storres | internalIpower += pIdegree |
555 | 105 | storres | # End for tPower |
556 | 105 | storres | # The i^iPower * p^pPower * N^(alpha-pPower) i-shift. |
557 | 105 | storres | if columnsWidth != 0: |
558 | 111 | storres | polExpStr = spo_expression_as_string(iPower, |
559 | 111 | storres | iBound, |
560 | 111 | storres | 0, |
561 | 111 | storres | tBound, |
562 | 111 | storres | pPower, |
563 | 105 | storres | alpha-pPower) |
564 | 105 | storres | print "->", polExpStr |
565 | 111 | storres | currentExpression = i^iPower * nAtPower * iBound^iPower |
566 | 105 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
567 | 105 | storres | polynomialsList.append(currentPolynomial) |
568 | 105 | storres | # End for iPower |
569 | 105 | storres | polynomialAtPower *= p |
570 | 105 | storres | nAtPower /= N |
571 | 105 | storres | # End for pPower loop |
572 | 109 | storres | return polynomialsList |
573 | 105 | storres | # End spo_polynomial_to_proto_matrix_2 |
574 | 105 | storres | |
575 | 111 | storres | def spo_polynomial_to_polynomials_list_3(p, alpha, N, iBound, tBound, |
576 | 109 | storres | columnsWidth=0): |
577 | 108 | storres | """ |
578 | 108 | storres | From p, alpha, N build a list of polynomials... |
579 | 108 | storres | TODO: more in depth rationale... |
580 | 108 | storres | |
581 | 108 | storres | Our goal is to introduce each monomial with the smallest coefficient. |
582 | 108 | storres | |
583 | 108 | storres | |
584 | 108 | storres | |
585 | 108 | storres | Parameters |
586 | 108 | storres | ---------- |
587 | 108 | storres | p: the (bivariate) polynomial; |
588 | 108 | storres | pRing: the ring over which p is defined; |
589 | 108 | storres | alpha: |
590 | 108 | storres | N: |
591 | 108 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
592 | 108 | storres | printed in colums of columnsWitdth width. |
593 | 108 | storres | """ |
594 | 108 | storres | pRing = p.parent() |
595 | 108 | storres | polynomialsList = [] |
596 | 108 | storres | pVariables = p.variables() |
597 | 108 | storres | iVariable = pVariables[0] |
598 | 108 | storres | tVariable = pVariables[1] |
599 | 108 | storres | polynomialAtPower = pRing(1) |
600 | 108 | storres | currentPolynomial = pRing(1) |
601 | 108 | storres | pIdegree = p.degree(iVariable) |
602 | 108 | storres | pTdegree = p.degree(tVariable) |
603 | 108 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
604 | 108 | storres | nAtAlpha = N^alpha |
605 | 108 | storres | nAtPower = nAtAlpha |
606 | 108 | storres | polExpStr = "" |
607 | 108 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
608 | 108 | storres | for pPower in xrange(0, alpha + 1): |
609 | 108 | storres | # pPower == 0 is a special case. We introduce all the monomials in i |
610 | 108 | storres | # up to i^pIdegree. |
611 | 108 | storres | if pPower == 0: |
612 | 108 | storres | # Notice who iPower runs up to i^pIdegree. |
613 | 108 | storres | for iPower in xrange(0, pIdegree + 1): |
614 | 108 | storres | # No t power is taken into account as we limit our selves to |
615 | 108 | storres | # degree 1 in t and make no "t-shifts". |
616 | 108 | storres | if columnsWidth != 0: |
617 | 108 | storres | polExpStr = spo_expression_as_string(iPower, |
618 | 111 | storres | iBound, |
619 | 108 | storres | 0, |
620 | 111 | storres | tBound, |
621 | 108 | storres | 0, |
622 | 108 | storres | alpha) |
623 | 108 | storres | print "->", polExpStr |
624 | 111 | storres | currentExpression = iVariable^iPower * nAtAlpha * iBound^iPower |
625 | 108 | storres | # polynomialAtPower == 1 here. Next line should be commented |
626 | 108 | storres | # out but it does not work! Some conversion problem? |
627 | 108 | storres | currentPolynomial = pRing(currentExpression) |
628 | 108 | storres | polynomialsList.append(currentPolynomial) |
629 | 108 | storres | # End for iPower. |
630 | 108 | storres | else: # pPower > 0: (p^1..p^alpha) |
631 | 108 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
632 | 108 | storres | # polynomial. This is also where the t^pPower monomials shows up for |
633 | 108 | storres | # the first time. It app |
634 | 108 | storres | if columnsWidth != 0: |
635 | 111 | storres | polExpStr = spo_expression_as_string(0, iBound, |
636 | 111 | storres | 0, tBound, |
637 | 111 | storres | pPower, alpha-pPower) |
638 | 108 | storres | print "->", polExpStr |
639 | 108 | storres | currentPolynomial = polynomialAtPower * nAtPower |
640 | 108 | storres | polynomialsList.append(currentPolynomial) |
641 | 108 | storres | # Exit when pPower == alpha |
642 | 108 | storres | if pPower == alpha: |
643 | 111 | storres | return polynomialsList |
644 | 108 | storres | # This is where the "i-shifts" take place. Mixed terms, i^k * t^l |
645 | 108 | storres | # (that were introduced at a previous |
646 | 108 | storres | # stage or are introduced now) are only shifted to already existing |
647 | 108 | storres | # ones with the notable exception of i^iPower * t^pPower, which |
648 | 108 | storres | # must be manually introduced. |
649 | 108 | storres | # p^pPower is "shifted" to higher degrees in i as far as |
650 | 108 | storres | # possible, up to of the degree in i of p^(pPower+1). |
651 | 108 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
652 | 108 | storres | for iPower in xrange(1, pIdegree + 1): |
653 | 108 | storres | # The i^iPower * t^pPower monomial. Notice the alpha exponent |
654 | 108 | storres | # for N. |
655 | 108 | storres | internalIpower = iPower |
656 | 108 | storres | for tPower in xrange(pPower,0,-1): |
657 | 108 | storres | if columnsWidth != 0: |
658 | 111 | storres | polExpStr = spo_expression_as_string(internalIpower, |
659 | 111 | storres | iBound, |
660 | 111 | storres | tPower, |
661 | 111 | storres | tBound, |
662 | 111 | storres | 0, |
663 | 108 | storres | alpha) |
664 | 108 | storres | print "->", polExpStr |
665 | 111 | storres | currentExpression = i^internalIpower * t^tPower * nAtAlpha * \ |
666 | 111 | storres | iBound^internalIpower * tBound^tPower |
667 | 108 | storres | currentPolynomial = pRing(currentExpression) |
668 | 108 | storres | polynomialsList.append(currentPolynomial) |
669 | 108 | storres | internalIpower += pIdegree |
670 | 108 | storres | # End for tPower |
671 | 109 | storres | # Here we have to choose between a |
672 | 109 | storres | # i^iPower * p^pPower * N^(alpha-pPower) i-shift and |
673 | 111 | storres | # i^iPower * i^(d_i(p) * pPower) * N^alpha, depending on which |
674 | 109 | storres | # coefficient is smallest. |
675 | 109 | storres | IcurrentExponent = iPower + \ |
676 | 111 | storres | (pPower * polynomialAtPower.degree(iVariable)) |
677 | 111 | storres | currentMonomial = pRing(iVariable^IcurrentExponent) |
678 | 111 | storres | currentPolynomial = pRing(iVariable^iPower * nAtPower * \ |
679 | 111 | storres | iBound^iPower) * \ |
680 | 111 | storres | polynomialAtPower |
681 | 109 | storres | currMonomials = currentPolynomial.monomials() |
682 | 109 | storres | currCoefficients = currentPolynomial.coefficients() |
683 | 109 | storres | currentCoefficient = spo_get_coefficient_for_monomial( \ |
684 | 109 | storres | currMonomials, |
685 | 109 | storres | currCoefficients, |
686 | 109 | storres | currentMonomial) |
687 | 111 | storres | print "Current coefficient:", currentCoefficient |
688 | 111 | storres | alterCoefficient = iBound^IcurrentExponent * nAtAlpha |
689 | 111 | storres | print "N^alpha * ibound^", IcurrentExponent, ":", \ |
690 | 111 | storres | alterCoefficient |
691 | 111 | storres | if currentCoefficient > alterCoefficient : |
692 | 109 | storres | if columnsWidth != 0: |
693 | 111 | storres | polExpStr = spo_expression_as_string(IcurrentExponent, |
694 | 111 | storres | iBound, |
695 | 111 | storres | 0, |
696 | 111 | storres | tBound, |
697 | 111 | storres | 0, |
698 | 109 | storres | alpha) |
699 | 111 | storres | print "->", polExpStr |
700 | 111 | storres | polynomialsList.append(currentMonomial * \ |
701 | 111 | storres | alterCoefficient) |
702 | 109 | storres | else: |
703 | 109 | storres | if columnsWidth != 0: |
704 | 111 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
705 | 111 | storres | 0, tBound, |
706 | 111 | storres | pPower, |
707 | 109 | storres | alpha-pPower) |
708 | 111 | storres | print "->", polExpStr |
709 | 109 | storres | polynomialsList.append(currentPolynomial) |
710 | 108 | storres | # End for iPower |
711 | 108 | storres | polynomialAtPower *= p |
712 | 108 | storres | nAtPower /= N |
713 | 108 | storres | # End for pPower loop |
714 | 109 | storres | return polynomialsList |
715 | 108 | storres | # End spo_polynomial_to_proto_matrix_3 |
716 | 108 | storres | |
717 | 111 | storres | def spo_polynomial_to_polynomials_list_4(p, alpha, N, iBound, tBound, |
718 | 111 | storres | columnsWidth=0): |
719 | 83 | storres | """ |
720 | 111 | storres | From p, alpha, N build a list of polynomials... |
721 | 111 | storres | TODO: more in depth rationale... |
722 | 83 | storres | |
723 | 111 | storres | Our goal is to introduce each monomial with the smallest coefficient. |
724 | 111 | storres | |
725 | 111 | storres | |
726 | 111 | storres | |
727 | 83 | storres | Parameters |
728 | 83 | storres | ---------- |
729 | 111 | storres | p: the (bivariate) polynomial; |
730 | 111 | storres | pRing: the ring over which p is defined; |
731 | 111 | storres | alpha: |
732 | 111 | storres | N: |
733 | 111 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
734 | 111 | storres | printed in colums of columnsWitdth width. |
735 | 111 | storres | """ |
736 | 111 | storres | pRing = p.parent() |
737 | 111 | storres | polynomialsList = [] |
738 | 111 | storres | pVariables = p.variables() |
739 | 111 | storres | iVariable = pVariables[0] |
740 | 111 | storres | tVariable = pVariables[1] |
741 | 111 | storres | polynomialAtPower = copy(p) |
742 | 111 | storres | currentPolynomial = pRing(1) |
743 | 111 | storres | pIdegree = p.degree(iVariable) |
744 | 111 | storres | pTdegree = p.degree(tVariable) |
745 | 111 | storres | maxIdegree = pIdegree * alpha |
746 | 111 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
747 | 111 | storres | nAtAlpha = N^alpha |
748 | 111 | storres | nAtPower = nAtAlpha |
749 | 111 | storres | polExpStr = "" |
750 | 111 | storres | # We first introduce all the monomials in i alone multiplied by N^alpha. |
751 | 111 | storres | for iPower in xrange(0, maxIdegree + 1): |
752 | 111 | storres | if columnsWidth !=0: |
753 | 111 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
754 | 111 | storres | 0, tBound, |
755 | 111 | storres | 0, alpha) |
756 | 111 | storres | print "->", polExpStr |
757 | 111 | storres | currentExpression = iVariable^iPower * nAtAlpha * iBound^iPower |
758 | 111 | storres | currentPolynomial = pRing(currentExpression) |
759 | 111 | storres | polynomialsList.append(currentPolynomial) |
760 | 111 | storres | # End for iPower |
761 | 111 | storres | # We work from p^1 * N^alpha-1 to p^alpha * N^0 |
762 | 111 | storres | for pPower in xrange(1, alpha + 1): |
763 | 111 | storres | # First of all the p^pPower * N^(alpha-pPower) polynomial. |
764 | 111 | storres | nAtPower /= N |
765 | 111 | storres | if columnsWidth !=0: |
766 | 111 | storres | polExpStr = spo_expression_as_string(0, iBound, |
767 | 111 | storres | 0, tBound, |
768 | 111 | storres | pPower, alpha-pPower) |
769 | 111 | storres | print "->", polExpStr |
770 | 111 | storres | currentPolynomial = polynomialAtPower * nAtPower |
771 | 111 | storres | polynomialsList.append(currentPolynomial) |
772 | 111 | storres | # Exit when pPower == alpha |
773 | 111 | storres | if pPower == alpha: |
774 | 111 | storres | return polynomialsList |
775 | 111 | storres | # We now introduce the mixed i^k * t^l monomials by i^m * p^n * N^(alpha-n) |
776 | 111 | storres | for iPower in xrange(1, pIdegree + 1): |
777 | 111 | storres | if columnsWidth != 0: |
778 | 111 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
779 | 111 | storres | 0, tBound, |
780 | 111 | storres | pPower, alpha-pPower) |
781 | 111 | storres | print "->", polExpStr |
782 | 111 | storres | currentExpression = i^iPower * iBound^iPower * nAtPower |
783 | 111 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
784 | 111 | storres | polynomialsList.append(currentPolynomial) |
785 | 111 | storres | # End for iPower |
786 | 111 | storres | polynomialAtPower *= p |
787 | 111 | storres | # End for pPower loop |
788 | 111 | storres | return polynomialsList |
789 | 111 | storres | # End spo_polynomial_to_proto_matrix_4 |
790 | 111 | storres | |
791 | 113 | storres | def spo_polynomial_to_polynomials_list_5(p, alpha, N, iBound, tBound, |
792 | 113 | storres | columnsWidth=0): |
793 | 113 | storres | """ |
794 | 113 | storres | From p, alpha, N build a list of polynomials use to create a base |
795 | 113 | storres | that will eventually be reduced with LLL. |
796 | 113 | storres | |
797 | 113 | storres | The bounds are computed for the coefficients that will be used to |
798 | 113 | storres | form the base. |
799 | 113 | storres | |
800 | 113 | storres | We try to introduce only one new monomial at a time, to obtain a |
801 | 113 | storres | triangular matrix (it is easy to compute the volume of the underlining |
802 | 113 | storres | latice if the matrix is triangular). |
803 | 113 | storres | |
804 | 113 | storres | There are many possibilities to introduce the monomials: our goal is also |
805 | 113 | storres | to introduce each of them on the diagonal with the smallest coefficient. |
806 | 113 | storres | |
807 | 113 | storres | The method depends on the structure of the polynomial. Here it is adapted |
808 | 113 | storres | to the a_n*i^n + ... + a_1 * i - t + b form. |
809 | 113 | storres | |
810 | 113 | storres | Parameters |
811 | 113 | storres | ---------- |
812 | 113 | storres | p: the (bivariate) polynomial; |
813 | 113 | storres | alpha: |
814 | 113 | storres | N: |
815 | 113 | storres | iBound: |
816 | 113 | storres | tBound: |
817 | 113 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
818 | 113 | storres | printed in colums of columnsWitdth width. |
819 | 113 | storres | """ |
820 | 113 | storres | pRing = p.parent() |
821 | 113 | storres | polynomialsList = [] |
822 | 113 | storres | pVariables = p.variables() |
823 | 113 | storres | iVariable = pVariables[0] |
824 | 113 | storres | tVariable = pVariables[1] |
825 | 113 | storres | polynomialAtPower = copy(p) |
826 | 113 | storres | currentPolynomial = pRing(1) |
827 | 113 | storres | pIdegree = p.degree(iVariable) |
828 | 113 | storres | pTdegree = p.degree(tVariable) |
829 | 113 | storres | maxIdegree = pIdegree * alpha |
830 | 113 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
831 | 113 | storres | nAtAlpha = N^alpha |
832 | 113 | storres | nAtPower = nAtAlpha |
833 | 113 | storres | polExpStr = "" |
834 | 113 | storres | # We first introduce all the monomials in i alone multiplied by N^alpha. |
835 | 113 | storres | for iPower in xrange(0, maxIdegree + 1): |
836 | 113 | storres | if columnsWidth !=0: |
837 | 113 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
838 | 113 | storres | 0, tBound, |
839 | 113 | storres | 0, alpha) |
840 | 113 | storres | print "->", polExpStr |
841 | 113 | storres | currentExpression = iVariable^iPower * nAtAlpha * iBound^iPower |
842 | 113 | storres | currentPolynomial = pRing(currentExpression) |
843 | 113 | storres | polynomialsList.append(currentPolynomial) |
844 | 113 | storres | # End for iPower |
845 | 113 | storres | # We work from p^1 * N^alpha-1 to p^alpha * N^0 |
846 | 113 | storres | for pPower in xrange(1, alpha + 1): |
847 | 113 | storres | # First of all the p^pPower * N^(alpha-pPower) polynomial. |
848 | 113 | storres | nAtPower /= N |
849 | 113 | storres | if columnsWidth !=0: |
850 | 113 | storres | polExpStr = spo_expression_as_string(0, iBound, |
851 | 113 | storres | 0, tBound, |
852 | 113 | storres | pPower, alpha-pPower) |
853 | 113 | storres | print "->", polExpStr |
854 | 113 | storres | currentPolynomial = polynomialAtPower * nAtPower |
855 | 113 | storres | polynomialsList.append(currentPolynomial) |
856 | 113 | storres | # Exit when pPower == alpha |
857 | 113 | storres | if pPower == alpha: |
858 | 113 | storres | return polynomialsList |
859 | 113 | storres | for iPower in xrange(1, pIdegree + 1): |
860 | 113 | storres | iCurrentPower = pIdegree + iPower |
861 | 113 | storres | for tPower in xrange(pPower-1, 0, -1): |
862 | 114 | storres | #print "tPower:", tPower |
863 | 113 | storres | if columnsWidth != 0: |
864 | 113 | storres | polExpStr = spo_expression_as_string(iCurrentPower, iBound, |
865 | 113 | storres | tPower, tBound, |
866 | 113 | storres | 0, alpha) |
867 | 113 | storres | print "->", polExpStr |
868 | 113 | storres | currentExpression = i^iCurrentPower * iBound^iCurrentPower * t^tPower * tBound^tPower *nAtAlpha |
869 | 113 | storres | currentPolynomial = pRing(currentExpression) |
870 | 113 | storres | polynomialsList.append(currentPolynomial) |
871 | 113 | storres | iCurrentPower += pIdegree |
872 | 113 | storres | # End for tPower |
873 | 113 | storres | # We now introduce the mixed i^k * t^l monomials by i^m * p^n * N^(alpha-n) |
874 | 113 | storres | if columnsWidth != 0: |
875 | 113 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
876 | 113 | storres | 0, tBound, |
877 | 113 | storres | pPower, alpha-pPower) |
878 | 113 | storres | print "->", polExpStr |
879 | 113 | storres | currentExpression = i^iPower * iBound^iPower * nAtPower |
880 | 113 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
881 | 113 | storres | polynomialsList.append(currentPolynomial) |
882 | 113 | storres | # End for iPower |
883 | 113 | storres | polynomialAtPower *= p |
884 | 113 | storres | # End for pPower loop |
885 | 113 | storres | return polynomialsList |
886 | 113 | storres | # End spo_polynomial_to_proto_matrix_5 |
887 | 113 | storres | |
888 | 155 | storres | def spo_polynomial_to_polynomials_list_6(p, alpha, N, iBound, tBound, |
889 | 155 | storres | columnsWidth=0): |
890 | 155 | storres | """ |
891 | 155 | storres | From p, alpha, N build a list of polynomials use to create a base |
892 | 155 | storres | that will eventually be reduced with LLL. |
893 | 155 | storres | |
894 | 155 | storres | The bounds are computed for the coefficients that will be used to |
895 | 155 | storres | form the base. |
896 | 155 | storres | |
897 | 155 | storres | We try to introduce only one new monomial at a time, whithout trying to |
898 | 155 | storres | obtain a triangular matrix. |
899 | 155 | storres | |
900 | 155 | storres | There are many possibilities to introduce the monomials: our goal is also |
901 | 155 | storres | to introduce each of them on the diagonal with the smallest coefficient. |
902 | 155 | storres | |
903 | 155 | storres | The method depends on the structure of the polynomial. Here it is adapted |
904 | 155 | storres | to the a_n*i^n + ... + a_1 * i - t + b form. |
905 | 155 | storres | |
906 | 155 | storres | Parameters |
907 | 155 | storres | ---------- |
908 | 155 | storres | p: the (bivariate) polynomial; |
909 | 155 | storres | alpha: |
910 | 155 | storres | N: |
911 | 155 | storres | iBound: |
912 | 155 | storres | tBound: |
913 | 155 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
914 | 155 | storres | printed in colums of columnsWitdth width. |
915 | 155 | storres | """ |
916 | 155 | storres | pRing = p.parent() |
917 | 155 | storres | polynomialsList = [] |
918 | 155 | storres | pVariables = p.variables() |
919 | 155 | storres | iVariable = pVariables[0] |
920 | 155 | storres | tVariable = pVariables[1] |
921 | 155 | storres | polynomialAtPower = copy(p) |
922 | 155 | storres | currentPolynomial = pRing(1) # Constant term. |
923 | 155 | storres | pIdegree = p.degree(iVariable) |
924 | 155 | storres | pTdegree = p.degree(tVariable) |
925 | 155 | storres | maxIdegree = pIdegree * alpha |
926 | 155 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
927 | 155 | storres | nAtAlpha = N^alpha |
928 | 155 | storres | nAtPower = nAtAlpha |
929 | 155 | storres | polExpStr = "" |
930 | 155 | storres | # |
931 | 157 | storres | """ |
932 | 157 | storres | ## Bound for iPower + pIdegree*tPower <= alpha*pIdegree |
933 | 157 | storres | print "degree in i:", pIdegree |
934 | 157 | storres | powersRangeUpperBound = alpha * pIdegree + 1 # +1 for the range. |
935 | 157 | storres | for iPower in xrange(0, powersRangeUpperBound): |
936 | 157 | storres | tPower = 0 |
937 | 157 | storres | while (iPower + tPower * pIdegree) < powersRangeUpperBound: |
938 | 155 | storres | print "iPower:", iPower, " tPower:", tPower |
939 | 155 | storres | q = pRing(iVariable * iBound)^iPower * ((p * N)^tPower) |
940 | 157 | storres | print "q monomials:", q.monomials() |
941 | 155 | storres | polynomialsList.append(q) |
942 | 157 | storres | tPower += 1 |
943 | 157 | storres | """ |
944 | 168 | storres | """ |
945 | 168 | storres | Start from iExp = 0 since starting from 1 does not allow for |
946 | 168 | storres | resultants != 0. |
947 | 168 | storres | """ |
948 | 157 | storres | for iExp in xrange(0, alpha+1): |
949 | 157 | storres | tExp = 0 |
950 | 157 | storres | while iExp + tExp <= alpha: |
951 | 157 | storres | q = pRing(iVariable * iBound)^iExp * ((p * N)^tExp) |
952 | 168 | storres | sys.stdout.write("q " + str(iExp) + "," + str(tExp) + ": ") |
953 | 168 | storres | print q |
954 | 157 | storres | polynomialsList.append(q) |
955 | 157 | storres | tExp += 1 |
956 | 155 | storres | return polynomialsList |
957 | 155 | storres | |
958 | 155 | storres | """ |
959 | 155 | storres | # We first introduce all the monomials in i alone multiplied by N^alpha. |
960 | 155 | storres | for iPower in xrange(0, maxIdegree + 1): |
961 | 155 | storres | if columnsWidth !=0: |
962 | 155 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
963 | 155 | storres | 0, tBound, |
964 | 155 | storres | 0, alpha) |
965 | 155 | storres | print "->", polExpStr |
966 | 155 | storres | currentExpression = iVariable^iPower * nAtAlpha * iBound^iPower |
967 | 155 | storres | currentPolynomial = pRing(currentExpression) |
968 | 155 | storres | polynomialsList.append(currentPolynomial) |
969 | 155 | storres | # End for iPower |
970 | 155 | storres | # We work from p^1 * N^alpha-1 to p^alpha * N^0 |
971 | 155 | storres | for pPower in xrange(1, alpha + 1): |
972 | 155 | storres | # First of all the p^pPower * N^(alpha-pPower) polynomial. |
973 | 155 | storres | nAtPower /= N |
974 | 155 | storres | if columnsWidth !=0: |
975 | 155 | storres | polExpStr = spo_expression_as_string(0, iBound, |
976 | 155 | storres | 0, tBound, |
977 | 155 | storres | pPower, alpha-pPower) |
978 | 155 | storres | print "->", polExpStr |
979 | 155 | storres | currentPolynomial = polynomialAtPower * nAtPower |
980 | 155 | storres | polynomialsList.append(currentPolynomial) |
981 | 155 | storres | # Exit when pPower == alpha |
982 | 155 | storres | if pPower == alpha: |
983 | 155 | storres | return polynomialsList |
984 | 155 | storres | for iPower in xrange(1, pIdegree + 1): |
985 | 155 | storres | iCurrentPower = pIdegree + iPower |
986 | 155 | storres | for tPower in xrange(pPower-1, 0, -1): |
987 | 155 | storres | #print "tPower:", tPower |
988 | 155 | storres | if columnsWidth != 0: |
989 | 155 | storres | polExpStr = spo_expression_as_string(iCurrentPower, iBound, |
990 | 155 | storres | tPower, tBound, |
991 | 155 | storres | 0, alpha) |
992 | 155 | storres | print "->", polExpStr |
993 | 155 | storres | currentExpression = i^iCurrentPower * iBound^iCurrentPower * t^tPower * tBound^tPower *nAtAlpha |
994 | 155 | storres | currentPolynomial = pRing(currentExpression) |
995 | 155 | storres | polynomialsList.append(currentPolynomial) |
996 | 155 | storres | iCurrentPower += pIdegree |
997 | 155 | storres | # End for tPower |
998 | 155 | storres | # We now introduce the mixed i^k * t^l monomials by i^m * p^n * N^(alpha-n) |
999 | 155 | storres | if columnsWidth != 0: |
1000 | 155 | storres | polExpStr = spo_expression_as_string(iPower, iBound, |
1001 | 155 | storres | 0, tBound, |
1002 | 155 | storres | pPower, alpha-pPower) |
1003 | 155 | storres | print "->", polExpStr |
1004 | 155 | storres | currentExpression = i^iPower * iBound^iPower * nAtPower |
1005 | 155 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
1006 | 155 | storres | polynomialsList.append(currentPolynomial) |
1007 | 155 | storres | # End for iPower |
1008 | 155 | storres | polynomialAtPower *= p |
1009 | 155 | storres | # End for pPower loop |
1010 | 155 | storres | """ |
1011 | 155 | storres | return polynomialsList |
1012 | 155 | storres | # End spo_polynomial_to_proto_matrix_6 |
1013 | 155 | storres | |
1014 | 168 | storres | def spo_polynomial_to_polynomials_list_7(p, alpha, N, iBound, tBound, |
1015 | 168 | storres | columnsWidth=0): |
1016 | 168 | storres | """ |
1017 | 171 | storres | As per Random Bits... direct loops nesting. |
1018 | 168 | storres | """ |
1019 | 168 | storres | pRing = p.parent() |
1020 | 168 | storres | polynomialsList = [] |
1021 | 168 | storres | pVariables = p.variables() |
1022 | 168 | storres | iVariable = pVariables[0] |
1023 | 168 | storres | tVariable = pVariables[1] |
1024 | 168 | storres | polynomialAtPower = copy(p) |
1025 | 168 | storres | currentPolynomial = pRing(1) # Constant term. |
1026 | 168 | storres | |
1027 | 168 | storres | for iExp in xrange(0, alpha+1): |
1028 | 168 | storres | pExp = 0 |
1029 | 168 | storres | while (iExp + pExp) <= alpha: |
1030 | 168 | storres | print "iExp:", iExp, \ |
1031 | 168 | storres | "- pExp:", pExp, \ |
1032 | 168 | storres | "- alpha-pExp:", alpha-pExp |
1033 | 168 | storres | q = pRing(iVariable * iBound)^iExp * p^pExp * N^(alpha-pExp) |
1034 | 169 | storres | print q.monomials() |
1035 | 168 | storres | polynomialsList.append(q) |
1036 | 168 | storres | pExp += 1 |
1037 | 168 | storres | return polynomialsList |
1038 | 168 | storres | # End spo_polynomial_to_polynomials_list_7 |
1039 | 168 | storres | |
1040 | 169 | storres | def spo_polynomial_to_polynomials_list_8(p, alpha, N, iBound, tBound, |
1041 | 169 | storres | columnsWidth=0): |
1042 | 169 | storres | """ |
1043 | 171 | storres | As per Random Bits... (reversed loop nesting) |
1044 | 169 | storres | """ |
1045 | 169 | storres | pRing = p.parent() |
1046 | 169 | storres | polynomialsList = [] |
1047 | 169 | storres | pVariables = p.variables() |
1048 | 169 | storres | iVariable = pVariables[0] |
1049 | 169 | storres | tVariable = pVariables[1] |
1050 | 169 | storres | polynomialAtPower = copy(p) |
1051 | 169 | storres | currentPolynomial = pRing(1) # Constant term. |
1052 | 169 | storres | |
1053 | 169 | storres | for pExp in xrange(0, alpha+1): |
1054 | 169 | storres | iExp = 0 |
1055 | 169 | storres | while (iExp + pExp) <= alpha: |
1056 | 179 | storres | #print "iExp:", iExp, \ |
1057 | 179 | storres | # "- pExp:", pExp, \ |
1058 | 179 | storres | # "- alpha-pExp:", alpha-pExp |
1059 | 169 | storres | q = pRing(iVariable * iBound)^iExp * p^pExp * N^(alpha-pExp) |
1060 | 179 | storres | #print q.monomials() |
1061 | 223 | storres | if columnsWidth != 0: |
1062 | 223 | storres | print " " * columnsWidth, q, |
1063 | 169 | storres | polynomialsList.append(q) |
1064 | 169 | storres | iExp += 1 |
1065 | 169 | storres | return polynomialsList |
1066 | 169 | storres | # End spo_polynomial_to_polynomials_list_8 |
1067 | 169 | storres | |
1068 | 111 | storres | def spo_proto_to_column_matrix(protoMatrixColumns): |
1069 | 111 | storres | """ |
1070 | 111 | storres | Create a column (each row holds the coefficients for one monomial) matrix. |
1071 | 111 | storres | |
1072 | 111 | storres | Parameters |
1073 | 111 | storres | ---------- |
1074 | 87 | storres | protoMatrixColumns: a list of coefficient lists. |
1075 | 83 | storres | """ |
1076 | 87 | storres | numColumns = len(protoMatrixColumns) |
1077 | 87 | storres | if numColumns == 0: |
1078 | 83 | storres | return None |
1079 | 87 | storres | # The last column holds has the maximum length. |
1080 | 87 | storres | numRows = len(protoMatrixColumns[numColumns-1]) |
1081 | 83 | storres | if numColumns == 0: |
1082 | 83 | storres | return None |
1083 | 83 | storres | baseMatrix = matrix(ZZ, numRows, numColumns) |
1084 | 87 | storres | for colIndex in xrange(0, numColumns): |
1085 | 87 | storres | for rowIndex in xrange(0, len(protoMatrixColumns[colIndex])): |
1086 | 90 | storres | if protoMatrixColumns[colIndex][rowIndex] != 0: |
1087 | 90 | storres | baseMatrix[rowIndex, colIndex] = \ |
1088 | 111 | storres | protoMatrixColumns[colIndex][rowIndex] |
1089 | 83 | storres | return baseMatrix |
1090 | 83 | storres | # End spo_proto_to_column_matrix. |
1091 | 83 | storres | # |
1092 | 111 | storres | def spo_proto_to_row_matrix(protoMatrixRows): |
1093 | 83 | storres | """ |
1094 | 111 | storres | Create a row (each column holds the coefficients corresponding to one |
1095 | 111 | storres | monomial) matrix from the protoMatrixRows list. |
1096 | 83 | storres | |
1097 | 83 | storres | Parameters |
1098 | 83 | storres | ---------- |
1099 | 83 | storres | protoMatrixRows: a list of coefficient lists. |
1100 | 83 | storres | """ |
1101 | 83 | storres | numRows = len(protoMatrixRows) |
1102 | 83 | storres | if numRows == 0: |
1103 | 83 | storres | return None |
1104 | 171 | storres | # Search for the longest row to get the number of columns. |
1105 | 171 | storres | numColumns = 0 |
1106 | 171 | storres | for row in protoMatrixRows: |
1107 | 171 | storres | rowLength = len(row) |
1108 | 171 | storres | if numColumns < rowLength: |
1109 | 171 | storres | numColumns = rowLength |
1110 | 83 | storres | if numColumns == 0: |
1111 | 83 | storres | return None |
1112 | 83 | storres | baseMatrix = matrix(ZZ, numRows, numColumns) |
1113 | 83 | storres | for rowIndex in xrange(0, numRows): |
1114 | 83 | storres | for colIndex in xrange(0, len(protoMatrixRows[rowIndex])): |
1115 | 89 | storres | if protoMatrixRows[rowIndex][colIndex] != 0: |
1116 | 89 | storres | baseMatrix[rowIndex, colIndex] = \ |
1117 | 111 | storres | protoMatrixRows[rowIndex][colIndex] |
1118 | 89 | storres | #print rowIndex, colIndex, |
1119 | 89 | storres | #print protoMatrixRows[rowIndex][colIndex], |
1120 | 89 | storres | #print knownMonomialsList[colIndex](boundVar1,boundVar2) |
1121 | 83 | storres | return baseMatrix |
1122 | 83 | storres | # End spo_proto_to_row_matrix. |
1123 | 83 | storres | # |
1124 | 246 | storres | sys.stderr.write("\t...sagePolynomialOperations loaded.\n") |