Chapter 1

Desolv Documentation

The Desolv package is implemented in the symbolic software, Maple. Desolv is simply a solver which attempts to find the solutions of a linear or nonlinear system of overdetermined differential equations of polynomial type. The implementation of Desolv is based on a heuristic procedure, which does not always promise to produce complete solutions of a system. If the heuristic procedure fails to simplify a system completely, it will transform the system of unsolved equations into a standard form. The heuristics is made up of seven modules which are the one-term module, the direct separation module, the solving unknown module, the indirect separation module, the integration module, the solving ODE module and the decoupling module; these have been described in Chapter Three.

The development of Desolv is based on the symbolic package LIE written by Alan Head [?] in MuMATH; LIE is a solver for systems of linear PDEs. Desolv has a simpler structure than LIE. Desolv partly adopts the methods of separation and integration in CRACK [?] (Thomas Wolf and Andre Brand) written in REDUCE and the method of decoupling in both standard form [?] (Gregory Reid and Allan Wittkopf) and diffgrob2 [?] (Elizabeth Mansfield), which are both written in Maple. Desolv is designed to handle large systems efficiently without running into memory problems.

Desolv includes other algorithms which help to analyse Lie point symmetries, non-classical symmetries and potential symmetries.

The Desolv package runs on Maple Version V, Release 3 or later releases. The descriptions of main functions in Desolv are given in this chapter. The

global parameters of Desolv which can be controlled by users are also described in details. Examples are used to demonstrate how to find point symmetries, non-classical symmetries and potential factors. A classification problem of point symmetries is illustrated in section (6.6). Section (6.7) shows that the function decouple() can work differently with and without changing the global parameters.

1.1 The Main Function pdesolv()

The function pdesolv() solves a system of overdetermined differential equations. It requires three arguments and outputs a list of four elements. The input and output of

pdesolv(leqns,lfn,lvar)

are as follows:

\underline{INPUT:}

- leqns is a list of differential expressions

- lfn is a list of unknown functions

- lvar is a list of variables

\underline{OUTPUT:} a list of 4 elements [reqns,linq,lsol,lunk]

- reqns is a list of differential expressions

- linq is a list of differential expressions

- lsol is a list of solutions of functions in lfn

- lunk is a list of unknowns occuring in reqns, linq and lsol

The equations in the list leqns are in linear homogeneous form. The functions and the constants in lfn are the unknowns in leqns. The variables in lvar are the variables on which the differential equations leqns depend. In the cases that pdesolv() does not solve the system leqns completely, a list of unsolved equations is returned as reqns. During any calculation, pdesolv() normally makes the assumptions that some expressions are non-zero. These expression are called inequations. If these inequations contain any unknown functions/constants, they are stored in the list linq. If the input list leqns is a linear system of differential equations, then linq is an empty list. The list lsol contains the solutions to the unknowns in lfn. The list lunk consists of unknown functions or constants appearing in reqns, linq and lsol. Any unknowns arising during the simplification process are denoted as F i for functions and C k for constants, where i and k are integers. The highest values of i and k are stored in the global variables desolv function and desolv constant respectively.

1.2 The Parameters in Desolv

Desolv has some parameters that allows the user to have control of the simplification process.

desolv module - The number of modules used in the simplification process. The function pdesolv() consists of seven modules which are in the following order: the solving one-term module, the direct separating module, the solving unknown module, the indirect separating module, the integrating module, the solving ordinary DE module and the decoupling module. By default, desolv module is set to be 7. For example, if desolv module is set to be 4, then only the first four modules are used to simplify a system of differential equations.

desolv modstat - A list of integers which tells how many times each module is successfully applied.

desolv Nterm - The number of terms in an equation. There are two different ways of scanning a system of differential equations through the modules. One way is to take the whole system and scan it through each module. This is called horizontal scan. The other way is take one equation at a time and scan it through each module. This is called vertical scan. pdesolv() splits a system into two lists, one list S contains equations having desolv Nterm terms or less and one list L contains equation having more than desolv Nterm terms. pdesolv() uses vertical scan for the list S and horizontal scan for the list L. The default value is 1.

desolv function - The number of new unknown functions which are introduced into the system during the call of pdesolv. desolv function is set to be 1 each time the function pdesolv() is called.

desolv constant - The number of new unknown constants which are introduced into the system during the call of pdesolv(). desolv constant is set to be 1 each time the function pdesolv() is called.

desolv reduce - The default value is true. When an integrability condition has been computed, it is immediately reduced with respect to the current system. If desolv reduce is set to be false, no reduction is done at this stage.

desolv sublimit - The default value is 10^{5}. During the calculation of the decoupling module, the coefficient of the highest derivative of an equation may get very large. If the length of the coefficient is greater than desolv sublimit, the coefficient is replaced by an arbitrary function which has the same dependencies as the coefficient. The arbitrary function replacement of the coefficient prevents the system from large expansion.

desolv efactor - The expansion factor. When the solution of an unknown function or constant is found in the simplification process, its solution is substituted into equations that contain the function. As a substitution is made to an equation, the function pdesolv() compares the length of the equation before and after the substitution. If the equation is expanded after the substitution by a factor greater than desolv efactor, then the substitution will not be done to the equation and delay until all the modules have unsuccessfully applied to the equation. desolv efactor is set to 1000 by default.

