Symmetry Analysis of Ordinary Differential Equations (ODEs)

1. Introduction

The subject of this section is the discussion of a general procedure allowing the solution of an ordinary differential equation (ODE) of arbitrary order. We discuss the solution of first order, second order  and higher order ODEs. The technique of first integrals is one of several methods to integrate an ODE provided  the symmetries of the equation are known.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_1.gif]

The basis of our discussion is a theory conceived  by S. Lie (1842-1899). Lie developed a general theory dealing with symmetries and group properties of differential equations. The theory of Lie is a valuable tool for solving ODEs and partial differential equations. In the current text we restrict our considerations to ordinary differential equations. The discussion of partial differential equations can be found in standard books by Olver [1986]; Bluman and Kumei [1989]; Ibragimov [1994,95,96]; and Baumann [1998]. The impact of Lie's theory is not restricted to differential equations but is omnipresent in the field of mathematics. Lie deals with fundamental areas of modern mathematics combining group theory and analysis in a nut shell. One century after Lie, continuous groups are today the central and most important part of modern mathematics.

Our aim is to show that Lie's theory is instrumental for solving an ordinary differential equation. We will show you a method capable of deriving solutions of  an ODE of arbitrary order. The basis  of these considerations are the symmetries of the equation.

2. Symmetries

2.1 Overview

The word symmetry is used in our everyday language in different meanings. In the one sense symmetric means something like well-proportioned, well-balanced, and symmetry denotes that sort of concordance of several parts by which they integrate into a whole. Beauty is bound up with symmetry. Ebenmass is a good German equivalent for the Greek symmetry.

Symmetry, as wide or as narrow as you may define its meaning, is one idea by which man through the ages has tried to comprehend and create order, beauty, and perfection.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_2.gif]

Current research topics very closely related to symmetry are theories concerned with elementary particles (quantum field theories, theory of relativity, string theories, T.O.E. etc.). Actually symmetries are not only important in current theories but also in ancient considerations undertaken by Galilei, Kepler, and Newton. For example, Kepler demonstrated the existence of symmetries in his work Harmonice Mundi by deriving the distance period relation of a planet. The famous law T^2/a^3 = const . was derived by dimensional analysis of Brahes observation data.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_4.gif]

Since the formulation of physical laws by using differential equations by Leibniz and Newton continuous symmetries were introduced in the mathematical and physical laws. Leibniz and Newton were certainly not aware of this fact but they actually  established the basis of today symmetry analysis with their calculi.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_5.gif]

However, the mathematical notion of a symmetry in connection with differential equations was introduced by Lie 200 years after Leibniz and Newton. Before we apply Lie´s method to ordinary differential equations, we will briefly discuss symmetries in connection with functions. This section introduces the main concepts of the theory.

One of the remarkable aspects of Lie's theory of symmetry transformations is the invariance of a function under the change of coordinates. When dealing with differential equations, oftentimes the equation is simplified by an appropriate change of variables. This transformation generally involves both the independent and the dependent variables. In Mathematica, we can represent such transformations by the following list of rules

x = X(x, u) <br /> u = U(x, u)

This kind of transformation involving the original independent and dependent variables (x,u)  is called a point transformation. A point (x, u) of the manifold M is transformed into another point (X, U). Point transformations only take into account changes of the coordinates. The point transformations considered by Lie are transformations depending on at least one parameter ε. As we will see, the parameter ε is the parameter of the related group. Thus we call ε the group parameter of the transformation. A one-parameter transformation is thus given  by

x = X(x, u ; ϵ) <br /> u = U(x, u ; ϵ)

Such transformations have the following properties:

They are invertible if a corresponding Jacobi determinant exists.

Their repeated applications yield  transformations of the same kind.

An identity of the transformation for ε=0 is given.

These three properties are the basis of the definition of a Lie group. The three properties are summarized in the definition of a symmetry transformation.

Definition:  Symmetry Transformation

A set G of transformations given, for example by (2), is a one-parameter group if it contains the identical transformation I = T_0 and includes the inverse T_ϵ^(-1) and the composition T_ϵ ⊗ T_β ∈ . By choosing the group parameter to be ε, the main group property T_ϵ ⊗ T_β ∈   can be written as

 T_ϵ ⊗ T_β = T_ (ϵ + β)

that is

X(X(x,u,ε),U(x,u,e),β)  =  X(x,u,ε+β)

U(X(x,u,ε),U(x,u,ε),β)  ==  U(x,u,ε+β).

In special applications, conditions (4) and (5) only hold for sufficiently small values of ε and β. In such cases, we arrive at what is called a local one parameter group G or an infinitesimal group. ●

These properties ensure that the transformations given in (2) form a one-parameter group of point transformations.  

2.2 Infinitesimal transformations

One of the most important aspects of Lie's theory is the discovery that a symmetry transformation as given above can be represented in a much simpler form. Lie called this simplified representation an infinitesimal transformation. The main idea behind this simplification is the existence of an identity transformation. An identity transformation denotes the lowest approximation of a transformation around the identity. Correspondingly, the finite transformation can be expanded into a Taylor series around the identity transformation. In Mathematica we have direct access to such an expansion when applying the function Series[] to the global representation. The following is an example for a general point transformation in the (x, u) plane. We derive the infinitesimal representation by expanding the one-parameter transformation around ϵ = 0 up to the first order in ε

infinitesimalTrafo = Series[{X(x, u, ϵ), U(x, u, ϵ)}, {ϵ, 0, 1}] ; infinitesimalTrafo//TableForm

X[x, u, 0] + X^(0, 0, 1)[x, u, 0] ϵ + O[ϵ]^2
U[x, u, 0] + U^(0, 0, 1)[x, u, 0] ϵ + O[ϵ]^2

If we use the group properties of the transformation, especially the identity X(x, u, ϵ = 0) = x and U(x, u, ϵ = 0) = u, we can reduce the expression to

infinitesimalTrafo = MatrixForm[infinitesimalTrafo/. ({∂X(x, u, ϵ)/∂ ... 1013;ϕ(x, u), X(x, u, 0) x, U(x, u, 0) u}/.ϵ0)]

