Implementing support for 2nd Order ODE and Hidden Separable Equations

About Me

Basic Information

Name : Rahil Hastu

E-mail : [email protected]

University : PESIT South Campus, Bangalore

Github : rahilhastu

Timezone : IST (UTC+5:30)

Personal Background

Who am I ?

I’m a 3rd year undergraduate pursuing B.E. in Computer Science Engineering from PESIT South Campus, Bangalore.My fields of interest include coding and mathematics. My first tryst with python came in my freshman year, where I tried rewriting the programming assignments given in C or C++, in python and from there on I have only used Python for my projects, some of them are listed below:

1. https://github.com/rahilhastu/subjectiveEvaluation :It was a project done to evaluate the answer scripts automatically using NLTK(Natural Language Toolkit).

2. https://github.com/rahilhastu/competitiveSitesData/blob/master/spojData.py : This is a script for scraping the users of SPOJ (This is a website for competitive programming that is built basically to help students and professionals to enhance their coding skills.)

3. https://github.com/rahilhastu/competitiveSitesData/blob/master/coderByteData.py : This is a script for scraping the users of CoderByte (CoderByte is just like Spoj)

A few of my other projects can be seen at https://github.com/rahilhastu. Also I have been consistently exposed to the field of Mathematics since my high school, and I have taken courses in Advanced Engineering Mathematics-I,II,III,IV, Data Structures and Algorithms Design in my undergrad.I have also taken a course in my undergrad for Differential Equations.

Course Details : 15MAT21

Programming Details

I work on Ubuntu 17.10 LTS machine with Sublime Text 3 as my primary text editor because of its user-friendliness. I am very much familiar with version control system and have been using Git(command line) and GitHub for quite some time(more than a year). I have been programming for over 2 years now and have been coding in Python and few of its frameworks. Python has relatively fewer lines of code, which makes it less prone to issues, easier to debug and is more maintainable.Python can scale to solve complex problems, as illustrated by the fact that it powers most of YouTube and DropBox, not to mention Reddit, Quora, Disqus and FriendFeed. Even the mighty Google has made Python one of its official programming languages.

One of my favourite features of python is the NLTK library which I have used in one of my projects listed above (1).

My favourite feature about SymPy is the dsolve and classify_ode feature. It was one of the first features of SymPy I used while trying to familiarize myself with the library. The ability of the software to recognize and evaluate a differential equation amazed me:

>>> eq = Derivative(f(x),x)-6*f(x)**2*x

>>> classify_ode(eq)

('separable', '1st_exact', 'Bernoulli', 'separable_reduced', '1st_power_series', 'lie_group', 'separable_Integral', '1st_exact_Integral', 'Bernoulli_Integral', 'separable_reduced_Integral')

>>> dsolve(eq,f(x))

Eq(f(x), -1/(C1 + 3*x**2))

Contributions to Sympy

I have been an active member for around 1½ months. The total number of patches may be less in number but I consider myself as a fast learner. Following is the list of contributions to Sympy in forms of Patches.

Patches

• (Merged) Moved arithematic operators to sympy/sets/sets.py #14372
• (Merged) kinematic_equation returns zero if speeds are zero #14280
• (Open) Extracting common codes from Latex.py #14222

I am also active on Gitter (@rahilhastu) and a conversation led to opening an issue related to integral solving while dealing with 1st order ODE.(issue was opened by asmeurer).

Issue

• (Open) integrate(exp(x**2)*cos(x), x) cannot compute the integral #14398

Project Idea

Abstract

The current ODE solver in sympy was the work started by Aaron Meurer as his GSoC Project in 2009. The ODE solver, though robust, has a good number of possible features that can be added. And I’m enthusiastic and excited to work on this module as further improvement has been lingering on since 2014. My project proposal would be dealing with the implementation of two of these features:

• Discovering a differential equation is separable or not when the nonlinearity is not explicitly factorized

The current ODE module of Sympy successfully classifies the first order differential equations which are classified into linear,exact,bernoulli,liouville but sometimes a hidden separable differential equation fails to get recognized . By the end of my project the ode solver would be able to recognize the hidden separable 1st order differential equation .

• Solving non classifiable second order ODE using Lie Groups

Sympy has a support for solving the first order differential equations that can be classified as separable, linear, Bernoulli, Riccati, exact, homogeneous, almost-linear, separable-reduced, and linear coefficients. By the end of my project, the ode solver would be able to solve 2nd order differential equation which cannot be classified into any of the above given categories.

Inspiration

My project proposal idea is partly inspired by a course that I had taken in my 3rd semester 15MAT21 which involved differential equations and partly by the fact that since I have come to know about Sympy, I mostly like to work with its ODE module.

Theory

• Separable Differential Equation