desolv nzassume - A list of expressions that are assumed to be non-zero during the simplification process.

desolv liassume - A list of lists of expressions that are assumed to be linearly independent of each other during the simplification process.

infolevel[pdesolv] - The setting of infolevel[pdesolv] causes the display of intermediate results during the calculation of the function pdesolv(). Higher values of infolevel[pdesolv] will cause more information to be displayed. Other options for infolevel are provided, such as, infolevel[gendef], infolevel[genvec], infolevel[comtab], infolevel[separation], infolevel[decouple].

1.3 Classical symmetries

As a prototypical example of how to find point symmetry with the help of the Desolv package, the linear heat equation

u_{xx} − u_{t} = 0  (1.1)

is considered. The defining equations of the point symmetry can be generated by the function gendef() which requires three input arguments, as follows:

gendef(deqns, dvar, ivar)

\underline{INPUT:}

- deqns is a list of differential equations

- dvar is a list of dependent variables

- ivar is a list of independent variables \underline{OUTPUT:} a list of 3 elements [leqns, lfn, lvar]

- leqns is a list of homogeneous differential equations

- lfn is a list of infinitesimal functions

- lvar is a list of variables

The differential equations in deqns usually have the highest derivatives on the left hand side (lhs). If the highest derivatives are not on the left hand sides of the equations, gendef() will automatically put the equations into the solved form for the highest derivatives. The system of defining equations leqns are in the homogeneous form. The system leqns has been reduced with respect to its one-term equations. The list lfn contains infinitesimal functions corresponding to the variables in lvar. The infinitesimal function of an independent variable x is denoted by xi[x] in maple-syntax whereas the infinitesimal function of a dependent variable u is denoted by eta[u]. In the MapleV Release 3 and later versions, the name xi[x] is displayed as ξ_{x}. The symmetry generator of the heat equation is expressed as

v = xi[x]\frac{∂}{∂x} + xi[t]\frac{∂}{∂t} + eta[u]\frac{∂}{∂u}  (1.2)

\underline{Example}: Generation of the defining equations for the linear heat equation:

> read desolv:
> de := [ diff(u(x,t),x,x) - diff(u(x,t),t) = 0 ];

de := [(\frac{∂^{2}}{∂x^{2}} u(x, t)) − (\frac{∂}{∂t}u(x, t)) = 0]

> le := gendef(de,[u],[x,t]);

le := [[\frac{∂}{∂u}ξ_{x}(x, t, u), \frac{∂}{∂u}ξ_{t}(x, t, u), \frac{∂}{∂x}ξ_{t}(x, t, u), \frac{∂^{2}}{∂u^{2}} η_{u}(x, t, u),

−2 (\frac{∂}{∂x}ξ_{x}(x, t, u)) + (\frac{∂}{∂t}ξ_{t}(x, t, u)) , − (\frac{∂}{∂t}η_{u}(x, t, u)) + (\frac{∂^{2}}{∂x^{2}} η_{u}(x, t, u)) ,

2 (\frac{∂^{2}}{∂u∂x}η_{u}(x, t, u)) + (\frac{∂}{∂t}ξ_{x}(x, t, u)) − (\frac{∂^{2}}{∂x^{2}} ξ_{x}(x, t, u)) ],

[ξ_{x}(x, t, u), ξ_{t}(x, t, u), η_{u}(x, t, u)] , [x, t, u] ]

This system of overdetermined linear PDEs can be simplified by the function pdesolv(). The function pdesolv() returns one unsolved equation. There is no assumption made on any unknown functions or constants in the solution. The solution contains one unknown function F 4(x, t) and six unknown constants, C 1 to C_{6}.

\underline{Example}: (Heat equation)

> sol := pdesolv(le[1],le[2],le[3]);

sol := [[ − (\frac{∂}{∂t}F 4(x, t)) + (\frac{∂^{2}}{∂x^{2}} F 4(x, t)) ], [ ],

[ξ_{x}(x, t, u) = \frac{1}{2}xC 4 + xtC 5 + \frac{1}{2}C 1 + \frac{1}{2}tC 2,

ξ_{t}(x, t, u) = C 3 + tC 4 + t^{2}C 5,

η_{u}(x, t, u) = F 4(x, t) − \frac{1}{4}uC 5x^{2} − \frac{1}{4}uC 2x − \frac{1}{2}utC 5 + \frac{1}{8}uC 6],

[F 4(x, t), C 1, C 2, C 3, C 4, C 5, C 6] ]

¿From the solutions of the infinitesimal functions, Lie vectors can be formed by using the function genvec(). The function genvec() has three arguments and returns a list of vectors, as follows:

genvec(lsol, lpar, lvar)

\underline{INPUT:}

- lsol is a list of solutions

- lpar is a list of parameters

- lvar is a list of variables \underline{OUTPUT:}

- a list of vectors.

The left hand side of a solution in lsol must be a function. Each function on the left hand side in lsol corresponds to a variable in lvar. The corresponding pair of function and variable must have the same position in lsol and lvar. For example,

