Improvements to how series are represented

237 views
Skip to first unread message

Ondřej Čertík

unread,
Nov 22, 2013, 1:03:20 PM11/22/13
to sympy
Hi,

Currently the series are represented as a regular symbolic expression,
with the Order() instance added to it.

I think we should rather represent series as a Series class, roughly
equivalent to SeriesData in Mathematica:

http://reference.wolfram.com/mathematica/ref/SeriesData.html

The Series class would interact with SymPy similar to how Poly works in SymPy.
The Series would represent a series expansion in one variable (I don't
have opinions about multiple variables yet), so it would remember the
terms in some suitable format (see below, just like for Polys, there
are more options on the internal representation) and the order of the
series. Then, when you print it, it will print it just like we
currently do, with the "O(x^2)" term, though that's just how it is
printed.

The Series would have methods like differentiate, invert, etc.

If you look at how GiNaC does it:

https://github.com/certik/ginac/blob/master/ginac/pseries.h
https://github.com/certik/ginac/blob/master/ginac/pseries.cpp

They have a class pseries for it. Their representation is a list of
(coefficient, power) pairs (plus series variable and expansion point):

https://github.com/certik/ginac/blob/master/ginac/pseries.h#L112

In Mathematica, they only have a list of coefficients for all powers
from nmin/den to nmax/den (so some coefficients could be zero), plus
series variable and expansion point. In Ginac, when they convert the
series to an expression:

https://github.com/certik/ginac/blob/master/ginac/pseries.cpp#L593

they simply attach the Order instance.

Conclusion: Mathematica is using dense representation of series, while
GiNaC is using a sparse representation (instead of a list, one can
also use a dictionary). I currently don't know which one is better,
maybe we should have several implementations, just like Polys.



Open ideas:

* I don't like that Order requires special handling in Add:

https://github.com/sympy/sympy/blob/master/sympy/core/add.py#L237

So maybe we can just append Order to an expression in Add, but the
simplification will be done by the Series class, i.e. no special
handling in Add. So there will be a method Series.to_expr(), that
would return an Add instance, with Order term appended. Unfortunately
if you want to do some stuff like O(x)+O(x^2), then it wouldn't work
if O(x) is just an Order. However, one idea could be that if we
completely remove Order, and only have Series, then O(x) would simply
be an instance of Series. And then just like Polys work, O(x)+O(x^2)
would by dispatched to Series, that knows how to simplify it.
So Add internally doesn't have to handle simplification of this, only
"holding" the terms.

* I am thinking how to implement this in CSymPy as well. Thinking
about it helps with the design, because in CSymPy, it needs to be as
fast as possible. So handling of Order in Add is a nono. In general,
the best way is to always have a special algorithm+efficient
representation for the given thing. So Add should only handle the
usual x+x or x*x thing. Series handling should be done by a special
class Series, that knows everything about series expansion and can use
efficient representation (both the dense and sparse representation
above is more efficient than just the general expression with Order
term appended, for things like O(x)+O(x^2) or multiplication or
inversion of series). At the same time, what's best for CSymPy might
not necessarily be the best for SymPy. But so far it seems to me that
we always want to have the best algorithm+datastructures in SymPy, so
that it's fast, algorithmically. C++ then only provides a constant
speedup, as it doesn't have the overhead of Python.

* See also this:

https://github.com/sympy/sympy/wiki/UD-Under-Development#series

Ondrej

Ondřej Čertík

unread,
Nov 22, 2013, 1:05:48 PM11/22/13
to sympy
P.S. This again touches the multiple dispatch problem in SymPy.
We simply have to tackle it, I don't think there is any way around it.

Ondrej

Sergey Kirpichev

unread,
Nov 23, 2013, 5:36:53 AM11/23/13
to sy...@googlegroups.com
 
I think we should rather represent series as a Series class, roughly
equivalent to SeriesData in Mathematica:

http://reference.wolfram.com/mathematica/ref/SeriesData.html

This is a bit undocumented (as all in Mathematica, usually).  For example,
it's not so clear what kind of coefficients are allowed.
 
The Series class would interact with SymPy similar to how Poly works in SymPy.
The Series would represent a series expansion in one variable (I don't
have opinions about multiple variables yet), so it would remember the
terms in some suitable format (see below, just like for Polys, there
are more options on the internal representation) and the order of the
series. Then, when you print it, it will print it just like we
currently do, with the "O(x^2)" term, though that's just how it is
printed.

Not sure if this is a Series class job.  I mean, to remember terms.

The Series should be build on (or to be itself) the ground of
some stream-like entity, that can do actual arithmetic.
Memoization/caching - another issue.
 
If you look at how GiNaC does it

It looks, like it can only formal power series and laurent (just like
sage?).  We want more?

Open ideas:

* I don't like that Order requires special handling in Add:

https://github.com/sympy/sympy/blob/master/sympy/core/add.py#L237

That's another issue.  But arithmetics (not differential
calculus!) with O makes sense.

Aaron Meurer

unread,
Nov 23, 2013, 12:53:32 PM11/23/13
to sy...@googlegroups.com
> On Nov 22, 2013, at 11:03 AM, "Ondřej Čertík" <ondrej...@gmail.com> wrote:
>
> Hi,
>
> Currently the series are represented as a regular symbolic expression,
> with the Order() instance added to it.
>
> I think we should rather represent series as a Series class, roughly
> equivalent to SeriesData in Mathematica:
>
> http://reference.wolfram.com/mathematica/ref/SeriesData.html
>
> The Series class would interact with SymPy similar to how Poly works in SymPy.
> The Series would represent a series expansion in one variable (I don't
> have opinions about multiple variables yet), so it would remember the
> terms in some suitable format (see below, just like for Polys, there
> are more options on the internal representation) and the order of the
> series. Then, when you print it, it will print it just like we
> currently do, with the "O(x^2)" term, though that's just how it is
> printed.

And also remember the original expression!

>
> The Series would have methods like differentiate, invert, etc.
>
> If you look at how GiNaC does it:
>
> https://github.com/certik/ginac/blob/master/ginac/pseries.h
> https://github.com/certik/ginac/blob/master/ginac/pseries.cpp
>
> They have a class pseries for it. Their representation is a list of
> (coefficient, power) pairs (plus series variable and expansion point):
>
> https://github.com/certik/ginac/blob/master/ginac/pseries.h#L112
>
> In Mathematica, they only have a list of coefficients for all powers
> from nmin/den to nmax/den (so some coefficients could be zero), plus
> series variable and expansion point. In Ginac, when they convert the
> series to an expression:
>
> https://github.com/certik/ginac/blob/master/ginac/pseries.cpp#L593
>
> they simply attach the Order instance.
>
> Conclusion: Mathematica is using dense representation of series, while
> GiNaC is using a sparse representation (instead of a list, one can
> also use a dictionary). I currently don't know which one is better,
> maybe we should have several implementations, just like Polys.

I don't see the point of sparse, since most series are dense. But see
what I wrote below.
Also take a look at Mario's work (I think it may be buried in some
pull request). His implementation will probably be the fastest.

Aaron Meurer


>
> Ondrej
>
> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> For more options, visit https://groups.google.com/groups/opt_out.

Aaron Meurer

unread,
Nov 23, 2013, 12:58:43 PM11/23/13
to sy...@googlegroups.com

On Nov 23, 2013, at 3:36 AM, Sergey Kirpichev <skirp...@gmail.com> wrote:

 
I think we should rather represent series as a Series class, roughly
equivalent to SeriesData in Mathematica:

http://reference.wolfram.com/mathematica/ref/SeriesData.html

This is a bit undocumented (as all in Mathematica, usually).  For example,
it's not so clear what kind of coefficients are allowed.

