Methods in Mathematica for Solving Ordinary Differential Equations

An overview of the solution methods for ordinary differential equations in the Mathematica function DSolve is presented.


INTRODUCTION
As they arise in the mathematical formulations of real-life problems, differential equations play a central role in displaying the interrelations between mathematics, physical and biological sciences or engineering [1]. Hence studying the solution methods of differential equationbs has always been an important research area [2,3,4]. With the invention of computer algebra systems, ongoing efforts for finding new methods for computation of solutions of differential equations led to exciting developments. An understanding of the scope and built-in algorithms in such systems is very useful while applying them in practice as they typically allow for a variety of approaches (symbolic, numerical, and graphical) for solving differential equations.
In this paper, we give an overview of available methods for solving ordinary differential equations (ODEs) in closed form and give examples of these methods in action as they are being used in DSolve, the function for solving differential equations in Mathematica [5], a major computer algebra system. In section 2, we give a list of methods for solving first-order ODEs. Section 3 contains the methods for solving second or higher order linear ODEs. In section 4, we deal with second or higher order ODEs which are nonlinear. Section 5 contains the methods for handling systems of ODEs. The examples have been chosen to illustrate the structure of typical members of each class of ODEs and we have tried to give insight into the key ideas and algorithms which are used to solve them. We have also included applications of differential equations to population dynamics and differential geometry. The final section offers suggestions for pursuing the study of differential equations in greater depth using Mathematica.

FIRST-ORDER ODEs
As they also become useful when solving higher order equations and systems of ODEs, studying the solution methods of first-order ODEs is really important. For some classes of first-order ODEs solution methods are known and well-studied [3,4]. Among these classes of equations, we can list: linear, separable, Bernoulli, homogeneous, inverse linear, exact and Clairaut type first-order ODEs. Riccati, Abel [6] and Chini type first-order ODEs are also well-studied and many sub-classes of these equations for which a solution can be found are identified. For first-order ODEs which do not fit into one of these classes, one can try: recognizing the defining differential equations of the special functions of physics, finding integrating factors or computing Lie symmetries of the ODE [7,8,9].
In the well-known books by Kamke [3] and Murphy [4], the standard methods are covered. When dealing with a first-order ODE, based on the statistics of equations listed in [3], one can choose the sequence of the methods to apply. A great portion of the 576 first-order ODEs listed in [3] can be solved with the standard methods for solving linear, separable, Riccati and Abel type equations. For example, 103 equations fit into the separable class.

