diff --git a/src/sage/data_structures/stream.py b/src/sage/data_structures/stream.py index 87b5b1e7ecb..673c78e1476 100644 --- a/src/sage/data_structures/stream.py +++ b/src/sage/data_structures/stream.py @@ -3258,3 +3258,323 @@ def is_nonzero(self): True """ return self._series.is_nonzero() + + +class Stream_infinite_operator(Stream): + r""" + Stream defined by applying an operator an infinite number of times. + + The ``iterator`` returns elements `s_i` to compute an infinite operator. + The valuation of `s_i` is weakly increasing as we iterate over `I` and + there are only finitely many terms with any fixed valuation. + In particular, this *assumes* the result is nonzero. + + .. WARNING:: + + This does not check that the input is valid. + + INPUT: + + - ``iterator`` -- the iterator for the factors + """ + def __init__(self, iterator): + r""" + Initialize ``self``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + """ + self._op_iter = iterator + self._cur = None + self._cur_order = -infinity + super().__init__(False) + + @lazy_attribute + def _approximate_order(self): + r""" + Compute and return the approximate order of ``self``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum, Stream_infinite_product + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f._approximate_order + 1 + sage: it = (t^n for n in PositiveIntegers()) + sage: f = Stream_infinite_product(it) + sage: f._approximate_order + 0 + sage: it = (t^(n-10) for n in PositiveIntegers()) + sage: f = Stream_infinite_product(it) + sage: f._approximate_order + -45 + """ + if self._cur is None: + self._advance() + while self._cur_order <= 0: + self._advance() + return self._cur._coeff_stream._approximate_order + + def _advance(self): + r""" + Advance the iterator so that the approximate order increases + by at least one. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (n * t^n for n in range(10)) + sage: f = Stream_infinite_sum(it) + sage: f._cur is None + True + sage: f._advance() + sage: f._cur + t + 2*t^2 + sage: f._cur_order + 2 + sage: for _ in range(20): + ....: f._advance() + sage: f._cur + t + 2*t^2 + 3*t^3 + 4*t^4 + 5*t^5 + 6*t^6 + 7*t^7 + 8*t^8 + 9*t^9 + sage: f._cur_order + +Infinity + + sage: it = (t^(n//3) / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f._advance() + sage: f._cur + 2 + 3*t + 3*t^2 + 3*t^3 + O(t^4) + sage: f._advance() + sage: f._cur + 2 + 5*t + 6*t^2 + 6*t^3 + 6*t^4 + O(t^5) + """ + if self._cur is None: + temp = next(self._op_iter) + if isinstance(temp._coeff_stream, Stream_zero): + self._advance() + return + self.initial(temp) + self._cur_order = temp._coeff_stream._approximate_order + + order = self._cur_order + while order == self._cur_order: + try: + next_factor = next(self._op_iter) + except StopIteration: + self._cur_order = infinity + return + if isinstance(next_factor._coeff_stream, Stream_zero): + continue + coeff_stream = next_factor._coeff_stream + while coeff_stream._approximate_order < order: + # This check also updates the next_factor._approximate_order + if coeff_stream[coeff_stream._approximate_order]: + order = coeff_stream._approximate_order + raise ValueError(f"invalid product computation with invalid order {order} < {self._cur_order}") + self.apply_operator(next_factor) + order = coeff_stream._approximate_order + # We check to see if we need to increment the order + if order == self._cur_order and not coeff_stream[order]: + order += 1 + self._cur_order = order + + def __getitem__(self, n): + r""" + Return the ``n``-th coefficient of ``self``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f[2] + 2 + sage: f[5] + 5 + """ + while n >= self._cur_order: + self._advance() + return self._cur[n] + + def order(self): + r""" + Return the order of ``self``, which is the minimum index ``n`` such + that ``self[n]`` is nonzero. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^(5+n) / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f.order() + 6 + """ + return self._approximate_order + + def __hash__(self): + r""" + Return the hash of ``self``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: hash(f) == hash((type(f), f._op_iter)) + True + """ + return hash((type(self), self._op_iter)) + + def __ne__(self, other): + r""" + Return whether ``self`` and ``other`` are known to be equal. + + INPUT: + + - ``other`` -- a stream + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: itf = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(itf) + sage: itg = (t^(2*n-1) / (1 - t) for n in PositiveIntegers()) + sage: g = Stream_infinite_sum(itg) + sage: f != g + False + sage: f[10] + 10 + sage: g[10] + 5 + sage: f != g + True + """ + if not isinstance(other, type(self)): + return True + ao = min(self._approximate_order, other._approximate_order) + if any(self[i] != other[i] for i in range(ao, min(self._cur_order, other._cur_order))): + return True + return False + + def is_nonzero(self): + r""" + Return ``True`` if and only if this stream is known + to be nonzero. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f.is_nonzero() + True + """ + return True + + +class Stream_infinite_sum(Stream_infinite_operator): + r""" + Stream defined by an infinite sum. + + The ``iterator`` returns elements `s_i` to compute the product + `\sum_{i \in I} s_i`. See :class:`Stream_infinite_operator` + for restrictions on the `s_i`. + + INPUT: + + - ``iterator`` -- the iterator for the factors + """ + def initial(self, obj): + r""" + Set the initial data. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f._cur is None + True + sage: f._advance() # indirect doctest + sage: f._cur + t + 2*t^2 + 2*t^3 + 2*t^4 + O(t^5) + """ + self._cur = obj + + def apply_operator(self, next_obj): + r""" + Apply the operator to ``next_obj``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_sum + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^(n//2) / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_sum(it) + sage: f._advance() + sage: f._advance() # indirect doctest + sage: f._cur + 1 + 3*t + 4*t^2 + 4*t^3 + 4*t^4 + O(t^5) + """ + self._cur += next_obj + + +class Stream_infinite_product(Stream_infinite_operator): + r""" + Stream defined by an infinite product. + + The ``iterator`` returns elements `p_i` to compute the product + `\prod_{i \in I} (1 + p_i)`. See :class:`Stream_infinite_operator` + for restrictions on the `p_i`. + + INPUT: + + - ``iterator`` -- the iterator for the factors + """ + def initial(self, obj): + r""" + Set the initial data. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_product + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_product(it) + sage: f._cur is None + True + sage: f._advance() # indirect doctest + sage: f._cur + 1 + t + 2*t^2 + 3*t^3 + 4*t^4 + 5*t^5 + 6*t^6 + O(t^7) + """ + self._cur = obj + 1 + + def apply_operator(self, next_obj): + r""" + Apply the operator to ``next_obj``. + + EXAMPLES:: + + sage: from sage.data_structures.stream import Stream_infinite_product + sage: L. = LazyLaurentSeriesRing(QQ) + sage: it = (t^n / (1 - t) for n in PositiveIntegers()) + sage: f = Stream_infinite_product(it) + sage: f._advance() + sage: f._advance() # indirect doctest + sage: f._cur + 1 + t + 2*t^2 + 4*t^3 + 6*t^4 + 9*t^5 + 13*t^6 + O(t^7) + """ + self._cur = self._cur + self._cur * next_obj diff --git a/src/sage/rings/lazy_series_ring.py b/src/sage/rings/lazy_series_ring.py index 32545f65574..a4fb39436f3 100644 --- a/src/sage/rings/lazy_series_ring.py +++ b/src/sage/rings/lazy_series_ring.py @@ -55,6 +55,7 @@ from sage.misc.lazy_attribute import lazy_attribute from sage.rings.integer_ring import ZZ +from sage.rings.infinity import infinity from sage.rings.polynomial.laurent_polynomial_ring import LaurentPolynomialRing from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing from sage.rings.lazy_series import (LazyModuleElement, @@ -890,6 +891,193 @@ def is_exact(self): """ return self.base_ring().is_exact() + def prod(self, f, a=None, b=infinity, add_one=False): + r""" + The product of elements of ``self``. + + INPUT: + + - ``f`` -- a list (or iterable) of elements of ``self`` + - ``a``, ``b`` -- optional arguments + - ``add_one`` -- (default: ``False``); if ``True``, then converts a + lazy series `p_i` from ``args`` into `1 + p_i` for the product + + If ``a`` and ``b`` are both integers, then this returns the product + `\prod_{i=a}^b f(i)`, where `f(i) = p_i` if ``add_one=False`` or + `f(i) = 1 + p_i` otherwise. If ``b`` is not specified, then we consider + `b = \infty`. Note this corresponds to the Python ``range(a, b+1)``. + + If `a` is any other iterable, then this returns the product + `\prod_{i \in a} f(i)`, where `f(i) = p_i` if ``add_one=False`` or + `f(i) = 1 + p_i`. + + .. NOTE:: + + For infinite products, it is faster to use ``add_one=True`` since + the implementation is based on `p_i` in `\prod_i (1 + p_i)`. + + .. WARNING:: + + When ``f`` is an infinite generator, then the first argument + ``a`` must be ``True``. Otherwise this will loop forever. + + .. WARNING:: + + For an *infinite* product of the form `\prod_i (1 + p_i)`, + if `p_i = 0`, then this will loop forever. + + EXAMPLES:: + + sage: L. = LazyLaurentSeriesRing(QQ) + sage: euler = L.prod(lambda n: 1 - t^n, PositiveIntegers()) + sage: euler + 1 - t - t^2 + t^5 + O(t^7) + sage: 1 / euler + 1 + t + 2*t^2 + 3*t^3 + 5*t^4 + 7*t^5 + 11*t^6 + O(t^7) + sage: euler - L.euler() + O(t^7) + sage: L.prod(lambda n: -t^n, 1, add_one=True) + 1 - t - t^2 + t^5 + O(t^7) + + sage: L.prod((1 - t^n for n in PositiveIntegers()), True) + 1 - t - t^2 + t^5 + O(t^7) + sage: L.prod((-t^n for n in PositiveIntegers()), True, add_one=True) + 1 - t - t^2 + t^5 + O(t^7) + + sage: L.prod((1 + t^(n-3) for n in PositiveIntegers()), True) + 2*t^-3 + 4*t^-2 + 4*t^-1 + 4 + 6*t + 10*t^2 + 16*t^3 + O(t^4) + + sage: L.prod(lambda n: 2 + t^n, -3, 5) + 96*t^-6 + 240*t^-5 + 336*t^-4 + 840*t^-3 + 984*t^-2 + 1248*t^-1 + + 1980 + 1668*t + 1824*t^2 + 1872*t^3 + 1782*t^4 + 1710*t^5 + + 1314*t^6 + 1122*t^7 + 858*t^8 + 711*t^9 + 438*t^10 + 282*t^11 + + 210*t^12 + 84*t^13 + 60*t^14 + 24*t^15 + sage: L.prod(lambda n: t^n / (1 + abs(n)), -2, 2, add_one=True) + 1/3*t^-3 + 5/6*t^-2 + 13/9*t^-1 + 25/9 + 13/9*t + 5/6*t^2 + 1/3*t^3 + sage: L.prod(lambda n: t^-2 + t^n / n, -4, -2) + 1/24*t^-9 - 1/8*t^-8 - 1/6*t^-7 + 1/2*t^-6 + + sage: D = LazyDirichletSeriesRing(QQ, "s") + sage: D.prod(lambda p: (1+D(1, valuation=p)).inverse(), Primes()) + 1 - 1/(2^s) - 1/(3^s) + 1/(4^s) - 1/(5^s) + 1/(6^s) - 1/(7^s) + O(1/(8^s)) + + sage: D.prod(lambda p: D(1, valuation=p), Primes(), add_one=True) + 1 + 1/(2^s) + 1/(3^s) + 1/(5^s) + 1/(6^s) + 1/(7^s) + O(1/(8^s)) + """ + if a is None: + if add_one: + return super().prod(self.one() + g for g in f) + return super().prod(f) + + if a is True: + it = f + elif a in ZZ: + if b != infinity: + if add_one: + return super().prod(self.one() + f(i) for i in range(a, b+1)) + return super().prod(f(i) for i in range(a, b+1)) + from sage.sets.non_negative_integers import NonNegativeIntegers + it = (f(i+a) for i in NonNegativeIntegers()) + else: + it = (f(i) for i in a) + + # NOTE: We must have a new variable name for each new iterator + if not add_one: + data = (g - self.one() for g in it) + else: + data = it + + from sage.data_structures.stream import Stream_infinite_product + coeff_stream = Stream_infinite_product(data) + return self.element_class(self, coeff_stream) + + def sum(self, f, a=None, b=infinity): + r""" + The sum of elements of ``self``. + + INPUT: + + - ``f`` -- a list (or iterable or function) of elements of ``self`` + - ``a``, ``b`` -- optional arguments + + If ``a`` and ``b`` are both integers, then this returns the sum + `\sum_{i=a}^b f(i)`. If ``b`` is not specified, then we consider + `b = \infty`. Note this corresponds to the Python ``range(a, b+1)``. + + If `a` is any other iterable, then this returns the sum + `\sum{i \in a} f(i)`. + + .. WARNING:: + + When ``f`` is an infinite generator, then the first argument + ``a`` must be ``True``. Otherwise this will loop forever. + + .. WARNING:: + + For an *infinite* sum of the form `\sum_i s_i`, + if `s_i = 0`, then this will loop forever. + + EXAMPLES:: + + sage: L. = LazyLaurentSeriesRing(QQ) + sage: L.sum(lambda n: t^n / (n+1), PositiveIntegers()) + 1/2*t + 1/3*t^2 + 1/4*t^3 + 1/5*t^4 + 1/6*t^5 + 1/7*t^6 + 1/8*t^7 + O(t^8) + + sage: L. = LazyPowerSeriesRing(QQ) + sage: T = L.undefined(1) + sage: D = L.undefined(0) + sage: H = L.sum(lambda k: T(z^k)/k, 2) + sage: T.define(z*exp(T)*D) + sage: D.define(exp(H)) + sage: T + z + z^2 + 2*z^3 + 4*z^4 + 9*z^5 + 20*z^6 + 48*z^7 + O(z^8) + sage: D + 1 + 1/2*z^2 + 1/3*z^3 + 7/8*z^4 + 11/30*z^5 + 281/144*z^6 + O(z^7) + + We verify the Rogers-Ramanujan identities up to degree 100:: + + sage: L. = LazyPowerSeriesRing(QQ) + sage: Gpi = L.prod(lambda k: -q^(1+5*k), 0, oo, add_one=True) + sage: Gpi *= L.prod(lambda k: -q^(4+5*k), 0, oo, add_one=True) + sage: Gp = 1 / Gpi + sage: G = L.sum(lambda n: q^(n^2) / prod(1 - q^(k+1) for k in range(n)), 0, oo) + sage: G - Gp + O(q^7) + sage: all(G[k] == Gp[k] for k in range(100)) + True + + sage: Hpi = L.prod(lambda k: -q^(2+5*k), 0, oo, add_one=True) + sage: Hpi *= L.prod(lambda k: -q^(3+5*k), 0, oo, add_one=True) + sage: Hp = 1 / Hpi + sage: H = L.sum(lambda n: q^(n^2+n) / prod(1 - q^(k+1) for k in range(n)), 0, oo) + sage: H - Hp + O(q^7) + sage: all(H[k] == Hp[k] for k in range(100)) + True + + :: + + sage: D = LazyDirichletSeriesRing(QQ, "s") + sage: D.sum(lambda p: D(1, valuation=p), Primes()) + 1/(2^s) + 1/(3^s) + 1/(5^s) + 1/(7^s) + O(1/(9^s)) + """ + if a is None: + return super().sum(f) + + if a is True: + it = f + elif a in ZZ: + if b != infinity: + return super().sum(f(i) for i in range(a, b+1)) + from sage.sets.non_negative_integers import NonNegativeIntegers + it = (f(i+a) for i in NonNegativeIntegers()) + else: + it = (f(i) for i in a) + + from sage.data_structures.stream import Stream_infinite_sum + coeff_stream = Stream_infinite_sum(it) + return self.element_class(self, coeff_stream) + def _test_invert(self, **options): """ Test multiplicative inversion of elements of ``self``.