My impression was that Mathematica had excellent documentation. In fact, it seems to be the best in the field. But I've never used the actual functions outside of wolfram alpha, so I guess I missed things like under specification. 

But definitely if you consider Mathematica docs, plus the functions site, plus mathworld, plus the more or less interactive documentation that is wolfram alpha, the Mathematica documentation is the best in the field. 

Aaron Meurer

 
The Series class would interact with SymPy similar to how Poly works in SymPy.
The Series would represent a series expansion in one variable (I don't
have opinions about multiple variables yet), so it would remember the
terms in some suitable format (see below, just like for Polys, there
are more options on the internal representation) and the order of the
series. Then, when you print it, it will print it just like we
currently do, with the "O(x^2)" term, though that's just how it is
printed.

Not sure if this is a Series class job.  I mean, to remember terms.

The Series should be build on (or to be itself) the ground of
some stream-like entity, that can do actual arithmetic.
Memoization/caching - another issue.
 
If you look at how GiNaC does it

It looks, like it can only formal power series and laurent (just like
sage?).  We want more?

Open ideas:

* I don't like that Order requires special handling in Add:

https://github.com/sympy/sympy/blob/master/sympy/core/add.py#L237

That's another issue.  But arithmetics (not differential
calculus!) with O makes sense.

--

Ondřej Čertík

unread,
Nov 24, 2013, 1:23:13 AM11/24/13
to sympy
On Sat, Nov 23, 2013 at 10:53 AM, Aaron Meurer <asme...@gmail.com> wrote:
>> On Nov 22, 2013, at 11:03 AM, "Ondřej Čertík" <ondrej...@gmail.com> wrote:
>>
>> Hi,
>>
>> Currently the series are represented as a regular symbolic expression,
>> with the Order() instance added to it.
>>
>> I think we should rather represent series as a Series class, roughly
>> equivalent to SeriesData in Mathematica:
>>
>> http://reference.wolfram.com/mathematica/ref/SeriesData.html
>>
>> The Series class would interact with SymPy similar to how Poly works in SymPy.
>> The Series would represent a series expansion in one variable (I don't
>> have opinions about multiple variables yet), so it would remember the
>> terms in some suitable format (see below, just like for Polys, there
>> are more options on the internal representation) and the order of the
>> series. Then, when you print it, it will print it just like we
>> currently do, with the "O(x^2)" term, though that's just how it is
>> printed.
>
> And also remember the original expression!

If it is available (i.e. when you sum or multiply two Series, there
is no original expression).

>
>>
>> The Series would have methods like differentiate, invert, etc.
>>
>> If you look at how GiNaC does it:
>>
>> https://github.com/certik/ginac/blob/master/ginac/pseries.h
>> https://github.com/certik/ginac/blob/master/ginac/pseries.cpp
>>
>> They have a class pseries for it. Their representation is a list of
>> (coefficient, power) pairs (plus series variable and expansion point):
>>
>> https://github.com/certik/ginac/blob/master/ginac/pseries.h#L112
>>
>> In Mathematica, they only have a list of coefficients for all powers
>> from nmin/den to nmax/den (so some coefficients could be zero), plus
>> series variable and expansion point. In Ginac, when they convert the
>> series to an expression:
>>
>> https://github.com/certik/ginac/blob/master/ginac/pseries.cpp#L593
>>
>> they simply attach the Order instance.
>>
>> Conclusion: Mathematica is using dense representation of series, while
>> GiNaC is using a sparse representation (instead of a list, one can
>> also use a dictionary). I currently don't know which one is better,
>> maybe we should have several implementations, just like Polys.
>
> I don't see the point of sparse, since most series are dense. But see
> what I wrote below.

While most series are probably dense, e.g. even this one is dense:

In [1]: (sin(x)**10).series(x, 0, 20)
x**10 - 5*x**12/3 + 4*x**14/3 - 43*x**16/63 + 713*x**18/2835 + O(x**20)

