GSoC idea for a project: holonomic functions

[Ondřej Čertík]

Hi,

I've added a project idea here:

https://github.com/sympy/sympy/wiki/GSoC-2016-Ideas#implement-holonomic-functions

I would really like to have this capability in SymPy.

I think it should allow a robust symbolic integration of a large class
of functions, as well as analytic series expansion (having a
recurrence relation for the Taylor coefficients, and there are
algorithms to get a closed formula from it).

Ondrej

[Fredrik Johansson]

I will probably be available to give advice on such a project (not sure if I will have time to be a mentor... possibly co-mentor).

Fredrik

[Ondřej Čertík]

Thanks Fredrik. I put your name in (an indicated that you might only
advice). If no student picks this up, I want to implement this myself,
so you can advice me. :)

Ondrej

[meghana.m...@gmail.com]

Hi,
I went through the ideas page and https://gist.github.com/certik/847cafaf3111730c5b3f and this project seems very interesting to me. I am thinking about applying to GSoC 16 with this idea. I have gone through the links and thought it would be a good idea to familiarize myself with the hypergeometric module and work from there. Is there anything else that I should know (apart from what is mentioned in the outline given in the ideas page and the link) ?

-Meghana

[Ondřej Čertík]

Hi Meghana,

On Sat, Feb 20, 2016 at 11:10 AM, <meghana.m...@gmail.com> wrote:
> Hi,
> I went through the ideas page and
> https://gist.github.com/certik/847cafaf3111730c5b3f and this project seems
> very interesting to me. I am thinking about applying to GSoC 16 with this
> idea. I have gone through the links and thought it would be a good idea to
> familiarize myself with the hypergeometric module and work from there. Is
> there anything else that I should know (apart from what is mentioned in the
> outline given in the ideas page and the link) ?

You can also play with the ore_algebra package in Sage, as well as see
how Mathematica does things with the DifferentialRoot. You can use
Mathematica online at:
https://lab.open.wolframcloud.com/objects/wpl/GetStarted.nb, the free
version is limited, but it is enough to get an idea what
DifferentialRoot can and cannot do.

Then for the application, you should try to create a prototype, i.e. a
class that can represent a holonomic function (and investigate how
ore_algebra stores it inside and do the same or better), and implement
some operations that are easy to implement, like integration.

If you have any questions, please ask, Fredrik or I will be happy to help.

Ondrej

Ondrej

[meghana.m...@gmail.com]

Thank you Ondřej for the information.

-Meghana

[shubham tibra]

Hi,

I am interested in this project and went through the references mentioned in the ideas page.

From What I understood, We are going to represent the holonomic function in the form of it's differential equation and will use numerical integration methods to find values of the function at a required points. Can all the operations listed on Ideas page be implemented without getting the symbolic form of the function when it is solvable?

And regarding the conversion of Holonomic to Hypergeometric function, we need to expand the holonomic function to a power series. How we can do that without having the actual symbolic representation of the function?

How we will handle the special cases when the Holonomic function represents an elementary function or a polynomial?

Thanks

Shubham

[Ralf Stephan]

On Wednesday, February 24, 2016 at 3:59:08 PM UTC+1, shubham tibra wrote:

And regarding the conversion of Holonomic to Hypergeometric function, we need to expand the holonomic function to a power series. How we can do that without having the actual symbolic representation of the function?

Series coefficients of holonomic functions satisfy linear recurrences (with constant

and polynomial coefficients). It follows that the series coeffs can be represented

as formal sums, and equalities of different formal sums can be decided by comparing

the underlying holonomic function (so you can prove hypergeometric sum equalities).

See e.g. papers by Zeilberger for applications of this fascinating subject. 

 

How we will handle the special cases when the Holonomic function represents an elementary function or a polynomial?

They have hypergeometric representations as well as far as I recall.

Best,

[shubham tibra]

Correct me if I am wrong, then the coefficient ( constant / polynomial) of the i'th series coefficient in the recurrence relation must be Sum_{ j }  A(i, j)*x^j

in the context of Ondrej's write up https://gist.github.com/certik/847cafaf3111730c5b3f .

[Ralf Stephan]

Yes, I should add that the recurrence formula in some way is directly

equivalent to the diff.eq. only in discrete space so to say. And then if

you just have the (rational) coeff values you can via their Ore algebra

infer both the diff.eq. and the recurrence, i.e. "guess" the formula if

there is one.

[Ondřej Čertík]

Hi Shubham,