Linear equations
Linear first-order ODEs are identified as ) and the solution is given after integration: In the above solution, ] 1 [ K and ] 2 [ K denote the dummy integration variables. To suppress possible messages generated by DSolve, we initially ran the following command:

Separable equations
First-order ODEs which can be written as are called separable equations and the solution is again given after integration: The logistic equation Methods in Mathematica for Solving Ordinary Differential Equations

Bernoulli type equations
Equations of the form are called the Bernoulli type equations and the solution is found after integration. Following example is the equation 1.34 from [3]:

Homogeneous equations
A first-order ODE of the form

Inverse linear equations
ODEs of the form

Riccati type equations
First-order ODEs of the form are called Riccati type equations. The methods for solving the Riccati type equations are studied in [3,4]: The following Riccati type ODE is first transformed into a linear second-order ODE and once the second-order ODE is solved, solution to the original first-order Riccati type ODE is found by inverse transformation: Note that, in the above general solution if an initial condition such as 1 ) 0 (  y has to be satisfied, this could be handled by using limit:

Chini type equations
Equations of the form are called Chini type firstorder ODEs and various cases when a solution is directly formed are given in [3,4]:

Abel type equations
Equations of the form are called Abel equations of the first kind. Abel equations of the second kind are of the form . The cases where direct solution methods are available are studied in [3,4] and [6]: To verify this implicit solution, we can do:

Exact equations (and computing integrating factors)
A first-order ODE is exact if is an integrating factor. This means that you make the ODE exact if you can find an integrating factor [3,4]:

Lie symmetry methods
The knowledge of a symmetry in the form of an infinitesimal generator reduces the order of an ODE which typically simplifies the problem. DSolve function checks for standard types of symmetries in the given ODE and uses them to return the solution [7,8,9].
For the following ODE 1.357 of [3], DSolve uses the symmetry methods to find the solution: For the ODE 1.188 of [3], DSolve returns an implicit solution:

Clairaut equations
ODEs of the form are called Clairaut equations [3]:

Recognizing defining equations of special functions
Several elliptic functions are defined by first-order ODEs as illustrated below for the WeierstrassPPrime function:

LINEAR SECOND OR HIGHER ORDER ODEs
There are various standard methods for solving linear second or higher order ODEs. Equations with constant coefficients are solved using the roots of the characteristic equation. The Euler-Legendre type equations can be transformed into equations with constant coefficients. Some equations are exact and so can be directly integrated whereas some can be directly integrated after multiplying with an integrating factor. For second-order equations with rational coefficients there is the well-known Kovacic algorithm for finding the Liouvillian solutions [10]. For second or higher order equations with rational coefficients, in [11,12] and [13] the algorithms for finding rational and exponential solutions are given. The hypergeometric PFQ type solutions are found using the algorithms in [14] and [15]. In [16] an algorithm using the concept of symmetric powers, and in [17] an algorithm using the concept of symmetric products are given. There are also various factorization algorithms for higher order ODEs in [18] and [19].

ODEs with constant coefficients
The general solution of linear ODEs with constant coefficients are found by using the roots of the characteristic equation for the ODE:

Euler-Legendre type equations
The following equation is an Euler-Legendre ODE which can be solved by transforming it to a linear ODE with constant coefficients:

Exact equations (and integrating factors)
The following is an example of an exact ODE since the left-hand side can be integrated to a first-order expression whose solution gives one element of the general solution of the second-order ODE:

Recognizing the defining equations of special functions
There are many physical processes which are modeled by linear second or higher order ODEs. DSolve can recognize the defining ODEs for many special functions:

Kovacic algorithm
This is a standard algorithm for solving second-order linear homogeneous ODEs with rational function coefficients [10]: In the above solution, the piece with C [1] is found by the Kovacic algorithm, and the second independent solution is found by reduction of order.

Rational and exponential solutions
For finding rational solutions of linear ODEs with rational coefficients the algorithms in [11] and [12], and for finding exponential solutions the algorithm in [13] can be used. Here is a nice example where we can see a nice harmony of the methods for rational solutions, reduction of order and recognizing special functions:

Hypergeometric PFQ type solutions
Hypergeometric PFQ type solutions are found by using the algorithms in [14] and [15] . This fourth order ODE is called the symmetric product of the secondorder equations [17]: Here is the solution of the symmetric product of these ODEs:

Factorization
DSolve has the implementations of factorization algorithms in [18] and [19]:

Equations solved after a transformation
DSolve uses a number of transformation rules in order to solve ODEs whose coefficients are rational in either trigonometric, hyperbolic or exponential functions. The basic idea is to transform the given ODE into one whose coefficients are rational in the new independent variable.
In the following example, we use DSolve to find quantum eigenfunctions for a modified Pöschl-Teller potential, which requires the solution of a linear second-order ODE with hyperbolic coefficients:

NONLINEAR SECOND OR HIGHER ORDER ODEs
DSolve has special methods for solving important classes of nonlinear ODEs which arise in applications. This includes the standard methods for handling the equations with missing variables, exact equations, and homogeneous equations [3,4], computing integrating factors due to the algorithms in [20], and computing the Lie symmetries [9] and applying transformation rules that reduce the equation to one of the standard types.

Equations with missing variables
In the following example, the differential equation does not depend explicitly on the independent variable x , and a first integral can be found after making the substitution

Transformations
It is sometimes possible to remove the explicit dependence of an ODE on the independent variable by using a transformation. For instance, the following ODE is solved by making the transformation , which leads to a differential equation for ) (x u with missing independent variable x :

Homogeneous equations
The left-hand side of the following nonlinear ODE is a homogeneous function of degree 2 in the variables

Exact equations (and integrating factors)
In [20] methods for finding integrating factors of the form x y x   is used as an integrating factor for solving the following ODE:

SYSTEMS OF ODEs
DSolve has a variety of methods for solving systems of ODEs with constant or variable coefficients. Systems with higher-order derivatives are internally reduced to first-order systems and, wherever possible, the system is decoupled to reduce the problem to solving a set of independent single ODEs. We will now give a few examples for solving systems of ODEs with increasing levels of complexity to illustrate the techniques used for solving them.

Systems with constant coefficients
Linear systems with constant coefficients are solved by applying matrix exponentiation to the coefficient matrix:

Systems of decoupled equations
The following example illustrates a system composed of decoupled equations in which each equation involves a single dependent variable only. In such cases, the equations are solved independently using the available methods for single ODEs.

Systems solved by one at a time integration
The following system of ODEs is solved by integrating the first equation which involves only )

Systems with patterns
The following linear system of ODEs is solved in closed form as the coefficient matrix has a special structure:

Rational solutions
For linear systems with rational coefficients, the algorithms in [11,12] for finding rational solutions are used:

Application: Solving the equations for a Cornu spiral
The following equations describe a Cornu spiral. This is a plane curve through the origin whose curvature is equal to the parameter value s at every point.

The solution of the ODEs contains Fresnel functions.
A graph of the curve clearly shows that the curvature becomes large as s approaches infinity:

CONCLUSION
Further information on DSolve is available on the documentation page for this function [21]. The tutorial [22] gives a comprehensive overview of the functionality available in DSolve for symbolic solutions of differential equations, along with further references for this topic. Numerical solutions of differential equations can be computed using the NDSolve function [23] in Mathematica.
We hope that the examples and ideas outlined in this paper will be useful for elementary and advanced courses on ordinary differential equations, as well as for solving differential equations which occur in research and design problems in practice.