RowBox[{(, , GridBox[{{InterpretationBox[RowBox[{System`Convert`CommonDump`ConvertText ... {}]], SeriesData[õ, 0, {}, 0, 2, 1]]}], SeriesData[õ, 0, {u, Õ[x, u]}, 0, 2, 1]]}}], , )}]

The result shows that the transformation is represented by the original variables plus linear terms in ε. The coefficients of ε are called the infinitesimals of the transformation and are denoted by ξ and φ. Both expressions are the expansion coefficients of the Taylor expansion ξ = ∂_ϵX_ (ϵ = 0), and ϕ = ∂_ϵU_ (ϵ = 0), respectively.

The above is summarized in Lie's first theorem. The theorem considers the inverse problem of an infinitesimal representation. Since, in this case the infinitesimals ξ and φ are known, we are looking for a global transformation.

Theorem:  Lie´s first theorem

Let us assume that there exists a parameterization of a transformation such that the global transformation like rotation is equivalent to the solution of the initial value problem for the system of first order differential equations

dX/dϵ = ξ(X(ϵ), U(ϵ))

dU/dϵ = ϕ(X(ϵ), U(ϵ))

with the initial conditions

X(0) = x and U(0) = u.●

Lie's theorem states that the finite or global transformation represented by X and U is derived from an initial value problem. By integrating these equations with respect to ε and eliminating the group parameter, we end up with a global transformation.

As an example let us consider the transformations which are known as rotations:

Overscript[x, _] = x cos(ϵ) + y sin(ϵ)

Overscript[y, _] = y cos(ϵ) - x sin(ϵ)

The action of this global transformation is shown in the following figure.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_54.gif]

So far we have demonstrated that a global transformation is equivalent to its infinitesimal representation. However, the question of how symmetries and transformations are connected remains. Their link is the subject of the following subsection.

2.3 Group invariants

We know that a one-parameter group of transformations is represented by infinitesimal transformations. As the transformation is applied to a point of the manifold M it moves along a path in such a way that the path maps into itself. This basic concept of invariance needs a mathematical formulation. For the moment we restrict our considerations to curves in the plane.

Group invariants are the basic quantities of a symmetry analysis. A group invariant is defined as follows

Definition: Group invariant

A function F(x, u) is an invariant of the group of transformations, if

F(X(x,u,ε),U(x,u,ε)) = F(x,u)

is satisfied identically in x, u for all values of the group parameter ε.●

This definition contains the basic information for an invariant in symmetry analysis in a nut shell. The definition states that a function is an invariant if it has the same representation in the original and in the new coordinates. The definition also shows us a way how the invariant F may be calculated.

We know that any curve in the (x, u)-plane can be represented by a function F(x, u) = const. Following our definition this curve is considered to be invariant under a transformation if the curve remains the same in both coordinate systems. As already mentioned, global transformations need only an infinitesimal representation like

infini = {XFunction[{x, u, ϵ}, x + ϵ ξ(x, u)], UFunction[{x, u, ϵ}, u + ϵ ϕ(x, u)]}

{XFunction[{x, u, ϵ}, x + ϵ ξ[x, u]], UFunction[{x, u, ϵ}, u + ϵ ϕ[x, u]]}

where ξ and φ are the infinitesimals. The analytical expression for the invariance is simply stated by the relation

invar = F[x, u] == F[X[x, u, ϵ], U[x, u, ϵ]]

F[x, u] == F[X[x, u, ϵ], U[x, u, ϵ]]

showing that the form of the curve is the same in the old and new coordinates. Applying the transformation to this expression delivers the invariant condition in an implicit form

ivarI = invar/.infini

F[x, u] == F[x + ϵ ξ[x, u], u + ϵ ϕ[x, u]]

The explicit representation of the invariance condition is obtained by a Taylor expansion of the equation of invariance with respect to parameter ε.

Thread[Series[ivarI, {ϵ, 0, 1}], Equal]

F[x, u] == F[x, u] + (ϕ[x, u] F^(0, 1)[x, u] + ξ[x, u] F^(1, 0)[x, u]) ϵ + O[ϵ]^2

The invariance condition is satisfied once all terms containing the infinitesimal parameter ε have vanished. Introducing Lie's Symbol Overscript[v, ⇀] = ξ∂_x +ϕ∂_u, today known as tangent vector field or group generator, we can write the invariance condition as

F(x, u) = F(x, u) + Overscript[v, ⇀] F(x, u) ϵ + o(ϵ^2).

This calculation shows that a sufficient condition of invariance is Overscript[v, ⇀] F = 0. This condition is accessible by an algorithmic procedure and thus applicable to any function F.

Assuming that the tangent vector field applied to the curve F has to vanish, we can immediately conclude that all higher terms in  (12) containing multiple applications of the vector field Overscript[v , ⇀] also have to vanish. As a result, we have a sufficient condition of invariance. This condition is central to the determination of symmetries and can be generalized to a manifold of q independent variables  x = (x_1, x_2, …, x_q) and p dependent variables  u = (u^1, u^2, …, u^p). In the q × p-dimensional space, the invariance condition reads

Underoverscript[∑, i = 1, arg3] ξ_i (∂F (x, u))/∂x_i + Underoverscript[∑, α = 1, arg3] ϕ_α (∂F (x, u))/∂u^α = 0

where ξ_i and ϕ_α are the unknown infinitesimals of a point transformation.

The term of invariance boils down to the point  that each  one-parameter group of point transformations in the plane has one independent invariant which can be taken to be the left-hand side of any first integral J(x, u) = C of the characteristic equation (13) of the tangent vector field.

2.4 Symmetries of ordinary differential equations

So far we have introduced symmetry transformations of curves. The essential difference between curves and ODEs is that the manifold M is extended or prolonged by the derivatives. The manifold  of an ODE does not only contain the independent and dependent variable but also the derivative u ' = du/dx = p. Thus we have in the simplest case a three-dimensional space with coordinates (x, u, p). We are now looking for the symmetries of this kind of manifold. Again we can define the invariance of a differential equation in a similar way as for functions. Let us  state the definition for the most general case of an nth order ODE.

Definition: Symmetry of a differential equation

A group G of point transformations for variables in the manifold M is a symmetry group of a nth-order ordinary differential equation

Δ(x, u(x), du/dx, …, (d^nu)/dx^n) = 0

admitting G if the nth extended manifold ^n is invariant with respect to the nth prolongation _n of the group G.●

This definition contains the special case for first order ordinary differential equations.  A first-order differential equation

Δ(x, u(x), u ' (x)) = 0

admits a group G if the once extended manifold ^1, the surface in the space x, u, u ', is invariant with respect to the first prolongation _1 of G. This means that equation (15) is invariant under the coordinate transformations X = X(x, u, ϵ) and U = U(x, u, ϵ). The invariance condition can be formulated as

Δ(X, U, U ') = Δ(x, u, u ').

This kind of relation is also true for the general case of a nth order ODE. Consequently by knowing the transformations X and U, we also know the symmetries of the equation. Lie demonstrated that the knowledge of the infinitesimal transformations is sufficient to solve the equations. We will  state here only the fact that the knowledge of the infinitesimals is sufficient to solve the equations. We now show how the infinitesimals can be calculated with the help of MathLie Baumann [2000]. Once we know the symmetries, it is quite easy to derive the solutions  of the equations. The following sections will illustrate solutions for first-order, second-order, and higher-order equations.

3. First Order ODEs

The general form of a first order ordinary differential equation is given by


F(x, u(x), du/dx) = 0

where F is an arbitrary function combining  the independent and the dependent variables x, u, and the derivatives  du/dx in a general way. Our aim is to use point symmetries for the solution. Knowing the infinitesimal transformation of equation (17) allows two basic strategies for a solution. The first strategy is based on canonical variables extensively discussed by Ibragimov [1994,95,96]; Bluman and Kumei [1989]; as well as Olver [1986]. The second strategy originally invented by Lie in 1874 uses the well-known procedure of an integrating factor. A third procedure combines the two basic procedures and derives differential equations by an admitted group.

3.1 Integrating Factor Method

The common belief in literature is that the method involving an integrating factor is only useful in connection with a first order ordinary differential equation. In contrast Lie already developed a generalization of the procedure in the late 1870.

In the previous section, we demonstrated that a curve is invariant under a symmetry transformation if the tangent vector applied to the curve vanishes. Let us here assume that the curves discussed are integral curves or solutions of a first order differential equation. In this case, we can state the following theorem

Theorem: Integrating Factor

If the first order ordinary differential equation

du/dx = ω(x, u)

allows a one-parameter group given by the vector field Overscript[v, ⇀] = ξ ∂_x + ϕ∂_u , then and only then the following function

μ = 1/(ξ ω - ϕ)

is an integrating factor with ξω - ϕ≠0. If this is the case, the original equation is solved by a quadrature

∫ (ωdx - du) μ = const .

The theorem states that the knowledge of the infinitesimals ξ and φ allows the derivation of a solution just by integration. According to Lie, relation (20) can be symbolically simplified if we introduce two determinants of the form

dΔ = det (dx             du          )                 1              ω(x, u)

and

Δ = det(ξ         ϕ      ) = 1/μ               1              ω(x, u).

Then relation (20) reduces to the simple form

∫ (dΔ/Δ) = const .

Equation (23) combines all information necessary for our solution in a nut shell. All we need to know are the infinitesimals and the right hand side of  equation (18).  Integration over the manifold provides us with the solution.

As a first example, we discuss the Riccati equation used by Ibragimov [1994] to demonstrate the procedure of canonical transformations. The equation under consideration is

riccati = u(x)^2 + ∂u(x)/∂x - 2/x^2 == 0 ; LTF(riccati)

u^2 - 2/x^2 + u_x == 0

The infinitesimal symmetries of this equation are scaling symmetries for the coordinates given by ξ = x and ϕ = -u. These two expressions and the right hand side of the equation are the main components for deriving the Δ matrix of equation (21). With ,  we are able to create the Lie-matrix by using the independent and the dependent variable, the right hand side of the equation ω, and the infinitesimals ξ and φ.

Δmatrix = DeltaMatrix(x, u, 2/x^2 - u(x)^2, 1, ( {x}       {-u(x)} )) ; TableForm[Δmatrix]

x -u[x]
1 2/x^2 - u[x]^2

The above matrix is the basis for deriving a first integral via (23). The integral is the result of a series of integrations along a path in the (x, u)-plane. FirstIntegral[] of MathLie carries out the calculation

int1 = FirstIntegral(x, u, Δmatrix, 1) == c1

Log[x] + 1/3 Log[-2 + x u[x]] - 1/3 Log[1 + x u[x]] == c1

The input parameters for function FirstIntegral[] are the independent and dependent variables, the Lie matrix, and the number of the first integral. The number of first integrals is given by the order of the ODE. In case of first-order ODEs there exists only one integral. For second-order ODEs there exist two integrals and for nth-order ODEs there may exist n first integrals. If we solve the above expression with respect to the unknown function u, we end up with the solution of the problem.

sriccati = Solve[int1, u(x)]

{{u[x]  (^(3 c1) + 2 x^3)/(-^(3 c1) x + x^4)}}

Here c1 is an arbitrary constant of integration. The following section gives a generalization of the integrating factor method for a second order ODE.

3.2 Canonical Variables Method

Before we are dealing with the solution of first-order ODEs let us introduce two general notions of the theory of differential equations. These terms are the skeleton and the class of solution. We will show that both properties are instrumental in the geometric examination of an ODE and allow a straight forward simplification by means of canonical variables. Introducing these basic concepts it is easy to understand the geometrical meaning of a canonical transformation. Let us discuss these two properties in case of first order ODEs given in general by

F(x, u(x), u ') = 0

where u ' = du/dx = p is the first derivative.  In the following, we will discuss the two terms and give their definitions.

Definition: The Skeleton

The skeleton of a first order ordinary differential equation is defined as the surface

F(x, u, p) = 0

in the space of the three independent variables x, u, p. The corresponding differential equation is gained from the skeleton with the replacement p = u'.●

This definition can be generalized to higher-order ODEs if the corresponding derivatives are replaced by new variables independent of each other. The definition of the skeleton introduces nothing more than an extension of the manifold  = {x, u} by the variable p. This once extended manifold is very useful in the discussion of first order ordinary differential equations. Another term we will shortly define is the class of solutions.

Definition: The Class of Solutions

A class of solutions is a smooth solution which is a continuously differentiable function h(x) such that the curve u = h(x), u ' = dh(x)/dx belongs to the skeleton, that is, F(x, h(x), dh(x)/dx) = 0 identically in x for some interval.●

Thus the class of solutions is defined in accordance with certain natural mathematical assumptions or chosen from a physical significance of differential equations under discussion.

However, the crucial step in integrating differential equations is the simplification of the skeleton by a suitable change of the variables x and u. The change of variables is governed by the symmetries of the ODE. Meaning that the skeleton is invariant under the group transformation and its prolongations. The term invariance of an ODE is discussed in [7]. Provided that a symmetry group is known, a simplification of the skeleton is accomplished by introducing the canonical variables discussed above. Integration of the ODE in a canonical representation is simplified in such a way that the solution follows by a quadrature. To demonstrate the contents of these remarks let us examine three examples.

Example: 1 (Riccati equation)

As a first example let us inspect the Riccati equation discussed by Ibragimov [9]

riccati = u(x)^2 + ∂u(x)/∂x - 2/x^2 == 0

-2/x^2 + u[x]^2 + u^′[x] == 0

The skeleton of this Riccati equation is defined by the algebraic equation

u^2 + p - 2/x^2 == 0

p + u^2 - 2/x^2 == 0

The surface of this relation represents a so-called hyperbolic paraboloid. The three-dimensional visualization of this geometrical object is simplified if we define the function  

f(x_, u_) := 2/x^2 - u^2

The corresponding surface of the skeleton in three-dimensions  is shown in Fig. 1.3.1.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_148.gif]

Skeleton of the Riccati equation in the original (left) and the target  (right) coordinates.

Figure 1.3.1 shows that the skeleton in original coordinates x, u, and p has a singularity if x0. Also the figure clearly shows the parabolic shape of the surface for larger x-values. Thus the surface is bent in two directions which complicates the uncovering of the solution. Our goal is to find a transformation which reduces the  parabolic shape to a simpler and smoother representation. For the Riccati equation, a one-parameter symmetry group is provided by the following scaling transformations known as non-homogeneous dilations.

transformation = {xt ^(-a), uFunction[x, w(x ^a) ^a]}

{x^(-a) t, uFunction[x, w[x ^a] ^a]}

The invariance of the Riccati equation is checked by taking into account that the derivatives also need to be transformed by the rule

dtrafo = v_^(n_)(a_. x_) a^(-n) v^(n)(a x)

v_^(n_)[a_. x_] a^(-n) v^(n)[a x]

The rule dtrafo reflects the fact that the nth derivative of a function under a scaling is replaced by the nth derivative divided by the scaling factor a^n. The application of  transformation and dtrafo to the Riccati equation gives us

triccati = riccati/.transformation/.dtrafo ; LTF(triccati)

-(2 ^(2 a))/t^2 + ^(2 a) w^2 + ^(2 a) w_t == 0

We note that the original equation riccati is reproduced up to a common factor E^(2 a) containing the scaling parameter a. In the definition of the transformation in Mathematica it is essential that we use the target variables in the representation  for the original variables. x is simply replaced by the target r multiplied by the scaling factor E^(-a). The dependent variable u is replaced by w E^a. Since w depends on the target r, we have to represent r in Mathematica by x E^a. We also have to take into account that the derivative needs a special treatment which is realized in the rule dtrafo. Applying both rules to the original equation, we end up with the equation given in triccati. The discussed replacements of variables are actually the steps needed in a pencil and paper calculation. The application of  transformation to a first derivative shows that first order derivatives transform like

∂u(x)/∂x/.transformation/.dtrafo

^(2 a) w^′[t]

which in conventional notation reads  w '  u ' e^(-2a). Thus we observe that the equation's skeleton is invariant under the inhomogeneous scaling tx e^a, wue^(-a), w ' u ' e^(-2a) obtained by extending the transformations of the group to the first derivative u '. We also can check the invariance of the skeleton by applying the extended vector field to the skeleton.  The related vector field of the scaling transformation is defined by

Vect(f_) := -2 p ∂f/∂p - u ∂f/∂u + x ∂f/∂x

The skeleton of the Riccati equation is

skeleton = riccati/. {∂u(x)/∂xp, u(x) u}

p + u^2 - 2/x^2 == 0

The application of the  vector field to the skeleton reveals

Vect/@skeleton

-2 p - 2 u^2 + 4/x^2 == 0

the invariance. If we compare both expressions in more detail we recognize that the original and the transformed skeleton is reproduced up to a factor -2. Thus, if the skeleton vanishes the application of the vector field to the skeleton also vanishes. This result demonstrates that the skeleton of the Riccati equation is invariant under a scaling transformation.

Example: 2

Another example of a first order ordinary differential equation is the equation

example2 = (u(x)^3/x + u(x)^2/x^2) - ∂u(x)/∂x/x^2 == 0

u[x]^2/x^2 + u[x]^3/x - u^′[x]/x^2 == 0

This example is also invariant with respect to an inhomogeneous scaling transformation. We define this sort of transformation by a transformation rule like

scalingtrafo = {x^(-a) t, uFunction[x, w (x ^a) ^a]} ; ... ^(n_)(a_. x_) a^(-n) v^(n)(a x) ;) scaling(x_) := x/.scalingtrafo/.dtrafo

The application of the function scaling[] to the second example demonstrates again the invariance of the equation up to a common factor E^(4 a)

scaling(example2)

-(^(4 a) w)/t^2 + ^(4 a) w^2 + ^(4 a) t^2 w^3 == 0

The graphical representation of the  skeleton for this equation is given in Fig. 2.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_191.gif]

On the left the skeleton of the equation -(∂_xu(x))/x + u^2/x^2 + u^3/x = 0. The surface has some resemblance of a stingray missing the tail. The right part shows the canonical transformation of the left.

The structure of  this surface is simplified if we apply the following canonical transformation to the skeleton.

canonical = {x^t, uFunction[x, w(log(x))/x]} ; canonicaltransform(x_) := Simplify[PowerExpand[x/.canonical]]

The transformation is carried out by

canonicaltransform(Thread[example2 ^(4 t), Equal])

w[t] (1 + w[t] + w[t]^2) == w^′[t]

The related skeleton simplifies to a third order parabola translated along the x-axis, see Fig. 2 right. This example demonstrates that canonical coordinates dramatically simplify the shape of the skeleton and thus allows access to the solution in a simpler way.

Example: 3

Another example to demonstrate the concept of the skeleton is given by the first order equation

example3 = x ∂u(x)/∂x - u(x) + u(x)/x^(1/2) == 0

-u[x] + u[x]/x^(1/2) + x u^′[x] == 0

This type of equation is invariant with respect to the transformation Ξ = x/(1 - ϵ x) and Φ = y/(1 - ϵ x) the related canonical variables are w = u/x and t = -1/x. The equation in canonical variables reads

cexample3 = ∂w(t)/∂t + w(t)^(1/2) == 0

w[t]^(1/2) + w^′[t] == 0

The skeletons of the original and target equation is given in Fig. 3.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_204.gif]

Left: the skeleton of the equation ∂_xu(x) - u + (u/x)^(1/2) = 0  in the original representation. Right: the skeleton in canonical variable. The related equation reads ∂_tw(t) + w^(1/2) = 0.

It is obvious that the complicated skeleton of the original equation is reduced to a much simpler geometrical form. The convex shape of the original skeleton reduces to a concave surface invariant with respect to translations along the t-axis.

The three examples demonstrate that canonical transformations simplify the skeleton. Knowing this fact, we approach the second term, the solution. The question is how can we use the information of the symmetries to find solutions of any first order equation possessing symmetries. The idea is to use canonical variables allowing a transformation to a simpler representation of the equation.

Let us start our examination with the Riccati equation of the first example. We know that the canonical variables for a nonhomogeneous dilation is given by

trafo1 = {x^t, uFunction[x, w(log(x))/x]}

{x^t, uFunction[x, w[Log[x]]/x]}

This type of transformations was derived in section 2 by solving the defining equations for t and w. The application of this transformation to the Riccati equation gives us the following result

riccati1 = Simplify[Thread[PowerExpand[riccati/.trafo1] ^(2 t), Equal]]

w[t]^2 + w^′[t] == 2 + w[t]

Thus it straightens out the skeleton of the original Riccati equation, taking it to a parabolic cylinder, see Fig.1 right. This geometrical simplification is the reason why the Riccati equation takes on an integrable form when written in canonical variables.

Actually the scaling of the original Riccati equation is replaced by a group of translations  Overscript[t, _] = t + ϵ, w = u, and w ' = u '. The calculations so far executed by hand can be collected in a Mathematica function. The information we need for this kind of calculation are the original equation, and the dependent and independent variables. We also need the canonical transformations which are derived by our function CanonicalVariables[]. This function also has to know the names of the target coordinates.

We demonstrate the use of this function by applying it to the Riccati equation. We know that the Riccati equation is invariant with respect to an inhomogeneous scaling transformation. The related canonical variables are given by

ccoord = CanonicalVariables({u}, {x}, {x}, {-u}, {w}, {t})

{wFunction[{x, u}, u x], tFunction[{x, u}, Log[x]]}

Using this transformations, we can reduce the Riccati equation to

cricc = CanonicalRepresentation(riccati, u, x, {wx u, tlog(x)}, w, t)

-2 - w[t] + w[t]^2 + w^′[t] == 0

The resulting equation has the same structure as the equation derived by hand.

This equation can be solved by quadrature or by separating variables. We use here the Mathematica  function DSolve[] to integrate the equation.

sricc = DSolve[cricc, w, t]

{{wFunction[{t}, (2 ^(3 t) + ^(3 C[1]))/(^(3 t) - ^(3 C[1]))]}}

The solution in the original variables follows by applying the canonical transformation a second time

solricc = Simplify[w(x, u) == (w(t)/.sricc) 〚1〛/.tt(x, u)/.ccoord]

u x == (^(3 C[1]) + 2 x^3)/(-^(3 C[1]) + x^3)

Solving this implicit solution with respect to u delivers

Solve[solricc, u]

{{u (^(3 C[1]) + 2 x^3)/(x (-^(3 C[1]) + x^3))}}

which is just the solution found by Ibragimov. The second example discussed above allows the same scaling symmetries as the Riccati equation. The target equation follows by

cex2 = CanonicalRepresentation(example2, u, x, {wx u, tlog(x)}, w, t)

w[t] + w[t]^2 + w[t]^3 - w^′[t] == 0

The solution w follows by separation of variables and by integrating the left hand side and the right hand side

sex2 = Simplify[∫1/(w^3 + w^2 + w) w == c + ∫1t]

Log[w] == c + t + ArcTan[(1 + 2 w)/3^(1/2)]/3^(1/2) + 1/2 Log[1 + w + w^2]

Since the solution contains transcendental functions it is hard to get an explicit solution for the stingray equation. The inversion of the canonical transformation does not resolve this problem.

iex2 = sex2/. {ww(x, u), tt(x, u)}/.ccoord

Log[u x] == c + ArcTan[(1 + 2 u x)/3^(1/2)]/3^(1/2) + Log[x] + 1/2 Log[1 + u x + u^2 x^2]

However, we can solve this problem by a graphical representation of the solution in a two dimensional plot. We create a contour plot using the implicit  function for certain parameters

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_234.gif]

Contour plot of the solution for the stingray equation.

We observe that u increases if x increases. For small values of u there exist a linear relation between u and x.

The third example discussed  above is represented in canonical variables by

cex3 = CanonicalRepresentation(example3, u, x, {wu/x, t -1/x}, w, t)

w[t]^(1/2) + w^′[t] == 0

The solution of this equation for w follows by

sex3 = DSolve[cex3, w, t]

{{wFunction[{t}, 1/4 (t^2 - 2 t C[1] + C[1]^2)]}}

Inverting the transformation and solving with respect to u gives us the explicit solution

solex3 = u/x == (w(t)/.sex3) 〚1〛/.t -1/x

u/x == 1/4 (1/x^2 + (2 C[1])/x + C[1]^2)

sol3 = Solve[solex3, u]

{{u -1/4 x (-1/x^2 - (2 C[1])/x - C[1]^2)}}

So far we discussed three examples of first-order ODEs to show how the use of canonical variables can simplify the integration. The method of canonical variables is useful not only in the integration process of first-order ODEs but also in the integration of higher order ODEs.

3.3 Differential Equations Admitting a Group (Pattern Matching Method)

The following table contains some patterns for first order ODEs. These patterns follow by applying the methods of integrating factors and the method of canonical variables to a certain type of infinitesimals.

Infinitesimals First Order ODE
ξ = 1φ = 0 y ' = F (y)
ξ = 0φ = 1 y ' = F (x)
ξ = 0φ = G (y) y ' = H (x) G (y)
ξ = 0φ = F (y) G (x) y ' = (g (x) f (y) + θ (x))/f_y
ξ = 0φ = F (y) + G (x) y ' = (∫ (G_x (x))/(G (x) + F (y))^2y - θ (x)) (G (x) + F (y))
ξ = 0ϕ = 1/(F (y) + G (x)) y ' = (K (x) - yG ' (x))/(F (y) + G (x))
ξ = 0ϕ = ^(∫F (x) x) (y - H (x))^n y ' = (F (x) (H (x) - y))/(n - 1) + ((y - H (x))^nU (x) + H_x (x))
ξ = F (x) φ = 0 y ' = (K (y))/(F (x))
ξ = ^(∫g (x) x)/(f_y (y)) ϕ = 0 y ' = (g (x) f (y) + θ (x))/f_y
ξ = F (x) G (y) ϕ = 0 orξ = ^(∫g (y) y)/(f_x (x)) ϕ = 0 y ' = f_x/(g (y) f (x) + θ (y))
ξ = F (x) + G (y) φ = 0 y ' = 1/((∫ (G_y (y))/(G (y) + F (x))^2x - θ (y)) (G (y) + F (x)))
ξ = 1/(F (x) + G (y)) ϕ = 0 y ' = (K (x) - yG ' (x))/(F (y) + G (x))
ξ = x^2ϕ = xy y ' = (y + F (y/x))/x
ξ = xyϕ = y^2 y ' = y/(x + F (y/x))
ξ = F (x) φ = H (x) y ' = 1/(F (x)) (H (x) + K (y - ∫ (H (x))/(F (x)) x))
ξ = H (y) φ = F (y) y ' = (F (y))/(H (y) + K (x - ∫ (H (y))/(F (y)) y))
ξ = F (x) φ = G (y) y ' = (G (y))/(F (x)) K (∫1/(F (x)) x - ∫1/(G (y)) y)
ξ = x/yϕ = a y ' = (a^2x)/(ay + K (ax^2y^2))
ξ = aϕ = y/x y ' = (ax + K (ay^2 - x^2))/(a^2y)
ξ = y/xϕ = a y ' = (-a^2yK (F (y/x^a)))/(x (x^a - a K (F (y/x^a))))
ξ = aϕ = x/y y ' = y/ax - (y^(1 + a) K (x/y^a))/ax

The classification given by Ibragimov and Lie are special cases contained in the table above. In a future version of MathLie these procedure will be available.

4. Second Order ODEs

Second order ordinary differential equations are very important for applications in physics and engineering. All equations based on Newton's second law are second order equations. Thus mechanics, for example, is mainly based on second order ordinary differential equations. The integration theory of second-order ODEs was developed by Lie during the years 1871-1874 and published by Lie and Scheffers in 1891. Scheffers a student of Lie,  assembled all of Lie's work on ordinary differential equations in a beginners book. In this book Lie described the problem  of integrating a differential equation as follows:

I noticed that the majority of ordinary differential equations which were integrable by the old methods were left invariant under certain transformations, and that these integration methods consisted in using that property. Once I had thus represented many old integration methods from a common viewpoint, I set myself the natural problem: to develop a general theory of integration for all ordinary differential equations admitting finite or infinitesimal transformations ( ).

We will exemplify Lie's line of thought below. The most general form of a second-order ODE  is given by

FormBox[RowBox[{F(x, u, u ', u ' '),  , =,  , 0.}], TraditionalForm]

For our purposes, we assume that equation (26) is solvable with respect to the second-order derivative. Thus we consider equations in the form

u ' ' = ω(x, u, u ')

where ω is a given function of x, u, and u '.

According to Lie we can now state: If a second-order equation admits a finite symmetry of dimension r≥2, it can be integrated by a group theoretic quadrature method. This can be done in various ways, one of which is given by the following Algorithm:

1. compute the Lie algebra L_r. A basis of  L_r is the set  Overscript[v, ⇀] _1,  Overscript[v, ⇀] _2, …, Overscript[v, ⇀] _r. The tangent vector fields Overscript[v, ⇀] _i follow by appropriately specifying the group constants.

2. If r = 2, go to the next step;
     if r>2,  distinguish any two-dimensional subalgebra L_2 of L_r.
     If r = 1, the order of the equation may be lowered;
     if r = 0, the group method is not useful.

3. Calculate the Lie determinants dΔ_i and Δ and determine the two first integrals by integration. The Lie determinants are defined by

dΔ_1 = det(                           ),          ...                        1           u '         ω                   1         u '       ω

and

Δ = det (                          ')               ξ    ϕ    ϕ            ... 1;    ϕ                     2         2         2                  1         u '       ω

where ξ_i and ϕ_i are the infinitesimals corresponding to the vector field Overscript[v, ⇀] _i, and ϕ_i^' denote the first extensions of the infinitesimals ϕ_i. The first integrals ψ_i related to the Lie determinants dΔ_i are given by

FormBox[RowBox[{ψ_i = ∫dΔ_i/Δ = c_i, ,,          , i = 1, ,, 2.}], TraditionalForm]

The constants c_i denote the integration constants of the ODE.

4. Solve one of the two integrals with respect to the first order derivative and substitute the result into the remaining relation.

5. If we can solve the resulting relation with respect to the unknown function u, we end up with an explicit solution. Otherwise, we get the solution in an implicit form.

These five steps are completely algorithmic and are implemented among others in the package MathLie. The functions of MathLie carry out the needed calculations completely automatically. To show you how the solution procedure works, we demonstrate the application of the algorithm to an example.

Let us consider a nonlinear second order equation in which the nonlinearity is given by the square of the first derivative. This nonlinear term is multiplied by a real constant α. Our example equation is a generalization of the equations No. 6.1 and 6.2 by Kamke [1977] with α = 1 and α = 6, respectively. The generalized equation reads

firstExample = ∂^2u(x)/∂x^2 - α (∂u(x)/∂x)^2 ; LTF(firstExample)

-α u_x^2 + u_ (x, x) == 0

The infinitesimal symmetries are derived by the function Infinitesimals[] which is part of the package MathLie.    

infi = Infinitesimals(firstExample, u, x, {α}) ; LTF(infi)

ϕ_1 == (-^(-u α) k7 + ^(u α) (k6 + k1 x) + (k8 + k2 x) α)/α
ξ_1 == k4 + (k5 + ^(-u α) k7) x - (^(-u α) k3)/α - k2 x^2 α

The result is a symmetry group with 8 group parameters ki, i=1,2,…8. The group constants ki determine the symmetries of the equation. We note that this equation allows a transformation which reduces the original equation to the simple form FormBox[RowBox[{y ' ',  , =, 0.}], TraditionalForm] This reduction is always possible if a second order equation allows a symmetry group of order eight ( Lie [1891]).

To solve the nonlinear second order ODE, we select a solvable subgroup of order two from the eight-dimensional group. The last argument of the function DeltaMatrix[] contains the groups related to the solvable algebra. We choose  subgroups k8=β and k2=1.

inf1 = ( xi(1)  (x, u)  )/.infi/. {k10, k2> ... ;0, k70, k8β}/.uu(x)                   phi(1)  (x, u)

{{0}, {β}}

inf2 = ( xi(1)  (x, u)  )/.infi/. {k10, k2> ... 62754;0, k70, k80}/.uu(x)                   phi(1)  (x, u)

{{-x^2 α}, {x}}

Lie's matrix replacing the integrating factor is calculated by using the function DeltaMatrix[]

Δmatrix = DeltaMatrix(x, u, α (∂u(x)/∂x)^2, 2, {inf1, inf2}) ; TableForm[Δmatrix]

0 β 0
-x^2 α x 1 + 2 x α u^′[x]
1 u^′[x] α u^′[x]^2

The result is a 3×3 matrix containing the infinitesimals and the first prolongation of the two subgroups and the right hand side of the equation. The determinant of this matrix

Det[Δmatrix]

β + 2 x α β u^′[x] + x^2 α^2 β u^′[x]^2

is a non-vanishing expression containing a second order polynomial in the first derivatives. The coefficients of this polynomial depends on the independent variable x and the parameters α and β. Inserting the Lie matrix  into equation (26), we are able to calculate the first integrals of the nonlinear ODE

integrals = Thread[FirstIntegral(x, u, Δmatrix, {1, 2}) == {c1, c2}]

{-Log[1 + x α u^′[x]]/(α β) + u[x]/β == c1, 1/(x α) - 1/(x α (1 + x α u^′[x])) == c2}

The result contains two expressions for the integrals combining the variable u and its first order derivative in an algebraic way. These two first integrals define two manifolds in the space (x, u, u '). The projection of the intersection onto the (x, u)-plane of these two manifolds defines the solution we are looking for. The following figure represents a case with fixed values c1 and c2.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_336.gif]

Figure x shows the intersection of  the two integrals  for c1=c2=1. The parameters α and β  take the values α = 1/10 and β=1. The intersecting line represents the solution of the equation if we project the intersection to the (x, u)-plane.

Since we know two first integrals, we are able to eliminate the derivatives from the two integrals. Let us solve the first relation with respect to u '

sol1 = Simplify[Solve[integrals〚1〛, u^′(x)]]

{{u^′[x]  (-1 + ^(α (-c1 β + u[x])))/(x α)}}

Substituting this result into both expressions into integrals gives us

int21 = Simplify[PowerExpand[integrals/.sol1〚1, 1〛]]

{True, (1 - ^(α (c1 β - u[x])))/(x α) == c2}

The resulting list contains True and an implicit representation of the solution. The explicit solution follows from the second relation if we solve it with respect to u.

solution = Solve[int21〚2〛, u(x)]

{{u[x]  (c1 α β - Log[x (-c2 + 1/(x α)) α])/α}}

The constants c1 and c2 denote constants of integration. α and β are the model constant and the group parameter, respectively. The occurrence of the group parameter β as a multiplier  reminds us that the additive integration constant is related to the group of translations.

The example above demonstrates that the technique of integrating factors can be generalized to second order ODEs. The following section will discuss the extension of the integrating factor technique to a general nth order ODE. Lie called this procedure the method of generalized multipliers.

5. Higher Order ODEs

Differential equations of higher order arise naturally in physics. For example, third order ODEs come up in fluid dynamics and fourth order ODEs in elasticity. For general higher order equations there exist hardly any techniques for obtaining explicit symbolic solutions. Higher order ODEs are thus not well studied in the literature. In the following, we will present a symbolic technique for producing explicit solutions independent of the order of the equation. The method described in the preceding sections can, without essential changes, be applied to the solution of differential equations of higher order. The essential change in the method is the extension of the Lie matrix to higher prolongations. Extending the Lie matrix to higher prolongations is the key step for generalizing the procedure of integrating factors. Lie called this extension the determination of the multiplier of the differential equation.

The five steps of integration discussed for second order equations remain the same for higher order equations. However, the dimension of the Lie matrix changes from a 3 × 3 matrix to a (n + 1) × (n + 1) matrix. Before we describe the algorithm to solve higher order ODEs let us first state the general settings for higher order equations.

The most general form of a nth-order ODE is given by

FormBox[RowBox[{F(x, u, u ', u ' ', …, u^(n)),  , =,  , 0.}], TraditionalForm]

where u^(n) = d^nu/dx^n denotes the nth derivative of u with respect to x. In the following, we assume that equation (24) can be solved with respect to the nth derivative. The actual equation under consideration is thus

u^(n) = ω(x, u, u ', …, u^(n - 1))

where ω is a given function of x, u, u ', …, u^(n - 1).

If an nth-order equation admits a finite symmetry of dimension r≥n, then the equation can be integrated by the group theoretic quadrature method. This again can be done in various ways. One of the group theoretic algorithms is based on first integrals. The algorithm is summarized as follows:

1. compute the Lie algebra L_r. A basis for  L_r is the set  Overscript[v, ⇀] _1,  Overscript[v, ⇀] _2, …, Overscript[v, ⇀] _r. The tangent vector fields Overscript[v, ⇀] _i follow from appropriately specifying the group constants.

2. If r = n, go to the next step;
     if r>n,  distinguish any n-dimensional subalgebra L_n of L_r.
     If r = n - 1, the order of the equation may be lowered;
     if r = 0, the group method is not useful.

3. Calculate the Lie determinants dΔ_i and Δ and if  possible determine the n first integrals by integration. The Lie determinants are defined by

dΔ_1 = det(                                                                             ) ...  ⋮                    1               u '             u ' '           …         ω

and

dΔ_i = det(                                                                      (n - 1)) ...  ⋮                    1               u '             u ' '           ⋯         ω

where i denotes the row of the Lie matrix in which infinitesimals are replaced by differentials. The n × n Lie determinant  reads

Δ = det (                                      '                               (n - 1))   ...    ⋮                  1               u '             u ' '           …         ω

where ξ_i and ϕ_i are the infinitesimals of the vector field Overscript[v, ⇀] _i and ϕ_i^' denote the first and ϕ_i^(n - 1) the (n-1)th extension of the infinitesimals ϕ_i. The corresponding first integrals ψ_i are given by

ψ_i = ∫dΔ_i/Δ = c_i,          i = 1, 2, …, n .

The constants c_i denote the integration constants of the ODE.

4. Solve one of the n integrals with respect to the (n-1)th order derivative and substitute the result into the remaining relations. Repeat this procedure until no derivative remains in the relations.

5. If we can solve the resulting relation with respect to the unknown function u, we have found an explicit solution. Otherwise, our solution is implicit.

This procedure becomes very cumbersome with increasing orders of the equation if done by hand. In principle, the procedure can be  applied to any kind of linear or nonlinear ODE. How the algorithm works in particular is demonstrated below.

An example to test the algorithm from above is a third-order equation listed by Kamke as No. 7.13.

thirdOrderExample = ∂^2u(x)/(∂x ∂x) ∂^3u(x)/∂x^3 - α (β^2 (∂^2u(x)/(∂x ∂x))^2 + 1)^(1/2) ; LTF(thirdOrderExample)

-α (1 + β^2 u_ (x, x)^2)^(1/2) + u_ (x, x) u_ (x, x, x) == 0

The parameters α and β are real constants. According to Kamke this third-order ODE is solvable and the solution can only be represented in parametric form. We will show that an explicit solution of the equation is possible.

The question now is can we derive the necessary number of symmetries in order to integrate the equation? The derivation of the symmetries is the first step in the general algorithm. The calculation of symmetries is carried out by the MathLie function Infinitesimals[]

infi = Infinitesimals(thirdOrderExample, u, x, {α, β}) ; LTF(infi)

ϕ_1 == k2 + k3 x
ξ_1 == k1

The result is a symmetry group of order three. The specific symmetries are denoted by the group constants k1, k2, and k3. Each of these parameters is related to a vector field Overscript[v, ⇀] _i, i = 1, 2, 3. Since the number of vector fields is equal to the order of  equation, we can go to step three of the algorithm. In the third step we determine the Lie matrix by inserting the prolongations of the infinitesimals and the equation itself

inf1 = ( xi(1)  (x, u)  )/.infi/. {k11, k20, k30}/.uu(x)                   phi(1)  (x, u)

{{1}, {0}}

inf2 = ( xi(1)  (x, u)  )/.infi/. {k10, k21, k30}/.uu(x)                   phi(1)  (x, u)

{{0}, {1}}

inf3 = ( xi(1)  (x, u)  )/.infi/. {k10, k20, k31}/.uu(x)                   phi(1)  (x, u)

{{0}, {x}}

The Lie matrix is derived by

Δmatrix = DeltaMatrix(x, u, (α (β^2 (∂^2u(x)/(∂x ∂x))^2 ... )^(1/2))/∂^2u(x)/(∂x ∂x), 3, {inf1, inf2, inf3}) ; TableForm[Δmatrix]

1 0 0 0
0 1 0 0
0 x 1 0
1 u^′[x] u^′′[x] (α (1 + β^2 u^′′[x]^2)^(1/2))/u^′′[x]

One of the three first integrals is

integ1 = Simplify[FirstIntegral(x, u, Δmatrix, 1) == c1]

x == c1 + (1 + β^2 u^′′[x]^2)^(1/2)/(α β^2)

The result depends on u ' ' and now allows us to rewrite all terms containing u ' ' in the Lie matrices. Next we solve the first integral with respect to u ' '. Since the integral depends quadratically on u ' ', we get two solutions

sol1 = Solve[integ1, ∂^2u(x)/(∂x ∂x)]

{{u^′′[x]  -(-1/β^2 + c1^2 α^2 β^2 - 2 c1 x α^2 β ...  (-1/β^2 + c1^2 α^2 β^2 - 2 c1 x α^2 β^2 + x^2 α^2 β^2)^(1/2)}}

The first of the two solutions is used to replace u ' ' in the Lie matrix

dmat = Δmatrix/.sol1〚1〛

{{1, 0, 0, 0}, {0, 1, 0, 0}, {0, x, 1, 0}, {1, u^′[x], -(-1/β^2 + c1^2 α^2 	 ... /(-1/β^2 + c1^2 α^2 β^2 - 2 c1 x α^2 β^2 + x^2 α^2 β^2)^(1/2)}}

The simplified Lie matrix again is used to calculate the second first integral of the third-order equation.

integ2 = Simplify[FirstIntegral(x, u, dmat, 3) == c2]

((-1 + c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 β^4)^2 + α & ... + c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 β^4)/β^2^(1/2)) == c2

As expected we find the integral depending only on first derivatives of u. Since this integral is linear in u ', it is uniquely solvable in u '

sol2 = Simplify[Solve[integ2, u^′(x)]]

{{u^′[x] c2 - ((-1 + c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 β^4)/β^2)^(3/2)/(α ((c1 - x)^2 α^2 β^4)^(1/2))}}

The resulting expression contains radicals of quadratic polynomials in x. Inserting this result into the Lie matrix, we can eliminate the derevatives u '. We  find

dmat1 = Simplify[dmat/.sol2〚1, 1〛]

{{1, 0, 0, 0}, {0, 1, 0, 0}, {0, x, 1, 0}, {1, c2 - ((-1 + c1^2 α^2 β^4 - 2 c1 x 	 ... /(-1 + c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 β^4)/β^2^(1/2)}}

In a last step, we insert the Lie matrices into the third first integral depending on u and x

integ3 = Simplify[FirstIntegral(x, u, dmat1, 2) == c3]

-1/(α ((c1 - x)^2 α^2 β^4)^(1/2)) ((-1 + c1^2 α^2 β^4 - 2 c1 x α ...  c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 β^4)/β^2^(1/2))) == c3

The derived integral contains u linearly. In turn, we end up with a unique solution for the Kamke equation 7.13 which is

solution = Simplify[Solve[integ3, u(x)]]

{{u[x] c3 + (x ((-1 + c1^2 α^2 β^4 - 2 c1 x α^2 β^4 + x^2 α^2 ... ^2 α^2 β^4)/β^2^(1/2))))/(α^2 β^4 ((c1 - x)^2 α^2 β^4)^(3/2))}}

The solution depends on three constants c1, c2, and c3, all of which are constants of integration. The parameters α  and β are the parameters of  the original equation. To get a feeling how the solution evolves, we plot the solution for different parameter sets c1, c2, c3 at fixed α and β.

[Graphics:HTMLFiles/GBSymmetryAnalysisODES_432.gif]

Figure x represents the real valued solutions of the equation  ∂_ (x, x) u[x] ∂_ {x, 3} u[x] - α (1 + β^2 (∂_ (x, x) u[x])^2)^(1/2) = 0. The different curves represent the solutions for values of  c3 ∈ {0, 1, 2, 3, 4, 5} and fixed values for c2=1/2,c1=1. The parameters of the equation are α=1 and β=1/2.

We note  that the solution is explicitly represented by a complicated expression containing radicals and polynomials. This result is new since Kamke only offers a parametric representation of the solution. The example shows that we are able to arrive at a solution for higher order ODEs with Lie's procedure. Mathematica is not yet able to tread such types of equations.

6. References

W.F. Ames, Nonlinear Ordinary Differential Equations in Transport Processes, Academic Press, New York, 1968.

G. Baumann, MathLie, University of Ulm, 1997.

G.Baumann, Symmetry Analysis of Differential Equations with Mathematica, Springer, New York, 2000.

G.W. Bluman and S. Kumei, Symmetries and Differential Equations, Springer, New York, 1989.

N.H. Ibragimov, CRC Handbook of Lie Group Analysis of Differential Equations, Vol. 1-3, CRC-Press, Boca Raton, 1994, 1995, 1996.

N.H. Ibragimov, Sophus Lie and Harmony in Mathematical Physics, on the 150th Anniversary of His Birth, The Math. Intel. 16, 20-28, (1994).

E. Kamke, Differentialgleichungen, G.B. Teubner, Stuttgart, 1977.

S. Lie, Verh. Gesell. d. Wissenschaften zu Christiania, 1874

S. Lie, Theorie der Transformationsgruppen, Vol. I-II, Leipzig, (1888, 1890, 1893), reprinted by Chelsea Publ. Com., New York , (1970)

S. Lie und G. Scheffers, Vorlesungen über Differentialgleichungen mit bekannten infinitesimalen Transformationen, B.G. Teubner, Leipzig, 1891.

P.J. Olver, Applications of Lie groups to differential equations, Springer, Berlin, 1986.


Created by Mathematica  (September 15, 2003)