On Wed, Feb 24, 2016 at 7:59 AM, shubham tibra <shubh...@gmail.com> wrote:
> Hi,
> I am interested in this project and went through the references mentioned in
> the ideas page.

Excellent. Ralf already answered most of your questions, so let me
just answer the rest.

>
> From What I understood, We are going to represent the holonomic function in
> the form of it's differential equation and will use numerical integration
> methods to find values of the function at a required points. Can all the

The numerical integration will only be used if the user wants a
numerical value at a point, i.e. when the user calls `.n()` on it.
Otherwise we want to manipulate the holonomic function symbolically.

> operations listed on Ideas page be implemented without getting the symbolic
> form of the function when it is solvable?

Yes, all these:

* addition
* multiplication
* power
* application of holonomic function to an algebraic one
* differentiation
* integration
* numerical evaluation for any complex value of x

Can be done for any holonomic function, at least that is my
understanding. I provided links on the wiki page to papers with the
actual algorithms, check them out and let me know what you think.

> And regarding the conversion of Holonomic to Hypergeometric function, we
> need to expand the holonomic function to a power series. How we can do that
> without having the actual symbolic representation of the function?

As Ralf said, you can get a recurrence relation of the Taylor
coefficients. The ore_algebra package in Sage does it, for example, or
lookup the other papers I link.
So this recurrence relation allows you to evaluate any Taylor
coefficient. I think there should be a way to obtain a ratio of those,
in order to recover a hypergeometric function.

Also in some cases, one can solve the recurrence relation to obtain a
closed form formula for the Taylor series, I think that's possible in
simple cases like sin(x) or cos(x). That would also be useful to have.
But mainly, we need to figure out how to robustly convert to a
hypergeometric function, because then we can reuse the SymPy's
hyperexpand.

>
> How we will handle the special cases when the Holonomic function represents
> an elementary function or a polynomial?

As Ralf said, most of these are hypergeometric functions, so those
cases are already covered. For holonomic functions that are not
hypergeometric, but that still have a representation in terms of more
elementary functions, we'll have to do some kind of pattern matching,
or a database. Or possibly our differential equation solver dsolve()
can do some of those cases.

Another application that I have in mind is definite integration. I
think you can literally subtract the two endpoints of a holonomic
function (the antiderivative), so then we just need to be able to
symbolically evaluate the holonomic function at a point. I don't know
if there are algorithms for that. Another perhaps a better option is
to convert the integrand, which is a holonomic function, to a
hypergeometric or MeijerG (or a product of two MeijerG) functions. And
then reuse the robust definite integration of those from sympy. In
other words, my hope is that holonomic functions will greatly help in
rewriting any expression into hypergeometric or MeijerG functions,
which currently is very heuristic and sometimes doesn't work.


Let me know if you have more questions.

Ondrej

>
> Thanks
> Shubham
>
> On Friday, February 19, 2016 at 2:01:16 PM UTC+5:30, Ondřej Čertík wrote:
>>
>> Hi,
>>
>> I've added a project idea here:
>>
>>
>> https://github.com/sympy/sympy/wiki/GSoC-2016-Ideas#implement-holonomic-functions
>>
>> I would really like to have this capability in SymPy.
>>
>> I think it should allow a robust symbolic integration of a large class
>> of functions, as well as analytic series expansion (having a
>> recurrence relation for the Taylor coefficients, and there are
>> algorithms to get a closed formula from it).
>>
>> 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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/f8add7b7-e38d-4093-80a2-e571a162918b%40googlegroups.com.
>
> For more options, visit https://groups.google.com/d/optout.

[shubham tibra]

Thanks Ondrej and Ralf.

Is it possible to get the source code of ore_algebra package of sage's? Seeing an example of a already implemented algorithm can help in a lot of ways.

[Ralf Stephan]

Via Google:

http://arxiv.org/abs/1306.4263

http://www.kauers.de/software.html

On Thu, Feb 25, 2016 at 3:50 PM shubham tibra <shubh...@gmail.com> wrote:

Thanks Ondrej and Ralf.

Is it possible to get the source code of ore_algebra package of sage's? Seeing an example of a already implemented algorithm can help in a lot of ways.

--
You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/JkCkgltXeJM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.

To post to this group, send email to sy...@googlegroups.com.
Visit this group at https://groups.google.com/group/sympy.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/ea31688e-3768-4f21-b97e-411eaadffbb1%40googlegroups.com.

[shubham tibra]

Thanks Ralph. 

I had already downloaded the spkg file but didn't know that we can see the source code from it. 

