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