leqns = [xi[x](⋯) = ⋯, xi[t](⋯) = ⋯, eta[u](⋯) = ⋯]
lvar = [x, t, u]

A vector is expressed in the form of a sum of a list of a coefficient and a list of variables. As an example,

x\frac{∂}{∂t} + sin(u)\frac{∂}{∂x} − t^{2} \frac{∂}{∂u} ⟶ [x, [t]] + [sin(u), [x]] + [−t^{2}, [u]]

> lv := genvec(sol[3],sol[4],le[3]);

lv := [ [F 4(x, t), [u]] , [1, [x]] , [\frac{1}{2}t, [x]] + [−\frac{1}{4}ux, [u]] ,

[1, [t]] , [1, [u]] , [\frac{1}{2}x, [x]] + [t, [t]] ,

[xt, [x]] + [t^{2}, [t]] + [−\frac{1}{4}ux^{2} − \frac{1}{2}ut, [u]] ]

Thus in our example the Lie-vectors of the heat equation are

\frac{∂}{∂x} , \frac{∂}{∂t} , \frac{∂}{∂u} ,
\frac{t}{2} \frac{∂}{∂x} − \frac{xu}{4} \frac{∂}{∂u}
\frac{x}{2} \frac{∂}{∂x} + t\frac{∂}{∂t} ,
xt\frac{∂}{∂x} + t^{2} \frac{∂}{∂t} − \frac{u}{4}(x^{2} + 2t)\frac{∂}{∂u} , F 4(x, t)\frac{∂}{∂u}

The function comtab() calculates the table of commutators of a list of vectors in the form that generates from genvec(). The function comtab() requires two arguments and returns a square matrix or a square matrix with a list of equations, as follows:

comtab(lvec,lvar)

\underline{INPUT:}

- lvec is a list of vectors

- lvar is a list of variables

\underline{OUTPUT:}

- a square matrix or a square matrix and a list of equations.

Each element of the output matrix is expressed linearly in terms of psi[i] which denotes the i^{th}-vector in lvec. The following Maple session shows the calculation of the commutator table of the vectors lv without the vector, [F 4(x, t), [u]], because it contains the unknown function F 4(x, t).

\underline{Example}: (Heat equation)

> comtab(subsop(1=NULL,lv),le[3]);

[ \begin{matrix}0 & −\frac{1}{4}ψ_{6} & 0 & \frac{1}{2}ψ_{1} & 2ψ_{2} & 0 \\ \frac{1}{4}ψ_{6} & 0 & −\frac{1}{2}ψ_{1} & −\frac{1}{2}ψ_{2} & 0 & 0 \\ 0 & \frac{1}{2}ψ_{1} & 0 & ψ_{3} & 2ψ_{4} − \frac{1}{2}ψ_{6} & 0 \\ −\frac{1}{2}ψ_{1} & \frac{1}{2}ψ_{2} & ψ_{3} & 0 & ψ_{5} & 0 \\ −2ψ_{2} & 0 & −2ψ_{4} + \frac{1}{2}ψ_{6} & −ψ_{5} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 \end{matrix} ]

When a commutator of two vectors cannot be expressed as a linear combination of psi[i], it is replaced by a variable, called gamma[k] (γ_{k}) in the commutator table where k is an integer. For this case, comtab() outputs is a matrix of commutator and a list of expressions of gamma[k].

\underline{Example}: (Heat equation)

> ct := comtab(subsop(1=NULL,lv),le[3]):
> print(ct[1]);

[ \begin{matrix}0 & γ_{1} & γ_{2} & γ_{3} & γ_{4} & γ_{5} & ψ_{1} \\ −γ_{1} & 0 & −\frac{1}{4}ψ_{6} & 0 & \frac{1}{2}ψ_{1} & 2ψ_{2} & 0 \\ −γ_{2} & \frac{1}{4}ψ_{6} & 0 & −\frac{1}{2}ψ_{1} & −\frac{1}{2}ψ_{2} & 0 & 0 \\ −γ_{3} & 0 & \frac{1}{2}ψ_{1} & 0 & ψ_{3} & 2ψ_{4} − \frac{1}{2}ψ_{6} & 0 \\ −γ_{4} & −\frac{1}{2}ψ_{1} & \frac{1}{2}ψ_{2} & ψ_{3} & 0 & ψ_{5} & 0 \\ −γ_{5} & −2ψ_{2} & 0 & −2ψ_{4} + \frac{1}{2}ψ_{6} & −ψ_{5} & 0 & 0 \\ −ψ_{1} & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix} ]

> ct[2];

[γ_{1} = [− (\frac{∂}{∂x}F 4(x, t)) , [u]] , γ_{2} = [−\frac{1}{4}F 4(x, t)x − \frac{1}{2}t (\frac{∂}{∂x}F 4(x, t)) , [u]] ,

γ_{3} = [− (\frac{∂}{∂t}F 4(x, t)) , [u]] , γ_{4} = [−\frac{1}{2}x (\frac{∂}{∂x}F 4(x, t)) − t (\frac{∂}{∂t}F 4(x, t)) , [u]] ,

γ_{5} = [−\frac{1}{4}F 4(x, t)x^{2} − \frac{1}{2}F 4(x, t)t − xt (\frac{∂}{∂x}F 4(x, t)) − t^{2} (\frac{∂}{∂t}F 4(x, t)) , [u]] ]