A 1st order Differential Equation, y'= f(x,y) is separable if the right hand side of the equation can be factorized i.e . After factoring the RHS, we just need to integrate both sides of the equation.

But sometimes the given differential equation is mistaken to be any other kind of Differential Equation, but, instead it is a hidden separable differential equation .

For ex.

Eq is a kind of separable differential equation but on verifying the same in ODE module we get

>>>classify_ode(Eq)
>>>('1st_exact', '1st_power_series', 'lie_group', '1st_exact_Integral')

• Lie group methods for 2nd order ODEs

Take the case of most general non constant linear 2nd order differential equation

and a case for constant 2nd order differential equation

To solve these type of 2nd order differential equations we will implement Lie Symmetry methods(SM). The same was achieved for the 1st order Differential Equation in GSoC 2013 by Manoj Kumar. The key difference here will be the implementation of SM for high order ODEs and solving the determining equations for the infinitesimals. The trickiest thing to solve the 2nd order ODEs will be that that the knowledge of n symmetries does not directly reduce the problem to a line integral as in the 1st order case.

The key point of Lie’s method for solving ODEs is that the knowledge of a (Lie) group of transformations which leaves a given ODE invariant may help in reducing the problem of finding its solution to quadratures.

Considering the most general second order differential equation of the form

Eq.(1)

`

The basic result in this case is that if we know a one-parameter symmetry group of this equation,then we can reduce the order of the equation by one.

The algorithm for finding the solution using Lie Groups is as follows:

• Finding pair of infinitesimals of a one-parameter Symmetry Group which leaves Eq.(1) invariant. This is done by solving the PDE

Eq.(2)

(2*y1*(∂^2η/(∂y*∂y1)) + 2*(∂^2η/(∂x*∂y1)) - 3*y1*(∂ξ/∂y) - 2*(∂ξ/∂x) - 2*y1*(∂^2ξ/(∂x*∂y1)) + (∂η/∂y) -
 2*y1^2*(∂^2ξ/(∂y*∂y1)))φ + (-y1*(∂η/∂y) + y1^2*(∂ξ/∂y) + y1(∂ξ/∂x) - (∂η/∂x))*(∂φ/∂y1) + (∂^2η/∂x^2) + 
(∂^2η/∂y^2)*y1^2 + ( y1(∂η/∂y1) + y1^2*(∂ξ/∂y1) - η )*(∂φ/∂y) + ( (∂^2η/∂y1^2) - y1(∂^2ξ/∂y1^2) - 2(∂ξ/∂y1) 
)*φ^2 + ( (∂η/∂y1) - y1(∂ξ/∂y1) -ξ )*(∂φ/∂x) - (∂^2ξ/∂y^2)*y1^3 - 2*(∂^2ξ/(∂y*∂x))*y1^2 + 2(∂^2η/∂y*∂x)y1 - 
y1(∂^2ξ/∂x^2) = 0

where are infinitesimals and .

Integrating these infinitesimals by any of the following three methods :

• Determining the differential invariants of order 0 and 1 of the symmetry group
• Determining first integral (from the knowledge of )
• Determining the canonical coordinates {r,s(r)} of the associated Lie group.

Ex

Let there be an n-order differential equation involving the single dependent variable u (in our case that n has to be 2):

Where

We can get a solution to this differential equation if we have one-parameter symmetry group so that it can help in reducing the order of the equation by 1. So if we have a 2nd order differential equation then its order, on finding the symmetry group is reduced to one, of which we have a way to find a solution. We choose coordinates such that the group transforms into a group of translations with infinitesimals generator ∂/∂w.We express the derivatives of u with respect to x in terms of y, w and the derivatives of w with respect to y,

For some certain function . Substituting these expressions into our Eq to get an equivalent equation:

`-(3)

in terms of new coordinates, y and w and the infinitesimal has trivial prolongation v = ∂/∂w

The infinitesimal criterion of invariance implies

This means that there is an equivalent equation

Which is independent of w,i.e, if and only if .Now on setting z =w we have an (n-1) order equation for z,

-(4)

The solutions of (4) provide a general solution to our original equation.If z = h(y) is a solution of (4) ,then is a solution of (3). Therefore on replacing w and y by their expressions in terms of x and u ,implicitly defines a solution of the original equation.

Implementation

• Separable Differential Equation

Let D be a subset of and . We know from a theorem[1] that

when .

For proving it the other way around ,assume f(x,y) is not zero on D.

Let be such that and let there be a function defined by and , which helps in reaching the fact that

This is applicable to wider class of functions as there are no assumptions on f.

This method is not only meant for 2nd order but is generalized for functions having n variables (so it can be used in the near future when Sympy will be dealing with higher order differential equations).

So the algorithm that will help in deriving the above mentioned separable condition for differential equation will be :

