root / pobysoPythonSage / src / sageSLZ / sagePolynomialOperations.sage @ 108
Historique | Voir | Annoter | Télécharger (32,54 ko)
1 | 74 | storres | load "/home/storres/recherche/arithmetique/pobysoPythonSage/src/sageSLZ/sageMatrixOperations.sage" |
---|---|---|---|
2 | 87 | storres | print "sagePolynomialOperations loading..." |
3 | 106 | storres | def spo_add_polynomial_coeffs_to_matrix_row(poly, |
4 | 83 | storres | knownMonomials, |
5 | 83 | storres | protoMatrixRows, |
6 | 83 | storres | columnsWidth=0): |
7 | 80 | storres | """ |
8 | 106 | storres | For a given polynomial , |
9 | 80 | storres | add the coefficients of the protoMatrix (a list of proto matrix rows). |
10 | 80 | storres | Coefficients are added to the protoMatrix row in the order imposed by the |
11 | 80 | storres | monomials discovery list (the knownMonomials list) built as construction |
12 | 80 | storres | goes on. |
13 | 83 | storres | As a bonus, data can be printed out for a visual check. |
14 | 106 | storres | poly : the polynomial; in argument; |
15 | 106 | storres | knownMonomials : the list of the already known monomials; will determine |
16 | 106 | storres | the order of the coefficients appending to a row; in-out |
17 | 106 | storres | argument (new monomials may be discovered and then |
18 | 106 | storres | appended the the knowMonomials list); |
19 | 80 | storres | protoMatrixRows: a list of lists, each one holding the coefficients of the |
20 | 106 | storres | monomials of a polynomial; in-out argument: a new row is |
21 | 106 | storres | added at each call; |
22 | 80 | storres | columnWith : the width, in characters, of the displayed column ; if 0, |
23 | 106 | storres | do not display anything; in argument. |
24 | 80 | storres | """ |
25 | 106 | storres | pMonomials = poly.monomials() |
26 | 106 | storres | pCoefficients = poly.coefficients() |
27 | 80 | storres | # We have started with the smaller degrees in the first variable. |
28 | 80 | storres | pMonomials.reverse() |
29 | 80 | storres | pCoefficients.reverse() |
30 | 80 | storres | # New empty proto matrix row. |
31 | 80 | storres | protoMatrixRowCoefficients = [] |
32 | 80 | storres | # We work according to the order of the already known monomials |
33 | 80 | storres | # No known monomials yet: add the pMonomials to knownMonomials |
34 | 80 | storres | # and add the coefficients to the proto matrix row. |
35 | 80 | storres | if len(knownMonomials) == 0: |
36 | 80 | storres | for pmIdx in xrange(0, len(pMonomials)): |
37 | 80 | storres | knownMonomials.append(pMonomials[pmIdx]) |
38 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[pmIdx]) |
39 | 80 | storres | if columnsWidth != 0: |
40 | 80 | storres | monomialAsString = str(pCoefficients[pmIdx]) + " " + \ |
41 | 80 | storres | str(pMonomials[pmIdx]) |
42 | 80 | storres | print monomialAsString, " " * \ |
43 | 80 | storres | (columnsWidth - len(monomialAsString)), |
44 | 80 | storres | # There are some known monomials. We search for them in pMonomials and |
45 | 80 | storres | # add their coefficients to the proto matrix row. |
46 | 80 | storres | else: |
47 | 80 | storres | for knownMonomialIndex in xrange(0,len(knownMonomials)): |
48 | 80 | storres | # We lazily use an exception here since pMonomials.index() function |
49 | 80 | storres | # may fail throwing the ValueError exception. |
50 | 80 | storres | try: |
51 | 80 | storres | indexInPmonomials = \ |
52 | 80 | storres | pMonomials.index(knownMonomials[knownMonomialIndex]) |
53 | 80 | storres | if columnsWidth != 0: |
54 | 80 | storres | monomialAsString = str(pCoefficients[indexInPmonomials]) + \ |
55 | 80 | storres | " " + str(knownMonomials[knownMonomialIndex]) |
56 | 80 | storres | print monomialAsString, " " * \ |
57 | 80 | storres | (columnsWidth - len(monomialAsString)), |
58 | 80 | storres | # Add the coefficient to the proto matrix row and delete the \ |
59 | 80 | storres | # known monomial from the current pMonomial list |
60 | 80 | storres | #(and the corresponding coefficient as well). |
61 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[indexInPmonomials]) |
62 | 80 | storres | del pMonomials[indexInPmonomials] |
63 | 80 | storres | del pCoefficients[indexInPmonomials] |
64 | 80 | storres | # The knownMonomials element is not in pMonomials |
65 | 80 | storres | except ValueError: |
66 | 80 | storres | protoMatrixRowCoefficients.append(0) |
67 | 80 | storres | if columnsWidth != 0: |
68 | 80 | storres | monomialAsString = "0" + " "+ \ |
69 | 80 | storres | str(knownMonomials[knownMonomialIndex]) |
70 | 80 | storres | print monomialAsString, " " * \ |
71 | 80 | storres | (columnsWidth - len(monomialAsString)), |
72 | 80 | storres | # End for knownMonomialKey loop. |
73 | 80 | storres | # We now append the remaining monomials of pMonomials to knownMonomials |
74 | 80 | storres | # and the corresponding coefficients to proto matrix row. |
75 | 80 | storres | for pmIdx in xrange(0, len(pMonomials)): |
76 | 80 | storres | knownMonomials.append(pMonomials[pmIdx]) |
77 | 80 | storres | protoMatrixRowCoefficients.append(pCoefficients[pmIdx]) |
78 | 80 | storres | if columnsWidth != 0: |
79 | 80 | storres | monomialAsString = str(pCoefficients[pmIdx]) + " " \ |
80 | 80 | storres | + str(pMonomials[pmIdx]) |
81 | 80 | storres | print monomialAsString, " " * \ |
82 | 80 | storres | (columnsWidth - len(monomialAsString)), |
83 | 80 | storres | # End for pmIdx loop. |
84 | 80 | storres | # Add the new list row elements to the proto matrix. |
85 | 80 | storres | protoMatrixRows.append(protoMatrixRowCoefficients) |
86 | 80 | storres | if columnsWidth != 0: |
87 | 80 | storres | |
88 | 83 | storres | # End spo_add_polynomial_coeffs_to_matrix_row |
89 | 80 | storres | |
90 | 105 | storres | def spo_expression_as_string(powI, powT, powP, powN): |
91 | 80 | storres | """ |
92 | 80 | storres | Computes a string version of the i^k + t^l + p^m + N^n expression for |
93 | 80 | storres | output. |
94 | 80 | storres | """ |
95 | 80 | storres | expressionAsString ="" |
96 | 80 | storres | if powI != 0: |
97 | 80 | storres | expressionAsString += "i^" + str(powI) |
98 | 80 | storres | if powT != 0: |
99 | 80 | storres | if len(expressionAsString) != 0: |
100 | 80 | storres | expressionAsString += " * " |
101 | 80 | storres | expressionAsString += "t^" + str(powT) |
102 | 80 | storres | if powP != 0: |
103 | 80 | storres | if len(expressionAsString) != 0: |
104 | 80 | storres | expressionAsString += " * " |
105 | 80 | storres | expressionAsString += "p^" + str(powP) |
106 | 105 | storres | if (powN) != 0 : |
107 | 80 | storres | if len(expressionAsString) != 0: |
108 | 80 | storres | expressionAsString += " * " |
109 | 105 | storres | expressionAsString += "N^" + str(powN) |
110 | 80 | storres | return(expressionAsString) |
111 | 80 | storres | # End spo_expression_as_string. |
112 | 80 | storres | |
113 | 87 | storres | def spo_norm(poly, p=2): |
114 | 81 | storres | """ |
115 | 81 | storres | Behaves more or less (no infinity defined) as the norm for the |
116 | 81 | storres | univariate polynomials. |
117 | 107 | storres | Quoting Sage documentation: |
118 | 107 | storres | "Definition: For integer p, the p-norm of a polynomial is the pth root of |
119 | 81 | storres | the sum of the pth powers of the absolute values of the coefficients of |
120 | 107 | storres | the polynomial." |
121 | 87 | storres | |
122 | 81 | storres | """ |
123 | 87 | storres | # TODO: check the arguments (for p see below).. |
124 | 81 | storres | norm = 0 |
125 | 87 | storres | # For infinity norm. |
126 | 87 | storres | if p == Infinity: |
127 | 87 | storres | for coefficient in poly.coefficients(): |
128 | 87 | storres | coefficientAbs = coefficient.abs() |
129 | 87 | storres | if coefficientAbs > norm: |
130 | 87 | storres | norm = coefficientAbs |
131 | 87 | storres | return norm |
132 | 87 | storres | # TODO: check here the value of p |
133 | 107 | storres | # p must be a positive integer >= 1. |
134 | 107 | storres | if p < 1 or (not p in ZZ): |
135 | 94 | storres | return None |
136 | 87 | storres | # For 1 norm. |
137 | 87 | storres | if p == 1: |
138 | 87 | storres | for coefficient in poly.coefficients(): |
139 | 87 | storres | norm += coefficient.abs() |
140 | 87 | storres | return norm |
141 | 87 | storres | # For other norms |
142 | 81 | storres | for coefficient in poly.coefficients(): |
143 | 103 | storres | norm += coefficient.abs()^p |
144 | 87 | storres | return pow(norm, 1/p) |
145 | 81 | storres | # end spo_norm |
146 | 81 | storres | |
147 | 100 | storres | def spo_polynomial_to_proto_matrix(p, alpha, N, columnsWidth=0): |
148 | 74 | storres | """ |
149 | 83 | storres | From a (bivariate) polynomial and some other parameters build a proto |
150 | 87 | storres | matrix (an array of "rows") to be converted into a "true" matrix and |
151 | 83 | storres | eventually by reduced by fpLLL. |
152 | 102 | storres | The matrix is such as those found in Boneh-Durphee and Stehlé. |
153 | 74 | storres | |
154 | 83 | storres | Parameters |
155 | 83 | storres | ---------- |
156 | 87 | storres | p: the (bivariate) polynomial; |
157 | 87 | storres | pRing: the ring over which p is defined; |
158 | 74 | storres | alpha: |
159 | 74 | storres | N: |
160 | 83 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
161 | 83 | storres | printed in colums of columnsWitdth width. |
162 | 74 | storres | """ |
163 | 100 | storres | pRing = p.parent() |
164 | 77 | storres | knownMonomials = [] |
165 | 77 | storres | protoMatrixRows = [] |
166 | 92 | storres | polynomialsList = [] |
167 | 74 | storres | pVariables = p.variables() |
168 | 74 | storres | iVariable = pVariables[0] |
169 | 76 | storres | tVariable = pVariables[1] |
170 | 87 | storres | polynomialAtPower = pRing(1) |
171 | 87 | storres | currentPolynomial = pRing(1) |
172 | 74 | storres | pIdegree = p.degree(pVariables[0]) |
173 | 74 | storres | pTdegree = p.degree(pVariables[1]) |
174 | 87 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
175 | 105 | storres | nAtAlpha = N^alpha |
176 | 105 | storres | nAtPower = nAtAlpha |
177 | 92 | storres | polExpStr = "" |
178 | 74 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
179 | 74 | storres | for pPower in xrange(0, alpha + 1): |
180 | 76 | storres | # pPower == 0 is a special case. We introduce all the monomials but one |
181 | 78 | storres | # in i and those in t necessary to be able to introduce |
182 | 76 | storres | # p. We arbitrary choose to introduce the highest degree monomial in i |
183 | 76 | storres | # with p. We also introduce all the mixed i^k * t^l monomials with |
184 | 77 | storres | # k < p.degree(i) and l <= p.degree(t). |
185 | 78 | storres | # Mixed terms introduction is necessary here before we start "i shifts" |
186 | 78 | storres | # in the next iteration. |
187 | 74 | storres | if pPower == 0: |
188 | 78 | storres | # Notice that i^pIdegree is excluded as the bound of the xrange is |
189 | 78 | storres | # pIdegree |
190 | 74 | storres | for iPower in xrange(0, pIdegree): |
191 | 74 | storres | for tPower in xrange(0, pTdegree + 1): |
192 | 77 | storres | if columnsWidth != 0: |
193 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
194 | 76 | storres | tPower, |
195 | 76 | storres | pPower, |
196 | 105 | storres | alpha-pPower) |
197 | 92 | storres | print "->", polExpStr |
198 | 74 | storres | currentExpression = iVariable^iPower * \ |
199 | 91 | storres | tVariable^tPower * nAtAlpha |
200 | 78 | storres | # polynomialAtPower == 1 here. Next line should be commented |
201 | 78 | storres | # out but it does not work! Some conversion problem? |
202 | 91 | storres | currentPolynomial = pRing(currentExpression) |
203 | 106 | storres | polynomialsList.append(currentPolynomial) |
204 | 74 | storres | pMonomials = currentPolynomial.monomials() |
205 | 74 | storres | pCoefficients = currentPolynomial.coefficients() |
206 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
207 | 83 | storres | pCoefficients, |
208 | 83 | storres | knownMonomials, |
209 | 83 | storres | protoMatrixRows, |
210 | 83 | storres | columnsWidth) |
211 | 78 | storres | # End tPower. |
212 | 78 | storres | # End for iPower. |
213 | 77 | storres | else: # pPower > 0: (p^1..p^alpha) |
214 | 78 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
215 | 77 | storres | # polynomial. |
216 | 77 | storres | # This step could technically be fused as the first iteration |
217 | 77 | storres | # of the next loop (with iPower starting at 0). |
218 | 77 | storres | # We set it apart for clarity. |
219 | 77 | storres | if columnsWidth != 0: |
220 | 105 | storres | polExpStr = spo_expression_as_string(0, 0, pPower, alpha-pPower) |
221 | 92 | storres | print "->", polExpStr |
222 | 77 | storres | currentPolynomial = polynomialAtPower * nAtPower |
223 | 106 | storres | polynomialsList.append(currentPolynomial) |
224 | 77 | storres | pMonomials = currentPolynomial.monomials() |
225 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
226 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
227 | 83 | storres | pCoefficients, |
228 | 83 | storres | knownMonomials, |
229 | 83 | storres | protoMatrixRows, |
230 | 83 | storres | columnsWidth) |
231 | 77 | storres | |
232 | 77 | storres | # The i^iPower * p^pPower polynomials: they add i^k monomials to |
233 | 77 | storres | # p^pPower up to k < pIdegree * pPower. This only introduces i^k |
234 | 77 | storres | # monomials since mixed terms (that were introduced at a previous |
235 | 77 | storres | # stage) are only shifted to already existing |
236 | 77 | storres | # ones. p^pPower is "shifted" to higher degrees in i as far as |
237 | 77 | storres | # possible, one step short of the degree in i of p^(pPower+1) . |
238 | 77 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
239 | 77 | storres | for iPower in xrange(1, pIdegree): |
240 | 87 | storres | if columnsWidth != 0: |
241 | 92 | storres | polExpStr = spo_expression_as_string(iPower, \ |
242 | 87 | storres | 0, \ |
243 | 87 | storres | pPower, \ |
244 | 87 | storres | alpha) |
245 | 92 | storres | print "->", polExpStr |
246 | 77 | storres | currentExpression = i^iPower * nAtPower |
247 | 87 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
248 | 106 | storres | polynomialsList.append(currentPolynomial) |
249 | 77 | storres | pMonomials = currentPolynomial.monomials() |
250 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
251 | 87 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, \ |
252 | 87 | storres | pCoefficients, \ |
253 | 87 | storres | knownMonomials, \ |
254 | 87 | storres | protoMatrixRows, \ |
255 | 83 | storres | columnsWidth) |
256 | 77 | storres | # End for iPower |
257 | 77 | storres | # We want now to introduce a t * p^pPower polynomial. But before |
258 | 77 | storres | # that we must introduce some mixed monomials. |
259 | 77 | storres | # This loop is no triggered before pPower == 2. |
260 | 78 | storres | # It introduces a first set of high i degree mixed monomials. |
261 | 77 | storres | for iPower in xrange(1, pPower): |
262 | 77 | storres | tPower = pPower - iPower + 1 |
263 | 77 | storres | if columnsWidth != 0: |
264 | 92 | storres | polExpStr = spo_expression_as_string(iPower * pIdegree, |
265 | 77 | storres | tPower, |
266 | 77 | storres | 0, |
267 | 77 | storres | alpha) |
268 | 92 | storres | print "->", polExpStr |
269 | 91 | storres | currentExpression = i^(iPower * pIdegree) * t^tPower * nAtAlpha |
270 | 87 | storres | currentPolynomial = pRing(currentExpression) |
271 | 106 | storres | polynomialsList.append(currentPolynomial) |
272 | 77 | storres | pMonomials = currentPolynomial.monomials() |
273 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
274 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
275 | 83 | storres | pCoefficients, |
276 | 83 | storres | knownMonomials, |
277 | 83 | storres | protoMatrixRows, |
278 | 83 | storres | columnsWidth) |
279 | 77 | storres | # End for iPower |
280 | 78 | storres | # |
281 | 78 | storres | # This is the mixed monomials main loop. It introduces: |
282 | 77 | storres | # - the missing mixed monomials needed before the |
283 | 78 | storres | # t^l * p^pPower * N^(alpha-pPower) polynomial; |
284 | 78 | storres | # - the t^l * p^pPower * N^(alpha-pPower) itself; |
285 | 78 | storres | # - for each of i^k * t^l * p^pPower * N^(alpha-pPower) polynomials: |
286 | 78 | storres | # - the the missing mixed monomials needed polynomials, |
287 | 78 | storres | # - the i^k * t^l * p^pPower * N^(alpha-pPower) itself. |
288 | 78 | storres | # The t^l * p^pPower * N^(alpha-pPower) is introduced when |
289 | 78 | storres | # |
290 | 77 | storres | for iShift in xrange(0, pIdegree): |
291 | 77 | storres | # When pTdegree == 1, the following loop only introduces |
292 | 77 | storres | # a single new monomial. |
293 | 77 | storres | #print "++++++++++" |
294 | 77 | storres | for outerTpower in xrange(1, pTdegree + 1): |
295 | 77 | storres | # First one high i degree mixed monomial. |
296 | 77 | storres | iPower = iShift + pPower * pIdegree |
297 | 77 | storres | if columnsWidth != 0: |
298 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
299 | 77 | storres | outerTpower, |
300 | 77 | storres | 0, |
301 | 77 | storres | alpha) |
302 | 92 | storres | print "->", polExpStr |
303 | 91 | storres | currentExpression = i^iPower * t^outerTpower * nAtAlpha |
304 | 87 | storres | currentPolynomial = pRing(currentExpression) |
305 | 106 | storres | polynomialsList.append(currentPolynomial) |
306 | 77 | storres | pMonomials = currentPolynomial.monomials() |
307 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
308 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
309 | 83 | storres | pCoefficients, |
310 | 83 | storres | knownMonomials, |
311 | 83 | storres | protoMatrixRows, |
312 | 83 | storres | columnsWidth) |
313 | 77 | storres | #print "+++++" |
314 | 78 | storres | # At iShift == 0, the following innerTpower loop adds |
315 | 78 | storres | # duplicate monomials, since no extra i^l * t^k is needed |
316 | 78 | storres | # before introducing the |
317 | 77 | storres | # i^iShift * t^outerPpower * p^pPower * N^(alpha-pPower) |
318 | 77 | storres | # polynomial. |
319 | 77 | storres | # It introduces smaller i degree monomials than the |
320 | 77 | storres | # one(s) added previously (no pPower multiplication). |
321 | 77 | storres | # Here the exponent of t decreases as that of i increases. |
322 | 78 | storres | # This conditional is not entered before pPower == 1. |
323 | 78 | storres | # The innerTpower loop does not produce anything before |
324 | 78 | storres | # pPower == 2. We keep it anyway for other configuration of |
325 | 78 | storres | # p. |
326 | 77 | storres | if iShift > 0: |
327 | 77 | storres | iPower = pIdegree + iShift |
328 | 77 | storres | for innerTpower in xrange(pPower, 1, -1): |
329 | 77 | storres | if columnsWidth != 0: |
330 | 92 | storres | polExpStr = spo_expression_as_string(iPower, |
331 | 77 | storres | innerTpower, |
332 | 77 | storres | 0, |
333 | 77 | storres | alpha) |
334 | 77 | storres | currentExpression = \ |
335 | 91 | storres | i^(iPower) * t^(innerTpower) * nAtAlpha |
336 | 87 | storres | currentPolynomial = pRing(currentExpression) |
337 | 106 | storres | polynomialsList.append(currentPolynomial) |
338 | 77 | storres | pMonomials = currentPolynomial.monomials() |
339 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
340 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
341 | 77 | storres | pCoefficients, |
342 | 77 | storres | knownMonomials, |
343 | 77 | storres | protoMatrixRows, |
344 | 77 | storres | columnsWidth) |
345 | 77 | storres | iPower += pIdegree |
346 | 77 | storres | # End for innerTpower |
347 | 77 | storres | # End of if iShift > 0 |
348 | 78 | storres | # When iShift == 0, just after each of the |
349 | 78 | storres | # p^pPower * N^(alpha-pPower) polynomials has |
350 | 78 | storres | # been introduced (followed by a string of |
351 | 78 | storres | # i^k * p^pPower * N^(alpha-pPower) polynomials) a |
352 | 78 | storres | # t^l * p^pPower * N^(alpha-pPower) is introduced here. |
353 | 78 | storres | # |
354 | 77 | storres | # Eventually, the following section introduces the |
355 | 105 | storres | # i^iShift * t^outerTpower * p^iPower * N^(alpha-pPower) |
356 | 77 | storres | # polynomials. |
357 | 77 | storres | if columnsWidth != 0: |
358 | 92 | storres | polExpStr = spo_expression_as_string(iShift, |
359 | 77 | storres | outerTpower, |
360 | 77 | storres | pPower, |
361 | 105 | storres | alpha-pPower) |
362 | 92 | storres | print "->", polExpStr |
363 | 77 | storres | currentExpression = i^iShift * t^outerTpower * nAtPower |
364 | 105 | storres | currentPolynomial = pRing(currentExpression) * \ |
365 | 105 | storres | polynomialAtPower |
366 | 106 | storres | polynomialsList.append(currentPolynomial) |
367 | 77 | storres | pMonomials = currentPolynomial.monomials() |
368 | 77 | storres | pCoefficients = currentPolynomial.coefficients() |
369 | 83 | storres | spo_add_polynomial_coeffs_to_matrix_row(pMonomials, |
370 | 83 | storres | pCoefficients, |
371 | 83 | storres | knownMonomials, |
372 | 83 | storres | protoMatrixRows, |
373 | 83 | storres | columnsWidth) |
374 | 77 | storres | # End for outerTpower |
375 | 77 | storres | #print "++++++++++" |
376 | 77 | storres | # End for iShift |
377 | 77 | storres | polynomialAtPower *= p |
378 | 77 | storres | nAtPower /= N |
379 | 77 | storres | # End for pPower loop |
380 | 92 | storres | return ((protoMatrixRows, knownMonomials, polynomialsList)) |
381 | 83 | storres | # End spo_polynomial_to_proto_matrix |
382 | 81 | storres | |
383 | 106 | storres | def spo_polynomial_to_polynomials_list_2(p, alpha, N, columnsWidth=0): |
384 | 105 | storres | """ |
385 | 106 | storres | From p, alpha, N build a list of polynomials... |
386 | 106 | storres | TODO: clean up the comments below! |
387 | 106 | storres | |
388 | 105 | storres | From a (bivariate) polynomial and some other parameters build a proto |
389 | 105 | storres | matrix (an array of "rows") to be converted into a "true" matrix and |
390 | 105 | storres | eventually by reduced by fpLLL. |
391 | 105 | storres | The matrix is based on a list of polynomials that are built in a way |
392 | 105 | storres | that one and only monomial is added at each new polynomial. Among the many |
393 | 105 | storres | possible ways to build this list we pick one strongly dependent on the |
394 | 105 | storres | structure of the polynomial and of the problem. |
395 | 105 | storres | We consider here the polynomials of the form: |
396 | 105 | storres | a_k*i^k + a_(k-1)*i^(k-1) + ... + a_1*i + a_0 - t |
397 | 105 | storres | The values of i and t are bounded and we eventually look for (i_0,t_0) |
398 | 105 | storres | pairs such that: |
399 | 105 | storres | a_k*i_0^k + a_(k-1)*i_0^(k-1) + ... + a_1*i_0 + a_0 = t_0 |
400 | 105 | storres | Hence, departing from the procedure in described in Boneh-Durfee, we will |
401 | 105 | storres | not use "t-shifts" but only "i-shifts". |
402 | 105 | storres | |
403 | 105 | storres | Parameters |
404 | 105 | storres | ---------- |
405 | 105 | storres | p: the (bivariate) polynomial; |
406 | 105 | storres | pRing: the ring over which p is defined; |
407 | 105 | storres | alpha: |
408 | 105 | storres | N: |
409 | 105 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
410 | 105 | storres | printed in colums of columnsWitdth width. |
411 | 105 | storres | """ |
412 | 105 | storres | pRing = p.parent() |
413 | 105 | storres | knownMonomials = [] |
414 | 105 | storres | polynomialsList = [] |
415 | 105 | storres | pVariables = p.variables() |
416 | 105 | storres | iVariable = pVariables[0] |
417 | 105 | storres | tVariable = pVariables[1] |
418 | 105 | storres | polynomialAtPower = pRing(1) |
419 | 105 | storres | currentPolynomial = pRing(1) |
420 | 105 | storres | pIdegree = p.degree(iVariable) |
421 | 105 | storres | pTdegree = p.degree(tVariable) |
422 | 105 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
423 | 105 | storres | nAtAlpha = N^alpha |
424 | 105 | storres | nAtPower = nAtAlpha |
425 | 105 | storres | polExpStr = "" |
426 | 105 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
427 | 105 | storres | for pPower in xrange(0, alpha + 1): |
428 | 105 | storres | # pPower == 0 is a special case. We introduce all the monomials in i |
429 | 105 | storres | # up to i^pIdegree. |
430 | 105 | storres | if pPower == 0: |
431 | 105 | storres | # Notice who iPower runs up to i^pIdegree. |
432 | 105 | storres | for iPower in xrange(0, pIdegree + 1): |
433 | 105 | storres | # No t power is taken into account as we limit our selves to |
434 | 105 | storres | # degree 1 in t and make no "t-shifts". |
435 | 105 | storres | if columnsWidth != 0: |
436 | 105 | storres | polExpStr = spo_expression_as_string(iPower, |
437 | 105 | storres | 0, |
438 | 105 | storres | 0, |
439 | 105 | storres | alpha) |
440 | 105 | storres | print "->", polExpStr |
441 | 105 | storres | currentExpression = iVariable^iPower * nAtAlpha |
442 | 105 | storres | # polynomialAtPower == 1 here. Next line should be commented |
443 | 105 | storres | # out but it does not work! Some conversion problem? |
444 | 105 | storres | currentPolynomial = pRing(currentExpression) |
445 | 105 | storres | polynomialsList.append(currentPolynomial) |
446 | 105 | storres | # End for iPower. |
447 | 105 | storres | else: # pPower > 0: (p^1..p^alpha) |
448 | 105 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
449 | 105 | storres | # polynomial. This is also where the t^pPower monomials shows up for |
450 | 105 | storres | # the first time. |
451 | 105 | storres | if columnsWidth != 0: |
452 | 105 | storres | polExpStr = spo_expression_as_string(0, 0, pPower, alpha-pPower) |
453 | 105 | storres | print "->", polExpStr |
454 | 105 | storres | currentPolynomial = polynomialAtPower * nAtPower |
455 | 105 | storres | polynomialsList.append(currentPolynomial) |
456 | 106 | storres | # Exit when pPower == alpha |
457 | 106 | storres | if pPower == alpha: |
458 | 106 | storres | return((knownMonomials, polynomialsList)) |
459 | 105 | storres | # This is where the "i-shifts" take place. Mixed terms, i^k * t^l |
460 | 105 | storres | # (that were introduced at a previous |
461 | 105 | storres | # stage or are introduced now) are only shifted to already existing |
462 | 105 | storres | # ones with the notable exception of i^iPower * t^pPower, which |
463 | 105 | storres | # must be manually introduced. |
464 | 105 | storres | # p^pPower is "shifted" to higher degrees in i as far as |
465 | 105 | storres | # possible, up to of the degree in i of p^(pPower+1). |
466 | 105 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
467 | 105 | storres | for iPower in xrange(1, pIdegree + 1): |
468 | 105 | storres | # The i^iPower * t^pPower monomial. Notice the alpha exponent |
469 | 105 | storres | # for N. |
470 | 105 | storres | internalIpower = iPower |
471 | 105 | storres | for tPower in xrange(pPower,0,-1): |
472 | 105 | storres | if columnsWidth != 0: |
473 | 105 | storres | polExpStr = spo_expression_as_string(internalIpower, \ |
474 | 105 | storres | tPower, \ |
475 | 105 | storres | 0, \ |
476 | 105 | storres | alpha) |
477 | 105 | storres | print "->", polExpStr |
478 | 105 | storres | currentExpression = i^internalIpower * t^tPower * nAtAlpha |
479 | 105 | storres | currentPolynomial = pRing(currentExpression) |
480 | 105 | storres | polynomialsList.append(currentPolynomial) |
481 | 105 | storres | internalIpower += pIdegree |
482 | 105 | storres | # End for tPower |
483 | 105 | storres | # The i^iPower * p^pPower * N^(alpha-pPower) i-shift. |
484 | 105 | storres | if columnsWidth != 0: |
485 | 105 | storres | polExpStr = spo_expression_as_string(iPower, \ |
486 | 105 | storres | 0, \ |
487 | 105 | storres | pPower, \ |
488 | 105 | storres | alpha-pPower) |
489 | 105 | storres | print "->", polExpStr |
490 | 105 | storres | currentExpression = i^iPower * nAtPower |
491 | 105 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
492 | 105 | storres | polynomialsList.append(currentPolynomial) |
493 | 105 | storres | # End for iPower |
494 | 105 | storres | polynomialAtPower *= p |
495 | 105 | storres | nAtPower /= N |
496 | 105 | storres | # End for pPower loop |
497 | 106 | storres | return((knownMonomials, polynomialsList)) |
498 | 105 | storres | # End spo_polynomial_to_proto_matrix_2 |
499 | 105 | storres | |
500 | 108 | storres | def spo_polynomial_to_polynomials_list_3(p, alpha, N, columnsWidth=0): |
501 | 108 | storres | """ |
502 | 108 | storres | From p, alpha, N build a list of polynomials... |
503 | 108 | storres | TODO: more in depth rationale... |
504 | 108 | storres | |
505 | 108 | storres | Our goal is to introduce each monomial with the smallest coefficient. |
506 | 108 | storres | |
507 | 108 | storres | |
508 | 108 | storres | |
509 | 108 | storres | Parameters |
510 | 108 | storres | ---------- |
511 | 108 | storres | p: the (bivariate) polynomial; |
512 | 108 | storres | pRing: the ring over which p is defined; |
513 | 108 | storres | alpha: |
514 | 108 | storres | N: |
515 | 108 | storres | columsWidth: if == 0, no information is displayed, otherwise data is |
516 | 108 | storres | printed in colums of columnsWitdth width. |
517 | 108 | storres | """ |
518 | 108 | storres | pRing = p.parent() |
519 | 108 | storres | knownMonomials = [] |
520 | 108 | storres | polynomialsList = [] |
521 | 108 | storres | pVariables = p.variables() |
522 | 108 | storres | iVariable = pVariables[0] |
523 | 108 | storres | tVariable = pVariables[1] |
524 | 108 | storres | polynomialAtPower = pRing(1) |
525 | 108 | storres | currentPolynomial = pRing(1) |
526 | 108 | storres | pIdegree = p.degree(iVariable) |
527 | 108 | storres | pTdegree = p.degree(tVariable) |
528 | 108 | storres | currentIdegree = currentPolynomial.degree(iVariable) |
529 | 108 | storres | nAtAlpha = N^alpha |
530 | 108 | storres | nAtPower = nAtAlpha |
531 | 108 | storres | polExpStr = "" |
532 | 108 | storres | # We work from p^0 * N^alpha to p^alpha * N^0 |
533 | 108 | storres | for pPower in xrange(0, alpha + 1): |
534 | 108 | storres | # pPower == 0 is a special case. We introduce all the monomials in i |
535 | 108 | storres | # up to i^pIdegree. |
536 | 108 | storres | if pPower == 0: |
537 | 108 | storres | # Notice who iPower runs up to i^pIdegree. |
538 | 108 | storres | for iPower in xrange(0, pIdegree + 1): |
539 | 108 | storres | # No t power is taken into account as we limit our selves to |
540 | 108 | storres | # degree 1 in t and make no "t-shifts". |
541 | 108 | storres | if columnsWidth != 0: |
542 | 108 | storres | polExpStr = spo_expression_as_string(iPower, |
543 | 108 | storres | 0, |
544 | 108 | storres | 0, |
545 | 108 | storres | alpha) |
546 | 108 | storres | print "->", polExpStr |
547 | 108 | storres | currentExpression = iVariable^iPower * nAtAlpha |
548 | 108 | storres | # polynomialAtPower == 1 here. Next line should be commented |
549 | 108 | storres | # out but it does not work! Some conversion problem? |
550 | 108 | storres | currentPolynomial = pRing(currentExpression) |
551 | 108 | storres | polynomialsList.append(currentPolynomial) |
552 | 108 | storres | # End for iPower. |
553 | 108 | storres | else: # pPower > 0: (p^1..p^alpha) |
554 | 108 | storres | # This where we introduce the p^pPower * N^(alpha-pPower) |
555 | 108 | storres | # polynomial. This is also where the t^pPower monomials shows up for |
556 | 108 | storres | # the first time. It app |
557 | 108 | storres | if columnsWidth != 0: |
558 | 108 | storres | polExpStr = spo_expression_as_string(0, 0, pPower, alpha-pPower) |
559 | 108 | storres | print "->", polExpStr |
560 | 108 | storres | currentPolynomial = polynomialAtPower * nAtPower |
561 | 108 | storres | polynomialsList.append(currentPolynomial) |
562 | 108 | storres | # Exit when pPower == alpha |
563 | 108 | storres | if pPower == alpha: |
564 | 108 | storres | return((knownMonomials, polynomialsList)) |
565 | 108 | storres | # This is where the "i-shifts" take place. Mixed terms, i^k * t^l |
566 | 108 | storres | # (that were introduced at a previous |
567 | 108 | storres | # stage or are introduced now) are only shifted to already existing |
568 | 108 | storres | # ones with the notable exception of i^iPower * t^pPower, which |
569 | 108 | storres | # must be manually introduced. |
570 | 108 | storres | # p^pPower is "shifted" to higher degrees in i as far as |
571 | 108 | storres | # possible, up to of the degree in i of p^(pPower+1). |
572 | 108 | storres | # These "pure" i^k monomials can only show up with i multiplications. |
573 | 108 | storres | for iPower in xrange(1, pIdegree + 1): |
574 | 108 | storres | # The i^iPower * t^pPower monomial. Notice the alpha exponent |
575 | 108 | storres | # for N. |
576 | 108 | storres | internalIpower = iPower |
577 | 108 | storres | for tPower in xrange(pPower,0,-1): |
578 | 108 | storres | if columnsWidth != 0: |
579 | 108 | storres | polExpStr = spo_expression_as_string(internalIpower, \ |
580 | 108 | storres | tPower, \ |
581 | 108 | storres | 0, \ |
582 | 108 | storres | alpha) |
583 | 108 | storres | print "->", polExpStr |
584 | 108 | storres | currentExpression = i^internalIpower * t^tPower * nAtAlpha |
585 | 108 | storres | currentPolynomial = pRing(currentExpression) |
586 | 108 | storres | polynomialsList.append(currentPolynomial) |
587 | 108 | storres | internalIpower += pIdegree |
588 | 108 | storres | # End for tPower |
589 | 108 | storres | # The i^iPower * p^pPower * N^(alpha-pPower) i-shift. |
590 | 108 | storres | if columnsWidth != 0: |
591 | 108 | storres | polExpStr = spo_expression_as_string(iPower, \ |
592 | 108 | storres | 0, \ |
593 | 108 | storres | pPower, \ |
594 | 108 | storres | alpha-pPower) |
595 | 108 | storres | print "->", polExpStr |
596 | 108 | storres | currentExpression = i^iPower * nAtPower |
597 | 108 | storres | currentPolynomial = pRing(currentExpression) * polynomialAtPower |
598 | 108 | storres | polynomialsList.append(currentPolynomial) |
599 | 108 | storres | # End for iPower |
600 | 108 | storres | polynomialAtPower *= p |
601 | 108 | storres | nAtPower /= N |
602 | 108 | storres | # End for pPower loop |
603 | 108 | storres | return((knownMonomials, polynomialsList)) |
604 | 108 | storres | # End spo_polynomial_to_proto_matrix_3 |
605 | 108 | storres | |
606 | 90 | storres | def spo_proto_to_column_matrix(protoMatrixColumns, \ |
607 | 90 | storres | knownMonomialsList, \ |
608 | 90 | storres | boundVar1, \ |
609 | 90 | storres | boundVar2): |
610 | 83 | storres | """ |
611 | 87 | storres | Create a column (each row holds the coefficients of one monomial) matrix. |
612 | 83 | storres | |
613 | 83 | storres | Parameters |
614 | 83 | storres | ---------- |
615 | 87 | storres | protoMatrixColumns: a list of coefficient lists. |
616 | 83 | storres | """ |
617 | 87 | storres | numColumns = len(protoMatrixColumns) |
618 | 87 | storres | if numColumns == 0: |
619 | 83 | storres | return None |
620 | 87 | storres | # The last column holds has the maximum length. |
621 | 87 | storres | numRows = len(protoMatrixColumns[numColumns-1]) |
622 | 83 | storres | if numColumns == 0: |
623 | 83 | storres | return None |
624 | 83 | storres | baseMatrix = matrix(ZZ, numRows, numColumns) |
625 | 87 | storres | for colIndex in xrange(0, numColumns): |
626 | 87 | storres | for rowIndex in xrange(0, len(protoMatrixColumns[colIndex])): |
627 | 90 | storres | if protoMatrixColumns[colIndex][rowIndex] != 0: |
628 | 90 | storres | baseMatrix[rowIndex, colIndex] = \ |
629 | 90 | storres | protoMatrixColumns[colIndex][rowIndex] * \ |
630 | 90 | storres | knownMonomialsList[rowIndex](boundVar1, boundVar2) |
631 | 83 | storres | return baseMatrix |
632 | 83 | storres | # End spo_proto_to_column_matrix. |
633 | 83 | storres | # |
634 | 88 | storres | def spo_proto_to_row_matrix(protoMatrixRows, \ |
635 | 88 | storres | knownMonomialsList, \ |
636 | 88 | storres | boundVar1, \ |
637 | 88 | storres | boundVar2): |
638 | 83 | storres | """ |
639 | 90 | storres | Create a row (each column holds the evaluation one monomial at boundVar1 and |
640 | 90 | storres | boundVar2 values) matrix. |
641 | 83 | storres | |
642 | 83 | storres | Parameters |
643 | 83 | storres | ---------- |
644 | 83 | storres | protoMatrixRows: a list of coefficient lists. |
645 | 83 | storres | """ |
646 | 83 | storres | numRows = len(protoMatrixRows) |
647 | 83 | storres | if numRows == 0: |
648 | 83 | storres | return None |
649 | 91 | storres | # The last row is the longest one. |
650 | 83 | storres | numColumns = len(protoMatrixRows[numRows-1]) |
651 | 83 | storres | if numColumns == 0: |
652 | 83 | storres | return None |
653 | 83 | storres | baseMatrix = matrix(ZZ, numRows, numColumns) |
654 | 83 | storres | for rowIndex in xrange(0, numRows): |
655 | 83 | storres | for colIndex in xrange(0, len(protoMatrixRows[rowIndex])): |
656 | 89 | storres | if protoMatrixRows[rowIndex][colIndex] != 0: |
657 | 89 | storres | baseMatrix[rowIndex, colIndex] = \ |
658 | 89 | storres | protoMatrixRows[rowIndex][colIndex] * \ |
659 | 89 | storres | knownMonomialsList[colIndex](boundVar1,boundVar2) |
660 | 89 | storres | #print rowIndex, colIndex, |
661 | 89 | storres | #print protoMatrixRows[rowIndex][colIndex], |
662 | 89 | storres | #print knownMonomialsList[colIndex](boundVar1,boundVar2) |
663 | 83 | storres | return baseMatrix |
664 | 83 | storres | # End spo_proto_to_row_matrix. |
665 | 83 | storres | # |
666 | 87 | storres | print "\t...sagePolynomialOperations loaded" |