Found out after googling, just rename the ending .spkg to .tar.bz2 and then we can extract it.

[shubham tibra]

I know it's a newbie question but I am stuck at it. 

In the paper mentioned in the ideas page http://www.risc.jku.at/publications/download/risc_2244/DIPLFORM.pdf

the algorithm for addition of holonomic functions is described in Example 1.4.2, page 20. How do we get the linear system as mentioned in Equation (1.4.13)

from the derivatives of functions `f` and `g` ?

[Ralf Stephan]

The rows are the ansatz for the holonomic diff.eq. of order 3.

We got this order from step 1. The linear system solves for the p_i

--

You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/JkCkgltXeJM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.
To post to this group, send email to sy...@googlegroups.com.
Visit this group at https://groups.google.com/group/sympy.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/84abd1b9-d6da-4c22-8951-159f7e7943f9%40googlegroups.com.

[Aaron Meurer]

Does manipulating holonomic functions require solving linear systems
with rational function coefficients, or are there more direct methods?
I ask because if it does a GSoC project for this for SymPy may
require improving SymPy's matrices to support working over rational
function domains.

Aaron Meurer

> 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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit

> https://groups.google.com/d/msgid/sympy/CAMNHCDuZchi5CNMbuyDAii8yhKMFWbBQ__Ng2g8iFeJgm0%2B6Gg%40mail.gmail.com.

[Fredrik Johansson]

On Friday, February 26, 2016 at 6:00:10 PM UTC+1, Aaron Meurer wrote:

Does manipulating holonomic functions require solving linear systems
with rational function coefficients, or are there more direct methods?
 I ask because if it does a GSoC project for this for SymPy may
require improving SymPy's matrices to support working over rational
function domains.

