|
The operations described in this section are relevant for polynomials
over series rings. There are two main algorithms involved.
Set the verbose printing to level v for PuiseuxExpansion,
ExpandToPrecision, DuvalPuiseuxExpansion, Roots
and ImplicitFunction. A level of 1 will mean that any partial
solutions that could not be expanded to the precision requested will be
printed before an error is returned except for Roots.
The polynomials used in forming
extensions will also be printed before the extension is computed.
In Roots, the algorithm used to compute the expansions will be printed.
When Walker's algorithm is being used the current value of the denominator
will be printed.
For ImplicitFunction a warning about a potentially bad value of d
will be printed if the value of d given is not divisible by the
exponent denominator of some coefficient of f.
A level of 2
will print the last polynomials calculated during the newton polygon part
of PuiseuxExpansion and DuvalPuiseuxExpansion and some
evaluated polynomials during ImplicitFunction.
PuiseuxExpansion(f, n) : RngUPolElt, RngIntElt -> SeqEnum[RngSerPuisElt]
PreciseRoot: BoolElt Default: false
TestSquarefree: BoolElt Default: true
NoExtensions: BoolElt Default: false
LowerFaces: BoolElt Default: true
OneRoot: BoolElt Default: false
SetVerbose("Newton", n): Maximum: 2
This function implements the algorithm described in [Wal78].
Return a sequence of partial expansions of the roots of the polynomial
f over a series ring as puiseux series. The roots are returned with relative
precision at least n/d where d is the least common multiple of exponent
denominator for the
series expansion and the exponent denominators of the coefficients of the f.
An input of n = 0 will return the expansions calculated
by the newton polygon part of the algorithm and these will be to the
precision of what is known. The coefficient ring of the series ring
containing the coefficients of f must always be a field
and unless extensions are not required it must be able to be extended.
If the coefficient ring of the series ring is a finite field whose
characteristic is less than or equal to the degree of the polynomial then the
denominators computed in the newton polygon part of the algorithm
may not be bounded and the function will return an error. However, it
is possible for some polynomials that the denominators will be bounded.
This is stated by [Gri95], pg 269 -- 272.
Care needs to be taken with polynomials whose coefficients have low precision.
The algorithm must extract
from f the squarefree part and in doing so lose even more precision.
The result is that the algorithm may not have enough precision to calculate the
expansions correctly.
A solution is to set TestSquarefree to false if the polynomial
is known to have no multiple roots. However this will not solve all precision
problems and the answer is only as good as the precision allows it to be.
If PreciseRoot is set to true then
the partial expansions to be returned are checked and if any are exact
roots of f they are returned with full precision.
If NoExtensions is set to true then expansions are found
within the puiseux series ring only. By default, the coefficient ring of the
series ring is extended to find all the partial expansions.
If LowerFaces is set to false then expansions with negative
valuations will not be found.
If OneRoot is set to true then representatives of conjugate roots
only will be found instead of each of the roots individually.
Note that it may useful to define the polynomial over an algebraically
closed field (via AlgebraicClosure), so that all roots may be
found.
PreciseRoot: BoolElt Default: false
TestSquarefree: BoolElt Default: true
SetVerbose("Newton", n): Maximum: 2
Given a polynomial f over a Puiseux series ring and a partial root c of
that polynomial (found by PuiseuxExpansion for example)
continue to expand that root until it has relative precision n/d
where d is the least common multiple of the exponent denominator of c
and the exponent denominators of the coefficients of f.
If c is given to greater
precision (or length greater than n if its precision is infinite), the
relative precision of c is reduced to n/d. An error results if
c is not a partial expansion to precision n/d.
If PreciseRoot is true then the partial expansion to be returned is
checked and if it is an exact root of f it is returned with full precision.
All input is checked for being an exact root regardless.
If TestSquarefree is false the polynomial will not be made
squarefree. This may avoid some loss of precision but may result in some
unique partial roots not being recognized as unique partial roots and
as such they cannot be expanded.
An error may result if c is a partial root of f but the exponent
denominator of the full expansion is greater than that of c.
Therefore for c to be expanded it must have the same denominator as
the expansion it is part of. This will rarely be an issue for those partial
expansions resulting from PusieuxExpansion which did not encounter
problems with precision in the newton polygon part since the algorithm
for this will find at least as
much of the expansion as necessary to compute the exponent denominator.
SetVerbose("Newton", n): Maximum: 2
Return a root of the polynomial f over a series ring.
The input d is the denominator (or a multiple of) the exponent denominator
of the root. The root is given to absolute precision n/d. The evaluation
of f at zero (polynomial) evaluated at zero (series) must be zero but that
of its derivative must be nonzero.
This example illustrates the joint use of PuiseuxExpansion and
ExpandToPrecision which can be used together to gain and improve
partial roots of a polynomial. The use of ExpandToPrecision following
PuiseuxExpansion avoids the recalculation of information already known.
> P<x> := PuiseuxSeriesRing(Rationals());
> R<y> := PolynomialRing(P);
> f := y^3 + 2*x^-1*y^2 + 1*x^-2*y + 2*x;
> c := PuiseuxExpansion(f, 0);
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> Q<q> := PolynomialRing(A);
> c;
[
-2*a^3 + O(a^4),
-a^-1 + n*a + O(a^2),
-a^-1 - n*a + O(a^2)
]
> [ExpandToPrecision(f, c[i], 10) : i in [1 .. #c]];
[
-2*a^3 - 8*a^7 - 56*a^11 + O(a^13),
-a^-1 + n*a + a^3 + 5/4*n*a^5 + 4*a^7 + O(a^9),
-a^-1 - n*a + a^3 - 5/4*n*a^5 + 4*a^7 + O(a^9)
]
The same results could have been gained using PuiseuxExpansion with the
required precision in the first place.
> c := PuiseuxExpansion(f, 10);
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> c;
[
-2*a^3 - 8*a^7 - 56*a^11 + O(a^13),
-a^-1 + n*a + a^3 + 5/4*n*a^5 + 4*a^7 + O(a^9),
-a^-1 - n*a + a^3 - 5/4*n*a^5 + 4*a^7 + O(a^9)
]
However, asking for more precision requires time so that if it is not
necessary the extra calculation can be avoided and if more precision
happens to be required then it can be gained without recalculation.
ExpandToPrecision is also called on only one root so that if
only one expansion is required using PusieuxExpansion and then
ExpandToPrecision will not calculate any unnecessary information.
> time c := PuiseuxExpansion(f, 100);
Time: 2.810
> time c := PuiseuxExpansion(f, 10);
Time: 0.060
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> time ExpandToPrecision(f, c[1], 100);
-2*a^3 - 8*a^7 - 56*a^11 - 480*a^15 - 4576*a^19 - 46592*a^23 -
496128*a^27 - 5457408*a^31 - 61529600*a^35 - 707266560*a^39 -
8257566720*a^43 - 97654702080*a^47 - 1167349284864*a^51 -
14082308833280*a^55 - 171221451538432*a^59 -
2096081963188224*a^63 - 25814314231136256*a^67 -
319605795242639360*a^71 - 3975750610806374400*a^75 -
49666299938073477120*a^79 - 622818862289639178240*a^83 -
7837247078959687925760*a^87 - 98931046460491133091840*a^91 -
1252424949872174982758400*a^95 - 15897106567806080658702336*a^99
+ O(a^103)
Time: 0.410
Return true if the series c can be expanded to at least one root of
the polynomial f.
TestSquarefree: BoolElt Default: true
Return true if the series
c can be expanded to exactly one distinct root of the polynomial f.
By default f will have multiple factors removed to allow partial expansions
of multiple roots to be recognized as being unique. If TestSquarefree
is set to false then f will be taken as given which may avoid errors
due to lost precision but may not pick partial expansions of multiple roots
as being unique and as such is best used when f is squarefree or the expansion
is known to be of a single root.
The above 2 functions can be used to reduce the occurrence of errors
from ExpandToPrecision by checking that the input can be expanded.
Errors resulting from a lack of precision which means that the expansion
cannot be calculated to the requested precision are the only errors that
cannot be removed. Only unique partial roots can be expanded.
If a partial root is not unique then calling PuiseuxExpansion will
provide several further partial expansions of the partial root that
will themselves be unique
and so can be used to calculate several expansions of the original.
> P<x> := PuiseuxSeriesRing(Rationals());
> R<y> := PolynomialRing(P);
> f := (y^2 - x^3)^2 - y*x^6;
> IsPartialRoot(f, x^(3/2));
true
> ExpandToPrecision(f, x^(3/2), 10);
>> ExpandToPrecision(f, x^(3/2), 10);
^
Runtime error in 'ExpandToPrecision': Element is not a unique partial
root of the polynomial
> IsUniquePartialRoot(f, x^(3/2));
false
> c := PuiseuxExpansion(f, 0);
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> Q<q> := PolynomialRing(A);
> c;
[
a^(3/2) + 1/2*a^(9/4) + O(a^(5/2)),
a^(3/2) - 1/2*a^(9/4) + O(a^(5/2)),
-a^(3/2) + 1/2*n*a^(9/4) + O(a^(5/2)),
-a^(3/2) - 1/2*n*a^(9/4) + O(a^(5/2))
]
> IsUniquePartialRoot(f, x^(3/2) + 1/2*x^(9/4));
true
> ExpandToPrecision(f, x^(3/2) + 1/2*x^(9/4), 10);
x^(3/2) + 1/2*x^(9/4) - 1/64*x^(15/4) + O(x^4)
> ExpandToPrecision(f, x^(3/2) + x^2, 30);
>> ExpandToPrecision(f, x^(3/2) + x^2, 30);
^
Runtime error in 'ExpandToPrecision': Element is not a partial root of
the polynomial
> IsPartialRoot(f, x^(3/2) + x^2);
false
So if IsPartialRoot returns false then no expansion can be made.
If IsUniquePartialRoot returns false (but IsPartialRoot
returns true) then several expansions can be made after calling
PuiseuxExpansion.
Given a series expansion return the sequence of exponents [a/b]
of the non zero terms of the series p up to and including the first one where
b is the global denominator for the series.
Given two series return the sequence of exponents [a/b] of the non
zero initial terms of the series p and q which are equal up to but not including the
first unequal terms.
This example illustrates how PuiseuxExponents and
PuiseuxExponentsCommon can be used on output from PusieuxExpansion.
(Similar can be done with related functions and general series).
> P<x> := PuiseuxSeriesRing(FiniteField(5, 3));
> R<y> := PolynomialRing(P);
> f := (1+x)*y^4 - x^(-1/3)*y^2 + y + x^(1/2);
> time c := PuiseuxExpansion(f, 5);
Time: 0.030
> c;
[
4*x^(1/2) + x^(2/3) + 3*x^(5/6) + x^(7/6) + O(x^(4/3)),
x^(1/3) + x^(1/2) + 4*x^(2/3) + 2*x^(5/6) + O(x^(7/6)),
4*x^(-1/6) + 2*x^(1/3) + O(x^(2/3)),
x^(-1/6) + 2*x^(1/3) + O(x^(2/3))
]
> PuiseuxExponents(c[1]);
[ 1/2, 2/3, 5/6 ]
> PuiseuxExponents(c[3]);
[ -1/6 ]
> P<x> := PuiseuxSeriesRing(FiniteField(5, 3));
> R<y> := PolynomialRing(P);
> f := ((y^2 - x^3)^2 - y*x^6)^2 - y*x^15;
> c := PuiseuxExpansion(f, 0);
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> Q<q> := PolynomialRing(A);
> c;
[
4*a^(3/2) + 4*a^(9/4) + 4*a^3 + O(a^(13/4)),
4*a^(3/2) + 4*a^(9/4) + a^3 + O(a^(13/4)),
4*a^(3/2) + a^(9/4) + 4*a^3 + O(a^(13/4)),
4*a^(3/2) + a^(9/4) + a^3 + O(a^(13/4)),
a^(3/2) + 3*a^(9/4) + 4*a^3 + O(a^(13/4)),
a^(3/2) + 3*a^(9/4) + a^3 + O(a^(13/4)),
a^(3/2) + 2*a^(9/4) + 4*a^3 + O(a^(13/4)),
a^(3/2) + 2*a^(9/4) + a^3 + O(a^(13/4))
]
> PuiseuxExponentsCommon(c[1], c[1]);
[ 3/2, 9/4, 3 ]
> PuiseuxExponentsCommon(c[1], c[2]);
[ 3/2, 9/4 ]
> PuiseuxExponentsCommon(c[1], c[3]);
[ 3/2 ]
> PuiseuxExponentsCommon(c[1], c[8]);
[]
The following functions have a similar use to those given above
but implement a different algorithm, namely that of [Duv89]
which is faster and can handle larger degree polynomials. However,
it can only be used with polynomials which are essentially over a laurent
series ring and the coefficient ring of that laurent series ring
has either characteristic zero or characteristic greater than the degrees
of the squarefree factors of the polynomial.
Version: MonStgElt Default: "Rational"
TestSquarefree: BoolElt Default: true
NoExtensions: BoolElt Default: false
LowerFaces: BoolElt Default: true
OneRoot: BoolElt Default: false
SetVerbose("Newton", n): Maximum: 2
A sequence of parametrizations of puiseux expansions of roots of f,
as puiseux series, where f is a polynomial over a series ring.
The expansions will have at least n
non zero terms (unless the expansion is finite and has less than n
non zero terms), with more than n occurring only if n is less than the
number of terms returned by the newton polygon part of the algorithm.
The coefficients of f must have exponent denominator 1.
This algorithm is faster than that given by Walker and implemented in
PuiseuxExpansion, since it doesn't calculate non zero terms explicitly
and doesn't make all necessary extensions during the algorithm leaving
some to be made when the series is computed from the parametrizations.
If f has coefficients with finite precision then the expansions can
only be computed to as many non zero terms as can be known for that expansion.
After this limit has been reached, an error results since the next non
zero term is not known. If f has roots with finite puiseux expansions
then if n is greater than the number of non zero terms in the expansion
the expansion is returned with infinite precision.
If Version is set to "Classical" then the (slower)
classical branch of the algorithm will be run which makes all extensions
necessary for the computation of the expansions. It is still faster than
PusieuxExpansion since it does not iterate through and calculate
zero terms but will encounter the same problems that PuiseuxExpansion
does with field extensions over the rationals. The classical version
will return as many parametrizations as there are expansions and some of
these parametrizations will give the same set of expansions.
If NoExtensions is set to true then only the expansions which lie
in the puiseux series ring corresponding to the coefficient ring of f
are calculated. Otherwise, all the expansions of roots of f are calculated
regardless of where they lie.
LowerFaces and OneRoot work as for PuiseuxExpansion.
This algorithm works with the squarefree part of f only. If any coefficient
of f has low precision then this step may make it impossible for any
information about the expansions to be gained due to a loss of further
precision. A way around this is to set TestSquarefree to false
if the polynomial is known to be squarefree. This may result in some
information being returned but such information is only as good
as the precision it was allowed.
The series that satisfy the parametrization T. These are found by
evaluating T[2] at t where T[1] = λ te.
A parametrization of the series S. It is the simplest one which takes the
denominator out of S and makes it the exponent of the first entry in the
parametrization.
This example illustrates the use of DuvalPuiseuxExpansion
and ParametrizationToPuiseux to gain the information given by
PuiseuxExpansion and also compares the performance of the two algorithms.
It also highlights some of the anomalies that may be encountered due
to precision concerns.
> P<x> := PuiseuxSeriesRing(Rationals());
> R<y> := PolynomialRing(P);
> f := (y^2 - x^3)^2 - y*x^6;
> time D := DuvalPuiseuxExpansion(f, 0);
Time: 0.000
> D;
[
<16*x^4, 64*x^6 + 256*x^9 + O(x^10)>
]
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.060
> A<a> := Parent(P[1]);
> N<n> := CoefficientRing(A);
> P;
[
a^(3/2) + 1/2*a^(9/4) + O(a^(5/2)),
a^(3/2) - 1/2*a^(9/4) + O(a^(5/2)),
-a^(3/2) + 1/2*n*a^(9/4) + O(a^(5/2)),
-a^(3/2) - 1/2*n*a^(9/4) + O(a^(5/2))
]
> time c := PuiseuxExpansion(f, 0);
Time: 0.050
Here it can be seen that the newton polygon part of the algorithm is
substantially faster using Duval's method, though the converting of
the parametrization to a series is not as fast. Asking for more terms
shows this more substantially.
> time D := DuvalPuiseuxExpansion(f, 10);
Time: 0.020
> D;
[
<16*x^4, 64*x^6 + 256*x^9 - 512*x^15 + 2048*x^18 - 4608*x^21 + 56320*x^27 -
294912*x^30 + 792064*x^33 - 12082176*x^39 + 68157440*x^42 + O(x^43)>
]
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.129
> time c := PuiseuxExpansion(f, 10);
Time: 0.100
> A<a> := Parent(c[1]);
> N<n> := CoefficientRing(A);
> c;
[
a^(3/2) + 1/2*a^(9/4) - 1/64*a^(15/4) + O(a^4),
a^(3/2) - 1/2*a^(9/4) + 1/64*a^(15/4) + O(a^4),
-a^(3/2) + 1/2*n*a^(9/4) + 1/64*n*a^(15/4) + O(a^4),
-a^(3/2) - 1/2*n*a^(9/4) - 1/64*n*a^(15/4) + O(a^4)
]
> A<a> := Parent(P[1]);
> N<n> := CoefficientRing(A);
> P;
[
a^(3/2) + 1/2*a^(9/4) - 1/64*a^(15/4) + 1/128*a^(9/2) - 9/4096*a^(21/4) +
55/131072*a^(27/4) - 9/32768*a^(15/2) + 1547/16777216*a^(33/4) -
11799/536870912*a^(39/4) + 65/4194304*a^(21/2) + O(a^(43/4)),
a^(3/2) - 1/2*a^(9/4) + 1/64*a^(15/4) + 1/128*a^(9/2) + 9/4096*a^(21/4) -
55/131072*a^(27/4) - 9/32768*a^(15/2) - 1547/16777216*a^(33/4) +
11799/536870912*a^(39/4) + 65/4194304*a^(21/2) + O(a^(43/4)),
-a^(3/2) + 1/2*n*a^(9/4) + 1/64*n*a^(15/4) - 1/128*a^(9/2) -
9/4096*n*a^(21/4) - 55/131072*n*a^(27/4) + 9/32768*a^(15/2) +
1547/16777216*n*a^(33/4) + 11799/536870912*n*a^(39/4) -
65/4194304*a^(21/2) + O(a^(43/4)),
-a^(3/2) - 1/2*n*a^(9/4) - 1/64*n*a^(15/4) - 1/128*a^(9/2) +
9/4096*n*a^(21/4) + 55/131072*n*a^(27/4) + 9/32768*a^(15/2) -
1547/16777216*n*a^(33/4) - 11799/536870912*n*a^(39/4) -
65/4194304*a^(21/2) + O(a^(43/4))
]
It can be seen that the computation of the information is a lot faster
using Duval's method. It is only the cosmetic of converting this
information into series that could make this algorithm seem slow. But
also note that there is much greater information given by Duval's algorithm.
The equivalent information is given below.
> time D := DuvalPuiseuxExpansion(f, 3);
Time: 0.009
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.049
> A<a> := Parent(P[1]);
> N<n> := CoefficientRing(A);
> P;
[
a^(3/2) + 1/2*a^(9/4) - 1/64*a^(15/4) + O(a^4),
a^(3/2) - 1/2*a^(9/4) + 1/64*a^(15/4) + O(a^4),
-a^(3/2) + 1/2*n*a^(9/4) + 1/64*n*a^(15/4) + O(a^4),
-a^(3/2) - 1/2*n*a^(9/4) - 1/64*n*a^(15/4) + O(a^4)
]
One thing that may be taken for granted from PuiseuxExpansion is
that all the expansions lie in the same puiseux series ring. However, for
DuvalPuiseuxExpansion this may not be the case. It will always be
true that each of the parametrizations will lie in the same puiseux series
ring but series resulting from different parametrizations may not. This occurs
since some extensions are left to the stage of calculating the series from
the parametrization to be made and for different parametrizations these
extensions may be different.
> f := (-x^3 + x^4) - 2*x^2*y - x*y^2 + 2*x*y^4 + y^5;
> time D := DuvalPuiseuxExpansion(f, 0);
Time: 0.010
> D;
[
<x^2, -x^2 - x^3 + O(x^4)>,
<x^3, x + O(x^2)>
]
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.000
> P;
[
-x - x^(3/2) + O(x^2),
-x + x^(3/2) + O(x^2)
]
> Parent(P[1]);
Puiseux series field in x over Rational Field
> time P := ParametrizationToPuiseux(D[2]);
Time: 0.030
> A<a> := Parent(P[1]);
> N<n> := CoefficientRing(A);
> P;
[
a^(1/3) + O(a^(2/3)),
n*a^(1/3) + O(a^(2/3)),
(-n - 1)*a^(1/3) + O(a^(2/3))
]
> Parent(P[1]);
Puiseux series field in a over N
> N;
Number Field with defining polynomial $.1^2 + $.1 + 1 over the Rational Field
DuvalPuiseuxExpansion reacts differently to PuiseuxExpansion
when given input which has finite expansions either due to finite precision
or exact roots. These differences are shown below and are due to the fact that
DuvalPuiseuxExpansion always looks for the next non zero term in
an expansion whereas PuiseuxExpansion will calculate zero terms.
> f := y - x^3 - x^7 - x^76 + O(x^200);
> D := DuvalPuiseuxExpansion(f, 0);
> D;
[
<x, x^3 + O(x^4)>
]
> D := DuvalPuiseuxExpansion(f, 3);
> D;
[
<x, x^3 + x^7 + x^76 + O(x^77)>
]
> D := DuvalPuiseuxExpansion(f, 4);
>> D := DuvalPuiseuxExpansion(f, 4);
^
Runtime error in 'DuvalPuiseuxExpansion': Insufficient precision to calculate to
requested precision
> c := PuiseuxExpansion(f, 197);
> c;
[
x^3 + x^7 + x^76 + O(x^200)
]
> c := PuiseuxExpansion(f, 200);
>> c := PuiseuxExpansion(f, 200);
^
Runtime error in 'PuiseuxExpansion': Insufficient precision to calculate to
requested precision
> f := y - x^3 - x^7 - x^76;
> D := DuvalPuiseuxExpansion(f, 0);
> D;
[
<x, x^3 + O(x^4)>
]
> D := DuvalPuiseuxExpansion(f, 3);
> D;
[
<x, x^3 + x^7 + x^76 + O(x^77)>
]
> D := DuvalPuiseuxExpansion(f, 4);
> D;
[
<x, x^3 + x^7 + x^76>
]
> c := PuiseuxExpansion(f, 10);
> c;
[
x^3 + x^7 + O(x^13)
]
> c := PuiseuxExpansion(f, 100);
> c;
[
x^3 + x^7 + x^76 + O(x^103)
]
> c := PuiseuxExpansion(f, 200);
> c;
[
x^3 + x^7 + x^76 + O(x^203)
]
> c := PuiseuxExpansion(f, 200 : PreciseRoot := true);
> c;
[
x^3 + x^7 + x^76
]
The two methods can be combined. Given a series there is no way that Duval's
REGULAR algorithm can be used to further the precision of an expansion.
But ExpandToPrecision can be used to gain the extra precision. Using
the REGULAR algorithm would be preferable since it is faster but this is not
possible for the type of input that is available. This is explored in the
case of our first example.
> f := (y^2 - x^3)^2 - y*x^6;
> time D := DuvalPuiseuxExpansion(f, 0);
Time: 0.009
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.050
> A<a> := Parent(P[1]);
> N<n> := CoefficientRing(A);
> P;
[
a^(3/2) + 1/2*a^(9/4) + O(a^(5/2)),
a^(3/2) - 1/2*a^(9/4) + O(a^(5/2)),
-a^(3/2) + 1/2*n*a^(9/4) + O(a^(5/2)),
-a^(3/2) - 1/2*n*a^(9/4) + O(a^(5/2))
]
> time ExpandToPrecision(f, P[1], 20);
a^(3/2) + 1/2*a^(9/4) - 1/64*a^(15/4) + 1/128*a^(9/2) - 9/4096*a^(21/4) +
O(a^(13/2))
Time: 0.070
> time D := DuvalPuiseuxExpansion(f, 5);
Time: 0.010
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.059
It can be seen that using ExpandToPrecision is slower even
than rerunning Duval's algorithm from the beginning. Even more so when it is
remembered that running Duval's algorithm with the extra precision will give
parametrizations of all the expansions of the roots of f and not just one
expansion. This seems to still be the case when larger examples are
considered so that if more terms of an expansion are required it is probably
best to start from the beginning asking for these extra terms.
> time p1 := ExpandToPrecision(f, P[1], 50);
Time: 0.439
> A<a> := Parent(p1);
> N<n> := CoefficientRing(A);
> p1;
a^(3/2) + 1/2*a^(9/4) - 1/64*a^(15/4) + 1/128*a^(9/2) - 9/4096*a^(21/4) +
55/131072*a^(27/4) - 9/32768*a^(15/2) + 1547/16777216*a^(33/4) -
11799/536870912*a^(39/4) + 65/4194304*a^(21/2) - 189805/34359738368*a^(45/4)
+ 1584999/1099511627776*a^(51/4) - 2261/2147483648*a^(27/2) + O(a^14)
> time D := DuvalPuiseuxExpansion(f, 13);
Time: 0.009
> time P := ParametrizationToPuiseux(D[1]);
Time: 0.200
This section describes two similar functions that can be used for finding
roots of polynomials over series rings in a similar way to finding roots
of polynomials over any other ring for which roots can be computed in Magma.
Roots(f, n) : RngUPolElt, RngIntElt -> [<RngSerElt, RngIntElt>]
SetVerbose("Newton", n): Maximum: 2
Find the roots of the polynomial f which lie in the coefficient ring of f. The first
form of this function can be used on any polynomial over any ring for which
Magma can compute roots. Since precision is an issue in series rings and
some roots may have infinite expansions, the second version of this function
which is specific to series rings allows a lower bound for
the precision to which these
roots will be known to be specified.
The first computes the roots to at least the default precision of the
ring if that ring
has infinite precision otherwise it computes them to at least the precision of the ring.
The first version will be enough when all the coefficients of the polynomial
have infinite precision. The second version may
be required when a precision other than that assumed by the first is
sought which may be due to the impossibility of computing the roots
to such a high precision as the default.
The precision is relative to the least common multiple of the exponent
denominators of the coefficients of f and the exponent denominator
of the root.
Roots which are known to be different but are identical to the precision
specified will be returned as two distinct roots.
Duval's algorithm as implemented in DuvalPuiseuxExpansion will usually
be used. Walker's algorithm as implemented in PuiseuxExpansion will
be used if the polynomial has coefficients involving fractional powers
or the characteristic of the coefficient ring of the series ring is
less than the degree of a squarefree factor.
If Walker's algorithm is used and the characteristic of the field is less
than the degree of the polynomial then the computation may not finish
(see remarks under PuiseuxExpansion) and control will return
to the user when interrupted.
If verbose printing of partial output that doesn't have enough precision is
required the functions PuiseuxExpansion and DuvalPuiseuxExpansion
should be used with the appropriate precision. Roots also requires that
there is enough precision in the roots so that the multiplicities can be
calculated correctly and the parts of the roots that are returned are distinct.
Therefore an error will be given if there is not enough precision to calculate
the part of the root that results from the newton polygon part of the
algorithm and the root is not known to be a single root since the
multiplicity may not be able to be calculated correctly.
Information at such a low precision can be gained correctly by using
the PuiseuxExpansion functions and determining the multiplicities
manually.
Return true if the polynomial f has a root in its coefficient ring and that root
can be found to the fixed or default precision of the ring as applicable.
A root is also returned in this case. If f is irreducible over its
coefficient ring then return false.
Below are some examples of the use of the Roots function.
> SetVerbose("Newton", 1);
> P<x> := PuiseuxSeriesRing(Rationals());
> R<y> := PolynomialRing(P);
> f := y^3 + 2*x^-1*y^2 + 1*x^-2*y + 2*x;
> Roots(f);
DUVAL :
[
<-2*x^3 - 8*x^7 - 56*x^11 - 480*x^15 - 4576*x^19 + O(x^23), 1>
]
> f := f^2;
> Roots(f);
DUVAL :
[
<-2*x^3 - 8*x^7 - 56*x^11 - 480*x^15 - 4576*x^19 + O(x^23), 2>
]
> f := y^3 + 2*x^-1*y^2 + 1*x^-2*y + 2*x;
> f +:= O(x^20)*(y^3 + y^2 + y + 1);
> f;
(1 + O(x^20))*y^3 + (2*x^-1 + O(x^20))*y^2 + (x^-2 + O(x^20))*y + 2*x + O(x^20)
> Roots(f);
DUVAL :
>> Roots(f);
^
Runtime error in 'Roots': Roots not calculable to default precision
> Roots(f, 10);
DUVAL :
[
<-2*x^3 - 8*x^7 - 56*x^11 + O(x^13), 1>
]
> f := f^2;
> f;
(1 + O(x^20))*y^6 + (4*x^-1 + O(x^19))*y^5 + (6*x^-2 + O(x^18))*y^4 + (4*x^-3 +
4*x + O(x^18))*y^3 + (x^-4 + 8 + O(x^18))*y^2 + (4*x^-1 + O(x^18))*y + 4*x^2
+ O(x^21)
> Roots(f, 10);
DUVAL :
[
<-2*x^3 - 8*x^7 - 56*x^11 + O(x^13), 2>
]
> f := (y - x^(1/4))*(y - x^(1/3));
> Roots(f);
WALKER :
[
<x^(1/3) + O(x^2), 1>,
<x^(1/4) + O(x^(23/12)), 1>
]
[Next][Prev] [Right] [Left] [Up] [Index] [Root]
|