You can easily construct series that are sparse:

In [2]: (sin(x**15)).series(x, 0, 100)
x**15 - x**45/6 + x**75/120 + O(x**100)

So I can definitely see that sparse representation will be faster for
some examples.

Ondrej

Aaron Meurer

unread,
Nov 25, 2013, 6:20:09 PM11/25/13
to sy...@googlegroups.com, Mario Pernici
On Sat, Nov 23, 2013 at 11:23 PM, Ondřej Čertík <ondrej...@gmail.com> wrote:
> On Sat, Nov 23, 2013 at 10:53 AM, Aaron Meurer <asme...@gmail.com> wrote:
>>> On Nov 22, 2013, at 11:03 AM, "Ondřej Čertík" <ondrej...@gmail.com> wrote:
>>>
>>> Hi,
>>>
>>> Currently the series are represented as a regular symbolic expression,
>>> with the Order() instance added to it.
>>>
>>> I think we should rather represent series as a Series class, roughly
>>> equivalent to SeriesData in Mathematica:
>>>
>>> http://reference.wolfram.com/mathematica/ref/SeriesData.html
>>>
>>> The Series class would interact with SymPy similar to how Poly works in SymPy.
>>> The Series would represent a series expansion in one variable (I don't
>>> have opinions about multiple variables yet), so it would remember the
>>> terms in some suitable format (see below, just like for Polys, there
>>> are more options on the internal representation) and the order of the
>>> series. Then, when you print it, it will print it just like we
>>> currently do, with the "O(x^2)" term, though that's just how it is
>>> printed.
>>
>> And also remember the original expression!
>
> If it is available (i.e. when you sum or multiply two Series, there
> is no original expression).

Wouldn't it just be the product of the original expressions? Keeping
track as much as possible is useful so that you can have some kind of
extend() method that just extends the series with more terms
(efficiently).
I would highly recommend looking at the work Mario did. It's going to
be more efficient than any representation that you come up with in
just a few minutes, and plus he already wrote a bunch of algorithms
around it.

The work was https://github.com/sympy/sympy/pull/609, but I have no
idea how much of that has been forward ported to ring() and friends.
It also looks like some more has been ported at
https://github.com/sympy/sympy/pull/2630. I guess Mario will need to
comment on the status of it all.

Aaron Meurer

Ondřej Čertík

unread,
Nov 25, 2013, 7:25:03 PM11/25/13
to sympy, Mario Pernici
Ah I see --- that's how lseries() works in SymPy, that it recursively
returns you the next element of the series, from the original expression.
I would call it the "top-down" approach. Asking for the next term
in the top expression and recursively propagating this to the inner
most expression (which might require couple terms, due to various
cancellations)

The nseries() in sympy is a "bottom-up" approach, asking for n-terms
in the inner most expression and multiplying things out.

I don't know if these two algorithms can be unified. I think that
each has advantages and disadvantages.
Thanks for the pointer.

Mario, if you could comment on these, that would be awesome.

Ondrej

mario

unread,
Nov 26, 2013, 7:38:10 AM11/26/13
to sy...@googlegroups.com, Mario Pernici

In PR 609 and PR 2630 power series are implemented with dictionaries, that is in a sparse representation.
The sparse representation is efficient when there are many variables. PR 2630 uses ``ring()``.

For power series in one variable the dense representation is usually faster.
In `compose_add` in PR 2630 series in one variable are used, so using the dense representation
would be faster, but in a negligible way, since `compose_add` takes roughly 5% of the time in `minpoly`.

In PR 609 there is lpoly.py, which deals with power series, and ltaylor.py.
lpoly.py can now be implemented using ``ring()``; this is partially done in ring_series.py in PR 2630,
where there are only the functions used in the computed sum. It seems to me that ring_series.py fits
well in SymPy; if one wants to compute efficiently power series,
one can use ``ring`` with ring_series.py functions (once all the series functions in lpoly.py are ported).