Fundamentally, manipulating linear differential operators with rational function coefficients is the same thing as doing linear algebra with rational function coefficients. It has been shown from a complexity point of view that multiplying linear differential operators is equivalent to multiplying matrices (http://arxiv.org/abs/0804.2181).

Both the direct and linear algebra point of views are valid. For example, to find an annihilator for f + g given annihilators A, B for f and g, one computes the LCLM (least common left multiple) of A and B. This can be done "directly" by a version of the Euclidean algorithm. It can also be done by solving a linear system, as above. There may be some operations that are more natural to do directly and vice versa; the advantage of linear algebra approach is that as soon as you have translated the problem to linear algebra, it is obvious how to do the actual computation and how to optimize that.

Fredrik

[shubham tibra]

I looked at the papers and would like to start with creating a class that can represent holonomic functions and implement algorithms like addition using LCLM( least common left multiple) of their annihilators (or should we use the linear system approach).

But I have no idea how to start, things like what are the functions to be defined, adding docstrings and all seems significant amount of work. Can someone guide me?

[Fredrik Johansson]

On Saturday, February 27, 2016 at 7:07:10 PM UTC+1, shubham tibra wrote:

I looked at the papers and would like to start with creating a class that can represent holonomic functions and implement algorithms like addition using LCLM( least common left multiple) of their annihilators (or should we use the linear system approach).

But I have no idea how to start, things like what are the functions to be defined, adding docstrings and all seems significant amount of work. Can someone guide me?

I suggest that you look at the documentation for the available packages (gfun, HolonomicFunctions, ore_algebra, etc). How do you define and add two functions with each of them? What are the pros and cons?

I think the most basic code you will need is classes to represent differential and difference operators. From there, it's mostly a matter of adding all the necessary methods.

Fredrik

[shubham tibra]

Which case of holonomic functions we need to implement, univariate or multivariate? Because the mathematica package HolonomicFunctions is about multivariate holonomic functions whereas gfun and ore_algebra are for univariate holonomic functions AFAIK.

[Ondřej Čertík]

We will start with univariate.

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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit

> https://groups.google.com/d/msgid/sympy/b833d9b7-49e6-48de-b4cb-82c6723f969c%40googlegroups.com.

[shubham tibra]

I looked at both the methods for addition algorithm( using Euclidean algorithm for computing LCLM, or making an ansatz and converting it into a linear system and then solving).

But the Ore algebra method will also require a class representing rings and some basic functions like getting parent ring, generators of the ring. Does Sympy have support for rings or Do we need to create a class for rings.

Also I think we should discuss on which approach to use, as pointed out by Aaron solving linear system with rational coefficients may require improving Sympy's support for matrices.

[Ondřej Čertík]

On Mon, Feb 29, 2016 at 5:48 AM, shubham tibra <shubh...@gmail.com> wrote:
> I looked at both the methods for addition algorithm( using Euclidean
> algorithm for computing LCLM, or making an ansatz and converting it into a
> linear system and then solving).
> But the Ore algebra method will also require a class representing rings and
> some basic functions like getting parent ring, generators of the ring. Does
> Sympy have support for rings or Do we need to create a class for rings.

Why don't you start a wiki page and put there the exact algorithm that
is needed. The polynomials in SymPy work in different rings, so in
this sense, it's already implemented. SymPy doesn't have the parent
ring structure as Sage. I would argue that it is not needed, at least
it wasn't so far. It's easy to get buried in this abstract design
stuff. Rather, let's figure out the actual technical algorithm,
understand it and write it down. Only then start thinking what the
right abstractions are.

Questions to figure out: what are the rings needed for in ore_algebra?
What kinds of rings are needed? How are they used in addition of two
holonomic functions?

Isn't it the case that you just choose polynomials with rational
coefficients and just stick to that? In fact, since the differential
equation can be multiplied by any integer, we probably only need
polynomials with integer coefficients to represent a holonomic
function (though rational coefficients might be needed for some
intermediate calculations).

> Also I think we should discuss on which approach to use, as pointed out by
> Aaron solving linear system with rational coefficients may require improving
> Sympy's support for matrices.

I think we can use our general solver for now, cannot we? There are
always things to improve and get faster or more robust results, but in
this holonomic functions project, the goal is to get everything up and
running (perhaps not optimally) and then iteratively improve upon
that.

Ondrej

[Fredrik Johansson]

On Monday, February 29, 2016 at 1:48:11 PM UTC+1, shubham tibra wrote:

I looked at both the methods for addition algorithm( using Euclidean algorithm for computing LCLM, or making an ansatz and converting it into a linear system and then solving).

But the Ore algebra method will also require a class representing rings and some basic functions like getting parent ring, generators of the ring. Does Sympy have support for rings or Do we need to create a class for rings.

This isn't required. You could just work with SymPy expressions with a distinguished symbol, e.g. x is the distinguished symbol of the differential operator Dx.

You might be able to use SymPy's polys as well, though I think you do need at least rational functions to do certain things conveniently (not sure how well supported this is in SymPy).

Fredrik

[shubham tibra]

Thanks for your response Ondrej and Fredrik.

I will soon be writing a wiki page for the algorithm. I agree with what Ondrej said, it's better to figure out the technical algorithms rather than abstracting the design.

I am creating a class which can represent a holonomic function. Since we need a Matrix and a Vector, should I inherit from already built classes of these in Sympy or 

create a completely new class a define things from scratch?

[Ondřej Čertík]

I would actually create it separate from sympy as a first step, and
implement the addition, multiplication and other operations on it.
Then in the second step we would figure out what the best way to hook
it up is. Most probably subclassing a Function.

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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit

> https://groups.google.com/d/msgid/sympy/1f53dfb6-0b13-4f22-ab19-cedf9a810e82%40googlegroups.com.

[meghana.m...@gmail.com]

Hi,

I have a few questions.
Would this project also involve representing holonomic sequences ? (in the case that one would want to obtain a relation for the Taylor series coefficients which essentially forms a holonomic sequence). In that case, the generators of the ore algebra can be either the shift operator or the differential operator right ? Also, in the ideas page, I didn't quite understand this point:
"- converting any symbolic expression tree bottom up (as long as elementary leafs can be converted and as long as the symbolic tree only uses +, *, ^ and composition)" Does this expression tree consist of the holonomic representations of the functions or just the symbolic form which must be converted ?

-Meghana

[Ondřej Čertík]

On Wed, Mar 2, 2016 at 8:45 AM, <meghana.m...@gmail.com> wrote:
> Hi,
>
> I have a few questions.
> Would this project also involve representing holonomic sequences ? (in the
> case that one would want to obtain a relation for the Taylor series
> coefficients which essentially forms a holonomic sequence). In that case,
> the generators of the ore algebra can be either the shift operator or the
> differential operator right ?

The priority will be holonomic functions. But I think the holonomic
sequences are involved somehow when you calculate the series expansion
of a holonomic function.

> Also, in the ideas page, I didn't quite
> understand this point:
> "- converting any symbolic expression tree bottom up (as long as elementary
> leafs can be converted and as long as the symbolic tree only uses +, *, ^
> and composition)" Does this expression tree consist of the holonomic
> representations of the functions or just the symbolic form which must be
> converted ?