1.3.1 Arbitrary Functions of Derivatives

In classification problems, some equations may contain arbitrary functions of derivatives of dependent variables. When the defining equations are calculated, these arbitrary functions must be treated carefully. In most of the available symbolic packages, defining equation generators have problems in dealing with such equations, but it is not the case for the function gendef(). A demonstration of Desolv calculating the point symmetry of such an equation is shown. The quasilinear hyperbolic equation of the general form is

v_{tt} = f(x, v_{x}) v_{xx} + g(x, v_{x})  (1.3)

This equation admits a symmetry group of 3 parameters.

\frac{∂}{∂t} , \frac{∂}{∂v} , t\frac{∂}{∂v}

\underline{Example}: the PDE v_{tt} = f(x, v_{x}) v_{xx} + g(x, v_{x})

> V := v(x,t):
> de := [ diff(V,t,t) = f(x,diff(V,x))*diff(V,x,x) + g(x,diff(V,x)) ];

de := [(\frac{∂^{2}}{∂t^{2}} v(x, t)) = f (x, \frac{∂}{∂x}v(x, t)) (\frac{∂^{2}}{∂x^{2}} v(x, t)) + g (x, \frac{∂}{∂x}v(x, t))]

> le := gendef(de,[v],[x,t]);

le := [[ξ_{x}(x, t, v), \frac{∂}{∂t}ξ_{t}(x, t, v), \frac{∂}{∂v} η_{v}(x, t, v), \frac{∂}{∂x}η_{v}(x, t, v),

\frac{∂}{∂v} ξ_{t}(x, t, v), \frac{∂}{∂x}ξ_{t}(x, t, v), \frac{∂^{2}}{∂t^{2}} η_{v}(x, t, v)],

[ξ_{x}(x, t, v), ξ_{t}(x, t, v), η_{v}(x, t, v)] , [x, t, v]]

> sol := pdesolv(op(le));