• Choose any (x0,y0) such that .

• Construct the functions and .

• The equation is separable if and only if .

• Lie group methods for 2nd order ODEs

As mentioned in the theory the first and most difficult step would be solving the PDE ,Eq(2). To find the infinitesimals which would be tougher than solving the differential equation itself. Interestingly, Maple has a set of six robust algorithms, which consists of heuristics or guessing solutions for η and ξ , so what I would be trying to write is essentially a clone of that in sympy(the same has been accomplished for 1st order differential equations in GSoC’13 by Manoj Kumar). The six algorithms as described by Maple are

• The first set of the algorithm is going to formulate a linear PDE system for the infinitesimals , ,then triangularize the system using differential algebra techniques. This will return the complete set of point symmetries.Ex: For a 2nd order differential equation like

y''+(y+3f(x))y'-y3+y2f(x)+y(f'+2f(x)2)=0

The infinitesimal generated will be

and

• The 2nd set of algorithms, deals with one of and , to be zero and the other to be a function of x or y or y1. This will reduce the PDE in Eq.(2) to a simpler ODE.

Ex : Consider a 2nd order nonlinear ODE

The pair of infinitesimal’s we will get are:

• The 3rd and 4th set of algorithms would be dealing with the polynomials in x,y,y1solutions to Eq.2

Let us take an example given in [2], which is . Putting =x and =y we get the solution to the ODE equal to .

• The fifth algorithm looks for rational ansatze in (x,y,y1) solutions to Eq.(2) for the infinitesimals,, whose maximum degree can be determined by observing the differential equation.

Let us look into an example provide in [2] i.e. . On using this algo we get and (we can get more than one infinitesimal). And using these infinitesimals we get to the solution .

• The sixth algorithm uses a polynomial of degree two constructed from a basis of functions and algebraic objects, together with their derivatives, taken from the given ODE[2]. The symmetry found can be of a dynamic type. And the last algorithm which I’ll be implementing will look for the form where q1 and q2 are functions of [x,y,y1].

The problem of finding a solution to the 2nd order ODE once the infinitesimals are found isn’t solved yet. The methods discussed now will be dealing with integration of the ODE even when the symmetries found are of dynamic type.

• Using differential invariants

Once a pair of infinitesimals for a given 2nd order ODE is found, it is possible to integrate the ODE if one succeeds in using these infinitesimals to find the differential invariants of order 0 and 1 and in solving an intermediate 1st order ODE.What I will be trying to achieve is, first I need to transform the symmetry to evolutionary form, which will result in keeping x as invariant(order = 0).

The algorithm which I will be following will be:

• Put the symmetry in the evolutionary form( )

• Use u = x as invariant of order 0 and find the invariant of order 1, v = v(x,y,y1)

• Then, I’ll construct the invariant of order 2: dv/du

• The variables need to be changed in terms of u , v, dv/du which will help in obtaining the expected reduction of order

• solve this 1st order ODE, and reintroduce the original variables (x,y,y1), obtaining another 1st order ODE, which can be integrated straightforwardly since it has the same symmetry of the original 2nd order ODE.

All the occurences of y1 needs to be changed in the infinitesimal in last step when the symmetry is of dynamic type.

Example

Consider the 2nd order equation

The given equation is invariant under the scaling group .

On integrating the above equation using the method of differential invariants, we get

The new equation involving w and u is

The new equation has two families of solutions; either w=y or dw/du = -w. On integrating dw/du = -w, we get ,c being a constant. On changing back to the original variables we get two homogeneous first order equations

The first has solution ; the second has implicit solutions

where u=y/x

• Using First Integrals

There is one more way to integrate a 2nd order ODE,i.e., to use the pair of infinitesimals to build 2 first integrals from which a closed solution can be obtained.Using this method we need only 1 pair of infinitesimal by looking for 1st integrals.

The algorithm which I’ll be sticking with is similar to the one used in the previous method:

• Put symmetry in evolutionary form. There will be a first integral which will satisfy

and

Where and

• Then solve X( )=0 to obtain invariants of order 0 and 1 : c0(q) and c1(q) where q =f(x,y,y1).

Because of the evolutionary form c0=x, X( ) is reduced to a single 1st order ODE.

• The variables in A( ) needs to be changed from (x,y,y1) to (c0,y,c1). Doing this results in

, a PDE in 2 variables. This can be solved to obtain .

• Change the variables back to the original variables .This is a 1st order ODE which can be integrated easily.

Example

Let the ODE be:

This has point symmetries. The infinitesimals found are

On integrating it the solution we find is:

• Using Canonical coordinates

As discussed in the theory, a pair of infinitesimals is used to reduce the order of a 2nd order ODE by introducing the canonical variables.But the use of this method is overshadowed in case of dynamical symmetries by the methods discussed previously.