The application is just the symbolic form that must be converted. You
do that by converting a few leafs. So then you have a tree with a mix
of symbolic forms from SymPy and holonomic functions (though you don't
need to represent this as a tree --- you can just traverse the
symbolic expression and construct the corresponding holonomic function
on the fly).
So the answer is yes to both.

Ondrej

> https://groups.google.com/d/msgid/sympy/e56e0bb3-0c62-4360-b874-c4b27afd666a%40googlegroups.com.

[Tanu Hari Dixit]

Hey Shubham,

How can I add pretty printing and other methods for this class?

 You can take a look at sympy/printing/pretty/pretty.py or maybe at this PR I created for pretty printing of trace of matrices, for an example.

[shubham tibra]

It was helpful. Thanks, Tanu.

[shubham tibra]

Regarding the initial conditions in multiplication and addition when they are not at the same point:

I think the best is(already discussed at ideas page) first to check if both the given holonomic functions are elementary function and calculate the initial conditions at a same point symbolically, where both the functions doesn't have a pole. So for the resulting function the initial condition would be addition/multiplication of calculated conditions at that point.

In case when one or both can't be converted to elementary we will calculate numerical value at a point and add/multiply them.

How are we gonna do this in INTEGRATION and DIFFERENTIATION?

In application with a algebraic function, let's say we have initial condition at x0, then if it's possible we can get the point where the algebraic function is equals to x0, let it be x1 ,and now we have the initial condition of resulting function at the point x1. For derivatives we will multiply the given initial condition with the corresponding derivative of algebraic function at x1. Is this method suitable enough?

I'd appreciate to know your views on it.

Thanks

[Ondřej Čertík]

On Sun, Mar 6, 2016 at 12:42 AM, shubham tibra <shubh...@gmail.com> wrote:
> Regarding the initial conditions in multiplication and addition when they
> are not at the same point:
>
> I think the best is(already discussed at ideas page) first to check if both
> the given holonomic functions are elementary function and calculate the
> initial conditions at a same point symbolically, where both the functions
> doesn't have a pole. So for the resulting function the initial condition
> would be addition/multiplication of calculated conditions at that point.
>
> In case when one or both can't be converted to elementary we will calculate
> numerical value at a point and add/multiply them.

Right. The numerical initial condition is the last resort, as at that
point, we can't really convert it back to a symbolic form. But I think
it's still useful, since if you know the holonomic function, but not
the initial condition (or only numerically), this still restricts what
the function looks like. For example the differential equation f''(x)
+ f(x) = 0 always has a solution of the form sin(a*x+b), so you know
that the function is a sinusoid, you just don't know the exact
parameters "a" and "b" without the initial condition, or only know
them numerically, so you can't convert it to, say, sin(x), but the
result is still useful.

>
> How are we gonna do this in INTEGRATION and DIFFERENTIATION?
>
> In application with a algebraic function, let's say we have initial
> condition at x0, then if it's possible we can get the point where the
> algebraic function is equals to x0, let it be x1 ,and now we have the
> initial condition of resulting function at the point x1. For derivatives we
> will multiply the given initial condition with the corresponding derivative
> of algebraic function at x1. Is this method suitable enough?

I don't understand your question. Can you formulate it on a particular example?

Ondrej

>
> I'd appreciate to know your views on it.
>
> Thanks
>

> --
> 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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit

> https://groups.google.com/d/msgid/sympy/21859ab1-f780-4d91-9e22-7f7303a47ad4%40googlegroups.com.

[Subham Tibra]

Hi Ondrej, thank you for your response.

ore_algebra in Sage has implemented `annihilator_of_composition` which finds annihilator of a holonomic function after application with an algebraic function but doesn't support initial conditions.

Do we want to get initial conditions of the resulting function when we will implement this operation or just the annihilator?

For instance, we have a holonomic function `sin(x)` defined by diff. eq. and we want to apply this with the function  z(x) = x**2 - 4