sol := [[ ], [ ], [ξ_{x}(x, t, v) = 0, ξ_{t}(x, t, v) = C 1, η_{v}(x, t, v) = C 2 + tC 3] ,

[C 1, C 2, C 3]

> genvec(sol[3],sol[4],le[3]);

[[1, [t]], [1, [v]], [t, [v]]]

1.4 Nonclassical Symmetries

For the Burgers’ equation [?]

u_{t} = u_{xx} + 2uu_{x}  (1.4)

the infinitesimal generator is

v = xi[x]\frac{∂}{∂x} + xi[t]\frac{∂}{∂t} + eta[u]\frac{∂}{∂u}  (1.5)

The invariant surface condition (ISC) is

eta[u] = xi[t]u_{t} + xi[x]u_{x}  (1.6)

For the case 1 with xi[t] = 1, the ISC becomes

eta[u] = u_{t} + xi[x]u_{x}  (1.7)

and for the case 2 with xi[x] = 1, xi[t] = 0, the ISC is

eta[u] = u_{x}  (1.8)

The function gennc() generates the defining equations of the nonclassical symmetry. gennc() requires four arguments with the fourth argument to be optional and returns a list of elements, as follows:

gennc(deqns, dvar, ivar, lcase)

\underline{INPUT:}

- deqns is a list of differential equations

- dvar is a list of dependent variables

- ivar is a list of independent variables

- lcase is a list of case-number (Optional) \underline{OUTPUT:} a list of 3 elements [leqns,lfn,lvar]

- leqns is a list of lists of differential expressions

- lfn is a list of infinitesimal functions

- lvar is a list of variables

gennc() will generate all the cases for the nonclassical method when the argument lcase is omitted. The case-numbers in lcase represent the cases that the function gennc() will output. The following Maple session shows the calculation of defining equations of nonclassical symmetry for the Burgers’ equation with the case 1.

\underline{Example:} Burger’s equation

> read desolv:
> U := u(x,t):
> de := [ diff(U,t) = diff(U,x,x) + 2*U*diff(U,x) ];

de := [\frac{∂}{∂t}u(x, t) = \frac{∂^{2}}{∂x^{2}} u(x, t) + 2u(x, t)\frac{∂}{∂x}u(x, t)]

> le := gennc(de,[u],[x,t],[1]):
> le[1][1];

[\frac{∂^{2}}{∂u^{2}} ξ_{x}(x, t, u) , ξ_{t}(x, t, u) − 1, −2 (\frac{∂}{∂x}ξ_{x}(x, t, u)) η_{u}(x, t, u) + (\frac{∂^{2}}{∂x^{2}} η_{u}(x, t, u))

− (\frac{∂}{∂t}η_{u}(x, t, u)) + 2u (\frac{∂}{∂x}η_{u}(x, t, u)) 4u (\frac{∂}{∂u}ξ_{x}(x, t, u)) + (\frac{∂^{2}}{∂u^{2}} η_{u}(x, t, u))

+2ξ_{x}(x, t, u) (\frac{∂}{∂u}ξ_{x}(x, t, u)) − 2 (\frac{∂^{2}}{∂u∂x}xi_{x}(x, t, u)) , 2u (\frac{∂}{∂x}ξ_{x}(x, t, u))

+2 (\frac{∂^{2}}{∂u∂x}η_{u}(x, t, u)) − (\frac{∂^{2}}{∂x^{2}} ξ_{x}(x, t, u)) + 2ξ_{x}(x, t, u) (\frac{∂}{∂x}ξ_{x}(x, t, u))

−2η_{u}(x, t, u) (\frac{∂}{∂u}ξ_{x}(x, t, u)) + 2η_{u}(x, t, u) + (\frac{∂}{∂t}ξ_{x}(x, t, u)) ]

This system of defining equations is nonlinear in xi[x], xi[t] and eta[u]. When the function pdesolv() is used to simplify an system of nonlinear differential equations, it is better to exclude the decoupling module from the simplification process. The decoupling module very often runs out of memory space when dealing with nonlinear systems because of the complexity and expansion. The global parameter desolv module which is set to 6, will force the function pdesolv() to use the first 6 modules, without the decoupling module.

\underline{Example}: Burger’s equation

> desolv_module := 6:
> sol := pdesolv(le[1][1],le[2],le[3]):

This Maple session yields three sets of solutions for the nonclassical symmetry of the Burgers’ equation with the case 1. The first solution contains 2 unsolved equations.

4F 5(t)^{3} + (\frac{∂^{2}}{∂t^{2}} F 5(t)) + 6F 5(t) (\frac{∂}{∂t}F 5(t)) = 0

− (\frac{∂^{2}}{∂t^{2}} F 6(t))−4F 5(t)^{2}F 6(t)−4F 5(t) (\frac{∂}{∂t}F 6(t))−2F 6(t) (\frac{∂}{∂t}F 5(t)) = 0

ξ_{t}(x, t, u) = 1
ξ_{x}(x, t, u) = F 5(t)x − F 6(t)
η_{u}(x, t, u) = −\frac{1}{2}x (\frac{∂}{∂t}F 5(t)) + \frac{1}{2} (\frac{∂}{∂t}F 6(t))
−F 5(t)^{2}x + F 5(t)F 6(t) − uF 5(t)

subject to F 5(t)x − F 6(t) ≠ 0.

The second solution is

ξ_{t}(x, t, u) = 1
ξ_{x}(x, t, u) = −2u
η_{u}(x, t, u) = 0

The third solution contains three unsolved equations subject to one condition, F 1(x, t) ≠ 0.

− (\frac{∂^{2}}{∂x^{2}} F 3(x, t)) + (\frac{∂}{∂t}F 3(x, t)) + 2 (\frac{∂}{∂x}F 1(x, t)) F 3(x, t) = 0

(\frac{∂}{∂t}F 4(x, t)) − (\frac{∂^{2}}{∂x^{2}} F 4(x, t)) − 2 (\frac{∂}{∂x}F 3(x, t)) + 2 (\frac{∂}{∂x}F 1(x, t)) F 4(x, t) = 0 (\frac{∂}{∂t}F 1(x, t)) − 2 (\frac{∂}{∂x}F 4(x, t)) − (\frac{∂^{2}}{∂x^{2}} F 1(x, t)) + 2F 1(x, t) (\frac{∂}{∂x}F 1(x, t)) = 0

ξ_{t}(x, t, u) = 1
ξ_{x}(x, t, u) = F 1(x, t) + u
η_{u}(x, t, u) = −u^{3} − F 1(x, t)u^{2} − F 3(x, t) − uF 4(x, t)

1.5 Potential symmetries

An important calculation in the case of potential symmetry is to find potential factors. Potential factors are needed to put an equation into a potential form, from which a potential system can be derived. In this section, an example of calculating potential factors using Desolv is demonstrated.

\underline{Example}: Consider the nonlinear diffusion equation S_{0}\{u\} [?] given by

u_{t} = (K(u)u_{x})_{x}  (1.9)

The Fréchet derivative of S_{0} can be computed by the function frechet():

frechet(leqns, dvar, ivar)

\underline{INPUT:}

- leqns is a list of differential equations

- dvar is a list of dependent variables

- ivar is a list of independent variables \underline{OUTPUT:}

- a square matrix

The output is a matrix of vectors which are expressed in the form of a sum of lists.

> de := [ diff(u(x,t),t) = diff( K(u(x,t))*diff(u(x,t),x),x ) ];

de := [\frac{∂}{∂t}u(x, t) = D(K)(u(x, t)) (\frac{∂}{∂x}u(x, t))^{2} K(u(x, t)) (\frac{∂^{2}}{∂x^{2}} u(x, t))]

> F := frechet(de,[u],[x,t]);

F := [[−K(u(x, t)), [x, x]] + [−2D(K)(u(x, t)) (\frac{∂}{∂x}u(x, t)) , [x]] + [1, [t]]

+[ − D^{(2)}(K)(u(x, t)) (\frac{∂}{∂x}u(x, t))^{2} − D(K)(u(x, t)) (\frac{∂^{2}}{∂x^{2}} u(x, t)) , 1]]

The adjoint of the Fréchet derivative F can be calculated by calling the function adjointD():

adjointD(M, lvar)

\underline{INPUT:}

- M is a square matrix of vectors

- lvar is a list of variables

\underline{OUTPUT:}

- a square matrix

> adjointD(F,[x,t]);

[[−K(u(x, t)), [x, x]] + [−1, [t]]]

The function genfac() generates the system of defining equations for potential factors of a system of differential equations. The function genfac() requires three arguments, a list of equations, a list of dependent variables and a list of independent variables:

genfac(deqns, dvar, ivar)

\underline{INPUT:}

- deqns is a list of differential equations

- dvar is a list of dependent variables

- ivar is a list of independent variables \underline{OUTPUT:} a list of 3 elements [leqns,lfn,lvar]

- leqns is a list of differential expressions

- lfn is a list of potential factor-functions

- lvar is a list of variables

The potential factor-functions are denoted as lambda[i] by default where i are integers.

> le := genfac(de,[u],[x,t]);

le := [[\frac{∂}{∂u}λ_{1}(x, t, u), −K(u) (\frac{∂^{2}}{∂x^{2}} λ_{1}(x, t, u)) − (\frac{∂}{∂t}λ_{1}(x, t, u))] , [λ_{1}(x, t, u)] , [x, t, u]]

> sol := pdesolv(op(le));

sol := [[ ], [ ], [λ_{1}(x, t, u) = C 1 + xC 2] , [C 1, C 2]]

Therefore the potential factors of the system S_{0} are {1, x}. The potential factor λ_{1} = x yields the potential system S_{1}{u,v},

v_{x} = xu  (1.10)
v_{t} = x(L(u))_{x} − L(u)  (1.11)

where K(u) = L′(u). By repeating the above steps, we find that the solutions for the potential factors of S_{1}{u,v} are (λ_{1}, λ_{2}) = (0, x^{−2}).

> U := u(x,t): V := v(x,t):
> de := [ diff(V,x) = x*U, diff(V,t) = x*diff( L(U),x) - L(U) ];

de := [\frac{∂}{∂x}v(x, t) = xu(x, t), \frac{∂}{∂t}v(x, t) = xD(L)(u(x, t)) (\frac{∂}{∂x}u(x, t)) − L(u(x, t))]

> F := frechet(de,u,v],[x,t]);