Examples(Prototype)

>>> from sympy import *
>>> from sympy.abc import x
>>> f = symbols('f ', cls=Function)
#Separable Equation
>>> eq = Derivative(f(x),x)-(E**(x**2+f(x)**2)*(cos(x+f(x))+cos(x-f(x))))
>>> classif_ode(eq)
Output : (‘Separable’ ,'1st_exact', '1st_power_series', 'lie_group', '1st_exact_Integral')

#2nd Order ODE
>>>dsolve(Derivative(f(x),x,x)-(1/(x+f(x))*(2*Derivative(f(x),x)+1)Derivative(f(x),x)),hint=’lie_group’)
Output : f(x) = -(x**2+c1)/(2*x+2*c2)
>>>infinitesimals_second(Derivative(f(x),x,x)-(1/(x+f(x))*(2*Derivative(f(x),x)+1)Derivative(f(x),x)))
Output : {[-1/(x+y),0], [-y/(x+y),0], [x*(x+2*y)/(x+y),0],[-1,1]}

#2nd non-linear ODE
>>> dsolve(Derivative(f(x),x,x)*((x+y)*x**2)-(x*Derivative(f(x),x)-y)**2)
Output : f(x)  = x*e**((-c1/x)+c2)-x
>>>infinitesimals_second(Derivative(f(x),x,x)*((x+y)*x**2)-(x*Derivative(f(x),x)-y)**2)
Output : {[x,y], [-x,x], [x**2,x*y]}

#Higher Order Differential Equation(>2)
>>> dsolve(Derivative(f(x),x,x,x)-(x**2*E**f(x)))
Output : NotImplementedError: solve: Cannot solve -x**2*exp(f(x)) + Derivative(f(x), (x, 3))

Timeline

Community bonding period(April 23, 2018 - May 14, 2018)

Goal : Community Bonding

• The principal focus in this period would be studying in detail ,the implementation of the algorithms and methods of integration for Lie group implementation of 2nd order ODE and ask guidance from my mentor, because I feel that would be the most challenging part of my project.

• If possible, I would also start coding in this phase only, so that I get a head-start.

Week 1 - Week 4 (May 14 - June 11)[Midterm Evaluation]

Goal : Implement the Separable Differential Equation and 2 algorithms for Lie Groups

• I plan to implement the Algorithm required to determine the separability of a differential equation.
• Implementation will correspond to a pull request.
• After I’m done with separable differential equation, I’ll move on to implement the first 2 algorithm’s for Lie Groups 2nd order ODE.
• The algorithm will correspond to a pull request.

Week 5 - Week 8 (June 13 - July 11)[Midterm Evaluation]

Goal : Implement the algorithms for Lie Groups 2nd Order ODE

• Finish coding the rest of the algorithms, if possible all the remaining of them will be implemented
• Each algorithm will correspond to a pull request.
• Fix bugs if any.

Week 9 - Week 11 (July 13 - Aug 2)

Goal : Implement Integration techniques

• If all the algorithms have been implemented, I’ll move towards the implementation of the integration techniques as mentioned in the Implementation Section, else one more week may be needed to complete the implementation of the algorithms for finding infinitesimals.

• Each implementation of the integration technique will correspond to a pull request.

Week 12 (Aug 2- Aug 11)

Goal : Final Evaluation

• Write improved tests, documentation and fix bugs that may arise.

Future work

All the methods will be structured in such a way that if anyone wants to add support for high order differential equations it can be done easily.

Note

From the start,i.e community bonding period, I’ll be devoting 6-7 hours daily as I don’t have any other big project devoted for the summer but I’ll be having my exams in June stretched for 2 weeks or less(dates not finalised yet, when finalised, the strategy for the Second evaluation may have to be discussed with the mentor). For those 2 weeks I’ll be able to spend no more than 4-5 hrs.

References

1. "A simple method to find out when an ordinary differential equation is separable" by José ́Ángel Cid

2. "Computer Algebra Solving of Second Order ODEs Using Symmetry Methods" by E.S. Cheb-Terrab, L.G.S. Duarte, L.A.C.P. da Mota

3. "Computer Algebra Solving of First Order ODEs Using Symmetry Methods" by E.S. Cheb-Terrab, L.G.S. Duarte and L.A.C.P. da Mota.

4. P.J. Olver,Applications of Lie Groups to Differential Equations, Springer-Verlag

5. G.W. Bluman and S. Kumei,Symmetries and Differential Equations, Applied Mathematical Sciences 81, Springer-Verlag, (1989)

6. Lie Groups and Differential Equations,Gilmore_LG

7. Maple Programming help.

8. GSoC 2013 Application Manoj Kumar: Improved ODE Solver