f(x) = holonomic(diff(diff(f)) + f, f(0) = 0, f'(0) = 1, x)
g(x) = f(z(x)) = f(x**2 - 4)
g'(x) = f'(z(x))*z'(x)                                            -(1)

Now, as we have value of f(x) at the point 0, we find the point where z(x) is 0 which is x = 2 or -2, if we take x=2 we have g(2) = 0, for finding g'(2) we put x = 2 in equation(1) so g'(2) = (1)*(2*2) = 4, and the annihilator of g will be obtained by the algorithm which we will be implemented.

Is this method decent enough for implementing, if we need the initial conditions?

I had some other problems, if you can please take a look at my previous post in this same thread that'd be very helpful.

I thought about the API and liked the way gfun does it in maple. I thought of something like this, initial condition will be given as a list [ f(x0), f'(x0), ... ]. 

In []: holonomic(diffeq, x, initial_condition, x0 = 0)
Out[]:holonomic(diffeq, x, f(x0) = f1, .... )

 
Is it possible to provide initial conditions like this way:

In []:holonomic(diffeq, x, f(x0) = f1, f(x1) = f2 .... )
Out[]:holonomic(diffeq, x, f(x0) = f1, .... )

where we can provide initial conditions as f(x0) = x1. We can also provide annihilators instead of differential equation as in sage. Which way will be more preferred for it's use in SymPy?

[Fredrik Johansson]

 
In general, one needs to allow generalized (singular) initial conditions. You could allow specifying a solution by a truncated series expansion, allowing singular terms, e.g. 1 + O(x^2) and x + O(x^2) for cos and sin, and as something like log(x) + x^(-1/3) + 5*x + O(x^2) for some function with a singularity.

I'm not sure if SymPy's way to represent series expansions is perfect for this (see the class for generalized series expansions in ore_algebra). Of course, you could just do it with a symbolic expression without explicitly giving the O term if that is a problem. However, ordinary initial conditions should be implemented first, and generalized conditions could be done later, so you could start by simply taking a list of initial values as input.

Fredrik

[Fredrik Johansson]

On Monday, March 7, 2016 at 7:12:44 PM UTC+1, Subham Tibra wrote:

Hi Ondrej, thank you for your response.

ore_algebra in Sage has implemented `annihilator_of_composition` which finds annihilator of a holonomic function after application with an algebraic function but doesn't support initial conditions.

Do we want to get initial conditions of the resulting function when we will implement this operation or just the annihilator?

For instance, we have a holonomic function `sin(x)` defined by diff. eq. and we want to apply this with the function  z(x) = x**2 - 4

f(x) = holonomic(diff(diff(f)) + f, f(0) = 0, f'(0) = 1, x)
g(x) = f(z(x)) = f(x**2 - 4)
g'(x) = f'(z(x))*z'(x)                                            -(1)

Now, as we have value of f(x) at the point 0, we find the point where z(x) is 0 which is x = 2 or -2, if we take x=2 we have g(2) = 0, for finding g'(2) we put x = 2 in equation(1) so g'(2) = (1)*(2*2) = 4, and the annihilator of g will be obtained by the algorithm which we will be implemented.

Is this method decent enough for implementing, if we need the initial conditions?

I had some other problems, if you can please take a look at my previous post in this same thread that'd be very helpful.

I thought about the API and liked the way gfun does it in maple. I thought of something like this, initial condition will be given as a list [ f(x0), f'(x0), ... ]. 

In []: holonomic(diffeq, x, initial_condition, x0 = 0)
Out[]:holonomic(diffeq, x, f(x0) = f1, .... )

 
Is it possible to provide initial conditions like this way:

In []:holonomic(diffeq, x, f(x0) = f1, f(x1) = f2 .... )
Out[]:holonomic(diffeq, x, f(x0) = f1, .... )

where we can provide initial conditions as f(x0) = x1. We can also provide annihilators instead of differential equation as in sage. Which way will be more preferred for it's use in SymPy?

The internal data structure should probably be the annihilator. The input could be given either as a differential equation or an annihilator. If the given differential equation is inhomogeneous, you need to convert it to homogeneous form, so a bit of preprocessing is needed anyway.

Fredrik

[Subham Tibra]

Thank you Fredrik for the clarification. I agree to use annihilators to store the function inside.

Regarding the numerical computation at a point, SymPy doesn't support numerical methods as of now AFAIK. As mpmath is included in SymPy, we can use  `odefun`  for numerical computation, but it doesn't support complex points. Is this good for now or do we need to implement a method when the point is complex?

[Fredrik Johansson]

On Tuesday, March 8, 2016 at 10:58:33 AM UTC+1, Subham Tibra wrote:

Thank you Fredrik for the clarification. I agree to use annihilators to store the function inside.

Regarding the numerical computation at a point, SymPy doesn't support numerical methods as of now AFAIK. As mpmath is included in SymPy, we can use  `odefun`  for numerical computation, but it doesn't support complex points. Is this good for now or do we need to implement a method when the point is complex?

mpmath.odefun does support complex points. You shouldn't use a general-purpose ODE solver, though. The holonomic ODE can be converted to a direct recurrence relation for the Taylor series coefficients, so you can analytically continue a solution very efficiently using Taylor polynomials of any order.

Fredrik

[Fredrik Johansson]

On Tuesday, March 8, 2016 at 2:20:14 PM UTC+1, Fredrik Johansson wrote:

On Tuesday, March 8, 2016 at 10:58:33 AM UTC+1, Subham Tibra wrote:

Thank you Fredrik for the clarification. I agree to use annihilators to store the function inside.

Regarding the numerical computation at a point, SymPy doesn't support numerical methods as of now AFAIK. As mpmath is included in SymPy, we can use  `odefun`  for numerical computation, but it doesn't support complex points. Is this good for now or do we need to implement a method when the point is complex?

mpmath.odefun does support complex points.

Oops, sorry, it doesn't do this directly. However, this is just a superficial restriction. It does support complex values, so you can always make a linear change of variables to move a complex path to a real interval.

Fredrik

[Subham Tibra]

On Tuesday, March 8, 2016 at 1:59:43 PM UTC+5:30, Fredrik Johansson wrote:

If the given differential equation is inhomogeneous, you need to convert it to homogeneous form, so a bit of preprocessing is needed anyway.

Does converting an inhomogeneous diff. eq. to homogeneous one, here means just removing the inhomogeneous part?

Inhomogeneous: a(x)*f(x) + a1(x)*f'(x) + ........ = g(x)

Homogeneous:   a(x)*f(x) + a1(x)*f'(x) + ........ = 0

[Fredrik Johansson]

No. Assuming that g(x) is holonomic (otherwise this won't work), you must use its annihilator to convert the inhomogeneous differential equation for f(x) to a homogeneous one.

By linearity, I think this just means computing an annihilator for h + g where h is a solution of the homogeneous part of the original equation. But don't quote me on that. I think the details are covered in some of the references.

Fredrik

[Tabot Kevin]

Hello Ondrej, I am really interested in this projects. Please can you point me to steps i can take to get started familiarizing myself with the SymPy environment for this project?

[Ondřej Čertík]

Hi Tabot,

The best is to try to figure out how the holonomic functions should be
implemented, as discussed in this thread.

Everyone --- if you are interested in this project, definitely
consider applying. We accept the best proposals, and if two or more
proposals are accepted for similar work, GSoC allows us to amend the
work, so that you work complements each other. We've done it in the
past, e.g. with regards to fast series expansion/polynomials.

Ondrej

[Subham Tibra]

Hi, I have created a pull request regarding this. Please give your suggestions and ideas here.

[Subham Tibra]

I have few questions.

When finding a recurrence relation for the series expansion, is it required just about the origin or at any arbitrary point?

In cases, when the ratio of terms can be found out directly by the recurrence relation or we can find the closed form of the recurrence, we can convert holonomic to hypergeometric. If these are not possible, then what would be our approach to convert to hypergeometric?

[Subham Tibra]

Hi Ondrej, 

Regarding the conversion of holonomic to hypergeometric, approaches I have in mind:

1. If the ratio of terms can be found out directly using the recurrence relation then one can get the hypergeometric representation
2. Getting a closed form of the recurrence and thus computing the ratio to convert to hypergeometric
3. Guessing the function from the coefficients in series expansions

Are there any more methods which I am missing? As you wrote here 

"To convert a holonomic function back to a hypergeometric one is also quite easy, as one can expand the holonomic function into a (formal) series, and then compute the ratio of terms and get the hypergeometric function out of it."

 

Does the ratio here means the numerical value? How one can get the hypergeometric function from the numerical values of the ratio?

[Ondřej Čertík]

On Tue, Mar 15, 2016 at 12:14 AM, Subham Tibra <shubh...@gmail.com> wrote:
> Hi Ondrej,
> Regarding the conversion of holonomic to hypergeometric, approaches I have
> in mind:
>
> If the ratio of terms can be found out directly using the recurrence
> relation then one can get the hypergeometric representation
> Getting a closed form of the recurrence and thus computing the ratio to
> convert to hypergeometric
> Guessing the function from the coefficients in series expansions
>
> Are there any more methods which I am missing? As you wrote here

I think those are the main once, at least as far as I know.

>>
>> "To convert a holonomic function back to a hypergeometric one is also
>> quite easy, as one can expand the holonomic function into a (formal) series,
>> and then compute the ratio of terms and get the hypergeometric function out
>> of it."
>
>
> Does the ratio here means the numerical value? How one can get the
> hypergeometric function from the numerical values of the ratio?

Not the numerical value, but a symbolic value of the ratio. Once you
have the ratio, then you factor the polynomial numerator and
denominator and read off the hypergeometric coefficients. See e.g.
here for details:

http://www.theoretical-physics.net/dev/math/hyper.html

> On Sunday, March 13, 2016 at 12:41:41 AM UTC+5:30, Subham Tibra wrote:
>>
>> I have few questions.
>>
>> When finding a recurrence relation for the series expansion, is it
>> required just about the origin or at any arbitrary point?

I am not sure about that. That's a good question. I would assume an
arbitrary point, but perhaps not.

>>
>> In cases, when the ratio of terms can be found out directly by the
>> recurrence relation or we can find the closed form of the recurrence, we can
>> convert holonomic to hypergeometric. If these are not possible, then what
>> would be our approach to convert to hypergeometric?

If it's not possible, then we have to use other means, perhaps some
kind of a pattern matching. A lot of functions can be written as a
linear combination of hypergeometric functions, so I would concentrate
on that.

Ondrej

>>
>> On Friday, March 11, 2016 at 9:48:14 PM UTC+5:30, Subham Tibra wrote:
>>>
>>> Hi, I have created a pull request regarding this. Please give your
>>> suggestions and ideas here.
>>>
>>> On Tuesday, March 8, 2016 at 10:17:24 PM UTC+5:30, Ondřej Čertík wrote:
>>>>
>>>> Hi Tabot,
>>>>
>>>> On Tue, Mar 8, 2016 at 4:59 AM, Tabot Kevin <tabot...@gmail.com> wrote:
>>>> > Hello Ondrej, I am really interested in this projects. Please can you
>>>> > point
>>>> > me to steps i can take to get started familiarizing myself with the
>>>> > SymPy
>>>> > environment for this project?
>>>>
>>>> The best is to try to figure out how the holonomic functions should be
>>>> implemented, as discussed in this thread.
>>>>
>>>> Everyone --- if you are interested in this project, definitely
>>>> consider applying. We accept the best proposals, and if two or more
>>>> proposals are accepted for similar work, GSoC allows us to amend the
>>>> work, so that you work complements each other. We've done it in the
>>>> past, e.g. with regards to fast series expansion/polynomials.
>>>>
>>>> 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 https://groups.google.com/group/sympy.
> To view this discussion on the web visit

> https://groups.google.com/d/msgid/sympy/31238a9b-1830-4e47-86a6-5b5f1f44ede1%40googlegroups.com.

[Subham Tibra]

I have a question regarding converting an expression tree to holonomic

For instance, we have an expression expr = sin(x)*exp(x), So traversing this gives sin(x) and exp(x). Once we have the holonomic representation of sin and exp we can use our multiplication algorithm to convert the whole expression to a holonomic one.

So how is one going to get the holonomic representation of sin and exp? Do we have create some kind of table to store all the holonomic representation of elementary functions or should we develop a function to find holonomic representations?

[Ondřej Čertík]

Hi Subham,

On Fri, Mar 18, 2016 at 9:36 AM, Subham Tibra <shubh...@gmail.com> wrote:
> I have a question regarding converting an expression tree to holonomic
>
> For instance, we have an expression expr = sin(x)*exp(x), So traversing this
> gives sin(x) and exp(x). Once we have the holonomic representation of sin
> and exp we can use our multiplication algorithm to convert the whole
> expression to a holonomic one.
> So how is one going to get the holonomic representation of sin and exp? Do
> we have create some kind of table to store all the holonomic representation
> of elementary functions or should we develop a function to find holonomic
> representations?

There are multiple ways, but the easiest is to use what's already in
SymPy via the hypergeometric functions, convert sin(x) and exp(x) to
hypergeometric functions, and then convert hypergeometric functions to
holonomic ones via the simple code I posted before.

Ondrej

[Subham Tibra]

Hi !

I have shared my draft with SymPy, the Proposal's Title is "Implementation of Holonomic Function in Sympy".

Any reviews and suggestions would be truly helpful for me.

Thanks

Subham

[Henrique Gasparini Fiuza do Nascimento]

I have just submited my proposal for this project and I would be very interesting to work on it. Please excuse-me for not discussing anything about it with you previously. I was told that the deadline was today yesterday evening and I don't wanna miss this opportunity.
Please read my proposal, I hope you appreciate it.
Thank you,
Henrique