F := [[−x, 1], [1, [x]]]

[[−xD(L)(u(x, t)), [x]] + [−xD^{(2)}(L)(u(x, t)) (\frac{∂}{∂x}u(x, t)) + D(L)(u(x, t)), 1] ,

[1, [t]]]

> adjointD(F,[x,t]);

[ \begin{matrix}[−x, 1] & [xD(L)(u(x, t)), [x]] + [2D(L)(u(x, t)), 1] \\ [−1, [x]] & [−1, [t]] \end{matrix} ]

> le := genfac(de,[u,v],[x,t]);

le := [[\frac{∂}{∂u}λ_{2}(x, t, u, v), − (\frac{∂}{∂u}λ_{1}(x, t, u, v)) − (\frac{∂}{∂v} λ_{2}(x, t, u, v)) x (\frac{∂}{∂u}L(u)) ,

−xλ_{1}(x, t, u, v) + x (\frac{∂}{∂u}L(u)) (\frac{∂}{∂x}λ_{2}(x, t, u, v)) + x^{2} (\frac{∂}{∂u}L(u))

(\frac{∂}{∂v} λ_{2}(x, t, u, v)) u + 2 (\frac{∂}{∂u}L(u)) λ_{2}(x, t, u, v),

− (\frac{∂}{∂x}λ_{1}(x, t, u, v)) − (\frac{∂}{∂v} λ_{1}(x, t, u, v)) xu − (\frac{∂}{∂t}λ_{2}(x, t, u, v))

+ (\frac{∂}{∂v} λ_{2}(x, t, u, v)) L(u)], [λ_{1}(x, t, u, v), λ_{2}(x, t, u, v)] , [x, t, u, v]]

> sol := pdesolv(op(le));

sol := [[ ], [ ], [λ_{1}(x, t, u) = 0, λ_{2}(x, t, u) = \frac{C 1}{x^{2}}] , [C 1]]

1.6 Classification

\underline{Example}: The nonhomogenous Monge-Ampère equation contains an arbitrary function a(x, y),

u_{xx}u_{yy} − u_{xy}^{2} + a^{2} = 0  (1.12)

Different values of the arbitrary function give different symmetry groups (see, e.g., the section 9.4 in the handbook of Lie group analysis of differential equations Vol.1, N.H.Ibragimov [?]). The case where a = x^{β}y^{γ} is examined:

> U := u(x,y):
> de := [diff(U,x,x)*diff(U,y,y)-diff(U,x,y)^2+(x^beta*y^gamma)^2=0];

de := [(\frac{∂^{2}}{∂x^{2}} u(x, y)) (\frac{∂^{2}}{∂y^{2}} u(x, y)) − (\frac{∂^{2}}{∂y∂x}u(x, y))^{2} + (x^{β})^{2} (y^{γ})^{2} = 0]

> le := gendef(de,[u],[x,y]);

