Asymptotics of Multivariate Generating Series¶
Let \(F(x) = \sum_{\nu \in \NN^d} F_{\nu} x^\nu\) be a multivariate power series with complex coefficients that converges in a neighborhood of the origin. Assume that \(F = G/H\) for some functions \(G\) and \(H\) holomorphic in a neighborhood of the origin. Assume also that \(H\) is a polynomial.
This computes asymptotics for the coefficients \(F_{r \alpha}\) as \(r \to \infty\) with \(r \alpha \in \NN^d\) for \(\alpha\) in a permissible subset of \(d\)-tuples of positive reals. More specifically, it computes arbitrary terms of the asymptotic expansion for \(F_{r \alpha}\) when the asymptotics are controlled by a strictly minimal multiple point of the algebraic variety \(H = 0\).
The algorithms and formulas implemented here come from [RW2008] and [RW2012]. For a general reference take a look in the book [PW2013].
Introductory Examples¶
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
>>> from sage.all import *
>>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
A univariate smooth point example:
sage: R.<x> = PolynomialRing(QQ)
sage: FFPD = FractionWithFactoredDenominatorRing(R, SR)
sage: H = (x - 1/2)^3
sage: Hfac = H.factor()
sage: G = -1/(x + 3)/Hfac.unit()
sage: F = FFPD(G, Hfac)
sage: F
(-1/(x + 3), [(x - 1/2, 3)])
sage: alpha = [1]
sage: decomp = F.asymptotic_decomposition(alpha)
sage: decomp
(0, []) +
(-1/2*r^2*(x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 6*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 9/(x^5 + 9*x^4 + 27*x^3 + 27*x^2))
- 1/2*r*(5*x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 24*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 27/(x^5 + 9*x^4 + 27*x^3 + 27*x^2))
- 3*x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
- 9*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
- 9/(x^5 + 9*x^4 + 27*x^3 + 27*x^2),
[(x - 1/2, 1)])
sage: F1 = decomp[1]
sage: p = {x: 1/2}
sage: asy = F1.asymptotics(p, alpha, 3)
sage: asy
(8/343*(49*r^2 + 161*r + 114)*2^r, 2, 8/7*r^2 + 184/49*r + 912/343)
sage: F.relative_error(asy[0], alpha, [1, 2, 4, 8, 16], asy[1])
[((1,), 7.555555556, [7.556851312], [-0.0001714971672]),
((2,), 14.74074074, [14.74052478], [0.00001465051901]),
((4,), 35.96502058, [35.96501458], [1.667911934e-7]),
((8,), 105.8425656, [105.8425656], [4.399565380e-11]),
((16,), 355.3119534, [355.3119534], [0.0000000000])]
>>> from sage.all import *
>>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1)
>>> FFPD = FractionWithFactoredDenominatorRing(R, SR)
>>> H = (x - Integer(1)/Integer(2))**Integer(3)
>>> Hfac = H.factor()
>>> G = -Integer(1)/(x + Integer(3))/Hfac.unit()
>>> F = FFPD(G, Hfac)
>>> F
(-1/(x + 3), [(x - 1/2, 3)])
>>> alpha = [Integer(1)]
>>> decomp = F.asymptotic_decomposition(alpha)
>>> decomp
(0, []) +
(-1/2*r^2*(x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 6*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 9/(x^5 + 9*x^4 + 27*x^3 + 27*x^2))
- 1/2*r*(5*x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 24*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
+ 27/(x^5 + 9*x^4 + 27*x^3 + 27*x^2))
- 3*x^2/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
- 9*x/(x^5 + 9*x^4 + 27*x^3 + 27*x^2)
- 9/(x^5 + 9*x^4 + 27*x^3 + 27*x^2),
[(x - 1/2, 1)])
>>> F1 = decomp[Integer(1)]
>>> p = {x: Integer(1)/Integer(2)}
>>> asy = F1.asymptotics(p, alpha, Integer(3))
>>> asy
(8/343*(49*r^2 + 161*r + 114)*2^r, 2, 8/7*r^2 + 184/49*r + 912/343)
>>> F.relative_error(asy[Integer(0)], alpha, [Integer(1), Integer(2), Integer(4), Integer(8), Integer(16)], asy[Integer(1)])
[((1,), 7.555555556, [7.556851312], [-0.0001714971672]),
((2,), 14.74074074, [14.74052478], [0.00001465051901]),
((4,), 35.96502058, [35.96501458], [1.667911934e-7]),
((8,), 105.8425656, [105.8425656], [4.399565380e-11]),
((16,), 355.3119534, [355.3119534], [0.0000000000])]
Another smooth point example (Example 5.4 of [RW2008]):
sage: R.<x,y> = PolynomialRing(QQ)
sage: FFPD = FractionWithFactoredDenominatorRing(R)
sage: q = 1/2
sage: qq = q.denominator()
sage: H = 1 - q*x + q*x*y - x^2*y
sage: Hfac = H.factor()
sage: G = (1 - q*x)/Hfac.unit()
sage: F = FFPD(G, Hfac)
sage: alpha = list(qq*vector([2, 1 - q]))
sage: alpha
[4, 1]
sage: I = F.smooth_critical_ideal(alpha)
sage: I
Ideal (y^2 - 2*y + 1, x + 1/4*y - 5/4) of
Multivariate Polynomial Ring in x, y over Rational Field
sage: s = solve([SR(z) for z in I.gens()],
....: [SR(z) for z in R.gens()], solution_dict=true)
sage: s == [{SR(x): 1, SR(y): 1}]
True
sage: p = s[0]
sage: asy = F.asymptotics(p, alpha, 1, verbose=True)
Creating auxiliary functions...
Computing derivatives of auxiliary functions...
Computing derivatives of more auxiliary functions...
Computing second order differential operator actions...
sage: asy
(1/24*2^(2/3)*(sqrt(3) + 4/(sqrt(3) + I) + I)*gamma(1/3)/(pi*r^(1/3)),
1,
1/24*2^(2/3)*(sqrt(3) + 4/(sqrt(3) + I) + I)*gamma(1/3)/(pi*r^(1/3)))
sage: r = SR('r')
sage: tuple((a*r^(1/3)).full_simplify() / r^(1/3) for a in asy) # make nicer coefficients
(1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)),
1,
1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)))
sage: F.relative_error(asy[0], alpha, [1, 2, 4, 8, 16], asy[1])
[((4, 1), 0.1875000000, [0.1953794675...], [-0.042023826...]),
((8, 2), 0.1523437500, [0.1550727862...], [-0.017913673...]),
((16, 4), 0.1221771240, [0.1230813519...], [-0.0074009592...]),
((32, 8), 0.09739671811, [0.09768973377...], [-0.0030084757...]),
((64, 16), 0.07744253816, [0.07753639308...], [-0.0012119297...])]
>>> from sage.all import *
>>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2)
>>> FFPD = FractionWithFactoredDenominatorRing(R)
>>> q = Integer(1)/Integer(2)
>>> qq = q.denominator()
>>> H = Integer(1) - q*x + q*x*y - x**Integer(2)*y
>>> Hfac = H.factor()
>>> G = (Integer(1) - q*x)/Hfac.unit()
>>> F = FFPD(G, Hfac)
>>> alpha = list(qq*vector([Integer(2), Integer(1) - q]))
>>> alpha
[4, 1]
>>> I = F.smooth_critical_ideal(alpha)
>>> I
Ideal (y^2 - 2*y + 1, x + 1/4*y - 5/4) of
Multivariate Polynomial Ring in x, y over Rational Field
>>> s = solve([SR(z) for z in I.gens()],
... [SR(z) for z in R.gens()], solution_dict=true)
>>> s == [{SR(x): Integer(1), SR(y): Integer(1)}]
True
>>> p = s[Integer(0)]
>>> asy = F.asymptotics(p, alpha, Integer(1), verbose=True)
Creating auxiliary functions...
Computing derivatives of auxiliary functions...
Computing derivatives of more auxiliary functions...
Computing second order differential operator actions...
>>> asy
(1/24*2^(2/3)*(sqrt(3) + 4/(sqrt(3) + I) + I)*gamma(1/3)/(pi*r^(1/3)),
1,
1/24*2^(2/3)*(sqrt(3) + 4/(sqrt(3) + I) + I)*gamma(1/3)/(pi*r^(1/3)))
>>> r = SR('r')
>>> tuple((a*r**(Integer(1)/Integer(3))).full_simplify() / r**(Integer(1)/Integer(3)) for a in asy) # make nicer coefficients
(1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)),
1,
1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)))
>>> F.relative_error(asy[Integer(0)], alpha, [Integer(1), Integer(2), Integer(4), Integer(8), Integer(16)], asy[Integer(1)])
[((4, 1), 0.1875000000, [0.1953794675...], [-0.042023826...]),
((8, 2), 0.1523437500, [0.1550727862...], [-0.017913673...]),
((16, 4), 0.1221771240, [0.1230813519...], [-0.0074009592...]),
((32, 8), 0.09739671811, [0.09768973377...], [-0.0030084757...]),
((64, 16), 0.07744253816, [0.07753639308...], [-0.0012119297...])]
A multiple point example (Example 6.5 of [RW2012]):
sage: R.<x,y> = PolynomialRing(QQ)
sage: FFPD = FractionWithFactoredDenominatorRing(R, SR)
sage: H = (1 - 2*x - y)**2 * (1 - x - 2*y)**2
sage: Hfac = H.factor()
sage: G = 1/Hfac.unit()
sage: F = FFPD(G, Hfac)
sage: F
(1, [(x + 2*y - 1, 2), (2*x + y - 1, 2)])
sage: I = F.singular_ideal()
sage: I
Ideal (x - 1/3, y - 1/3) of
Multivariate Polynomial Ring in x, y over Rational Field
sage: p = {x: 1/3, y: 1/3}
sage: F.is_convenient_multiple_point(p)
(True, 'convenient in variables [x, y]')
sage: alpha = (var('a'), var('b'))
sage: decomp = F.asymptotic_decomposition(alpha); decomp
(0, []) +
(-1/9*r^2*(2*a^2/x^2 + 2*b^2/y^2 - 5*a*b/(x*y))
- 1/9*r*(6*a/x^2 + 6*b/y^2 - 5*a/(x*y) - 5*b/(x*y))
- 4/9/x^2 - 4/9/y^2 + 5/9/(x*y),
[(x + 2*y - 1, 1), (2*x + y - 1, 1)])
sage: F1 = decomp[1]
sage: F1.asymptotics(p, alpha, 2)
(-3*((2*a^2 - 5*a*b + 2*b^2)*r^2 + (a + b)*r + 3)*(1/((1/3)^a*(1/3)^b))^r,
1/((1/3)^a*(1/3)^b), -3*(2*a^2 - 5*a*b + 2*b^2)*r^2 - 3*(a + b)*r - 9)
sage: alpha = [4, 3]
sage: decomp = F.asymptotic_decomposition(alpha)
sage: F1 = decomp[1]
sage: asy = F1.asymptotics(p, alpha, 2)
sage: asy
(3*(10*r^2 - 7*r - 3)*2187^r, 2187, 30*r^2 - 21*r - 9)
sage: F.relative_error(asy[0], alpha, [1, 2, 4, 8], asy[1])
[((4, 3), 30.72702332, [0.0000000000], [1.000000000]),
((8, 6), 111.9315678, [69.00000000], [0.3835519207]),
((16, 12), 442.7813138, [387.0000000], [0.1259793763]),
((32, 24), 1799.879232, [1743.000000], [0.03160169385])]
>>> from sage.all import *
>>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2)
>>> FFPD = FractionWithFactoredDenominatorRing(R, SR)
>>> H = (Integer(1) - Integer(2)*x - y)**Integer(2) * (Integer(1) - x - Integer(2)*y)**Integer(2)
>>> Hfac = H.factor()
>>> G = Integer(1)/Hfac.unit()
>>> F = FFPD(G, Hfac)
>>> F
(1, [(x + 2*y - 1, 2), (2*x + y - 1, 2)])
>>> I = F.singular_ideal()
>>> I
Ideal (x - 1/3, y - 1/3) of
Multivariate Polynomial Ring in x, y over Rational Field
>>> p = {x: Integer(1)/Integer(3), y: Integer(1)/Integer(3)}
>>> F.is_convenient_multiple_point(p)
(True, 'convenient in variables [x, y]')
>>> alpha = (var('a'), var('b'))
>>> decomp = F.asymptotic_decomposition(alpha); decomp
(0, []) +
(-1/9*r^2*(2*a^2/x^2 + 2*b^2/y^2 - 5*a*b/(x*y))
- 1/9*r*(6*a/x^2 + 6*b/y^2 - 5*a/(x*y) - 5*b/(x*y))
- 4/9/x^2 - 4/9/y^2 + 5/9/(x*y),
[(x + 2*y - 1, 1), (2*x + y - 1, 1)])
>>> F1 = decomp[Integer(1)]
>>> F1.asymptotics(p, alpha, Integer(2))
(-3*((2*a^2 - 5*a*b + 2*b^2)*r^2 + (a + b)*r + 3)*(1/((1/3)^a*(1/3)^b))^r,
1/((1/3)^a*(1/3)^b), -3*(2*a^2 - 5*a*b + 2*b^2)*r^2 - 3*(a + b)*r - 9)
>>> alpha = [Integer(4), Integer(3)]
>>> decomp = F.asymptotic_decomposition(alpha)
>>> F1 = decomp[Integer(1)]
>>> asy = F1.asymptotics(p, alpha, Integer(2))
>>> asy
(3*(10*r^2 - 7*r - 3)*2187^r, 2187, 30*r^2 - 21*r - 9)
>>> F.relative_error(asy[Integer(0)], alpha, [Integer(1), Integer(2), Integer(4), Integer(8)], asy[Integer(1)])
[((4, 3), 30.72702332, [0.0000000000], [1.000000000]),
((8, 6), 111.9315678, [69.00000000], [0.3835519207]),
((16, 12), 442.7813138, [387.0000000], [0.1259793763]),
((32, 24), 1799.879232, [1743.000000], [0.03160169385])]
Various¶
AUTHORS:
Alexander Raichev (2008)
Daniel Krenn (2014, 2016)
Classes and Methods¶
- class sage.rings.asymptotic.asymptotics_multivariate_generating_functions.FractionWithFactoredDenominator(parent, numerator, denominator_factored, reduce=True)[source]¶
Bases:
RingElementThis element represents a fraction with a factored polynomial denominator. See also its parent
FractionWithFactoredDenominatorRingfor details.Represents a fraction with factored polynomial denominator (FFPD) \(p/(q_1^{e_1} \cdots q_n^{e_n})\) by storing the parts \(p\) and \([(q_1, e_1), \ldots, (q_n, e_n)]\). Here \(q_1, \ldots, q_n\) are elements of a 0- or multi-variate factorial polynomial ring \(R\) , \(q_1, \ldots, q_n\) are distinct irreducible elements of \(R\) , \(e_1, \ldots, e_n\) are positive integers, and \(p\) is a function of the indeterminates of \(R\) (e.g., a Sage symbolic expression). An element \(r\) with no polynomial denominator is represented as
(r, []).INPUT:
numerator– an element \(p\); this can be of any ring from which parent’s base has coercion indenominator_factored– list of the form \([(q_1, e_1), \ldots, (q_n, e_n)]\), where the \(q_1, \ldots, q_n\) are distinct irreducible elements of \(R\) and the \(e_i\) are positive integersreduce– (optional) ifTrue, then represent \(p/(q_1^{e_1} \cdots q_n^{e_n})\) in lowest terms, otherwise this won’t attempt to divide \(p\) by any of the \(q_i\)
OUTPUT:
An element representing the rational expression \(p/(q_1^{e_1} \cdots q_n^{e_n})\).
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: df = [x, 1], [y, 1], [x*y+1, 1] sage: f = FFPD(x, df) sage: f (1, [(y, 1), (x*y + 1, 1)]) sage: ff = FFPD(x, df, reduce=False) sage: ff (x, [(y, 1), (x, 1), (x*y + 1, 1)]) sage: f = FFPD(x + y, [(x + y, 1)]) sage: f (1, [])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> df = [x, Integer(1)], [y, Integer(1)], [x*y+Integer(1), Integer(1)] >>> f = FFPD(x, df) >>> f (1, [(y, 1), (x*y + 1, 1)]) >>> ff = FFPD(x, df, reduce=False) >>> ff (x, [(y, 1), (x, 1), (x*y + 1, 1)]) >>> f = FFPD(x + y, [(x + y, Integer(1))]) >>> f (1, [])
sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 5*x^3 + 1/x + 1/(x-1) + 1/(3*x^2 + 1) sage: FFPD(f) (5*x^7 - 5*x^6 + 5/3*x^5 - 5/3*x^4 + 2*x^3 - 2/3*x^2 + 1/3*x - 1/3, [(x - 1, 1), (x, 1), (x^2 + 1/3, 1)])
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(5)*x**Integer(3) + Integer(1)/x + Integer(1)/(x-Integer(1)) + Integer(1)/(Integer(3)*x**Integer(2) + Integer(1)) >>> FFPD(f) (5*x^7 - 5*x^6 + 5/3*x^5 - 5/3*x^4 + 2*x^3 - 2/3*x^2 + 1/3*x - 1/3, [(x - 1, 1), (x, 1), (x^2 + 1/3, 1)])
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: f = 2*y/(5*(x^3 - 1)*(y + 1)) sage: FFPD(f) (2/5*y, [(y + 1, 1), (x - 1, 1), (x^2 + x + 1, 1)]) sage: p = 1/x^2 sage: q = 3*x**2*y sage: qs = q.factor() sage: f = FFPD(p/qs.unit(), qs) sage: f (1/3/x^2, [(y, 1), (x, 2)]) sage: f = FFPD(cos(x)*x*y^2, [(x, 2), (y, 1)]) sage: f (x*y^2*cos(x), [(y, 1), (x, 2)]) sage: G = exp(x + y) sage: H = (1 - 2*x - y) * (1 - x - 2*y) sage: a = FFPD(G/H) sage: a (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) sage: a.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field sage: b = FFPD(G, H.factor()) sage: b (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) sage: b.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> f = Integer(2)*y/(Integer(5)*(x**Integer(3) - Integer(1))*(y + Integer(1))) >>> FFPD(f) (2/5*y, [(y + 1, 1), (x - 1, 1), (x^2 + x + 1, 1)]) >>> p = Integer(1)/x**Integer(2) >>> q = Integer(3)*x**Integer(2)*y >>> qs = q.factor() >>> f = FFPD(p/qs.unit(), qs) >>> f (1/3/x^2, [(y, 1), (x, 2)]) >>> f = FFPD(cos(x)*x*y**Integer(2), [(x, Integer(2)), (y, Integer(1))]) >>> f (x*y^2*cos(x), [(y, 1), (x, 2)]) >>> G = exp(x + y) >>> H = (Integer(1) - Integer(2)*x - y) * (Integer(1) - x - Integer(2)*y) >>> a = FFPD(G/H) >>> a (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) >>> a.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field >>> b = FFPD(G, H.factor()) >>> b (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) >>> b.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
Singular throws a ‘not implemented’ error when trying to factor in a multivariate polynomial ring over an inexact field:
sage: R.<x,y> = PolynomialRing(CC) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = (x + 1)/(x*y*(x*y + 1)^2) sage: FFPD(f) Traceback (most recent call last): ... TypeError: Singular error: ? not implemented ? error occurred in or before STDIN line ...: `def sage...=factorize(sage...);`
>>> from sage.all import * >>> R = PolynomialRing(CC, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = (x + Integer(1))/(x*y*(x*y + Integer(1))**Integer(2)) >>> FFPD(f) Traceback (most recent call last): ... TypeError: Singular error: ? not implemented ? error occurred in or before STDIN line ...: `def sage...=factorize(sage...);`
AUTHORS:
Alexander Raichev (2012-07-26)
Daniel Krenn (2014-12-01)
- algebraic_dependence_certificate()[source]¶
Return the algebraic dependence certificate of
self.The algebraic dependence certificate is the ideal \(J\) of annihilating polynomials for the set of polynomials
[q^e for (q, e) in self.denominator_factored()], which could be the zero ideal. The ideal \(J\) lies in a polynomial ring over the fieldself.denominator_ring.base_ring()that hasm = len(self.denominator_factored())indeterminates.OUTPUT: an ideal
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/(x^2 * (x*y + 1) * y^3) sage: ff = FFPD(f) sage: J = ff.algebraic_dependence_certificate(); J Ideal (1 - 6*T2 + 15*T2^2 - 20*T2^3 + 15*T2^4 - T0^2*T1^3 - 6*T2^5 + T2^6) of Multivariate Polynomial Ring in T0, T1, T2 over Rational Field sage: g = J.gens()[0] sage: df = ff.denominator_factored() sage: g(*(q**e for q, e in df)) == 0 True
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/(x**Integer(2) * (x*y + Integer(1)) * y**Integer(3)) >>> ff = FFPD(f) >>> J = ff.algebraic_dependence_certificate(); J Ideal (1 - 6*T2 + 15*T2^2 - 20*T2^3 + 15*T2^4 - T0^2*T1^3 - 6*T2^5 + T2^6) of Multivariate Polynomial Ring in T0, T1, T2 over Rational Field >>> g = J.gens()[Integer(0)] >>> df = ff.denominator_factored() >>> g(*(q**e for q, e in df)) == Integer(0) True
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: G = exp(x + y) sage: H = x^2 * (x*y + 1) * y^3 sage: ff = FFPD(G, H.factor()) sage: J = ff.algebraic_dependence_certificate(); J Ideal (1 - 6*T2 + 15*T2^2 - 20*T2^3 + 15*T2^4 - T0^2*T1^3 - 6*T2^5 + T2^6) of Multivariate Polynomial Ring in T0, T1, T2 over Rational Field sage: g = J.gens()[0] sage: df = ff.denominator_factored() sage: g(*(q**e for q, e in df)) == 0 True
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> G = exp(x + y) >>> H = x**Integer(2) * (x*y + Integer(1)) * y**Integer(3) >>> ff = FFPD(G, H.factor()) >>> J = ff.algebraic_dependence_certificate(); J Ideal (1 - 6*T2 + 15*T2^2 - 20*T2^3 + 15*T2^4 - T0^2*T1^3 - 6*T2^5 + T2^6) of Multivariate Polynomial Ring in T0, T1, T2 over Rational Field >>> g = J.gens()[Integer(0)] >>> df = ff.denominator_factored() >>> g(*(q**e for q, e in df)) == Integer(0) True
sage: f = 1/(x^3 * y^2) sage: J = FFPD(f).algebraic_dependence_certificate() sage: J Ideal (0) of Multivariate Polynomial Ring in T0, T1 over Rational Field
>>> from sage.all import * >>> f = Integer(1)/(x**Integer(3) * y**Integer(2)) >>> J = FFPD(f).algebraic_dependence_certificate() >>> J Ideal (0) of Multivariate Polynomial Ring in T0, T1 over Rational Field
sage: f = sin(1)/(x^3 * y^2) sage: J = FFPD(f).algebraic_dependence_certificate() sage: J Ideal (0) of Multivariate Polynomial Ring in T0, T1 over Rational Field
>>> from sage.all import * >>> f = sin(Integer(1))/(x**Integer(3) * y**Integer(2)) >>> J = FFPD(f).algebraic_dependence_certificate() >>> J Ideal (0) of Multivariate Polynomial Ring in T0, T1 over Rational Field
- algebraic_dependence_decomposition(whole_and_parts=True)[source]¶
Return an algebraic dependence decomposition of
self.Let \(f = p/q\) where \(q\) lies in a \(d\)-variate polynomial ring \(K[X]\) for some field \(K\). Let \(q_1^{e_1} \cdots q_n^{e_n}\) be the unique factorization of \(q\) in \(K[X]\) into irreducible factors and let \(V_i\) be the algebraic variety \(\{x \in L^d \mid q_i(x) = 0\}\) of \(q_i\) over the algebraic closure \(L\) of \(K\). By [Rai2012], \(f\) can be written as
\[(*) \quad \sum_A \frac{p_A}{\prod_{i \in A} q_i^{b_i}},\]where the \(b_i\) are positive integers, each \(p_A\) is a products of \(p\) and an element in \(K[X]\), and the sum is taken over all subsets \(A \subseteq \{1, \ldots, m\}\) such that \(|A| \leq d\) and \(\{q_i \mid i \in A\}\) is algebraically independent.
We call \((*)\) an algebraic dependence decomposition of \(f\). Algebraic dependence decompositions are not unique.
The algorithm used comes from [Rai2012].
OUTPUT: an instance of
FractionWithFactoredDenominatorSumEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/(x^2 * (x*y + 1) * y^3) sage: ff = FFPD(f) sage: decomp = ff.algebraic_dependence_decomposition() sage: decomp (0, []) + (-x, [(x*y + 1, 1)]) + (x^2*y^2 - x*y + 1, [(y, 3), (x, 2)]) sage: decomp.sum().quotient() == f True sage: for r in decomp: ....: J = r.algebraic_dependence_certificate() ....: J is None or J == J.ring().ideal() # The zero ideal True True True
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/(x**Integer(2) * (x*y + Integer(1)) * y**Integer(3)) >>> ff = FFPD(f) >>> decomp = ff.algebraic_dependence_decomposition() >>> decomp (0, []) + (-x, [(x*y + 1, 1)]) + (x^2*y^2 - x*y + 1, [(y, 3), (x, 2)]) >>> decomp.sum().quotient() == f True >>> for r in decomp: ... J = r.algebraic_dependence_certificate() ... J is None or J == J.ring().ideal() # The zero ideal True True True
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: G = sin(x) sage: H = x^2 * (x*y + 1) * y^3 sage: f = FFPD(G, H.factor()) sage: decomp = f.algebraic_dependence_decomposition() sage: decomp (0, []) + (x^4*y^3*sin(x), [(x*y + 1, 1)]) + (-(x^5*y^5 - x^4*y^4 + x^3*y^3 - x^2*y^2 + x*y - 1)*sin(x), [(y, 3), (x, 2)]) sage: bool(decomp.sum().quotient() == G/H) True sage: for r in decomp: ....: J = r.algebraic_dependence_certificate() ....: J is None or J == J.ring().ideal() True True True
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> G = sin(x) >>> H = x**Integer(2) * (x*y + Integer(1)) * y**Integer(3) >>> f = FFPD(G, H.factor()) >>> decomp = f.algebraic_dependence_decomposition() >>> decomp (0, []) + (x^4*y^3*sin(x), [(x*y + 1, 1)]) + (-(x^5*y^5 - x^4*y^4 + x^3*y^3 - x^2*y^2 + x*y - 1)*sin(x), [(y, 3), (x, 2)]) >>> bool(decomp.sum().quotient() == G/H) True >>> for r in decomp: ... J = r.algebraic_dependence_certificate() ... J is None or J == J.ring().ideal() True True True
- asymptotic_decomposition(alpha, asy_var=None)[source]¶
Return the asymptotic decomposition of
self.The asymptotic decomposition of \(F\) is a sum that has the same asymptotic expansion as \(f\) in the direction
alphabut each summand has a denominator factorization of the form \([(q_1, 1), \ldots, (q_n, 1)]\), where \(n\) is at most thedimension()of \(F\).INPUT:
alpha– a \(d\)-tuple of positive integers or symbolic variablesasy_var– (default:None) a symbolic variable with respect to which to compute asymptotics; ifNoneis given, we setasy_var = var('r')
OUTPUT: an instance of
FractionWithFactoredDenominatorSumThe output results from a Leinartas decomposition followed by a cohomology decomposition.
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: f = (x^2 + 1)/((x - 1)^3*(x + 2)) sage: F = FFPD(f) sage: alpha = [var('a')] sage: F.asymptotic_decomposition(alpha) (0, []) + (1/54*(5*a^2 + 2*a^2/x + 11*a^2/x^2)*r^2 - 1/54*(5*a - 2*a/x - 33*a/x^2)*r + 11/27/x^2, [(x - 1, 1)]) + (-5/27, [(x + 2, 1)])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> f = (x**Integer(2) + Integer(1))/((x - Integer(1))**Integer(3)*(x + Integer(2))) >>> F = FFPD(f) >>> alpha = [var('a')] >>> F.asymptotic_decomposition(alpha) (0, []) + (1/54*(5*a^2 + 2*a^2/x + 11*a^2/x^2)*r^2 - 1/54*(5*a - 2*a/x - 33*a/x^2)*r + 11/27/x^2, [(x - 1, 1)]) + (-5/27, [(x + 2, 1)])
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - 2*x -y)*(1 - x -2*y)**2 sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = var('a, b') sage: F.asymptotic_decomposition(alpha) (0, []) + (-1/3*r*(a/x - 2*b/y) - 1/3/x + 2/3/y, [(x + 2*y - 1, 1), (2*x + y - 1, 1)])
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - Integer(2)*x -y)*(Integer(1) - x -Integer(2)*y)**Integer(2) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = var('a, b') >>> F.asymptotic_decomposition(alpha) (0, []) + (-1/3*r*(a/x - 2*b/y) - 1/3/x + 2/3/y, [(x + 2*y - 1, 1), (2*x + y - 1, 1)])
- asymptotics(p, alpha, N, asy_var=None, numerical=0, verbose=False)[source]¶
Return the asymptotics in the given direction.
This function returns the first \(N\) terms (some of which could be zero) of the asymptotic expansion of the Maclaurin ray coefficients \(F_{r \alpha}\) of the function \(F\) represented by
selfas \(r \to \infty\), where \(r\) isasy_varandalphais a tuple of positive integers of length \(d\) which isself.dimension(). Assume that\(F\) is holomorphic in a neighborhood of the origin;
the unique factorization of the denominator \(H\) of \(F\) in the local algebraic ring at \(p\) equals its unique factorization in the local analytic ring at \(p\);
the unique factorization of \(H\) in the local algebraic ring at \(p\) has at most
dirreducible factors, none of which are repeated (one can reduce to this case viaasymptotic_decomposition());\(p\) is a convenient strictly minimal smooth or multiple point with all nonzero coordinates that is critical and nondegenerate for
alpha.
The algorithms used here come from [RW2008] and [RW2012].
INPUT:
p– dictionary with keys that can be coerced to equalself.denominator_ring.gens()alpha– tuple of lengthself.dimension()of positive integers or, if \(p\) is a smooth point, possibly of symbolic variablesN– positive integerasy_var– (default:None) a symbolic variable for the asymptotic expansion; ifnoneis given, thenvar('r')will be assignednumerical– (default: 0) a natural number; ifnumericalis greater than 0, then return a numerical approximation of \(F_{r \alpha}\) withnumericaldigits of precision; otherwise return exact valuesverbose– boolean (default:False); print the current state of the algorithm
OUTPUT:
The tuple
(asy, exp_scale, subexp_part). Hereasyis the sum of the first \(N\) terms (some of which might be 0) of the asymptotic expansion of \(F_{r\alpha}\) as \(r \to \infty\);exp_scale**ris the exponential factor ofasy;subexp_partis the subexponential factor ofasy.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
A smooth point example:
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2 sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac); print(F) (1, [(x*y + x + y - 1, 2)]) sage: alpha = [4, 3] sage: decomp = F.asymptotic_decomposition(alpha); decomp (0, []) + (... - 1/2, [(x*y + x + y - 1, 1)]) sage: F1 = decomp[1] sage: p = {y: 1/3, x: 1/2} sage: asy = F1.asymptotics(p, alpha, 2, verbose=True) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... sage: asy (1/6000*(3600*sqrt(5)*sqrt(3)*sqrt(2)*sqrt(r)/sqrt(pi) + 463*sqrt(5)*sqrt(3)*sqrt(2)/(sqrt(pi)*sqrt(r)))*432^r, 432, 3/5*sqrt(5)*sqrt(3)*sqrt(2)*sqrt(r)/sqrt(pi) + 463/6000*sqrt(5)*sqrt(3)*sqrt(2)/(sqrt(pi)*sqrt(r))) sage: F.relative_error(asy[0], alpha, [1, 2, 4, 8, 16], asy[1]) # abs tol 1e-10 # long time [((4, 3), 2.083333333, [2.092576110], [-0.004436533009]), ((8, 6), 2.787374614, [2.790732875], [-0.001204811281]), ((16, 12), 3.826259447, [3.827462310], [-0.0003143703383]), ((32, 24), 5.328112821, [5.328540787], [-0.00008032230388]), ((64, 48), 7.475927885, [7.476079664], [-0.00002030232879])]
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac); print(F) (1, [(x*y + x + y - 1, 2)]) >>> alpha = [Integer(4), Integer(3)] >>> decomp = F.asymptotic_decomposition(alpha); decomp (0, []) + (... - 1/2, [(x*y + x + y - 1, 1)]) >>> F1 = decomp[Integer(1)] >>> p = {y: Integer(1)/Integer(3), x: Integer(1)/Integer(2)} >>> asy = F1.asymptotics(p, alpha, Integer(2), verbose=True) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... >>> asy (1/6000*(3600*sqrt(5)*sqrt(3)*sqrt(2)*sqrt(r)/sqrt(pi) + 463*sqrt(5)*sqrt(3)*sqrt(2)/(sqrt(pi)*sqrt(r)))*432^r, 432, 3/5*sqrt(5)*sqrt(3)*sqrt(2)*sqrt(r)/sqrt(pi) + 463/6000*sqrt(5)*sqrt(3)*sqrt(2)/(sqrt(pi)*sqrt(r))) >>> F.relative_error(asy[Integer(0)], alpha, [Integer(1), Integer(2), Integer(4), Integer(8), Integer(16)], asy[Integer(1)]) # abs tol 1e-10 # long time [((4, 3), 2.083333333, [2.092576110], [-0.004436533009]), ((8, 6), 2.787374614, [2.790732875], [-0.001204811281]), ((16, 12), 3.826259447, [3.827462310], [-0.0003143703383]), ((32, 24), 5.328112821, [5.328540787], [-0.00008032230388]), ((64, 48), 7.475927885, [7.476079664], [-0.00002030232879])]
A multiple point example:
sage: R.<x,y,z>= PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (4 - 2*x - y - z)**2*(4 - x - 2*y - z) sage: Hfac = H.factor() sage: G = 16/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (-16, [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 2)]) sage: alpha = [3, 3, 2] sage: decomp = F.asymptotic_decomposition(alpha); decomp (0, []) + (..., [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 1)]) sage: F1 = decomp[1] sage: p = {x: 1, y: 1, z: 1} sage: asy = F1.asymptotics(p, alpha, 2, verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... sage: asy # long time (4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r)), 1, 4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r))) sage: F.relative_error(asy[0], alpha, [1, 2, 4, 8], asy[1]) # long time [((3, 3, 2), 0.9812164307, [1.515572606], [-0.54458543...]), ((6, 6, 4), 1.576181132, [1.992989399], [-0.26444185...]), ((12, 12, 8), 2.485286378, [2.712196351], [-0.091301338...]), ((24, 24, 16), 3.700576827, [3.760447895], [-0.016178847...])]
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(4) - Integer(2)*x - y - z)**Integer(2)*(Integer(4) - x - Integer(2)*y - z) >>> Hfac = H.factor() >>> G = Integer(16)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (-16, [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 2)]) >>> alpha = [Integer(3), Integer(3), Integer(2)] >>> decomp = F.asymptotic_decomposition(alpha); decomp (0, []) + (..., [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 1)]) >>> F1 = decomp[Integer(1)] >>> p = {x: Integer(1), y: Integer(1), z: Integer(1)} >>> asy = F1.asymptotics(p, alpha, Integer(2), verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... >>> asy # long time (4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r)), 1, 4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r))) >>> F.relative_error(asy[Integer(0)], alpha, [Integer(1), Integer(2), Integer(4), Integer(8)], asy[Integer(1)]) # long time [((3, 3, 2), 0.9812164307, [1.515572606], [-0.54458543...]), ((6, 6, 4), 1.576181132, [1.992989399], [-0.26444185...]), ((12, 12, 8), 2.485286378, [2.712196351], [-0.091301338...]), ((24, 24, 16), 3.700576827, [3.760447895], [-0.016178847...])]
- asymptotics_multiple(p, alpha, N, asy_var, coordinate=None, numerical=0, verbose=False)[source]¶
Return the asymptotics in the given direction of a multiple point nondegenerate for
alpha.This is the same as
asymptotics(), but only in the case of a convenient multiple point nondegenerate foralpha. Assume also thatself.dimension >= 2and that thep.values()are not symbolic variables.The formulas used for computing the asymptotic expansion are Theorem 3.4 and Theorem 3.7 of [RW2012].
INPUT:
p– dictionary with keys that can be coerced to equalself.denominator_ring.gens()alpha– tuple of lengthd = self.dimension()of positive integers or, if \(p\) is a smooth point, possibly of symbolic variablesN– positive integerasy_var– (default:None) a symbolic variable; the variable of the asymptotic expansion, if none is given,var('r')will be assignedcoordinate– (default:None) an integer in \(\{0, \ldots, d-1\}\) indicating a convenient coordinate to base the asymptotic calculations on; ifNoneis assigned, then choosecoordinate=d-1numerical– (default: 0) a natural number; if numerical is greater than 0, then return a numerical approximation of the Maclaurin ray coefficients ofselfwithnumericaldigits of precision. Otherwise return exact values.verbose– boolean (default:False); print the current state of the algorithm
OUTPUT: the asymptotic expansion
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y,z>= PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = (4 - 2*x - y - z)*(4 - x -2*y - z) sage: Hfac = H.factor() sage: G = 16/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (16, [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 1)]) sage: p = {x: 1, y: 1, z: 1} sage: alpha = [3, 3, 2] sage: F.asymptotics_multiple(p, alpha, 2, var('r'), verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... (4/3*sqrt(3)/(sqrt(pi)*sqrt(r)) - 25/216*sqrt(3)/(sqrt(pi)*r^(3/2)), 1, 4/3*sqrt(3)/(sqrt(pi)*sqrt(r)) - 25/216*sqrt(3)/(sqrt(pi)*r^(3/2))) sage: H = (1 - x*(1 + y))*(1 - z*x**2*(1 + 2*y)) sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (1, [(x*y + x - 1, 1), (2*x^2*y*z + x^2*z - 1, 1)]) sage: p = {x: 1/2, z: 4/3, y: 1} sage: alpha = [8, 3, 3] sage: F.asymptotics_multiple(p, alpha, 2, var('r'), coordinate=1, verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... (1/172872*108^r*(24696*sqrt(7)*sqrt(3)/(sqrt(pi)*sqrt(r)) - 1231*sqrt(7)*sqrt(3)/(sqrt(pi)*r^(3/2))), 108, 1/7*sqrt(7)*sqrt(3)/(sqrt(pi)*sqrt(r)) - 1231/172872*sqrt(7)*sqrt(3)/(sqrt(pi)*r^(3/2)))
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = (Integer(4) - Integer(2)*x - y - z)*(Integer(4) - x -Integer(2)*y - z) >>> Hfac = H.factor() >>> G = Integer(16)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (16, [(x + 2*y + z - 4, 1), (2*x + y + z - 4, 1)]) >>> p = {x: Integer(1), y: Integer(1), z: Integer(1)} >>> alpha = [Integer(3), Integer(3), Integer(2)] >>> F.asymptotics_multiple(p, alpha, Integer(2), var('r'), verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... (4/3*sqrt(3)/(sqrt(pi)*sqrt(r)) - 25/216*sqrt(3)/(sqrt(pi)*r^(3/2)), 1, 4/3*sqrt(3)/(sqrt(pi)*sqrt(r)) - 25/216*sqrt(3)/(sqrt(pi)*r^(3/2))) >>> H = (Integer(1) - x*(Integer(1) + y))*(Integer(1) - z*x**Integer(2)*(Integer(1) + Integer(2)*y)) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (1, [(x*y + x - 1, 1), (2*x^2*y*z + x^2*z - 1, 1)]) >>> p = {x: Integer(1)/Integer(2), z: Integer(4)/Integer(3), y: Integer(1)} >>> alpha = [Integer(8), Integer(3), Integer(3)] >>> F.asymptotics_multiple(p, alpha, Integer(2), var('r'), coordinate=Integer(1), verbose=True) # long time Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second-order differential operator actions... (1/172872*108^r*(24696*sqrt(7)*sqrt(3)/(sqrt(pi)*sqrt(r)) - 1231*sqrt(7)*sqrt(3)/(sqrt(pi)*r^(3/2))), 108, 1/7*sqrt(7)*sqrt(3)/(sqrt(pi)*sqrt(r)) - 1231/172872*sqrt(7)*sqrt(3)/(sqrt(pi)*r^(3/2)))
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - 2*x - y) * (1 - x - 2*y) sage: Hfac = H.factor() sage: G = exp(x + y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) sage: p = {x: 1/3, y: 1/3} sage: alpha = (var('a'), var('b')) sage: F.asymptotics_multiple(p, alpha, 2, var('r')) # long time (3*(1/((1/3)^a*(1/3)^b))^r*e^(2/3), 1/((1/3)^a*(1/3)^b), 3*e^(2/3))
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - Integer(2)*x - y) * (Integer(1) - x - Integer(2)*y) >>> Hfac = H.factor() >>> G = exp(x + y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (e^(x + y), [(x + 2*y - 1, 1), (2*x + y - 1, 1)]) >>> p = {x: Integer(1)/Integer(3), y: Integer(1)/Integer(3)} >>> alpha = (var('a'), var('b')) >>> F.asymptotics_multiple(p, alpha, Integer(2), var('r')) # long time (3*(1/((1/3)^a*(1/3)^b))^r*e^(2/3), 1/((1/3)^a*(1/3)^b), 3*e^(2/3))
- asymptotics_smooth(p, alpha, N, asy_var, coordinate=None, numerical=0, verbose=False)[source]¶
Return the asymptotics in the given direction of a smooth point.
This is the same as
asymptotics(), but only in the case of a convenient smooth point.The formulas used for computing the asymptotic expansions are Theorems 3.2 and 3.3 [RW2008] with the exponent of \(H\) equal to 1. Theorem 3.2 is a specialization of Theorem 3.4 of [RW2012] with \(n = 1\).
INPUT:
p– dictionary with keys that can be coerced to equalself.denominator_ring.gens()alpha– tuple of lengthd = self.dimension()of positive integers or, if \(p\) is a smooth point, possibly of symbolic variablesN– positive integerasy_var– (default:None) a symbolic variable; the variable of the asymptotic expansion, if none is given,var('r')will be assignedcoordinate– (default:None) an integer in \(\{0, \ldots, d-1\}\) indicating a convenient coordinate to base the asymptotic calculations on; ifNoneis assigned, then choosecoordinate=d-1numerical– (default: 0) a natural number; if numerical is greater than 0, then return a numerical approximation of the Maclaurin ray coefficients ofselfwithnumericaldigits of precision; otherwise return exact valuesverbose– boolean (default:False); print the current state of the algorithm
OUTPUT: the asymptotic expansion
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = 2 - 3*x sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (-1/3, [(x - 2/3, 1)]) sage: alpha = [2] sage: p = {x: 2/3} sage: asy = F.asymptotics_smooth(p, alpha, 3, asy_var=var('r')) sage: asy (1/2*(9/4)^r, 9/4, 1/2)
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = Integer(2) - Integer(3)*x >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (-1/3, [(x - 2/3, 1)]) >>> alpha = [Integer(2)] >>> p = {x: Integer(2)/Integer(3)} >>> asy = F.asymptotics_smooth(p, alpha, Integer(3), asy_var=var('r')) >>> asy (1/2*(9/4)^r, 9/4, 1/2)
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = 1-x-y-x*y sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = [3, 2] sage: p = {y: 1/2*sqrt(13) - 3/2, x: 1/3*sqrt(13) - 2/3} sage: F.asymptotics_smooth(p, alpha, 2, var('r'), numerical=3, verbose=True) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... (71.2^r*(0.369/sqrt(r) - 0.018.../r^(3/2)), 71.2, 0.369/sqrt(r) - 0.018.../r^(3/2)) sage: q = 1/2 sage: qq = q.denominator() sage: H = 1 - q*x + q*x*y - x^2*y sage: Hfac = H.factor() sage: G = (1 - q*x)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = list(qq*vector([2, 1 - q])) sage: alpha [4, 1] sage: p = {x: 1, y: 1} sage: F.asymptotics_smooth(p, alpha, 5, var('r'), verbose=True) # not tested (140 seconds) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... (1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)) - 1/96*sqrt(3)*2^(1/3)*gamma(2/3)/(pi*r^(5/3)), 1, 1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)) - 1/96*sqrt(3)*2^(1/3)*gamma(2/3)/(pi*r^(5/3)))
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = Integer(1)-x-y-x*y >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = [Integer(3), Integer(2)] >>> p = {y: Integer(1)/Integer(2)*sqrt(Integer(13)) - Integer(3)/Integer(2), x: Integer(1)/Integer(3)*sqrt(Integer(13)) - Integer(2)/Integer(3)} >>> F.asymptotics_smooth(p, alpha, Integer(2), var('r'), numerical=Integer(3), verbose=True) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... (71.2^r*(0.369/sqrt(r) - 0.018.../r^(3/2)), 71.2, 0.369/sqrt(r) - 0.018.../r^(3/2)) >>> q = Integer(1)/Integer(2) >>> qq = q.denominator() >>> H = Integer(1) - q*x + q*x*y - x**Integer(2)*y >>> Hfac = H.factor() >>> G = (Integer(1) - q*x)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = list(qq*vector([Integer(2), Integer(1) - q])) >>> alpha [4, 1] >>> p = {x: Integer(1), y: Integer(1)} >>> F.asymptotics_smooth(p, alpha, Integer(5), var('r'), verbose=True) # not tested (140 seconds) Creating auxiliary functions... Computing derivatives of auxiliary functions... Computing derivatives of more auxiliary functions... Computing second order differential operator actions... (1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)) - 1/96*sqrt(3)*2^(1/3)*gamma(2/3)/(pi*r^(5/3)), 1, 1/12*sqrt(3)*2^(2/3)*gamma(1/3)/(pi*r^(1/3)) - 1/96*sqrt(3)*2^(1/3)*gamma(2/3)/(pi*r^(5/3)))
- cohomology_decomposition()[source]¶
Return the cohomology decomposition of
self.Let \(p / (q_1^{e_1} \cdots q_n^{e_n})\) be the fraction represented by
selfand let \(K[x_1, \ldots, x_d]\) be the polynomial ring in which the \(q_i\) lie. Assume that \(n \leq d\) and that the gradients of the \(q_i\) are linearly independent at all points in the intersection \(V_1 \cap \ldots \cap V_n\) of the algebraic varieties \(V_i = \{x \in L^d \mid q_i(x) = 0 \}\), where \(L\) is the algebraic closure of the field \(K\). Return aFractionWithFactoredDenominatorSum\(f\) such that the differential form \(f dx_1 \wedge \cdots \wedge dx_d\) is de Rham cohomologous to the differential form \(p / (q_1^{e_1} \cdots q_n^{e_n}) dx_1 \wedge \cdots \wedge dx_d\) and such that the denominator of each summand of \(f\) contains no repeated irreducible factors.The algorithm used here comes from the proof of Theorem 17.4 of [AY1983].
OUTPUT: an instance of
FractionWithFactoredDenominatorSumEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/(x^2 + x + 1)^3 sage: decomp = FFPD(f).cohomology_decomposition() sage: decomp (0, []) + (2/3, [(x^2 + x + 1, 1)])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/(x**Integer(2) + x + Integer(1))**Integer(3) >>> decomp = FFPD(f).cohomology_decomposition() >>> decomp (0, []) + (2/3, [(x^2 + x + 1, 1)])
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: FFPD(1, [(x, 1), (y, 2)]).cohomology_decomposition() (0, [])
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> FFPD(Integer(1), [(x, Integer(1)), (y, Integer(2))]).cohomology_decomposition() (0, [])
The following example was fixed in Issue #29465:
sage: p = 1 sage: qs = [(x*y - 1, 1), (x**2 + y**2 - 1, 2)] sage: f = FFPD(p, qs) sage: f.cohomology_decomposition() (0, []) + (-4/3*x*y, [(x^2 + y^2 - 1, 1)]) + (1/3, [(x*y - 1, 1), (x^2 + y^2 - 1, 1)])
>>> from sage.all import * >>> p = Integer(1) >>> qs = [(x*y - Integer(1), Integer(1)), (x**Integer(2) + y**Integer(2) - Integer(1), Integer(2))] >>> f = FFPD(p, qs) >>> f.cohomology_decomposition() (0, []) + (-4/3*x*y, [(x^2 + y^2 - 1, 1)]) + (1/3, [(x*y - 1, 1), (x^2 + y^2 - 1, 1)])
- critical_cone(p, coordinate=None)[source]¶
Return the critical cone of the convenient multiple point
p.INPUT:
p– dictionary with keys that can be coerced to equalself.denominator_ring.gens()and values in a fieldcoordinate– (default:None) a natural number
OUTPUT: list of vectors
This list of vectors generate the critical cone of
pand the cone itself, which isNoneif the values ofpdon’t lie in \(\QQ\). Divide logarithmic gradients by their componentcoordinateentries. Ifcoordinate = None, then search from \(d-1\) down to 0 for the first indexjsuch that for alliwe haveself.log_grads()[i][j] != 0and setcoordinate = j.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y,z> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: G = 1 sage: H = (1 - x*(1 + y)) * (1 - z*x**2*(1 + 2*y)) sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: p = {x: 1/2, y: 1, z: 4/3} sage: F.critical_cone(p) ([(2, 1, 0), (3, 1, 3/2)], 2-d cone in 3-d lattice N)
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> G = Integer(1) >>> H = (Integer(1) - x*(Integer(1) + y)) * (Integer(1) - z*x**Integer(2)*(Integer(1) + Integer(2)*y)) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> p = {x: Integer(1)/Integer(2), y: Integer(1), z: Integer(4)/Integer(3)} >>> F.critical_cone(p) ([(2, 1, 0), (3, 1, 3/2)], 2-d cone in 3-d lattice N)
- denominator()[source]¶
Return the denominator of
self.OUTPUT:
The denominator (i.e., the product of the factored denominator).
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.denominator() x^3*y^2 + 2*x^3*y + x^2*y^2 + x^3 - 2*x^2*y - x*y^2 - 3*x^2 - 2*x*y - y^2 + 3*x + 2*y - 1
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.denominator() x^3*y^2 + 2*x^3*y + x^2*y^2 + x^3 - 2*x^2*y - x*y^2 - 3*x^2 - 2*x*y - y^2 + 3*x + 2*y - 1
- denominator_factored()[source]¶
Return the factorization in
self.denominator_ringof the denominator ofselfbut without the unit part.OUTPUT:
The factored denominator as a list of tuple
(f, m), where \(f\) is a factor and \(m\) its multiplicity.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.denominator_factored() [(x - 1, 1), (x*y + x + y - 1, 2)]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.denominator_factored() [(x - 1, 1), (x*y + x + y - 1, 2)]
- property denominator_ring¶
Return the ring of the denominator.
OUTPUT: a ring
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field sage: F = FFPD(G/H) sage: F (e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) sage: F.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field >>> F = FFPD(G/H) >>> F (e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) >>> F.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
- dimension()[source]¶
Return the number of indeterminates of
self.denominator_ring.OUTPUT: integer
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.dimension() 2
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.dimension() 2
- grads(p)[source]¶
Return a list of the gradients of the polynomials
[q for (q, e) in self.denominator_factored()]evaluated atp.INPUT:
p– (default:None) a dictionary whose keys are the generators ofself.denominator_ring
OUTPUT: list
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: p = exp(x) sage: df = [(x^3 + 3*y^2, 5), (x*y, 2), (y, 1)] sage: f = FFPD(p, df) sage: f (e^x, [(y, 1), (x*y, 2), (x^3 + 3*y^2, 5)]) sage: R.gens() (x, y) sage: p = None sage: f.grads(p) [(0, 1), (y, x), (3*x^2, 6*y)] sage: p = {x: sqrt(2), y: var('a')} sage: f.grads(p) [(0, 1), (a, sqrt(2)), (6, 6*a)]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> p = exp(x) >>> df = [(x**Integer(3) + Integer(3)*y**Integer(2), Integer(5)), (x*y, Integer(2)), (y, Integer(1))] >>> f = FFPD(p, df) >>> f (e^x, [(y, 1), (x*y, 2), (x^3 + 3*y^2, 5)]) >>> R.gens() (x, y) >>> p = None >>> f.grads(p) [(0, 1), (y, x), (3*x^2, 6*y)] >>> p = {x: sqrt(Integer(2)), y: var('a')} >>> f.grads(p) [(0, 1), (a, sqrt(2)), (6, 6*a)]
- is_convenient_multiple_point(p)[source]¶
Test if
pis a convenient multiple point ofself.In case
pis a convenient multiple point,verdict = Trueandcommentis a string stating which variables it’s convenient to use. In casepis not,verdict = Falseandcommentis a string explaining whypfails to be a convenient multiple point.See [RW2012] for more details.
INPUT:
p– dictionary with keys that can be coerced to equalself.denominator_ring.gens()
OUTPUT:
A pair
(verdict, comment).EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y,z> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = (1 - x*(1 + y)) * (1 - z*x**2*(1 + 2*y)) sage: df = H.factor() sage: G = 1 / df.unit() sage: F = FFPD(G, df) sage: p1 = {x: 1/2, y: 1, z: 4/3} sage: p2 = {x: 1, y: 2, z: 1/2} sage: F.is_convenient_multiple_point(p1) (True, 'convenient in variables [x, y]') sage: F.is_convenient_multiple_point(p2) (False, 'not a singular point')
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = (Integer(1) - x*(Integer(1) + y)) * (Integer(1) - z*x**Integer(2)*(Integer(1) + Integer(2)*y)) >>> df = H.factor() >>> G = Integer(1) / df.unit() >>> F = FFPD(G, df) >>> p1 = {x: Integer(1)/Integer(2), y: Integer(1), z: Integer(4)/Integer(3)} >>> p2 = {x: Integer(1), y: Integer(2), z: Integer(1)/Integer(2)} >>> F.is_convenient_multiple_point(p1) (True, 'convenient in variables [x, y]') >>> F.is_convenient_multiple_point(p2) (False, 'not a singular point')
- leinartas_decomposition()[source]¶
Return a Leinartas decomposition of
self.Let \(f = p/q\) where \(q\) lies in a \(d\) -variate polynomial ring \(K[X]\) for some field \(K\). Let \(q_1^{e_1} \cdots q_n^{e_n}\) be the unique factorization of \(q\) in \(K[X]\) into irreducible factors and let \(V_i\) be the algebraic variety \(\{x\in L^d \mid q_i(x) = 0\}\) of \(q_i\) over the algebraic closure \(L\) of \(K\). By [Rai2012], \(f\) can be written as
\[(*) \quad \sum_A \frac{p_A}{\prod_{i \in A} q_i^{b_i}},\]where the \(b_i\) are positive integers, each \(p_A\) is a product of \(p\) and an element of \(K[X]\), and the sum is taken over all subsets \(A \subseteq \{1, \ldots, m\}\) such that
\(|A| \le d\),
\(\bigcap_{i\in A} T_i \neq \emptyset\), and
\(\{q_i \mid i\in A\}\) is algebraically independent.
In particular, any rational expression in \(d\) variables can be represented as a sum of rational expressions whose denominators each contain at most \(d\) distinct irreducible factors.
We call \((*)\) a Leinartas decomposition of \(f\). Leinartas decompositions are not unique.
The algorithm used comes from [Rai2012].
OUTPUT: an instance of
FractionWithFactoredDenominatorSumEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = (x^2 + 1)/((x + 2)*(x - 1)*(x^2 + x + 1)) sage: decomp = FFPD(f).leinartas_decomposition() sage: decomp (0, []) + (2/9, [(x - 1, 1)]) + (-5/9, [(x + 2, 1)]) + (1/3*x, [(x^2 + x + 1, 1)]) sage: decomp.sum().quotient() == f True
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = (x**Integer(2) + Integer(1))/((x + Integer(2))*(x - Integer(1))*(x**Integer(2) + x + Integer(1))) >>> decomp = FFPD(f).leinartas_decomposition() >>> decomp (0, []) + (2/9, [(x - 1, 1)]) + (-5/9, [(x + 2, 1)]) + (1/3*x, [(x^2 + x + 1, 1)]) >>> decomp.sum().quotient() == f True
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/x + 1/y + 1/(x*y + 1) sage: decomp = FFPD(f).leinartas_decomposition() sage: decomp (0, []) + (1, [(x*y + 1, 1)]) + (x + y, [(y, 1), (x, 1)]) sage: decomp.sum().quotient() == f True sage: def check_decomp(r): ....: L = r.nullstellensatz_certificate() ....: J = r.algebraic_dependence_certificate() ....: return L is None and (J is None or J == J.ring().ideal()) sage: all(check_decomp(r) for r in decomp) True
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/x + Integer(1)/y + Integer(1)/(x*y + Integer(1)) >>> decomp = FFPD(f).leinartas_decomposition() >>> decomp (0, []) + (1, [(x*y + 1, 1)]) + (x + y, [(y, 1), (x, 1)]) >>> decomp.sum().quotient() == f True >>> def check_decomp(r): ... L = r.nullstellensatz_certificate() ... J = r.algebraic_dependence_certificate() ... return L is None and (J is None or J == J.ring().ideal()) >>> all(check_decomp(r) for r in decomp) True
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: f = sin(x)/x + 1/y + 1/(x*y + 1) sage: G = f.numerator() sage: H = R(f.denominator()) sage: ff = FFPD(G, H.factor()) sage: decomp = ff.leinartas_decomposition() sage: decomp # random - non canonical depends on singular version (0, []) + (-(x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x)*y, [(y, 1)]) + ((x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x)*x*y, [(x*y + 1, 1)]) + (x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x, [(y, 1), (x, 1)]) sage: bool(decomp.sum().quotient() == f) True sage: all(check_decomp(r) for r in decomp) True
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> f = sin(x)/x + Integer(1)/y + Integer(1)/(x*y + Integer(1)) >>> G = f.numerator() >>> H = R(f.denominator()) >>> ff = FFPD(G, H.factor()) >>> decomp = ff.leinartas_decomposition() >>> decomp # random - non canonical depends on singular version (0, []) + (-(x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x)*y, [(y, 1)]) + ((x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x)*x*y, [(x*y + 1, 1)]) + (x*y^2*sin(x) + x^2*y + x*y + y*sin(x) + x, [(y, 1), (x, 1)]) >>> bool(decomp.sum().quotient() == f) True >>> all(check_decomp(r) for r in decomp) True
sage: R.<x,y,z>= PolynomialRing(GF(2, 'a')) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/(x * y * z * (x*y + z)) sage: decomp = FFPD(f).leinartas_decomposition() sage: decomp (0, []) + (1, [(z, 2), (x*y + z, 1)]) + (1, [(z, 2), (y, 1), (x, 1)]) sage: decomp.sum().quotient() == f True
>>> from sage.all import * >>> R = PolynomialRing(GF(Integer(2), 'a'), names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/(x * y * z * (x*y + z)) >>> decomp = FFPD(f).leinartas_decomposition() >>> decomp (0, []) + (1, [(z, 2), (x*y + z, 1)]) + (1, [(z, 2), (y, 1), (x, 1)]) >>> decomp.sum().quotient() == f True
- log_grads(p)[source]¶
Return a list of the logarithmic gradients of the polynomials
[q for (q, e) in self.denominator_factored()]evaluated atp.The logarithmic gradient of a function \(f\) at point \(p\) is the vector \((x_1 \partial_1 f(x), \ldots, x_d \partial_d f(x) )\) evaluated at \(p\).
INPUT:
p– (default:None) a dictionary whose keys are the generators ofself.denominator_ring
OUTPUT: list
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: p = exp(x) sage: df = [(x^3 + 3*y^2, 5), (x*y, 2), (y, 1)] sage: f = FFPD(p, df) sage: f (e^x, [(y, 1), (x*y, 2), (x^3 + 3*y^2, 5)]) sage: R.gens() (x, y) sage: p = None sage: f.log_grads(p) [(0, y), (x*y, x*y), (3*x^3, 6*y^2)] sage: p = {x: sqrt(2), y: var('a')} sage: f.log_grads(p) [(0, a), (sqrt(2)*a, sqrt(2)*a), (6*sqrt(2), 6*a^2)]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> p = exp(x) >>> df = [(x**Integer(3) + Integer(3)*y**Integer(2), Integer(5)), (x*y, Integer(2)), (y, Integer(1))] >>> f = FFPD(p, df) >>> f (e^x, [(y, 1), (x*y, 2), (x^3 + 3*y^2, 5)]) >>> R.gens() (x, y) >>> p = None >>> f.log_grads(p) [(0, y), (x*y, x*y), (3*x^3, 6*y^2)] >>> p = {x: sqrt(Integer(2)), y: var('a')} >>> f.log_grads(p) [(0, a), (sqrt(2)*a, sqrt(2)*a), (6*sqrt(2), 6*a^2)]
- maclaurin_coefficients(multi_indices, numerical=0)[source]¶
Return the Maclaurin coefficients of
selfwith givenmulti_indices.INPUT:
multi_indices– list of tuples of positive integers, where each tuple has lengthself.dimension()numerical– (default: 0) a natural number; if positive, return numerical approximations of coefficients withnumericaldigits of accuracy
OUTPUT:
A dictionary whose value of the key
nuare the Maclaurin coefficient of indexnuofself.Note
Uses iterated univariate Maclaurin expansions. Slow.
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = 2 - 3*x sage: Hfac = H.factor() sage: G = 1 / Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (-1/3, [(x - 2/3, 1)]) sage: F.maclaurin_coefficients([(2*k,) for k in range(6)]) {(0,): 1/2, (2,): 9/8, (4,): 81/32, (6,): 729/128, (8,): 6561/512, (10,): 59049/2048}
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = Integer(2) - Integer(3)*x >>> Hfac = H.factor() >>> G = Integer(1) / Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (-1/3, [(x - 2/3, 1)]) >>> F.maclaurin_coefficients([(Integer(2)*k,) for k in range(Integer(6))]) {(0,): 1/2, (2,): 9/8, (4,): 81/32, (6,): 729/128, (8,): 6561/512, (10,): 59049/2048}
sage: R.<x,y,z> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = (4 - 2*x - y - z) * (4 - x - 2*y - z) sage: Hfac = H.factor() sage: G = 16 / Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = vector([3, 3, 2]) sage: interval = [1, 2, 4] sage: S = [r*alpha for r in interval] sage: F.maclaurin_coefficients(S, numerical=10) # long time {(3, 3, 2): 0.7849731445, (6, 6, 4): 0.7005249476, (12, 12, 8): 0.5847732654}
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = (Integer(4) - Integer(2)*x - y - z) * (Integer(4) - x - Integer(2)*y - z) >>> Hfac = H.factor() >>> G = Integer(16) / Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = vector([Integer(3), Integer(3), Integer(2)]) >>> interval = [Integer(1), Integer(2), Integer(4)] >>> S = [r*alpha for r in interval] >>> F.maclaurin_coefficients(S, numerical=Integer(10)) # long time {(3, 3, 2): 0.7849731445, (6, 6, 4): 0.7005249476, (12, 12, 8): 0.5847732654}
- nullstellensatz_certificate()[source]¶
Return a Nullstellensatz certificate of
selfif it exists.Let \([(q_1, e_1), \ldots, (q_n, e_n)]\) be the denominator factorization of
self. The Nullstellensatz certificate is a list of polynomials \(h_1, \ldots, h_m\) inself.denominator_ringthat satisfies \(h_1 q_1 + \cdots + h_m q_n = 1\) if it exists.Note
Only works for multivariate base rings.
OUTPUT:
A list of polynomials or
Noneif no Nullstellensatz certificate exists.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: G = sin(x) sage: H = x^2 * (x*y + 1) sage: f = FFPD(G, H.factor()) sage: L = f.nullstellensatz_certificate() sage: L [y^2, -x*y + 1] sage: df = f.denominator_factored() sage: sum(L[i]*df[i][0]**df[i][1] for i in range(len(df))) == 1 True
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> G = sin(x) >>> H = x**Integer(2) * (x*y + Integer(1)) >>> f = FFPD(G, H.factor()) >>> L = f.nullstellensatz_certificate() >>> L [y^2, -x*y + 1] >>> df = f.denominator_factored() >>> sum(L[i]*df[i][Integer(0)]**df[i][Integer(1)] for i in range(len(df))) == Integer(1) True
sage: f = 1/(x*y) sage: L = FFPD(f).nullstellensatz_certificate() sage: L is None True
>>> from sage.all import * >>> f = Integer(1)/(x*y) >>> L = FFPD(f).nullstellensatz_certificate() >>> L is None True
- nullstellensatz_decomposition()[source]¶
Return a Nullstellensatz decomposition of
self.Let \(f = p/q\) where \(q\) lies in a \(d\) -variate polynomial ring \(K[X]\) for some field \(K\) and \(d \geq 1\). Let \(q_1^{e_1} \cdots q_n^{e_n}\) be the unique factorization of \(q\) in \(K[X]\) into irreducible factors and let \(V_i\) be the algebraic variety \(\{x \in L^d \mid q_i(x) = 0\}\) of \(q_i\) over the algebraic closure \(L\) of \(K\). By [Rai2012], \(f\) can be written as
\[(*) \quad \sum_A \frac{p_A}{\prod_{i \in A} q_i^{e_i}},\]where the \(p_A\) are products of \(p\) and elements in \(K[X]\) and the sum is taken over all subsets \(A \subseteq \{1, \ldots, m\}\) such that \(\bigcap_{i\in A} T_i \neq \emptyset\).
We call \((*)\) a Nullstellensatz decomposition of \(f\). Nullstellensatz decompositions are not unique.
The algorithm used comes from [Rai2012].
Note
Recursive. Only works for multivariate
self.OUTPUT: an instance of
FractionWithFactoredDenominatorSumEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import * sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 1/(x*(x*y + 1)) sage: decomp = FFPD(f).nullstellensatz_decomposition() sage: decomp (0, []) + (1, [(x, 1)]) + (-y, [(x*y + 1, 1)]) sage: decomp.sum().quotient() == f True sage: [r.nullstellensatz_certificate() is None for r in decomp] [True, True, True]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(1)/(x*(x*y + Integer(1))) >>> decomp = FFPD(f).nullstellensatz_decomposition() >>> decomp (0, []) + (1, [(x, 1)]) + (-y, [(x*y + 1, 1)]) >>> decomp.sum().quotient() == f True >>> [r.nullstellensatz_certificate() is None for r in decomp] [True, True, True]
sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: G = sin(y) sage: H = x*(x*y + 1) sage: f = FFPD(G, H.factor()) sage: decomp = f.nullstellensatz_decomposition() sage: decomp (0, []) + (sin(y), [(x, 1)]) + (-y*sin(y), [(x*y + 1, 1)]) sage: bool(decomp.sum().quotient() == G/H) True sage: [r.nullstellensatz_certificate() is None for r in decomp] [True, True, True]
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> G = sin(y) >>> H = x*(x*y + Integer(1)) >>> f = FFPD(G, H.factor()) >>> decomp = f.nullstellensatz_decomposition() >>> decomp (0, []) + (sin(y), [(x, 1)]) + (-y*sin(y), [(x*y + 1, 1)]) >>> bool(decomp.sum().quotient() == G/H) True >>> [r.nullstellensatz_certificate() is None for r in decomp] [True, True, True]
- numerator()[source]¶
Return the numerator of
self.OUTPUT: the numerator
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.numerator() -e^y
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.numerator() -e^y
- property numerator_ring¶
Return the ring of the numerator.
OUTPUT: a ring
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F.numerator_ring Symbolic Ring sage: F = FFPD(G/H) sage: F (e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) sage: F.numerator_ring Symbolic Ring
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F.numerator_ring Symbolic Ring >>> F = FFPD(G/H) >>> F (e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) >>> F.numerator_ring Symbolic Ring
- quotient()[source]¶
Convert
selfinto a quotient.OUTPUT: an element
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: H = (1 - x - y - x*y)**2*(1-x) sage: Hfac = H.factor() sage: G = exp(y)/Hfac.unit() sage: F = FFPD(G, Hfac) sage: F (-e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) sage: F.quotient() -e^y/(x^3*y^2 + 2*x^3*y + x^2*y^2 + x^3 - 2*x^2*y - x*y^2 - 3*x^2 - 2*x*y - y^2 + 3*x + 2*y - 1)
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> H = (Integer(1) - x - y - x*y)**Integer(2)*(Integer(1)-x) >>> Hfac = H.factor() >>> G = exp(y)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> F (-e^y, [(x - 1, 1), (x*y + x + y - 1, 2)]) >>> F.quotient() -e^y/(x^3*y^2 + 2*x^3*y + x^2*y^2 + x^3 - 2*x^2*y - x*y^2 - 3*x^2 - 2*x*y - y^2 + 3*x + 2*y - 1)
- relative_error(approx, alpha, interval, exp_scale=1, digits=10)[source]¶
Return the relative error between the values of the Maclaurin coefficients of
selfwith multi-indicesr alphaforrinintervaland the values of the functions (of the variabler) inapprox.INPUT:
approx– an individual or list of symbolic expressions in one variablealpha– list of positive integers of lengthself.denominator_ring.ngens()interval– list of positive integersexp_scale– (default: 1) a number
OUTPUT: list of tuples with properties described below
This outputs a list whose entries are a tuple
(r*alpha, a_r, b_r, err_r)forrininterval. Herer*alphais a tuple;a_ris ther*alpha(multi-index) coefficient of the Maclaurin series forselfdivided byexp_scale**r;b_ris a list of the values of the functions inapproxevaluated atrand divided byexp_scale**m;err_ris the list of relative errors(a_r - f)/a_rforfinb_r. All outputs are decimal approximations.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = 1 - x - y - x*y sage: Hfac = H.factor() sage: G = 1 / Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = [1, 1] sage: r = var('r') sage: a1 = (0.573/sqrt(r))*5.83^r sage: a2 = (0.573/sqrt(r) - 0.0674/r^(3/2))*5.83^r sage: es = 5.83 sage: F.relative_error([a1, a2], alpha, [1, 2, 4, 8], es) # long time [((1, 1), 0.5145797599, [0.5730000000, 0.5056000000], [-0.1135300000, 0.01745066667]), ((2, 2), 0.3824778089, [0.4051721856, 0.3813426871], [-0.05933514614, 0.002967810973]), ((4, 4), 0.2778630595, [0.2865000000, 0.2780750000], [-0.03108344267, -0.0007627515584]), ((8, 8), 0.1991088276, [0.2025860928, 0.1996074055], [-0.01746414394, -0.002504047242])]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = Integer(1) - x - y - x*y >>> Hfac = H.factor() >>> G = Integer(1) / Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = [Integer(1), Integer(1)] >>> r = var('r') >>> a1 = (RealNumber('0.573')/sqrt(r))*RealNumber('5.83')**r >>> a2 = (RealNumber('0.573')/sqrt(r) - RealNumber('0.0674')/r**(Integer(3)/Integer(2)))*RealNumber('5.83')**r >>> es = RealNumber('5.83') >>> F.relative_error([a1, a2], alpha, [Integer(1), Integer(2), Integer(4), Integer(8)], es) # long time [((1, 1), 0.5145797599, [0.5730000000, 0.5056000000], [-0.1135300000, 0.01745066667]), ((2, 2), 0.3824778089, [0.4051721856, 0.3813426871], [-0.05933514614, 0.002967810973]), ((4, 4), 0.2778630595, [0.2865000000, 0.2780750000], [-0.03108344267, -0.0007627515584]), ((8, 8), 0.1991088276, [0.2025860928, 0.1996074055], [-0.01746414394, -0.002504047242])]
- singular_ideal()[source]¶
Return the singular ideal of
self.Let \(R\) be the ring of
selfand \(H\) its denominator. Let \(H_{red}\) be the reduction (square-free part) of \(H\). Return the ideal in \(R\) generated by \(H_{red}\) and its partial derivatives. If the coefficient field of \(R\) is algebraically closed, then the output is the ideal of the singular locus (which is a variety) of the variety of \(H\).OUTPUT: an ideal
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y,z> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = (1 - x*(1 + y))^3 * (1 - z*x**2*(1 + 2*y)) sage: df = H.factor() sage: G = 1 / df.unit() sage: F = FFPD(G, df) sage: F.singular_ideal() Ideal (x*y + x - 1, y^2 - 2*y*z + 2*y - z + 1, x*z + y - 2*z + 1) of Multivariate Polynomial Ring in x, y, z over Rational Field
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y', 'z',)); (x, y, z,) = R._first_ngens(3) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = (Integer(1) - x*(Integer(1) + y))**Integer(3) * (Integer(1) - z*x**Integer(2)*(Integer(1) + Integer(2)*y)) >>> df = H.factor() >>> G = Integer(1) / df.unit() >>> F = FFPD(G, df) >>> F.singular_ideal() Ideal (x*y + x - 1, y^2 - 2*y*z + 2*y - z + 1, x*z + y - 2*z + 1) of Multivariate Polynomial Ring in x, y, z over Rational Field
- smooth_critical_ideal(alpha)[source]¶
Return the smooth critical ideal of
self.Let \(R\) be the ring of
selfand \(H\) its denominator. Return the ideal in \(R\) of smooth critical points of the variety of \(H\) for the directionalpha. If the variety \(V\) of \(H\) has no smooth points, then return the ideal in \(R\) of \(V\).See [RW2012] for more details.
INPUT:
alpha– tuple of positive integers and/or symbolic entries of lengthself.denominator_ring.ngens()
OUTPUT: an ideal
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: H = (1 - x - y - x*y)^2 sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = var('a1, a2') sage: F.smooth_critical_ideal(alpha) Ideal (y^2 + (2*a1)/a2*y - 1, x + (-a2)/a1*y + (-a1 + a2)/a1) of Multivariate Polynomial Ring in x, y over Fraction Field of Multivariate Polynomial Ring in a1, a2 over Rational Field sage: H = (1-x-y-x*y)^2 sage: Hfac = H.factor() sage: G = 1/Hfac.unit() sage: F = FFPD(G, Hfac) sage: alpha = [7/3, var('a')] sage: F.smooth_critical_ideal(alpha) Ideal (y^2 + 14/(3*a)*y - 1, x + (-3*a)/7*y + (3*a - 7)/7) of Multivariate Polynomial Ring in x, y over Fraction Field of Univariate Polynomial Ring in a over Rational Field
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> H = (Integer(1) - x - y - x*y)**Integer(2) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = var('a1, a2') >>> F.smooth_critical_ideal(alpha) Ideal (y^2 + (2*a1)/a2*y - 1, x + (-a2)/a1*y + (-a1 + a2)/a1) of Multivariate Polynomial Ring in x, y over Fraction Field of Multivariate Polynomial Ring in a1, a2 over Rational Field >>> H = (Integer(1)-x-y-x*y)**Integer(2) >>> Hfac = H.factor() >>> G = Integer(1)/Hfac.unit() >>> F = FFPD(G, Hfac) >>> alpha = [Integer(7)/Integer(3), var('a')] >>> F.smooth_critical_ideal(alpha) Ideal (y^2 + 14/(3*a)*y - 1, x + (-3*a)/7*y + (3*a - 7)/7) of Multivariate Polynomial Ring in x, y over Fraction Field of Univariate Polynomial Ring in a over Rational Field
- univariate_decomposition()[source]¶
Return the usual univariate partial fraction decomposition of
self.Assume that the numerator of
selflies in the same univariate factorial polynomial ring as the factors of the denominator.Let \(f = p/q\) be a rational expression where \(p\) and \(q\) lie in a univariate factorial polynomial ring \(R\). Let \(q_1^{e_1} \cdots q_n^{e_n}\) be the unique factorization of \(q\) in \(R\) into irreducible factors. Then \(f\) can be written uniquely as:
\[(*) \quad p_0 + \sum_{i=1}^{m} \frac{p_i}{q_i^{e_i}},\]for some \(p_j \in R\). We call \((*)\) the usual partial fraction decomposition of \(f\).
Note
This partial fraction decomposition can be computed using
partial_fraction()orpartial_fraction_decomposition()as well. However, here we use the already obtained/cached factorization of the denominator. This gives a speed up for non-small instances.OUTPUT: an instance of
FractionWithFactoredDenominatorSumEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing
One variable:
sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 5*x^3 + 1/x + 1/(x-1) + 1/(3*x^2 + 1) sage: f (5*x^7 - 5*x^6 + 5/3*x^5 - 5/3*x^4 + 2*x^3 - 2/3*x^2 + 1/3*x - 1/3)/(x^4 - x^3 + 1/3*x^2 - 1/3*x) sage: decomp = FFPD(f).univariate_decomposition() sage: decomp (5*x^3, []) + (1, [(x - 1, 1)]) + (1, [(x, 1)]) + (1/3, [(x^2 + 1/3, 1)]) sage: decomp.sum().quotient() == f True
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(5)*x**Integer(3) + Integer(1)/x + Integer(1)/(x-Integer(1)) + Integer(1)/(Integer(3)*x**Integer(2) + Integer(1)) >>> f (5*x^7 - 5*x^6 + 5/3*x^5 - 5/3*x^4 + 2*x^3 - 2/3*x^2 + 1/3*x - 1/3)/(x^4 - x^3 + 1/3*x^2 - 1/3*x) >>> decomp = FFPD(f).univariate_decomposition() >>> decomp (5*x^3, []) + (1, [(x - 1, 1)]) + (1, [(x, 1)]) + (1/3, [(x^2 + 1/3, 1)]) >>> decomp.sum().quotient() == f True
One variable with numerator in symbolic ring:
sage: R.<x> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: f = 5*x^3 + 1/x + 1/(x-1) + exp(x)/(3*x^2 + 1) sage: f (5*x^5 - 5*x^4 + 2*x - 1)/(x^2 - x) + e^x/(3*x^2 + 1) sage: decomp = FFPD(f).univariate_decomposition() sage: decomp (0, []) + (15/4*x^7 - 15/4*x^6 + 5/4*x^5 - 5/4*x^4 + 3/2*x^3 + 1/4*x^2*e^x - 3/4*x^2 - 1/4*x*e^x + 1/2*x - 1/4, [(x - 1, 1)]) + (-15*x^7 + 15*x^6 - 5*x^5 + 5*x^4 - 6*x^3 - x^2*e^x + 3*x^2 + x*e^x - 2*x + 1, [(x, 1)]) + (1/4*(15*x^7 - 15*x^6 + 5*x^5 - 5*x^4 + 6*x^3 + x^2*e^x - 3*x^2 - x*e^x + 2*x - 1)*(3*x - 1), [(x^2 + 1/3, 1)])
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> f = Integer(5)*x**Integer(3) + Integer(1)/x + Integer(1)/(x-Integer(1)) + exp(x)/(Integer(3)*x**Integer(2) + Integer(1)) >>> f (5*x^5 - 5*x^4 + 2*x - 1)/(x^2 - x) + e^x/(3*x^2 + 1) >>> decomp = FFPD(f).univariate_decomposition() >>> decomp (0, []) + (15/4*x^7 - 15/4*x^6 + 5/4*x^5 - 5/4*x^4 + 3/2*x^3 + 1/4*x^2*e^x - 3/4*x^2 - 1/4*x*e^x + 1/2*x - 1/4, [(x - 1, 1)]) + (-15*x^7 + 15*x^6 - 5*x^5 + 5*x^4 - 6*x^3 - x^2*e^x + 3*x^2 + x*e^x - 2*x + 1, [(x, 1)]) + (1/4*(15*x^7 - 15*x^6 + 5*x^5 - 5*x^4 + 6*x^3 + x^2*e^x - 3*x^2 - x*e^x + 2*x - 1)*(3*x - 1), [(x^2 + 1/3, 1)])
One variable over a finite field:
sage: R.<x> = PolynomialRing(GF(2)) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 5*x^3 + 1/x + 1/(x-1) + 1/(3*x^2 + 1) sage: f (x^6 + x^4 + 1)/(x^3 + x) sage: decomp = FFPD(f).univariate_decomposition() sage: decomp (x^3, []) + (1, [(x, 1)]) + (x, [(x + 1, 2)]) sage: decomp.sum().quotient() == f True
>>> from sage.all import * >>> R = PolynomialRing(GF(Integer(2)), names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(5)*x**Integer(3) + Integer(1)/x + Integer(1)/(x-Integer(1)) + Integer(1)/(Integer(3)*x**Integer(2) + Integer(1)) >>> f (x^6 + x^4 + 1)/(x^3 + x) >>> decomp = FFPD(f).univariate_decomposition() >>> decomp (x^3, []) + (1, [(x, 1)]) + (x, [(x + 1, 2)]) >>> decomp.sum().quotient() == f True
One variable over an inexact field:
sage: R.<x> = PolynomialRing(CC) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = 5*x^3 + 1/x + 1/(x-1) + 1/(3*x^2 + 1) sage: f (5.00000000000000*x^7 - 5.00000000000000*x^6 + 1.66666666666667*x^5 - 1.66666666666667*x^4 + 2.00000000000000*x^3 - 0.666666666666667*x^2 + 0.333333333333333*x - 0.333333333333333)/(x^4 - x^3 + 0.333333333333333*x^2 - 0.333333333333333*x) sage: decomp = FFPD(f).univariate_decomposition() sage: decomp (5.00000000000000*x^3, []) + (1.00000000000000, [(x - 1.00000000000000, 1)]) + (-0.288675134594813*I, [(x - 0.577350269189626*I, 1)]) + (1.00000000000000, [(x, 1)]) + (0.288675134594813*I, [(x + 0.577350269189626*I, 1)]) sage: decomp.sum().quotient() == f # Rounding error coming False
>>> from sage.all import * >>> R = PolynomialRing(CC, names=('x',)); (x,) = R._first_ngens(1) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = Integer(5)*x**Integer(3) + Integer(1)/x + Integer(1)/(x-Integer(1)) + Integer(1)/(Integer(3)*x**Integer(2) + Integer(1)) >>> f (5.00000000000000*x^7 - 5.00000000000000*x^6 + 1.66666666666667*x^5 - 1.66666666666667*x^4 + 2.00000000000000*x^3 - 0.666666666666667*x^2 + 0.333333333333333*x - 0.333333333333333)/(x^4 - x^3 + 0.333333333333333*x^2 - 0.333333333333333*x) >>> decomp = FFPD(f).univariate_decomposition() >>> decomp (5.00000000000000*x^3, []) + (1.00000000000000, [(x - 1.00000000000000, 1)]) + (-0.288675134594813*I, [(x - 0.577350269189626*I, 1)]) + (1.00000000000000, [(x, 1)]) + (0.288675134594813*I, [(x + 0.577350269189626*I, 1)]) >>> decomp.sum().quotient() == f # Rounding error coming False
AUTHORS:
Robert Bradshaw (2007-05-31)
Alexander Raichev (2012-06-25)
Daniel Krenn (2014-12-01)
- class sage.rings.asymptotic.asymptotics_multivariate_generating_functions.FractionWithFactoredDenominatorRing(denominator_ring, numerator_ring=None, category=None)[source]¶
Bases:
UniqueRepresentation,ParentThis is the ring of fractions with factored denominator.
INPUT:
denominator_ring– the base ring (a polynomial ring)numerator_ring– (optional) the numerator ring; the default is thedenominator_ringcategory– (default:Rings) the category
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: df = [x, 1], [y, 1], [x*y+1, 1] sage: f = FFPD(x, df) # indirect doctest sage: f (1, [(y, 1), (x*y + 1, 1)])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> df = [x, Integer(1)], [y, Integer(1)], [x*y+Integer(1), Integer(1)] >>> f = FFPD(x, df) # indirect doctest >>> f (1, [(y, 1), (x*y + 1, 1)])
AUTHORS:
Daniel Krenn (2014-12-01)
- Element[source]¶
alias of
FractionWithFactoredDenominator
- base_ring()[source]¶
Return the base ring.
OUTPUT: a ring
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing sage: P.<X, Y> = ZZ[] sage: F = FractionWithFactoredDenominatorRing(P); F Ring of fractions with factored denominator over Multivariate Polynomial Ring in X, Y over Integer Ring sage: F.base_ring() Integer Ring sage: F.base() Multivariate Polynomial Ring in X, Y over Integer Ring
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing >>> P = ZZ['X, Y']; (X, Y,) = P._first_ngens(2) >>> F = FractionWithFactoredDenominatorRing(P); F Ring of fractions with factored denominator over Multivariate Polynomial Ring in X, Y over Integer Ring >>> F.base_ring() Integer Ring >>> F.base() Multivariate Polynomial Ring in X, Y over Integer Ring
- class sage.rings.asymptotic.asymptotics_multivariate_generating_functions.FractionWithFactoredDenominatorSum(iterable=(), /)[source]¶
Bases:
listA list representing the sum of
FractionWithFactoredDenominatorobjects with distinct denominator factorizations.AUTHORS:
Alexander Raichev (2012-06-25)
Daniel Krenn (2014-12-01)
- property denominator_ring¶
Return the polynomial ring of the denominators of
self.OUTPUT: a ring or
Noneif the list is emptyEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = FFPD(x + y, [(y, 1), (x, 1)]) sage: s = FractionWithFactoredDenominatorSum([f]) sage: s.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field sage: g = FFPD(x + y, []) sage: t = FractionWithFactoredDenominatorSum([g]) sage: t.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = FFPD(x + y, [(y, Integer(1)), (x, Integer(1))]) >>> s = FractionWithFactoredDenominatorSum([f]) >>> s.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field >>> g = FFPD(x + y, []) >>> t = FractionWithFactoredDenominatorSum([g]) >>> t.denominator_ring Multivariate Polynomial Ring in x, y over Rational Field
- sum()[source]¶
Return the sum of the elements in
self.OUTPUT: an instance of
FractionWithFactoredDenominatorEXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: df = (x, 1), (y, 1), (x*y + 1, 1) sage: f = FFPD(2, df) sage: g = FFPD(2*x*y, df) sage: FractionWithFactoredDenominatorSum([f, g]) (2, [(y, 1), (x, 1), (x*y + 1, 1)]) + (2, [(x*y + 1, 1)]) sage: FractionWithFactoredDenominatorSum([f, g]).sum() (2, [(y, 1), (x, 1)]) sage: f = FFPD(cos(x), [(x, 2)]) sage: g = FFPD(cos(y), [(x, 1), (y, 2)]) sage: FractionWithFactoredDenominatorSum([f, g]) (cos(x), [(x, 2)]) + (cos(y), [(y, 2), (x, 1)]) sage: FractionWithFactoredDenominatorSum([f, g]).sum() (y^2*cos(x) + x*cos(y), [(y, 2), (x, 2)])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> df = (x, Integer(1)), (y, Integer(1)), (x*y + Integer(1), Integer(1)) >>> f = FFPD(Integer(2), df) >>> g = FFPD(Integer(2)*x*y, df) >>> FractionWithFactoredDenominatorSum([f, g]) (2, [(y, 1), (x, 1), (x*y + 1, 1)]) + (2, [(x*y + 1, 1)]) >>> FractionWithFactoredDenominatorSum([f, g]).sum() (2, [(y, 1), (x, 1)]) >>> f = FFPD(cos(x), [(x, Integer(2))]) >>> g = FFPD(cos(y), [(x, Integer(1)), (y, Integer(2))]) >>> FractionWithFactoredDenominatorSum([f, g]) (cos(x), [(x, 2)]) + (cos(y), [(y, 2), (x, 1)]) >>> FractionWithFactoredDenominatorSum([f, g]).sum() (y^2*cos(x) + x*cos(y), [(y, 2), (x, 2)])
- whole_and_parts()[source]¶
Rewrite
selfas a sum of a (possibly zero) polynomial followed by reduced rational expressions.OUTPUT: an instance of
FractionWithFactoredDenominatorSumOnly useful for multivariate decompositions.
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R, SR) sage: f = x**2 + 3*y + 1/x + 1/y sage: f = FFPD(f); f (x^3*y + 3*x*y^2 + x + y, [(y, 1), (x, 1)]) sage: FractionWithFactoredDenominatorSum([f]).whole_and_parts() (x^2 + 3*y, []) + (x + y, [(y, 1), (x, 1)]) sage: f = cos(x)**2 + 3*y + 1/x + 1/y; f cos(x)^2 + 3*y + 1/x + 1/y sage: G = f.numerator() sage: H = R(f.denominator()) sage: f = FFPD(G, H.factor()); f (x*y*cos(x)^2 + 3*x*y^2 + x + y, [(y, 1), (x, 1)]) sage: FractionWithFactoredDenominatorSum([f]).whole_and_parts() (0, []) + (x*y*cos(x)^2 + 3*x*y^2 + x + y, [(y, 1), (x, 1)])
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, FractionWithFactoredDenominatorSum >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R, SR) >>> f = x**Integer(2) + Integer(3)*y + Integer(1)/x + Integer(1)/y >>> f = FFPD(f); f (x^3*y + 3*x*y^2 + x + y, [(y, 1), (x, 1)]) >>> FractionWithFactoredDenominatorSum([f]).whole_and_parts() (x^2 + 3*y, []) + (x + y, [(y, 1), (x, 1)]) >>> f = cos(x)**Integer(2) + Integer(3)*y + Integer(1)/x + Integer(1)/y; f cos(x)^2 + 3*y + 1/x + 1/y >>> G = f.numerator() >>> H = R(f.denominator()) >>> f = FFPD(G, H.factor()); f (x*y*cos(x)^2 + 3*x*y^2 + x + y, [(y, 1), (x, 1)]) >>> FractionWithFactoredDenominatorSum([f]).whole_and_parts() (0, []) + (x*y*cos(x)^2 + 3*x*y^2 + x + y, [(y, 1), (x, 1)])
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.coerce_point(R, p)[source]¶
Coerce the keys of the dictionary
pinto the ringR.Warning
This method assumes that it is possible.
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, coerce_point sage: R.<x,y> = PolynomialRing(QQ) sage: FFPD = FractionWithFactoredDenominatorRing(R) sage: f = FFPD() sage: p = {SR(x): 1, SR(y): 7/8} sage: for k in sorted(p, key=str): ....: print("{} {} {}".format(k, k.parent(), p[k])) x Symbolic Ring 1 y Symbolic Ring 7/8 sage: q = coerce_point(R, p) sage: for k in sorted(q, key=str): ....: print("{} {} {}".format(k, k.parent(), q[k])) x Multivariate Polynomial Ring in x, y over Rational Field 1 y Multivariate Polynomial Ring in x, y over Rational Field 7/8
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import FractionWithFactoredDenominatorRing, coerce_point >>> R = PolynomialRing(QQ, names=('x', 'y',)); (x, y,) = R._first_ngens(2) >>> FFPD = FractionWithFactoredDenominatorRing(R) >>> f = FFPD() >>> p = {SR(x): Integer(1), SR(y): Integer(7)/Integer(8)} >>> for k in sorted(p, key=str): ... print("{} {} {}".format(k, k.parent(), p[k])) x Symbolic Ring 1 y Symbolic Ring 7/8 >>> q = coerce_point(R, p) >>> for k in sorted(q, key=str): ... print("{} {} {}".format(k, k.parent(), q[k])) x Multivariate Polynomial Ring in x, y over Rational Field 1 y Multivariate Polynomial Ring in x, y over Rational Field 7/8
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.diff_all(f, V, n, ending=[], sub=None, sub_final=None, zero_order=0, rekey=None)[source]¶
Return a dictionary of representative mixed partial derivatives of \(f\) from order 1 up to order \(n\) with respect to the variables in \(V\).
The default is to key the dictionary by all nondecreasing sequences in \(V\) of length 1 up to length \(n\).
INPUT:
f– an individual or list of \(\mathcal{C}^{n+1}\) functionsV– list of variables occurring in \(f\)n– a natural numberending– list of variables in \(V\)sub– an individual or list of dictionariessub_final– an individual or list of dictionariesrekey– a callable symbolic function in \(V\) or list thereofzero_order– a natural number
OUTPUT: the dictionary
{s_1:deriv_1, ..., sr:deriv_r}Here
s_1, ..., s_ris a listing of all nondecreasing sequences of length 1 up to length \(n\) over the alphabet \(V\), where \(w > v\) in \(X\) if and only ifstr(w) > str(v), andderiv_jis the derivative of \(f\) with respect to the derivative sequences_jand simplified with respect to the substitutions insuband evaluated atsub_final. Moreover, all derivatives with respect to sequences of length less thanzero_order(derivatives of order less thanzero_order) will be made zero.If
rekeyis nonempty, thens_1, ..., s_rwill be replaced by the symbolic derivatives of the functions inrekey.If
endingis nonempty, then every derivative sequences_jwill be suffixed byending.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_all sage: f = function('f')(x) sage: dd = diff_all(f, [x], 3) sage: dd[(x, x, x)] diff(f(x), x, x, x) sage: d1 = {diff(f, x): 4*x^3} sage: dd = diff_all(f, [x], 3, sub=d1) sage: dd[(x, x, x)] 24*x sage: dd = diff_all(f, [x], 3, sub=d1, rekey=f) sage: dd[diff(f, x, 3)] 24*x sage: a = {x:1} sage: dd = diff_all(f, [x], 3, sub=d1, rekey=f, sub_final=a) sage: dd[diff(f, x, 3)] 24
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_all >>> f = function('f')(x) >>> dd = diff_all(f, [x], Integer(3)) >>> dd[(x, x, x)] diff(f(x), x, x, x) >>> d1 = {diff(f, x): Integer(4)*x**Integer(3)} >>> dd = diff_all(f, [x], Integer(3), sub=d1) >>> dd[(x, x, x)] 24*x >>> dd = diff_all(f, [x], Integer(3), sub=d1, rekey=f) >>> dd[diff(f, x, Integer(3))] 24*x >>> a = {x:Integer(1)} >>> dd = diff_all(f, [x], Integer(3), sub=d1, rekey=f, sub_final=a) >>> dd[diff(f, x, Integer(3))] 24
sage: X = var('x, y, z') sage: f = function('f')(*X) sage: dd = diff_all(f, X, 2, ending=[y, y, y]) sage: dd[(z, y, y, y)] diff(f(x, y, z), y, y, y, z)
>>> from sage.all import * >>> X = var('x, y, z') >>> f = function('f')(*X) >>> dd = diff_all(f, X, Integer(2), ending=[y, y, y]) >>> dd[(z, y, y, y)] diff(f(x, y, z), y, y, y, z)
sage: g = function('g')(*X) sage: dd = diff_all([f, g], X, 2) sage: dd[(0, y, z)] diff(f(x, y, z), y, z) sage: dd[(1, z, z)] diff(g(x, y, z), z, z) sage: f = exp(x*y*z) sage: ff = function('ff')(*X) sage: dd = diff_all(f, X, 2, rekey=ff) sage: dd[diff(ff, x, z)] x*y^2*z*e^(x*y*z) + y*e^(x*y*z)
>>> from sage.all import * >>> g = function('g')(*X) >>> dd = diff_all([f, g], X, Integer(2)) >>> dd[(Integer(0), y, z)] diff(f(x, y, z), y, z) >>> dd[(Integer(1), z, z)] diff(g(x, y, z), z, z) >>> f = exp(x*y*z) >>> ff = function('ff')(*X) >>> dd = diff_all(f, X, Integer(2), rekey=ff) >>> dd[diff(ff, x, z)] x*y^2*z*e^(x*y*z) + y*e^(x*y*z)
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.diff_op(A, B, AB_derivs, V, M, r, N)[source]¶
Return the derivatives \(DD^{(l+k)}(A[j] B^l)\) evaluated at a point \(p\) for various natural numbers \(j, k, l\) which depend on \(r\) and \(N\).
Here \(DD\) is a specific second-order linear differential operator that depends on \(M\) , \(A\) is a list of symbolic functions, \(B\) is symbolic function, and
AB_derivscontains all the derivatives of \(A\) and \(B\) evaluated at \(p\) that are necessary for the computation.INPUT:
A– a single or lengthrlist of symbolic functions in the variablesVB– a symbolic function in the variablesVAB_derivs– dictionary whose keys are the (symbolic) derivatives ofA[0], ..., A[r-1]up to order2 * N-2and the (symbolic) derivatives ofBup to order2 * N; the values of the dictionary are complex numbers that are the keys evaluated at a common point \(p\)V– the variables of theA[j]andBM– a symmetric \(l \times l\) matrix, where \(l\) is the length ofVr,N– natural numbers
OUTPUT: a dictionary
The output is a dictionary whose keys are natural number tuples of the form \((j, k, l)\), where \(l \leq 2k\), \(j \leq r-1\), and \(j+k \leq N-1\), and whose values are \(DD^(l+k)(A[j] B^l)\) evaluated at a point \(p\), where \(DD\) is the linear second-order differential operator \(-\sum_{i=0}^{l-1} \sum_{j=0}^{l-1} M[i][j] \partial^2 /(\partial V[j] \partial V[i])\).
Note
For internal use by
FractionWithFactoredDenominator.asymptotics_smooth()andFractionWithFactoredDenominator.asymptotics_multiple().EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_op sage: T = var('x, y') sage: A = function('A')(*tuple(T)) sage: B = function('B')(*tuple(T)) sage: AB_derivs = {} sage: M = matrix([[1, 2],[2, 1]]) sage: DD = diff_op(A, B, AB_derivs, T, M, 1, 2) # long time (see :issue:`35207`) sage: sorted(DD) # long time [(0, 0, 0), (0, 1, 0), (0, 1, 1), (0, 1, 2)] sage: DD[(0, 1, 2)].number_of_operands() # long time 246
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_op >>> T = var('x, y') >>> A = function('A')(*tuple(T)) >>> B = function('B')(*tuple(T)) >>> AB_derivs = {} >>> M = matrix([[Integer(1), Integer(2)],[Integer(2), Integer(1)]]) >>> DD = diff_op(A, B, AB_derivs, T, M, Integer(1), Integer(2)) # long time (see :issue:`35207`) >>> sorted(DD) # long time [(0, 0, 0), (0, 1, 0), (0, 1, 1), (0, 1, 2)] >>> DD[(Integer(0), Integer(1), Integer(2))].number_of_operands() # long time 246
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.diff_op_simple(A, B, AB_derivs, x, v, a, N)[source]¶
Return \(DD^(e k + v l)(A B^l)\) evaluated at a point \(p\) for various natural numbers \(e, k, l\) that depend on \(v\) and \(N\).
Here \(DD\) is a specific linear differential operator that depends on \(a\) and \(v\) , \(A\) and \(B\) are symbolic functions, and
AB_derivscontains all the derivatives of \(A\) and \(B\) evaluated at \(p\) that are necessary for the computation.Note
For internal use by the function
FractionWithFactoredDenominator.asymptotics_smooth().INPUT:
A,B– symbolic functions in the variablexAB_derivs– dictionary whose keys are the (symbolic) derivatives ofAup to order2 * Nifvis even orNifvis odd and the (symbolic) derivatives ofBup to order2 * N + vifvis even orN + vifvis odd; the values of the dictionary are complex numbers that are the keys evaluated at a common point \(p\)x– a symbolic variablea– a complex numberv,N– natural numbers
OUTPUT: a dictionary
The output is a dictionary whose keys are natural number pairs of the form \((k, l)\), where \(k < N\) and \(l \leq 2k\) and whose values are \(DD^(e k + v l)(A B^l)\) evaluated at a point \(p\). Here \(e=2\) if \(v\) is even, \(e=1\) if \(v\) is odd, and \(DD\) is the linear differential operator \((a^{-1/v} d/dt)\) if \(v\) is even and \((|a|^{-1/v} i \text{sgn}(a) d/dt)\) if \(v\) is odd.
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_op_simple sage: A = function('A')(x) sage: B = function('B')(x) sage: AB_derivs = {} sage: sorted(diff_op_simple(A, B, AB_derivs, x, 3, 2, 2).items()) [((0, 0), A(x)), ((1, 0), 1/2*I*2^(2/3)*diff(A(x), x)), ((1, 1), 1/4*2^(2/3)*(B(x)*diff(A(x), x, x, x, x) + 4*diff(A(x), x, x, x)*diff(B(x), x) + 6*diff(A(x), x, x)*diff(B(x), x, x) + 4*diff(A(x), x)*diff(B(x), x, x, x) + A(x)*diff(B(x), x, x, x, x)))]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_op_simple >>> A = function('A')(x) >>> B = function('B')(x) >>> AB_derivs = {} >>> sorted(diff_op_simple(A, B, AB_derivs, x, Integer(3), Integer(2), Integer(2)).items()) [((0, 0), A(x)), ((1, 0), 1/2*I*2^(2/3)*diff(A(x), x)), ((1, 1), 1/4*2^(2/3)*(B(x)*diff(A(x), x, x, x, x) + 4*diff(A(x), x, x, x)*diff(B(x), x) + 6*diff(A(x), x, x)*diff(B(x), x, x) + 4*diff(A(x), x)*diff(B(x), x, x, x) + A(x)*diff(B(x), x, x, x, x)))]
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.diff_prod(f_derivs, u, g, X, interval, end, uderivs, atc)[source]¶
Take various derivatives of the equation \(f = ug\), evaluate them at a point \(c\), and solve for the derivatives of \(u\).
INPUT:
f_derivs– dictionary whose keys are all tuples of the forms + end, wheresis a sequence of variables fromXwhose length lies ininterval, and whose values are the derivatives of a function \(f\) evaluated at \(c\)u– a callable symbolic functiong– an expression or callable symbolic functionX– list of symbolic variablesinterval– list of positive integers Call the first and last values \(n\) and \(nn\), respectivelyend– a possibly empty list of repetitions of the variablez, wherezis the last element ofXuderivs– dictionary whose keys are the symbolic derivatives of order 0 to order \(n-1\) ofuevaluated at \(c\) and whose values are the corresponding derivatives evaluated at \(c\)atc– dictionary whose keys are the keys of \(c\) and all the symbolic derivatives of order 0 to order \(nn\) ofgevaluated \(c\) and whose values are the corresponding derivatives evaluated at \(c\)
OUTPUT:
A dictionary whose keys are the derivatives of
uup to order \(nn\) and whose values are those derivatives evaluated at \(c\).This function works by differentiating the equation \(f = ug\) with respect to the variable sequence
s + end, for all tuplessofXof lengths ininterval, evaluating at the point \(c\) , and solving for the remaining derivatives ofu. This function assumes thatunever appears in the differentiations of \(f = ug\) after evaluating at \(c\).Note
For internal use by
FractionWithFactoredDenominator.asymptotics_multiple().EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_prod sage: u = function('u')(x) sage: g = function('g')(x) sage: fd = {(x,):1,(x, x):1} sage: ud = {u(x=2): 1} sage: atc = {x: 2, g(x=2): 3, diff(g, x)(x=2): 5} sage: atc[diff(g, x, x)(x=2)] = 7 sage: dd = diff_prod(fd, u, g, [x], [1, 2], [], ud, atc) sage: dd[diff(u, x, 2)(x=2)] 22/9
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_prod >>> u = function('u')(x) >>> g = function('g')(x) >>> fd = {(x,):Integer(1),(x, x):Integer(1)} >>> ud = {u(x=Integer(2)): Integer(1)} >>> atc = {x: Integer(2), g(x=Integer(2)): Integer(3), diff(g, x)(x=Integer(2)): Integer(5)} >>> atc[diff(g, x, x)(x=Integer(2))] = Integer(7) >>> dd = diff_prod(fd, u, g, [x], [Integer(1), Integer(2)], [], ud, atc) >>> dd[diff(u, x, Integer(2))(x=Integer(2))] 22/9
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.diff_seq(V, s)[source]¶
Given a list
sof tuples of natural numbers, return the list of elements ofVwith indices the elements of the elements ofs.INPUT:
V– lists– list of tuples of natural numbers in the intervalrange(len(V))
OUTPUT:
The tuple
tuple([V[tt] for tt in sorted(t)]), wheretis the list of elements of the elements ofs.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_seq sage: V = list(var('x, t, z')) sage: diff_seq(V,([0, 1],[0, 2, 1],[0, 0])) (x, x, x, x, t, t, z)
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import diff_seq >>> V = list(var('x, t, z')) >>> diff_seq(V,([Integer(0), Integer(1)],[Integer(0), Integer(2), Integer(1)],[Integer(0), Integer(0)])) (x, x, x, x, t, t, z)
Note
This function is for internal use by
diff_op().
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.direction(v, coordinate=None)[source]¶
Return
[vv/v[coordinate] for vv in v]wherecoordinateis the last index ofvif not specified otherwise.INPUT:
v– a vectorcoordinate– (default:None) an index forv
EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import direction sage: direction([2, 3, 5]) (2/5, 3/5, 1) sage: direction([2, 3, 5], 0) (1, 3/2, 5/2)
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import direction >>> direction([Integer(2), Integer(3), Integer(5)]) (2/5, 3/5, 1) >>> direction([Integer(2), Integer(3), Integer(5)], Integer(0)) (1, 3/2, 5/2)
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.permutation_sign(s, u)[source]¶
This function returns the sign of the permutation on
1, ..., len(u)that is induced by the sublistsofu.Note
This function was intended for internal use and is deprecated now (Issue #29465).
INPUT:
s– a sublist ofuu– list
OUTPUT:
The sign of the permutation obtained by taking indices within
uof the lists + sc, wherescisuwith the elements ofsremoved.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import permutation_sign sage: u = ['a', 'b', 'c', 'd', 'e'] sage: s = ['b', 'd'] sage: permutation_sign(s, u) doctest:...: DeprecationWarning: the function permutation_sign is deprecated See https://github.com/sagemath/sage/issues/29465 for details. -1 sage: s = ['d', 'b'] sage: permutation_sign(s, u) 1
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import permutation_sign >>> u = ['a', 'b', 'c', 'd', 'e'] >>> s = ['b', 'd'] >>> permutation_sign(s, u) doctest:...: DeprecationWarning: the function permutation_sign is deprecated See https://github.com/sagemath/sage/issues/29465 for details. -1 >>> s = ['d', 'b'] >>> permutation_sign(s, u) 1
- sage.rings.asymptotic.asymptotics_multivariate_generating_functions.subs_all(f, sub, simplify=False)[source]¶
Return the items of \(f\) substituted by the dictionaries of
subin order of their appearance insub.INPUT:
f– an individual or list of symbolic expressions or dictionariessub– an individual or list of dictionariessimplify– boolean (default:False); set toTrueto simplify the result
OUTPUT:
The items of
fsubstituted by the dictionaries ofsubin order of their appearance insub. Thesubs()command is used. If simplify isTrue, thensimplify()is used after substitution.EXAMPLES:
sage: from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import subs_all sage: var('x, y, z') (x, y, z) sage: a = {x:1} sage: b = {y:2} sage: c = {z:3} sage: subs_all(x + y + z, a) y + z + 1 sage: subs_all(x + y + z, [c, a]) y + 4 sage: subs_all([x + y + z, y^2], b) [x + z + 2, 4] sage: subs_all([x + y + z, y^2], [b, c]) [x + 5, 4]
>>> from sage.all import * >>> from sage.rings.asymptotic.asymptotics_multivariate_generating_functions import subs_all >>> var('x, y, z') (x, y, z) >>> a = {x:Integer(1)} >>> b = {y:Integer(2)} >>> c = {z:Integer(3)} >>> subs_all(x + y + z, a) y + z + 1 >>> subs_all(x + y + z, [c, a]) y + 4 >>> subs_all([x + y + z, y**Integer(2)], b) [x + z + 2, 4] >>> subs_all([x + y + z, y**Integer(2)], [b, c]) [x + 5, 4]
sage: var('x, y') (x, y) sage: a = {'foo': x**2 + y**2, 'bar': x - y} sage: b = {x: 1, y: 2} sage: subs_all(a, b) {'bar': -1, 'foo': 5}
>>> from sage.all import * >>> var('x, y') (x, y) >>> a = {'foo': x**Integer(2) + y**Integer(2), 'bar': x - y} >>> b = {x: Integer(1), y: Integer(2)} >>> subs_all(a, b) {'bar': -1, 'foo': 5}