To do a thorough treatment of power series, one should also implement
the dense representation; consider the case of Laurent series;
add functions to manipulate series at a higher level, possibly using a Series class with an Order function.

Mateusz Paprocki might want to comment on this.

PowerSeriesRing in Sage could be an example of higher level representation.

ltaylor.py uses lpoly.py to compute the taylor series of SymPy functions;
if the function is not smooth enough, the computation is passed to series.series.
So if the function is smooth enough, the taylor series is computed fast
(more or less the speed of `taylor` in Sage); otherwise it can be
slightly slower than series.series.
I do not know if it is a good approach; anyway, since it uses power series,
for the moment I will leave it aside.

Ondřej Čertík

unread,
Nov 27, 2013, 11:39:30 PM11/27/13
to sympy
On Tue, Nov 26, 2013 at 5:38 AM, mario <mario....@gmail.com> wrote:
>
> In PR 609 and PR 2630 power series are implemented with dictionaries, that
> is in a sparse representation.
> The sparse representation is efficient when there are many variables. PR
> 2630 uses ``ring()``.
>
> For power series in one variable the dense representation is usually faster.
> In `compose_add` in PR 2630 series in one variable are used, so using the
> dense representation
> would be faster, but in a negligible way, since `compose_add` takes roughly
> 5% of the time in `minpoly`.
>
> In PR 609 there is lpoly.py, which deals with power series, and ltaylor.py.
> lpoly.py can now be implemented using ``ring()``; this is partially done in
> ring_series.py in PR 2630,
> where there are only the functions used in the computed sum. It seems to me
> that ring_series.py fits
> well in SymPy; if one wants to compute efficiently power series,
> one can use ``ring`` with ring_series.py functions (once all the series
> functions in lpoly.py are ported).
>
> To do a thorough treatment of power series, one should also implement
> the dense representation; consider the case of Laurent series;
> add functions to manipulate series at a higher level, possibly using a
> Series class with an Order function.
>
> Mateusz Paprocki might want to comment on this.
>
> PowerSeriesRing in Sage could be an example of higher level representation.

Thanks a lot for the feedback. So I think we should try to merge your
code, or improve it, so that it can be merged.


> ltaylor.py uses lpoly.py to compute the taylor series of SymPy functions;
> if the function is not smooth enough, the computation is passed to
> series.series.
> So if the function is smooth enough, the taylor series is computed fast
> (more or less the speed of `taylor` in Sage); otherwise it can be
> slightly slower than series.series.
> I do not know if it is a good approach; anyway, since it uses power series,
> for the moment I will leave it aside.

Btw, just today I discovered a real life example where series() is really slow:

https://github.com/sympy/sympy/issues/2635

Ondrej

Pablo Puente

unread,
Nov 29, 2013, 1:07:38 PM11/29/13
to sy...@googlegroups.com
I definitely agree that Series infrastructure is important and SymPy should have it.

Maybe the base class can be just a pure interface (in Java terms, pure virtual functions in C++), not requiring actual storage sparse or dense. So that coefficients can be retrieved through a function. Then in Sparse or Dense subclass coefficients can be cached or pre-calculated in the appropriate manner.

Other nice to have series also supported in the new structure would be Fourier trigonometric series:
- SUM( an*sin(n*x)+bn*cos (n*x) )
- SUM( cn*exp(i*n*x) )

And thinking even more general any projection on a Hilbert space with a given base of functions (of which particular cases are Fourier series, Legendre polynomial series, etc.). 

Pablo.  

Ondřej Čertík

unread,
Dec 2, 2013, 3:48:49 PM12/2/13
to sympy
On Fri, Nov 29, 2013 at 11:07 AM, Pablo Puente <ppu...@googlemail.com> wrote:
> I definitely agree that Series infrastructure is important and SymPy should
> have it.
>
> Maybe the base class can be just a pure interface (in Java terms, pure
> virtual functions in C++), not requiring actual storage sparse or dense. So
> that coefficients can be retrieved through a function. Then in Sparse or
> Dense subclass coefficients can be cached or pre-calculated in the
> appropriate manner.