le := [[\frac{∂}{∂u}ξ_{y}(x, y, u), \frac{∂}{∂u}ξ_{x}(x, y, u), \frac{∂^{2}}{∂x^{2}} ξ_{y}(x, y, u), \frac{∂^{2}}{∂y∂x}η_{u}(x, y, u),

\frac{∂^{2}}{∂y^{2}} η_{u}(x, y, u), \frac{∂^{2}}{∂u^{2}} η_{u}(x, y, u), \frac{∂^{2}}{∂x^{2}} η_{u}(x, y, u), \frac{∂^{2}}{∂y^{2}} ξ_{x}(x, y, u),

2 (\frac{∂^{2}}{∂u∂y} η_{u}(x, y, u)) − (\frac{∂^{2}}{∂y^{2}} ξ_{y}(x, y, u)) , − (\frac{∂^{2}}{∂u∂x}η_{u}(x, y, u))

+ (\frac{∂^{2}}{∂y∂x}ξ_{y}(x, y, u)) , (\frac{∂^{2}}{∂u∂y} η_{u}(x, y, u)) + (\frac{∂^{2}}{∂y∂x}ξ_{x}(x, y, u)) ,

2 (\frac{∂^{2}}{∂u∂x}η_{u}(x, y, u)) − (\frac{∂^{2}}{∂x^{2}} ξ_{x}(x, y, u)) , − (\frac{∂}{∂y} ξ_{y}(x, y, u)) yx

− (\frac{∂}{∂x}ξ_{x}(x, y, u)) yx + (\frac{∂}{∂u}η_{u}(x, y, u)) yx − ξ_{y}(x, y, u)γx − ξ_{x}(x, y, u)βy],

[ξ_{x}(x, y, u), ξ_{y}(x, y, u), η_{u}(x, y, u)] , [x, y, u]]

> sol := pdesolv(op(le));

sol := [[ ], [ ], [ξ_{x}(x, y, u) = xC 10,

ξ_{y}(x, y, u) = −\frac{y(C 14 + C 10 + C 10β)}{γ + 1},

η_{u}(x, y, u) = −C 13 − uC 14 + yC 1 + xC 3],

[C 1, C 3, C 10, C 13, C 14]]

> lv := genvec(sol[3],sol[4],le[3]);

lv := [[y, [u]], [x, [u]], [x + xγ, [x]] + [−y − yβ, [y]], [1, [u]], [x, [x]] + [u + uβ, [u]]] The equation admits 5 symmetries,

y \frac{∂}{∂u} , x\frac{∂}{∂u} , \frac{∂}{∂u} ,

(1 + γ)x\frac{∂}{∂x} − (1 + β)y \frac{∂}{∂y} , x\frac{∂}{∂x} + (1 + β)u\frac{∂}{∂u}

During the simplification process, the function pdesolv() assumes some expressions to be non-zero such as denominators and coefficients of the leading derivatives. When pdesolv() applies the method of direct separation to an equation, the expressions which depend on the separating variable are assumed to be linearly independent of each other. The function classify() displays non-zero assumptions and linearly independent assumptions. The arguments of the function classify(), which are optional, are names of arbitrary functions or constants.

> classify(gamma,beta);

non-zero assumptions γ β

γ + 1
β + 1
β + 2 + γ

linearly independent assumptions

If these expressions of β and γ are equal to zero, the equation may admit different sets of symmetries. These expressions are stored in the global variable, desolv nzassume and desolv liassume. For γ = β = −1, the equation becomes

u_{xx}u_{yy} − u_{xy}^{2} + x^{−2}y^{−2} = 0  (1.13)

This equation admits 7 symmetries instead of 5 in the general case.

> U := u(x,y):
> de := [ diff(U,x,x)*diff(U,y,y) - diff(U,x,y)^2 + 1/(x*y)^2 = 0];
> le := gendef(de,[u],[x,y]);
> sol := pdesolv(op(le));

sol := [[ ], [ ], [ξ_{x}(x, y, u) = −x^{2}C 12 − xC 17 + xyC 2,
ξ_{y}(x, y, u) = y^{2}C 2 + yC 15 − xyC 12,
η_{u}(x, y, u) = C 5 + yC 9 + yuC 2 + xC 7 − xuC 12],
[C 2, C 5, C 7, C 9, C 12, C 15, C 17]]

> lv := genvec(sol[3],sol[4],le[3]);

lv := [[x, [u]], [1, [u]], [xy, [x]] + [y^{2}, y] + [yu, [u]], [y, [u]], [y, [y]],

[x^{2}, [x]] + [xy, [y]] + [xu, [u]], [x, [u]]]

1.7 Decoupling

A system of differential equations which contain arbitrary functions (or classified functions) often cause difficulties for simplification. A system of nonlinear telegraph equations [?] which contains two arbitrary functions f and g is

v_{t} = u_{x}  (1.14)
v_{x} = f(u)u_{t} + g(u)  (1.15)

The 8 defining equations for point symmetries are

−\frac{∂ξ_{x}}{∂u} + \frac{∂ξ_{t}}{∂v} = 0

−\frac{∂ξ_{t}}{∂u} + f \frac{∂ξ_{x}}{∂v} = 0

