3

I would like to have first several coefficients of formal power series defined implicitly by a differential equation.

Example.

import sympy as sp
sp.init_printing() # math as latex
from IPython.display import display
z = sp.Symbol('z')
F = sp.Function('F')(z)
F_ = sp.Derivative(F, z)
equation = sp.Eq(F**2 + 1 - F_, 0)
display(equation)

enter image description here

solution = sp.dsolve(equation)
display(solution)

enter image description here

sp.series(sp.tan(z), n = 8)

enter image description here

Question. How to compute formal power series solution of ODE without explicitly solving it in terms of elementary functions?

Sidenote. Some differential equations (even linear) have solutions in divergent power series only, for example an equation L = z^2 + z^2 L + z^4 L'. It is interesting to know whether sympy supports such equations along with usual ones.

Related.

The current question is a sequel of a more easy question.

Sympy: how to solve algebraic equation in formal power series?

Community
  • 1
  • 1
Sergey Dovgal
  • 614
  • 6
  • 21
  • I don't know why people prefer upvoting this particular question/answer, but actually, the second one has a better answer (which works faster) and more universal one. https://stackoverflow.com/questions/46422538/how-to-solve-an-algebraic-equation-in-formal-power-series I was going to adapt the answer to the current situation, but so far was lazy to do it, and post a separate solution. – Sergey Dovgal Oct 11 '17 at 06:07

1 Answers1

4

UPDATE: an answer to a similar question has been posted here. I use this second answer instead of the one presented below. The solution using standard functions from sympy works very slowly.

This seems to be (partially) possible in sympy-1.1.1. Be sure to update to the corresponding version first. I refer to official documentation part on ordinary differential equations

http://docs.sympy.org/latest/modules/solvers/ode.html

We can use the method dsolve with additional hints asking to represent the solution as formal power series. This is only possible for certain types of equations. For the above equation, we ask for possible hint types.

>>> sp.classify_ode(equation)

('separable',
'1st_exact',
'1st_power_series',
'lie_group',
'separable_Integral',
'1st_exact_Integral')

Continuing the example in the question, we specify the hint '1st_power_series' and initial conditions F(0) = 0:

solution = sp.dsolve(equation, hint='1st_power_series', ics={F.subs(z,0):0})
display(solution)

enter image description here

Issue 1. If we want more terms, even say n = 10, the function works for a very long time (normally, I expect 20 coefficients within several seconds). Combinatorially, one can write a very fast recurrence for linear ODE. I don't know whether it is implemented in sympy.

Issue 2. This technique doesn't seem to apply to Ricatti equations like the one mentioned in the original question L = z^2 + z^2 L + z^4 L'.

If someone knows a better solution, (s)he is warmly welcome!

Community
  • 1
  • 1
Sergey Dovgal
  • 614
  • 6
  • 21