fast implementation of exp

[mantepse]

Using the recursive definition

exp(f) = 1 + int(f' *exp(f))

we can make the computation of the exponential almost the same as a single multiplication.

The drawback is that we do not benefit from the error handling in compose.

Depends on #34552
Depends on #34636

CC: @tscrim

Component: combinatorics

Keywords: LazyPowerSeries

Author: Martin Rubey

Branch/Commit: 5020b9d

Reviewer: Travis Scrimshaw

Issue created by migration from https://trac.sagemath.org/ticket/34611

[mantepse]

Branch: u/mantepse/fast_implementation_of_exp

[mantepse]

Last 10 new commits:

cdea820	Merge branch 'public/rings/lazy_series_revert-34383' of https://github.com/sagemath/sagetrac-mirror into u/tscrim/derivatives_lazy_series-34413
b7f04ed	Merge branch 'develop' of trac.sagemath.org:sage into t/34413/derivatives_lazy_series-34413
2325826	Merge branch 'u/mantepse/derivatives_lazy_series-34413' of trac.sagemath.org:sage into t/34422/implement_functorial_composition_of_lazy_symmetric_functiosn
ebcf5c3	Merge branch 'u/mantepse/implement_functorial_composition_of_lazy_symmetric_functiosn' of trac.sagemath.org:sage into t/34423/implement_arithmetic_product_of_lazy_symmetric_functions
4172688	Merge branch 'u/mantepse/implement_arithmetic_product_of_lazy_symmetric_functions' of trac.sagemath.org:sage into t/32367/replace_lazy_power_series-32367
41ca99b	Merge branch 'u/mantepse/replace_lazy_power_series-32367' of trac.sagemath.org:sage into t/34470/categories_lazy_series-34470
fb3a5cd	Merge branch 'u/mantepse/categories_lazy_series-34470' of trac.sagemath.org:sage into t/34552/lazy_series_test_suite-34552
855d2bf	remove unused variable
75c275c	add documentation and doctests for _approximate_order
80d1e58	use that exp is holonomic

[mantepse]

Author: Martin Rubey

[mantepse]

Commit: 80d1e58

[mantepse]

Changed keywords from none to LazyPowerSeries

[mantepse]

comment:3

Here is a comparison:

sage: L.<z> = LazyLaurentSeriesRing(QQ)
sage: f = L(lambda n: randint(1, 100), valuation=1)
sage: g = exp(f)
sage: %time a = g[200]
CPU times: user 86.4 ms, sys: 29 µs, total: 86.4 ms
Wall time: 86.4 ms

sage: e = L(coefficients=lambda n: 1/factorial(ZZ(n)), valuation=0)
sage: h = e(f)
sage: %time b = h[200]
CPU times: user 1.7 s, sys: 20 µs, total: 1.7 s
Wall time: 1.7 s

sage: a == b
True

[sagetrac-git]

Changed commit from 80d1e58 to 08ca5a3

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

08ca5a3

provide fast logarithm

[mantepse]

Dependencies: #34552

[mantepse]

comment:6

There is a slight problem:

sage: L.<z> = LazyLaurentSeriesRing(QQ)
sage: f = L(lambda n: randint(1, 10), valuation=1)
sage: g = exp(f)
sage: %time _ = g[500]
...
RecursionError: maximum recursion depth exceeded while calling a Python object

[mantepse]

comment:7

This is not a problem if we use dense series (of course) - so there are at least three solutions:

• we could create a dense ring in exp,

• we could modify Stream_uninitialized.get_coefficient to compute coefficients in order,

• we could ignore the problem, and tell the user to use dense series.

[tscrim]

comment:8

There is also another solution: Tell the user to compute some earlier coefficients:

sage: L.<z> = LazyLaurentSeriesRing(QQ)
sage: f = L(lambda n: randint(1, 10), valuation=1)
sage: g = exp(f)
sage: %time _ = g[100]
CPU times: user 8.87 ms, sys: 18 µs, total: 8.88 ms
Wall time: 8.74 ms
sage: %time _ = g[200]
CPU times: user 38.1 ms, sys: 0 ns, total: 38.1 ms
Wall time: 37.8 ms
sage: %time _ = g[500]
CPU times: user 678 ms, sys: 0 ns, total: 678 ms
Wall time: 678 ms

However, I like option (2) the best because there is not a reasonable situation where we would not essentially end up computing all of the coefficients (ultimately in order):

sage: len(g._coeff_stream._cache)
501

[mantepse]

comment:9

A similar technique works for computing the Cauchy inverse.

[tscrim]

comment:10

We might be able to limit the number of recursion calls by being careful about the order of multiplication too.

[mantepse]

comment:11

You mean Stream_cauchy_mul(d_self, f._coeff_stream) vs. Stream_cauchy_mul(f._coeff_stream, d_self)?

All of this may also be affected by #34616.

[tscrim]

comment:12

Yes, that's right. Although since #34616 might only apply to the dense case, it might not be relevant to the issue we have with the sparse case.

[mantepse]

comment:13

What do you think about making Stream_uninitialized always dense. This is essentially the same as (2) above, just with less overhead.

This could be done either here or in #34636, which is ready for review. I could then make #34636 a dependency for this ticket. I could even do a doctest there, because the problem arises already with, for example, the implicit definition of the Catalan numbers:

sage: L.<z> = LazyPowerSeriesRing(ZZ); C = L.undefined(); C.define(1 + z*C^2)
sage: %time _ = C[500]
...
RecursionError: maximum recursion depth exceeded while calling a Python object

sage: L.<z> = LazyPowerSeriesRing(ZZ, sparse=True); C = L.undefined(); C.define(1 + z*C^2)
sage: %time _ = C[350]
CPU times: user 201 ms, sys: 0 ns, total: 201 ms
Wall time: 201 ms
sage: L.<z> = LazyPowerSeriesRing(ZZ, sparse=False); C = L.undefined(); C.define(1 + z*C^2)
sage: %time _ = C[350]
CPU times: user 160 ms, sys: 0 ns, total: 160 ms
Wall time: 160 ms

[tscrim]

comment:14

Replying to Martin Rubey:

What do you think about making Stream_uninitialized always dense. This is essentially the same as (2) above, just with less overhead.

+1

This could be done either here or in #34636, which is ready for review. I could then make #34636 a dependency for this ticket.

That would probably be better.

[sagetrac-git]

Changed commit from 08ca5a3 to 5020b9d

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. Last 10 new commits:

3e84c90	Merge branch 'u/mantepse/categories_lazy_series-34470' of trac.sagemath.org:sage into t/34552/lazy_series_test_suite-34552
bd3e3eb	make sparsity a decision of the user
814aa7c	remove Stream_cauchy_invert.get_coefficient, make sparse a mandatory argument, move _is_sparse an attribute of Stream_inexact
7043d1c	remove redundant `__init__` methods, remove finished TODOs
e6c4cae	make Stream_uninitialized always dense to avoid maximal recursion error
9577386	implement Polynomial_generic_sparse.__floordiv__
ec1eebd	Merge branch 'u/mantepse/__floordiv___for_sparse_polynomials' of trac.sagemath.org:sage into t/34636/make_sparsity_a_decision_of_the_user
de424bd	make internal rings sparse or dense if the lazy series ring is sparse or dense
66c1208	Merge branch 'u/mantepse/make_sparsity_a_decision_of_the_user' of trac.sagemath.org:sage into t/34611/fast_implementation_of_exp
5020b9d	adapt exp and log to new sparsity

[mantepse]

comment:16

A fast algorithm for inversion can be found in van der Hoeven's paper, but is not completely trivial to implement, so it should be left for another ticket.

[mantepse]

Changed dependencies from #34552 to #34552, #34636

[mantepse]

comment:17

As an aside, it should be very easy to get rid of the duplicated cache in Stream_uninitialized: currently we have a cache for self._target, and create another one for self. I have an experimental branch, which removes _target, and simply sets self._iter to the iterator of the target. However, I get a strange doctest failure when there are several interlinked uninitialized streams which, additionally, have positive valuation. I don't understand the problem yet. However, considering #34616, we probably want a slightly more advanced cache anyway, i.e., a cache whose size grows in chunks and not just one by one. So it may be better to consider these things together.

[tscrim]

Reviewer: Travis Scrimshaw

[tscrim]

comment:18

Or on a separate ticket altogether.

I think this is good to go for now as it greatly improves the performance (asymptotically; I didn't check for small values, but there is likely at worst a negligible difference).

[vbraun]

Changed branch from u/mantepse/fast_implementation_of_exp to 5020b9d