@isuruf You mean C++ wrappers right? So we want Flint polys to be optionally used or what i.e. will we have both symengine implementation of polys as well as flint poly wrappers?

any suggestions

Looks good

thanks

@isuruf pls have a look at my proposal I have made some changes.Am I missing any point?

@isuruf aren't FLINT polynomial representations dense? Won't they perform poorly for very sparse polynomials?

and hence, we use FLINT only after checking that the polynomial is dense enough

right>

is Taylor Series Expansion implemented in symengine? if not, when will be relased?

@mkarakoc it is, for a subset of the functions of SymPy

@rwst I guess you mean this (https://github.com/symengine/symengine/wiki/expression-series-expansion)? But when will this be fully completed and functional for python wrapper (symengine.py)?

@mkarakoc no I mean this symengine/symengine.py#28, i.e., it is already in symengine.py

@rwst "ImportError: cannot import name series" . How can I upgrade symengine? pip install says i have the latest version.

@mkarakoc, you'll need SymEngine master and recompile SymEngine with Piranha (SymEngine's implementation is still a WIP, so we are using piranha for polynomials for now. Once #868 is in, you won't need Piranha)

Then you'll need to clone symengine.py and install.

Then you'll need to clone symengine.py and install.

@srajangarg, it depends on how sparse your polynomial is

@isuruf Can you please tell what you meant by

you might want to look at wrapping flint polys as well, for very fast manipulation.

yes and yes

But flint already has C++ wrappers, doesn’t it?

Yes, I meant, making it usable from SymEngine

Yeah, I thought so. Cool.

@isuruf Did you see my proposal? Any views on the questions in the end?

@isuruf Thank you for explanation. But I have no idea what is Piranha and WIP. I guess than I have to wait not to deal with all those steps you told.

@isuruf is there (or wil there be) a poly.coeffs() and sym.roots() equivalent in symengine?

WIP - Work in progress

Piranha - https://github.com/bluescarni/piranha - It's a library for polynomial manipulations

Piranha - https://github.com/bluescarni/piranha - It's a library for polynomial manipulations

@mkarakoc, yes, we are hoping to add polynomial support this summer.

```
for n in range(1, nmax+2):
c += [[]]
d += [[]]
for i in range(nmax+1 - n):
c[n] += [(i+1)*c[n-1][i+1] + d[n-1][i] + sum([c[0][k]*c[n-1][i-k] for k in range(i+1)])]
d[n] += [(i+1)*d[n-1][i+1] + sum([d[0][k]*c[n-1][i-k] for k in range(i+1)])]
```

The c's and d's are lists which are contains polynomials. The time to calculate this loop becomes so long when nmax>=20. I have waited more than 8 hours but it did not finish when I tried last time. If i can make it symengine work for my code I am hoping it will be so fast.

BTW: Same process done in Mathematica in minutes.

c[0][0] is a polynomial?

it contains a polynomial

What kind of a polynomial? Integer coefficient one? a SymPy polynomial?

sorry, they are coefficents of a taylor series of a function. but then the loops creates polynomials. and yes in the end they are SymPy polinomials.

```
c[0] = [0, 18.0000000000000, -54.0000000000000, 162.000000000000, -486.000000000000, 1458.00000000000,
-4374.00000000000, 13122.0000000000, -39366.0000000000, 118098.000000000, -354294.000000000,
1062882.00000000, -3188646.00000000, 9565938.00000000, -28697814.0000000, 86093442.0000000,
-258280326.000000, 774840978.000000, -2324522934.00000, 6973568802.00000, -20920706406.0000]
d[0] = [-1.0*En - 2.22608382037941, -18.0765317747859, 54.0050736479500, -162.000252832090, 486.000010090681, -1458.00000033582, 4374.00000000958, -13122.0000000002, 39366.0000000000, -118098.000000000, 354294.000000000, -1062882.00000000, 3188646.00000000, -9565938.00000000, 28697814.0000000, -86093442.0000000, 258280326.000000, -774840978.000000, 2324522934.00000, -6973568802.00000, 20920706406.0000]
```

What's En?

think "En" as an "x" of polynomial.

`En = symbol("x")`

?
yes.

With nmax = 20, I get the result in around 15 seconds. With SymPy

because precision is not infinite

c[0]=[0, 18, -54, 162, -486, 1458, -4374, 13122, -39366, 118098, -354294, 1062882, -3188646, 9565938, -28697814, 86093442, -258280326, 774840978, -2324522934, 6973568802, -20920706406]

Ah, okay, this is just an example with floats

d[0]=[-(En*exp(1/15) - 9*exp(1/15) + 12)*exp(-1/15),*

-(-192 + 270exp(1/15))*exp(-1/15)/5,*

(-2886 + 4050exp(1/15))*exp(-1/15)/25,*

-(-43292 + 60750exp(1/15))*exp(-1/15)/125,*

(-1298761 + 1822500exp(1/15))*exp(-1/15)/1250,*

-(-48703538 + 68343750exp(1/15))*exp(-1/15)/15625,*

(-8766636841 + 12301875000exp(1/15))*exp(-1/15)/937500,*

-(-460248434153 + 645848437500exp(1/15))*exp(-1/15)/16406250,*

(-110459624196721 + 155003625000000exp(1/15))*exp(-1/15)/1312500000,*

-(-1864006158319667 + 2615686171875000exp(1/15))*exp(-1/15)/7382812500,*

(-2236807389983600401 + 3138823406250000000exp(1/15))*exp(-1/15)/2953125000000,*

-(-184536609673647033083 + 258952931015625000000exp(1/15))*exp(-1/15)/81210937500000,*

(-66433179482512931909881 + 93223055165625000000000exp(1/15))*exp(-1/15)/9745312500000000,*

-(-3238617499772505430606699 + 4544623939324218750000000exp(1/15))*exp(-1/15)/158361328125000000,*

(-2720438699808904561709627161 + 3817484109032343750000000000exp(1/15))*exp(-1/15)/44341171875000000000,*

-(-23542257979115520245564081201 + 33035920174318359375000000000exp(1/15))*exp(-1/15)/127907226562500000000,*

(-146903689789680846332319866694241 + 206144141887746562500000000000000exp(1/15))*exp(-1/15)/266047031250000000000000,*

-(-1170638778011519244210673937719733 + 1642711130667980419921875000000000exp(1/15))*exp(-1/15)/706687426757812500000000,*

(-778024541693809713075401755530653317 + 1091771089920873140625000000000000000exp(1/15))*exp(-1/15)/156558445312500000000000000,*

-(-1441290463487782493472181752120535269743 + 2022505944078417493007812500000000000000exp(1/15))*exp(-1/15)/96674839980468750000000000000,*

(-864774278092669496083309051272321161845801 + 1213503566447050495804687500000000000000000exp(1/15))*exp(-1/15)/19334967996093750000000000000000]

-(-192 + 270

(-2886 + 4050

-(-43292 + 60750

(-1298761 + 1822500

-(-48703538 + 68343750

(-8766636841 + 12301875000

-(-460248434153 + 645848437500

(-110459624196721 + 155003625000000

-(-1864006158319667 + 2615686171875000

(-2236807389983600401 + 3138823406250000000

-(-184536609673647033083 + 258952931015625000000

(-66433179482512931909881 + 93223055165625000000000

-(-3238617499772505430606699 + 4544623939324218750000000

(-2720438699808904561709627161 + 3817484109032343750000000000

-(-23542257979115520245564081201 + 33035920174318359375000000000

(-146903689789680846332319866694241 + 206144141887746562500000000000000

-(-1170638778011519244210673937719733 + 1642711130667980419921875000000000

(-778024541693809713075401755530653317 + 1091771089920873140625000000000000000

-(-1441290463487782493472181752120535269743 + 2022505944078417493007812500000000000000

(-864774278092669496083309051272321161845801 + 1213503566447050495804687500000000000000000

Is there a version where I can use the fact that c[0] is a polynomial?

With the code you gave

`c[0]`

is just a list. If we can use polynomial manipulation, considering `c[0]`

is a polynomial, then it could made much faster
@mkarakoc, how about this, https://gist.github.com/isuruf/03fe806581a23c61f9c8?

Runs in 3 seconds in my computer

yeah. same here but then i have this. :)

```
np.set_printoptions(precision=50)
n = nmax
sol = sym.poly(sym.N(d[n][0], 100)*sym.N(c[n-1][0], 100) - sym.N(d[n-1][0], 100)*sym.N(c[n][0], 100), En)
solc = sol.coeffs()
np.roots(solc)
```

Roots calculation is the step taking the longest time?

I guess, there is no sym.poly(), sol.coeffs() and np.roots() or sym.roots() equivalent in symengine.

yes, then that part takes much more time.

Yes,

`sym.poly()`

and `sol.coeffs()`

needs to be added
`np.roots()`

should be fast though, since numpy's backend is C
i made nmax=40 and i am still waiting.

it is done with Mathematica in 1.953125 seconds.

and the precision in Mathematica is set to 500, but it is way less (64 bit) precision in Python.

I get a index out of bounds error when nmax=40

SymEngine does support arbitrary precision floats when installed with MPFR library

I think this would make a good benchmark for SymEngine. We can compare Mathematica, SymEngine and other software and try to improve SymEngine

Btw, do you want to keep

`c`

and `d`

symbolically or with a high precision like 500?
c[0] = [0, 18.0000000000000, -54.0000000000000, 162.000000000000, -486.000000000000, 1458.00000000000, -4374.00000000000, 13122.0000000000, -39366.0000000000, 118098.000000000, -354294.000000000, 1062882.00000000, -3188646.00000000, 9565938.00000000, -28697814.0000000, 86093442.0000000, -258280326.000000, 774840978.000000, -2324522934.00000, 6973568802.00000, -20920706406.0000, 62762119218.0000, -188286357654.000, 564859072962.000, -1694577218886.00, 5083731656658.00, -15251194969974.0, 45753584909922.0, -137260754729766., 411782264189298., -1.23534679256789e+15, 3.70604037770368e+15, -1.11181211331110e+16, 3.33543633993331e+16, -1.00063090197999e+17, 3.00189270593998e+17, -9.00567811781995e+17, 2.70170343534598e+18, -8.10511030603795e+18, 2.43153309181139e+19, -7.29459927543416e+19]

d[0] = [-1.0*En - 2.22608382037941, -18.0765317747859, 54.0050736479500, -162.000252832090, 486.000010090681, -1458.00000033582, 4374.00000000958, -13122.0000000002, 39366.0000000000, -118098.000000000, 354294.000000000, -1062882.00000000, 3188646.00000000, -9565938.00000000, 28697814.0000000, -86093442.0000000, 258280326.000000, -774840978.000000, 2324522934.00000, -6973568802.00000, 20920706406.0000, -62762119218.0000, 188286357654.000, -564859072962.000, 1694577218886.00, -5083731656658.00, 15251194969974.0, -45753584909922.0, 137260754729766., -411782264189298., 1.23534679256789e+15, -3.70604037770368e+15, 1.11181211331110e+16, -3.33543633993331e+16, 1.00063090197999e+17, -3.00189270593998e+17, 9.00567811781994e+17, -2.70170343534598e+18, 8.10511030603795e+18, -2.43153309181139e+19, 7.29459927543416e+19]

500 is so good if i can make it work.

I've updated the gist

It's with precision 53. If you can give me the symbolic values for nmax=40, I can try with precision 500 as well

ok.

i have put them to comment section of the gist. But may be you should also add np.roots() part of the code to the gist.

Yes, will do

c[0] = [sympify(i) for i in c[0]]

d[0] = [sympify(i) for i in d[0]]

d[0] = [sympify(i) for i in d[0]]

what is this for?

Convert them into SymEngine expressions, so that I can call expand in the loop

do you think .expand in loop makes it faster?

For symbolic code, yes. Probably not for numerics

but why? expanding will consume some more time?

Yeah, but the expressions will be smaller

@isuruf how can i return from SymPy to SymEngine as you do with sympify()?

`sympify()`

can take in a SymPy expression
then why i am getting the following error when finding the roots.

`AttributeError: 'symengine.lib.symengine_wrapper.Symbol' object has no attribute 'is_commutative'`

That's a problem with SymPy poly

`En._sympy_()`

should be passed instead of `En`

as the 2nd argument
With nmax=40, prec=500, first step takes 15.3s

Full calculation takes 19s

nice.

I've updated the gist

with Mathematica first step is 0.859375 sec for first step with nmax=40.

Thanks. I'll try to figure out if there is room for improvement after playing with Mathematica

Thanks to you also. How can I install MPFR?

i cannot rerun your version without MPFR.

What's your OS?

@nishnik @CodeMaxx I have left a few comments in your proposals. Do let me know if there is any trouble.

@isuruf do we need sympify in

`SymEngine.jl`

?
Yes

@isuruf ubuntu. I have installed libmpfr and libmpc. but it still complains about mpfr. may be it needs some python/ipython bindings.