\frac{∂η_{u}}{∂u} + \frac{∂ξ_{t}}{∂t} − \frac{∂η_{v}}{∂v} + \frac{∂ξ_{x}}{∂x} = 0

\frac{∂η_{u}}{∂x} − \frac{∂η_{v}}{∂t} + g \frac{∂ξ_{x}}{∂t} + g \frac{∂η_{u}}{∂u} = 0

−\frac{∂ξ_{t}}{∂x} − g \frac{∂ξ_{t}}{∂v} − \frac{∂η_{v}}{∂u} + f \frac{∂ξ_{x}}{∂t} + g \frac{∂ξ_{x}}{∂u} + f \frac{∂η_{u}}{∂v} = 0

−g \frac{∂ξ_{x}}{∂u} − g \frac{∂ξ_{t}}{∂v} − f \frac{∂η_{u}}{∂v} + \frac{∂η_{v}}{∂u} + f \frac{∂ξ_{x}}{∂t} − \frac{∂ξ_{t}}{∂x} = 0

g \frac{∂η_{v}}{∂v} − g^{2} \frac{∂ξ_{x}}{∂v} + \frac{∂η_{v}}{∂x} − g′η_{v} − f \frac{∂η_{u}}{∂t} − g \frac{∂ξ_{x}}{∂x} = 0

None of these equations can be simplified by the first six modules of the function pdesolv(). Thus the defining system needs to be transformed into the standard form by the decoupling module. The function decouple() takes in three arguments, a list of equations, a list of unknowns and a list of variables, as follows,

decouple(leqns, lfn, lvar)

\underline{INPUT:}

- leqns is a list of differential equations

- lfn is a list of unknown functions

- lvar is a list of variables

\underline{OUTPUT:}

- a list of differential equations

As the function decouple() simplifies the above system, the calculation does not get very far before the computer is hung up due to a great expansion of the system. The parameter infolevel[decouple] can be set to be 3 to display intermediate computation of the function decouple().

> U := u(x,t): V := v(x,t):
> de := [ diff(V,t) = diff(U,x), diff(V,x) = f(U)*diff(U,t) + g(U)]:
> le := gendef(de,[u,v],[x,t]):
> infolevel[decouple] := 3:
> s  := decouple(le[1],le[2],le[3]):

The large expansion is caused by equations whose coefficients of highest derivatives depend on arbitrary functions f(u) and g(u). When these equations reduce another equation, they may be differentiated with respect to u a number of times. One of the techniques used by the function decouple() to overcome the expansion problem is to replace large coefficients of the highest derivatives with arbitrary functions, namely AAA k where k are integers. The functions AAA k have the same dependencies as the coefficients. A replacement is done depending on how large the coefficient is. The global variable desolv sublimit allows the user to control such replacements. If the length of a coefficient is larger than the value of desolv sublimit, a replacement will take place. By default, desolv sublimit is set to be at 10^{5}.

Another technique of avoiding the expansion is to delay the reduction of integrability conditions. When an integrability condition is computed, immediately it is reduced with respect to the current system. The reduction can be delayed until all integrability conditions are calculated by setting the parameter desolv reduce to be false. By default, the parameter desolv reduce is set to be true. The following Maple session shows how the function decouple() simplifies the above defining equations by setting desolv reduce = false and desolv sublimit = 100.

\underline{Example}:

> read desolv:
> de := [ diff(v(x,t),t) = diff(u(x,t),x),
diff(v(x,t),x) = f(u(x,t))*diff(u(x,t),t) + g(u(x,t)) ]:
> le := gendef(de,[u,v],[x,t]):
> desolv_reduce := false:
> desolv_sublimit := 100:
> s1 := decouple(le[1],le[2],le[3]);

s1 := [η_{u}(x, t, u, v), \frac{∂}{∂v} ξ_{x}(x, t, u, v), \frac{∂}{∂v} ξ_{t}(x, t, u, v), \frac{∂}{∂v} η_{v}(x, t, u, v),

\frac{∂}{∂u}ξ_{x}(x, t, u, v), \frac{∂}{∂u}ξ_{t}(x, t, u, v), \frac{∂}{∂t}ξ_{x}(x, t, u, v), \frac{∂}{∂u}η_{v}(x, t, u, v),

\frac{∂}{∂t}ξ_{t}(x, t, u, v), \frac{∂}{∂t}η_{v}(x, t, u, v), \frac{∂}{∂x}ξ_{t}(x, t, u, v),

\frac{∂}{∂x}ξ_{t}(x, t, u, v), \frac{∂}{∂x}η_{v}(x, t, u, v)]

> sol := pdesolv(s1,le[2],le[3]);

sol := [[ ], [ ], [ξ_{x}(x, t, u, v) = C 3,
ξ_{t}(x, t, u, v) = C 1, η_{u}(x, t, u, v) = 0, η_{v}(x, t, u, v) = C 2],
[C 1, C 2, C 3]]

bytes used=82878192, alloc=2489912, time=547.67

This Maple session took 547 seconds of CPU time to complete. The intermediated calculation of the function decouple() can be displayed by setting the variable infolevel[decouple]. The higher the value of infolevel[decouple], the more calculation is displayed.