I think you are right. And we can implement this in parallel to what
we have already in sympy.

>
> Other nice to have series also supported in the new structure would be
> Fourier trigonometric series:
> - SUM( an*sin(n*x)+bn*cos (n*x) )
> - SUM( cn*exp(i*n*x) )

Interesting --- so you would subclass Series as FourierSeries,
and let FourierSeries return a Fourier series of the given function?
So internally it would
calculate the cn coefficients (or an, bn) and internally it understands
that the series is not given using the x, x^2, x^3 basis, but using
exp(i*n*x) basis.

The Series would be abstract enough so that subclasses can choose
both how coefficients are calculated as well as how the series is represented,
i.e. linear combination of {x, x^2, x^3, ...} or {exp(i*x), exp(2*i*x), ...}.

So DenseSeries, SparseSeries and FourierSeries would be on the same level.

> And thinking even more general any projection on a Hilbert space with a
> given base of functions (of which particular cases are Fourier series,
> Legendre polynomial series, etc.).

Good idea.

Ondrej

Ondřej Čertík

unread,
Feb 11, 2015, 8:43:17 PM2/11/15
to sympy, Mario Pernici
Hi Mario,

On Tue, Nov 26, 2013 at 5:38 AM, mario <mario....@gmail.com> wrote:
>
I finally got to this. I checked out your branch:

https://github.com/sympy/sympy/pull/609

Well, actually the Mateusz' improved branch, and I updated it further
to run with the latest gmpy2:
https://github.com/certik/sympy/tree/lpoly2-improved

And tested an example on my computer with gmpy:

In [1]: %time s = series(sin(x)*cos(x), x, 0, 1000)
CPU times: user 4.52 s, sys: 14 ms, total: 4.54 s
Wall time: 4.5 s

In [2]: %time s = series(log(1+x)*exp(x), x, 0, 500)
CPU times: user 4.15 s, sys: 14 ms, total: 4.16 s
Wall time: 4.13 s


Then I ran the same benchmark using the latest sympy master, using
your code that got merged at https://github.com/sympy/sympy/pull/2630,
I put the scripts here:
https://github.com/certik/sympy/compare/formal_power_series, it's just
the a.py and b.py. Your merged code doesn't seem to support sin/cos
benchmark [1] above, so I just ran the benchmark [2]:

https://github.com/certik/sympy/blob/59492520b443ea5f0ef31fc018e9bc700b93b818/b.py

$ python b.py
import
done
initialize series
done
start
stop
1.85430693626

Seems to be a little bit faster. Finally, I then wrote my own implementation:

https://github.com/certik/sympy/blob/59492520b443ea5f0ef31fc018e9bc700b93b818/a.py

$ python a.py
import
done
initialize series
done
start
1.17188405991
stop
1.17299985886


It's using the same algorithm, the speed improvement is just caused by
sorting the dictionary first before multiplying things out in rs_mul:

items1 = list(p1._data.items())
items1.sort(key=lambda e: e[0])
items2 = list(p2._data.items())
items2.sort(key=lambda e: e[0])


You only sort the items2 in sympy, while I noticed a speedup if I sort
items1 as well. But this is minor.


All this is using sparse representation, as a dictionary. I also
implemented dense representation, e.g. you can try the dense
representation on the example [2] by applying the following patch:

diff --git a/a.py b/a.py
index 5fabdd9..f7caf4a 100644
--- a/a.py
+++ b/a.py
@@ -179,20 +179,20 @@ def div(a, b):
for i in range(1, n):
t = t/i
data.append(t)
-exp = FormalPowerSeriesSparse.from_dense(data)
+exp = FormalPowerSeries(data)
# log(1+x)
data = [QQ.dtype(0)]
t = QQ.dtype(1)
for i in range(1, n):
data.append(t/i)
t = -t
-log = FormalPowerSeriesSparse.from_dense(data)
+log = FormalPowerSeries(data)
#x = FormalPowerSeries([0, 1, 0, 0, 0, 0, 0, 0, 0])
#onemx = FormalPowerSeries([1, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])
print "done"
print "start"
t1 = clock()
-s = mul_sparse(log, exp, n)
+s = mul(log, exp)
#s = mul(sin, cos)
t2 = clock()
print "stop"


Then you get:


$ python a.py
import
done
initialize series
done
start
1.11554312706
stop
1.11559510231

So it's a little bit faster, but the difference is small.

Conclusion --- a.py now contains simple implementation of the taylor
series that is as fast as Mario's fast code.


Anyway, I think this is the right approach, we should support more
functions (sin, cos, ...) in sympy.polys.ring_series and then create a
class like the FormalPowerSeriesSparse() in a.py that would hold the
series and would know how to print it, and also had methods to
manipulate it (or those can be external functions like in a.py).
Finally, this should then be hooked into our series expansion in
sympy.

Now, this works great for a Taylor expansion. We then need to extend
this to support Laurent series and eventually Puiseux series (that's
what SymPy's series is doing). An example of a Puiseux series is:

In [1]: %time print (1/cos(x/log(x))).series(x, 0, 10)
1 + x**2/(2*log(x)**2) + 5*x**4/(24*log(x)**4) +
61*x**6/(720*log(x)**6) + 277*x**8/(8064*log(x)**8) + O(x**10)
CPU times: user 3 s, sys: 14 ms, total: 3.01 s
Wall time: 3.03 s

As you can see, currently it takes 3s, and it should be immediate in
my opinion. You can read the references that I quoted in the docstring
of FormalPowerSeries() in a.py. Essentially, I think these are all
series of the form a*X^n, where n = 0, 1, 2, ... for Taylor series, n
= n0, n0+1, n0+2, ... for negative n0 for Laurent series, and X is a
symbol. For Puiseux series, the X is then something more complicated,
either a fractional power like x^(1/2), or something like x/log(x) in
the last example. But the contents of X is stored separately, and the
rest is just like polynomials. I am currently not sure if things like
multiplication of two Puiseux series needs to be changed somehow.

Sparse/Dense: I would just go with Sparse, it seems the speed is very
good both for really sparse series, as well as dense series.

Overall, I think this is the way to go. Let me know what you think.

Ondrej

Richard Fateman

unread,
Feb 11, 2015, 11:18:32 PM2/11/15
to sy...@googlegroups.com, mario....@gmail.com
I haven't studied all the notes prior to this, but it may be helpful to
look at Macsyma/ Maxima.
Series can be extended to several variables in different ways, e.g.
series in x to order xn with coefs as series in y to order yn.... etc

or to total order   that is degree in(x) + degree in (y) + ....  <=N

Also the index set could be other than the integers 0,1,2, ...
including negative and fractional powers.

Also, a key procedure in some algorithms is to find the lowest non-zero
coefficient, a task which can be easy, hard, or non-computable.
So you might have to think about this.  Others have resolved it
by heuristics.

Apologies if these comments have already been taken into account or
if they are irrelevant in this context.
RJF

mario

unread,
Feb 14, 2015, 11:52:21 AM2/14/15
to sy...@googlegroups.com, mario....@gmail.com
I agree that `ring_series` should be completed.
A couple of observations:
`FormalPowerSeriesSparse` is for univariate series; `ring_series` can be used also
for multivariate series (in the former the keys are integers, in the latter a tuple).
In the former one can have more than one variable using rings of rings; in the
latter one can do the same, or keep the series terms distributed.
The series inversion algorithm in `a.py` is faster than the one in `ring_series`.
Reply all
Reply to author
Forward
0 new messages