Padé-Hermite Approximants

-Hermite approximants}

Contents

Introduction

A given rational function F(z)∈k(z) over a field k can be written as a power series f(z) in the completion k((z)) of k(z) at the place (z). It thus is a good approximation of f(z) in the sense that F(z) is equal to f(z) up to infinite order.

Padé-Hermite approximants deal with the converse. Given a (formal) power series f(z)∈k[[z]] and some non- negative integers nP and nQ, find polynomials P(z), Q(z)∈k[z] of degrees at most nP and nQ, respectively, such that P(z) * f(z) - Q(z) = O(znP + nQ + 1) holds. In other words, P/Q(z) is an approximation for f(z) for polynomials P and Q of limited degree in z. In this notation the pair [P, Q] is known as the Padé-Hermite approximant for the power series tuple (f, - 1).

The idea of approximating one power series by two polynomials can be extended to approximating several power series at the same time as follows. Consider the vector underline(f)T:= (f1, f2, ..., fm)T in the m-dimensional vector space k((z))m over a power series ring. Let underline(n)=(n1, n2, ..., nm) be an m-tuple of non-negative integers. A Padé-Hermite approximant of underline(f)T(z) of type underline(n) is a non-zero vector of polynomials underline(P)=(P1, P2, ..., Pm) in k[z]m such that underline(P).underline(f)T = O(zN) , N=n1 + n2 + ... + nm + m - 1, is satisfied. A non-trivial Padé-Hermite approximant always exists. The approximant is contained in the sub-space V_(underline(f), N)={ underline(Q)∈k[z]m : underline(Q).underline(f)T = O(zN)}.

The implementation of Padé-Hermite approximants is based on [Der94] and [BL94]. They are implemented in Magma as sequences rather than vectors.

Ordering of Sequences

A Padé-Hermite approximant to some sequence of power series does not have to be unique. They can be ordered according to their maximum degree and their type of sequence. A sequence underline(P)=[P1, P2, ..., Pm] is of degree i if the maximum of the weak degrees of P1, P2, ..., Pm is i. An extension of the definition maximum degree is when a distortion of non-negative integers underline(d)=[d1, d2, ..., dm] on the degrees is allowed. In this case the maximum degree is defined as the maximum of deg(Pi) - di. The type of underline(P) identifies the last Pi in underline(P) whose weak degree equals the maximum degree of underline(P). Again a distortion on the degrees is allowed.

The Padé-Hermite approximant of underline(f) can be seen as an element of the free module k[z]m of rank m, or as an element of the sub-space V_(underline(f), N).

A free sub-model V⊂k[z]m is generated by m polynomial vectors underline(Q)i, 1=1, 2, ..., m, such that underline(Q)1(z), underline(Q)2(z), ..., underline(Q)m(z) form a minimal vector sequence of V. Such a sequence is defined as a sequence S of m vectors in V, such that S[i] is a non-trivial polynomial vector in V of minimal degree of type i, for i=1, 2, ..., m. A minimal vector sequence is not unique.

Two variations on a minimal vector sequence are implemented. The first allows a distortion as an attribute. The sequence is the based on the distorted maximal degree. The second attribute sets the positive power p of z in each of the underline(Q)i as follows. Instead of considering underline(Q)i(z) for each i, one considers underline(Q)i(zp).

MaximumDegree(f) : SeqEnum -> RngIntElt
    Distortion: SeqEnum                 Default: []
MaximumDegree returns the degree of a sequence of polynomials or power series, defined as the maximum of the degrees of f[i] - d[i], where d is the distortion. The value -infinity is returned in the case that f is weakly equal to the zero-sequence.

Example FldFunRat_degree-of-sequence (H44E5)

> S<u> := PowerSeriesRing(Rationals());
> f := [u+u^2, 2+u^2+u^3,0];
> MaximumDegree(f);
3
> MaximumDegree(f:Distortion:=[]);
3
> MaximumDegree(f:Distortion:=[0,2,1]);
2
> MaximumDegree([S|0,0]);
-Infinity
> MaximumDegree([O(u)]);
-Infinity
TypeOfSequence(f) : SeqEnum -> RngIntElt, RngIntElt
    Distortion: SeqEnum                 Default: []
Returns the highest index i for those f[i] whose (distorted) degree is weakly equal to the maximum of the degrees of all entries. The second integer returned is the maximum degree of the sequence.

Example FldFunRat_type-of-sequence (H44E6)

> S<u> := PowerSeriesRing(Rationals());
> f := [u+u^2, 2+u^2+u^3,0];
> TypeOfSequence(f);
2 3
> TypeOfSequence(f:Distortion:=[]);
2 3
> TypeOfSequence(f:Distortion:=[0,2,1]);
1 2
> TypeOfSequence([S|0,0]);
2 -Infinity
MinimalVectorSequence(f,n) : SeqEnum, RngIntElt -> SeqEnum
    Distortion: SeqEnum                 Default: []
    Power: RngIntElt                    Default: 1
A minimal sequence of vectors underline(Q)1, underline(Q)2, ..., underline(Q)m with respect to the sequence f of length m whose entries are polynomials or power series. The order of underline(Q)i.f is at least m.

Example FldFunRat_minimal-vector-sequence (H44E7)

> S<u> := PowerSeriesRing(Rationals());
> f := [u+u^2, 2+u^2+u^3];
> seq := MinimalVectorSequence(f, 2);
> seq;
[
    (u 0),
    (   -1 1/2*u)
]
> sums := [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
> sums;
[
    u^2 + u^3,
    -u^2 + 1/2*u^3 + 1/2*u^4
]
> seq := MinimalVectorSequence(f,3);
> #seq eq 2,  seq[1], seq[2];
true (u^2   0)
(-1 + u  1/2*u)

Example FldFunRat_the-next_example (H44E8)

> L<x> := PolynomialRing(Rationals());
> f := [1+x, 3-x^2, 5+x+x^3-x^5];
> seq := MinimalVectorSequence(f, 10);
> seq[1];
(-2*x^2 + 6   -2*x - 2          0)
> sums := [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
> sums;
[
    0,
    0,
    x^10
]

Example FldFunRat_another-example (H44E9)

> S<u> := PowerSeriesRing(Rationals());
> f := [2*u^4,2+u^3+u^6];
> seq := MinimalVectorSequence(f, 10);
> seq;
[
    (1/2*u^6       0),
    (-1/2 - 1/4*u^3        1/2*u^4)
]
> sums := [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
> sums;
[
    u^10,
    1/2*u^10
]

Example FldFunRat_one-more (H44E10)

> S<u> := PowerSeriesRing(Rationals());
> f :=  [1+u-7*u^2, 6-3*u+1/2*u^2-u^3, 5-u+u^2];
> seq := MinimalVectorSequence(f, 5);
>  [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
[
    0,
    u^5,
    0
]
> seq := MinimalVectorSequence(f, 5:Distortion:=[2,0,1]);
> [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
[
    -7/2*u^5,
    u^5,
    0
]
> p:=2;
> seq := MinimalVectorSequence(f, 5:Distortion:=[2,0,1], Power:=p);
> sums := [&+([Q[i]*f[i]: i in [1..#f]]) : Q in seq];
> sums;
> mp:= map<S->S|  x :-> (IsWeaklyZero(x) select 0
>     else  &+([Coefficient(x,i)*(S.1)^(p*i) : i in Exponents(x)]))
>       + (ISA(Type(v),RngIntElt) select O((S.1)^(p*v))
>   else S!0 where v := AbsolutePrecision(x))>;
> sums := [&+([mp(Q[i])*f[i]: i in [1..#f]]) : Q in seq];
> sums;
[
    u^6 + u^7 - 7*u^8,
    u^5 - 23/3*u^6,
    -5/3*u^5 + 35/3*u^6
]

Approximants

The Padé-Hermite approximant of type underline(d)=[d1, d2, ..., dm] with respect to the tuple underline(f)T∈k((z))m is an element of the space V_(underline(f), N)={ underline(Q)∈k[z]m : underline(Q).underline(f)T = O(zN)} for N equal to d1 + d2 + ... + dm + m - 1. This space is generated by the vectors in the minimal vector sequence with respect to underline(f) with distortion underline(d). The routine PadeHermiteApproximant returns one that is smallest with respect to the degree on sequences. The input sequence underline(f) must be a sequence of polynomial ring elements, or be a power series sequence. While the Padé Hermite approximants theoretically are polynomials, Magma returns them as elements of the same ring the entries of underline(f) are contained in.

A variant of the Padé-Hermite approximant is when the exponent in the order term is set rather that the type of the sequence. It is also possible to let underline(f) be a sequence such that its entries themselves are vectors of polynomials or power series.

PadeHermiteApproximant(f,d) : SeqEnum, SeqEnum -> ModTupRngElt, SeqEnum, RngIntElt
    Power: RngIntElt                    Default: 1
Returns a Padé-Hermite form underline(P) of f with distortion d, smallest with respect to the degree on sequences, and the corresponding minimal vector sequence. The third argument returned is the order in the order term of underline(P).f.

Example FldFunRat_pade-hermite-approximants (H44E11)

This example can be found on page 813 in [BL94].
> S<u> := PowerSeriesRing(Rationals());
> f := [1,u,u/(1-u^4)+u^10+O(u^16),u/(1+u^4)+u^12+O(u^16)];
> pade, padebasis, ord := PadeHermiteApproximant(f,[2,2,2,2]);
> pade, ord;
( u -1  0  0)
11
> BaseRing(Parent(pade)) eq S;
true
> MinimalVectorSequence(f,10);
[
    ( u -1  0  0),
    (   0  u^4 -1/2  1/2),
    (         0    1 - u^4 -1/2 + u^4       -1/2),
    (    0    -u 1/2*u 1/2*u)
]
>
> p := 2;
> seq := MinimalVectorSequence(f,10: Distortion :=[2,2,2,2],Power := p);
> seq;
[
    (u^5   0   0   0),
    (   0  u^2 -1/2  1/2),
    (         0    1 - u^2 -1/2 + u^2       -1/2),
    (    0    -u 1/2*u 1/2*u)
]
> mp:= map<S->S|  x :-> (IsWeaklyZero(x) select 0
>     else  &+([Coefficient(x,i)*(S.1)^(p*i) : i in Exponents(x)]))
>       + (ISA(Type(v),RngIntElt) select O((S.1)^(p*v))
>   else S!0 where v := AbsolutePrecision(x))>;
> sums := [&+([mp(Q[i])*f[i]: i in [1..#f]]) : Q in seq];
> [Valuation(v) : v in sums];
[ 10, 10, 10, 11 ]
> sums;
[
    u^10,
    -1/2*u^10 + 1/2*u^12 - u^13 + O(u^16),
    -1/2*u^10 - 1/2*u^12 + u^13 + u^14 + O(u^16),
    u^11 + 1/2*u^12 + 1/2*u^14 + O(u^18)
]

Example FldFunRat_last-example (H44E12)

This example covers the example on page 815 in [BL94].
> S<u> := PowerSeriesRing(Rationals());
> f := [1,u, -1-u^4-2*u^8+u^10+u^11-u^12+O(u^16),-u-u^5-u^9-u^14-u^15+O(u^16)];
> dist:=[2,2,3,3];
> seq := MinimalVectorSequence(f,13:Distortion:=dist);
> pade, padebasis, ord :=  PadeHermiteApproximant(f,dist);
> pade, ord;
(-u  1  0  0)
13
> padebasis;
[
    (-u  1  0  0),
    (1/2*u -1/2*u^4 1/2*u + 1/2*u^3 + 1/2*u^4 -1/2*u^2 - 1/2*u^3 - u^4),
    (-1/2 1/2*u^3 -1/2 - 1/2*u^2 - 1/2*u^3 - u^4 1/2*u + 1/2*u^2 + 2*u^3),
    (      -u        0        0 -1 + u^4)
]
> padebasis eq seq;
true
> [[Valuation(w[i]): i in [1..Degree(w)]] : w in seq];
[
    [ 1, 0, Infinity, Infinity ],
    [ 1, 4, 1, 2 ],
    [ 0, 3, 0, 1 ],
    [ 1, Infinity, Infinity, 0 ]
]
> [[MaximumDegree([w[i]])-dist[i]: i in [1..Degree(w)] ] : w in seq];
[
    [ -1, -2, -Infinity, -Infinity ],
    [ -1, 2, 1, 1 ],
    [ -2, 1, 1, 0 ],
    [ -1, -Infinity, -Infinity, 1 ]
]
> p:=2;
> seq := MinimalVectorSequence(f,12:Distortion:=dist,Power:=p);
> seq;
[
    (  u - u^3         0 u - 2*u^3         0),
    (-1 - u + u^2 -u^3 -1 - u + 2*u^2 + u^3 -u^3),
    (u + u^2 - u^3 0 u + u^2 - 2*u^3 - u^4 0),
    (      0       1       0 1 - u^2)
]
> seq[1]-seq[3];
(      -u^2          0 -u^2 + u^4          0)
> mp:= map<S->S|  x :-> (IsWeaklyZero(x) select 0
>     else  &+([Coefficient(x,i)*(S.1)^(p*i) : i in Exponents(x)]))
>       + (ISA(Type(v),RngIntElt) select O((S.1)^(p*v))
>   else S!0 where v := AbsolutePrecision(x))>;
>  [Valuation(&+([mp(Q[i])*f[i]: i in [1..#f]])) : Q in seq];
[ 12, 12, 13, 13 ]

Example FldFunRat_ (H44E13)

This example considers the example on page 816 in [BL94].
> S<u> := PowerSeriesRing(Rationals());
> f := [1,u,-1-u^4-2*u^8+u^10+O(u^12),-u-u^5-u^9+u^10+O(u^12)];
> dist := [2,2,3,3];
> seq := MinimalVectorSequence(f,12: Distortion := dist);
> [[MaximumDegree([w[i]])-dist[i]: i in [1..Degree(w)] ] : w in seq];
[
    [ -1, -2, -Infinity, -Infinity ],
    [ -2, 1, 0, 0 ],
    [ -Infinity, -Infinity, 1, 0 ],
    [ -1, -Infinity, 0, 1 ]
]
> [Valuation(&+([(Q[i])*f[i]: i in [1..#f]]) ) : Q in seq];
[ Infinity, 12, 12, 12 ]
> [MaximumDegree([ &+([(Q[i])*f[i]: i in [1..#f]]) ]) : Q in seq];
[ -Infinity, -Infinity, 14, -Infinity ]
> PadeHermiteApproximant(f,[2,2,3,3]);
> p := 2;
> seq := MinimalVectorSequence(f,12:Distortion:=[2,2,3,3],Power:=p);
> seq;
[
    (  1 - u^2        -1 1 - 2*u^2  -1 + u^2),
    (   0 -u^4    0 -u^4),
    (      -u       -1 -u + u^3 -1 + u^2),
    (      0       u       0 u - u^3)
]
> [[MaximumDegree([w[i]])-dist[i]: i in [1..Degree(w)] ] : w in seq];
[
    [ 0, -2, -1, -1 ],
    [ -Infinity, 2, -Infinity, 1 ],
    [ -1, -2, 0, -1 ],
    [ -Infinity, -1, -Infinity, 0 ]
]
> mp:= map<S->S|  x :-> (IsWeaklyZero(x) select 0
>     else  &+([Coefficient(x,i)*(S.1)^(p*i) : i in Exponents(x)]))
>       + (ISA(Type(v),RngIntElt) select O((S.1)^(p*v))
>   else S!0 where v := AbsolutePrecision(x))>;
>
> [Valuation(&+([mp(Q[i])*f[i]: i in [1..#f]])) : Q in seq];
[ 12, 13, 12, 12 ]
PadeHermiteApproximant(f,m) : SeqEnum, RngIntElt -> ModTupRngElt, SeqEnum
    Power: RngIntElt                    Default: 1
Returns a Padé-Hermite form of minimal degree in the corresponding minimal vector sequence, such that its inproduct with f has order at least m. The second argument returned is the corresponding minimal vector sequence.

Example FldFunRat_pade-hermite-approximants-vectors (H44E14)

This example can be found on page 813 in [BL94].
> S<u> := PowerSeriesRing(Rationals());
> f := [Vector([1]), Vector([u])];
> pade, seq := PadeHermiteApproximant(f,3);
Calculating the Pade'-Hermite approximant for the sequence [
    1,
    u
]
with order term 3 and power 1 .
> pade;
( u -1)
> seq;
[
    ( u -1),
    (  0 u^2)
]
> mat := Matrix([Eltseq(v): v in f]);
> pade*mat;
(0)
> PadeHermiteApproximant([1,u],5);
( u -1)
[
    ( u -1),
    (  0 u^4)
]
> PadeHermiteApproximant(f,3:Power:=2);
> g:= [Vector([1,0,0]), Vector([0,1,0]), Vector([1+u,2+u^2,u^3])];
> pade := PadeHermiteApproximant(g,5);
Calculating the Pade'-Hermite approximant for the sequence [
    1,
    u,
    1 + 2*u + u^3 + u^7 + u^11
]
with order term 15 and power 3 .
> pade;
(-u^3 - u^4     -2*u^3        u^3)
> pade*Matrix([Eltseq(v): v in g]);
(  0 u^5 u^6)

Example FldFunRat_ (H44E15)

This example considers Padé-Hermite approximants for some series that have non-trivial power series expansions.
> S<u> := PowerSeriesRing(Rationals());
> f := [Sin(u), Cos(u), Exp(u)];
> [Valuation(f[i]) : i in [1..#f]], [Degree(f[i]) : i in [1..#f]];
[ 1, 0, 0 ]
[ 19, 20, 20 ]
> [AbsolutePrecision(f[i]) : i in [1..#f]];
[ 21, 22, 21 ]
> dist := [3,2,5];
> pade, seq, ord := PadeHermiteApproximant(f,dist);
> 1/420*pade;
(-1275 - 255*u + 45*u^2 + 5*u^3 120 + 495*u + 75*u^2 -120 + 900*u - 600*u^2 +
    160*u^3 - 20*u^4 + u^5)
> ord eq &+(dist)+#f-1, ord;
true 12
> [Degree(pade[i]) : i in [1..Degree(pade)]];
[ 3, 2, 5 ]
> g:= [Cos(2*u)*(u+1)+3,Cos(u)^2+u*Cos(u)+1,Cos(2*u)+1,Cos(u)];
> pade, basis := PadeHermiteApproximant(g,20);
> 131/75880*pade;
(          2    -4 + 2*u        -3*u 4*u - 2*u^2)
> h := [ 1+u^2-u^7+u^12, Sin(u), Exp(u) ];
> dist:=[3,1,2];
> seq := MinimalVectorSequence(h,8:Distortion := dist);
> sums := [&+([Q[i]*h[i]: i in [1..#f]]) : Q in seq];
> [Valuation(s) : s in sums];
[ 8, 8, 8 ]
> [[MaximumDegree([w[i]]): i in [1..Degree(w)] ] : w in seq];
[
    [ 4, 1, 2 ],
    [ 4, 2, 2 ],
    [ 3, 1, 2 ]
]
>  [[MaximumDegree([w[i]])-dist[i]: i in [1..Degree(w)] ] : w in seq];
[
    [ 1, 0, 0 ],
    [ 1, 1, 0 ],
    [ 0, 0, 0 ]
]
V2.28, 13 July 2023