13:57, 8 November 2005

On using materially-conserved quantities to construct solutions of differential equations

Simon Hood

Department of Mathematical Sciences, University of Liverpool, Liverpool, L69 3BX, United Kingdom. Tel.: 0151 794 4043; Fax.: 0151 794 4061

                           email: simonh@liv.ac.uk

Abstract

In this article we consider how one may systematically determine materially-conserved quantities (mcqs) from the equations describing the dynamics of a fluid and then how these can be used to construct quite general, exact, analytical solutions of these equations. Such solutions necessarily reflect underlying physical processes and are often general enough to satisfy strong boundary conditions. The method described for determining mcqs is systematic and essentially algorithmic, and is therefore a good candidate for implementation using a computer-algebra system.

We use this method to recover three mcqs of the equations describing large-scale fluid flow one the surface of a rotating sphere (geostrophic flow). We then prove that these geostophic equations admit no further mcqs. Next, we describe how to construct fairly general, exact solutions from these quantities.

Finally we discuss the application of these ideas to general systems of fluid-dynamical equations, in particular we consider under what conditions there exist mcqs which hold for a significant number of solutions of a system of governing equations and when these can be found without a priori determining these solutions.

1. Introduction

It is well known that the dynamics defined by the equations which describe large-scale flow of a thin layer of fluid on the surface of a rotating sphere (e.g. the world’s oceans), the geostrophic equations, admit certain materially conserved quantities (mcqs). These quantities have been used in various ad hoc ways to construct solutions of the equations. A systematic approach is needed. Moreover, this approach should be applicable to other fluid systems. In this paper we discuss how one can systematically determine the mcqs of a wide class of fluid-dynamical systems and how they may be used to construct surprisingly general analytical solutions.

Given n mcqs of a system of pdes in n independent variables the mcqs are necessarily functionally-related. These relations take the form of differential equations which are significantly simpler than the governing equations. Solution of these equations yields classes of solution ansatz which necessarily reflect fundamental underlying physics and which satisfy strong boundary conditions — stronger than those which may be satisfied by solutions

                                                                  src : publication_quality.tex

On using mcqs to construct solutions of differential equations

obtained from point symmetries (e.g., Bluman & Kumei, 1989; Olver, 1986; Stephani, 1990) for example (see Section 5 for further discussion on this).

We show how one can search, in a systematic, essentially algorithmic way, for such materially conserved quantities directly from the governing equations of a large class of fluiddynamical systems. The method is a good candidate for implementation using a computeralgebra system such as maple or macsyma. We apply this method to the geostrophic equations. We then discuss the number of mcqs a given system can have. We show that the geostrophic equations have only three mcqs (to the author’s knowledge no such proof has been previously given) and that in general a system of equations in n independent variables has n mcqs, though it is not obvious, in general, how to find these systematically without first solving the governing equations; we consider under what circumstances it is possible to do this.

The class of solutions which one finds from the mcqs of a system depends on the function relating them; this function is arbitrary. Hence the question of which such functions lead to physical solutions thus arises. We finally describe how one can do this by introducing additional physical processes into the system.

2. Background

The equations most often used to describe large-scale motion of a thin layer of fluid on the surface of a rotating sphere are, in non-dimensional form,

u = −\frac{p_{y}}{y} ,  v = \frac{p_{x}}{y} ,  ρ = p_{z},  (2.1a)
u_{x} + v_{y} + w_{z} = 0  (2.1b)

and

uρ_{x} + vρ_{y} + wρ_{z} = 0  (2.1c)

(three momentum equations, continuity and a thermodynamic equation, respectively), where u, v and w are the components of the velocity field in the x, y and z directions, ρ is density, and p is pressure. They may be written

M_{x}M_{zzz} + y(M_{xz}M_{yzz} − M_{yz}M_{xzz}) = 0,  (2.2i)
u = −\frac{M_{yz}}{y} ,  v = \frac{M_{xz}}{y} ,  w = \frac{M_{x}}{y^{2}} ,  p = M_{z},  ρ = M_{zz}  (2.2ii)

(Welander, 1959). These equations are well known to admit the three materially-conserved quantities (mcqs)

ρ = M_{zz},  B = M_{z} − zM_{zz}  and  q = yM_{zzz}.  (2.3)

Physically, this means that these three quantities take a constant value on the trajectory of every fluid-parcel. Each may be interpreted physically: ρ = M_{zz} is the density of the fluid; B = M_{z} − zM_{zz} is the Bernoulli function which is associated with the energy per unit mass of the fluid (e.g., Batchelor, 1967); q = yM_{zzz} is the potential vorticity, a quantity related to the usual vorticity of a fluid (e.g., Pedlosky, 1987).

Welander (1971) showed that these three mcqs must be functionally-related, i.e., that

yM_{zzz} = M(M_{zz}, M_{z} − zM_{zz}),  (2.4)

for some function M. Significantly, (2.4) is an ordinary differential equation for M(z; x, y) (x and y play only a parametric rôle). Hence, if this function can be determined then we have only to integrate an ode and ensure consistency with the governing equation, (2.2i). This is a huge simplification of the problem.

Motivated by these ideas, we would like some systematic method of determining mcqs of a system of governing equations and also a method of relating these which, ideally, reflects the physics of the situation and leads to classes of solutions which are relevant. Furthermore, we would like these methods to have some degree of generality, i.e, be applicable to other systems of partial differential equations. We address these issues in the following sections.

Notation The following notation is used: ∂_{x}u = u_{x}, u_{y}, u_{z} and M_{x} = M_{x}, M_{y}, M_{z}, i.e., boldface subscripts indicate all derivatives of the corresponding order. Hence M_{xx} = M_{xx}, M_{xy}, … etc.

3. Determining ���s

In this section we first introduce formally the concept of material conservation-laws of a differential equation. We then show how one can directly search for such laws and recover the three quantities, (2.3). Finally we prove that no more materially-conserved quantities of the geostrophic equations, (2.2), exist.

Motivated by the conservation of physical quantities such as energy and momentum, a dynamical quantity Q = Q(t, s(t), ṡ(t)) (where s = (s_{1}, s_{2}, …) are the trajectory coordinates) is conserved if the differential equation

\frac{dQ}{dt} = (\frac{∂}{∂t} + ṡ• \frac{∂}{∂s} + s̈• \frac{∂}{∂ṡ}) Q = 0  (3.1)

is satisfied. More generally, a quantity

Q = Q(t, s(t), ṡ(t), s̈(t), …)

is conserved if

D_{t}(Q) = 0  (3.2)

where

D_{t} = \frac{∂}{∂t} + ṡ• \frac{∂}{∂s} + s̈• \frac{∂}{∂ṡ} + …  (3.3)

For an arbitrary differential equation

Δ(x, u, ∂_{x}u, …) = 0  (3.4)

in which we do not distinguish between temporal and spatial independent variables we can generalise further. We consider conserved quantities of the form

Q = Q(x, u, ∂_{x}u, …),  (3.5)

and our conservation law becomes

∑_{i=0}^{n} D_{i}(Q^{i}) = 0,  (3.6i)

where

D_{i} = \frac{∂}{∂x^{i}} + s_{i}^{α} \frac{∂}{∂s^{α}} + s_{ij}^{α} \frac{∂}{∂s_{j}^{α}} + …  (3.6ii)

which must be satisfied for all solutions u(x) of (3.4) (Ibragimov, 1994, Chapter 6).

Definition A conserved quantity, Q, of order n, is a function Q(x, u, ∂_{x}u, …, ∂x^{n}u) satisfying (3.6).

In continuum mechanics we often consider materially-conserved quantities, that is quantities constant on each trajectory; in this case (3.6) is more usually written

(\frac{∂}{∂t} + ṡ•∇) Q(x, u, ∂_{x}u, …) = 0,  (3.7)

equality holding on solutions.

One can use Equations (3.1) to (3.7) directly to find materially conserved quantities of a system of equations. In doing this one must remember that as a conservation law, (3.7), must be satisfied for all solutions, u(x, y, z, t), of a differential equation and one must therefore take account of its frame.

Definition The frame of a differential equation consists of the equation itself together with its differential consequences.

It may be the case that Q is a function of derivatives of higher order than the highest order terms in the given differential equation, (3.4); we must then take into account not only the equation itself, but also its differential consequences.

3.1. Determining ���s Algorithmically

With the formal definition of a material-conservation law in mind we now describe how one may directly determine materially-conserved quantities of the system (2.2). There is a considerable amount of algebra involved for third-order and higher laws, so to illustrate the principles involved more clearly we first consider laws of up to second order, i.e., we look for quantities,

Q = Q(x, M, M_{x}, M_{xx})  (3.8)

for which (3.7) is satisfied.

Substituting (3.8) into (3.7), using (2.2i) – (2.2ii) and expanding yields

\begin{matrix}− M_{yz}\{Q_{x} + Q_{M}M_{x} + Q_{M_{x}}M_{xx} + Q_{M_{y}}M_{xy} + Q_{M_{z}}M_{xz} + Q_{M_{xx}}M_{xxx} \\ + Q_{M_{xy}}M_{xxy} + Q_{M_{xz}}M_{xxz} + Q_{M_{yy}}M_{xyy} + Q_{M_{yz}}M_{xyz} + Q_{M_{zz}}M_{xzz}\} \\ + M_{xz}\{Q_{y} + Q_{M}M_{y} + Q_{M_{x}}M_{xy} + Q_{M_{y}}M_{yy} + Q_{M_{z}}M_{yz} + Q_{M_{xx}}M_{xxy} \\ + Q_{M_{xy}}M_{xyy} + Q_{M_{xz}}M_{xyz} + Q_{M_{yy}}M_{yyy} + Q_{M_{yz}}M_{yzz} + Q_{M_{zz}}M_{yzz}\} \\ + \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + Q_{M_{x}}M_{xz} + Q_{M_{y}}M_{yz} + Q_{M_{z}}M_{zz} + Q_{M_{xx}}M_{xxz} \\ + Q_{M_{xy}}M_{xyz} + Q_{M_{xz}}M_{xzz} + Q_{M_{yy}}M_{yyz} + Q_{M_{yz}}M_{yzz} + Q_{M_{zz}}M_{zzz}\} = 0. \end{matrix} (3.9)

This equation must hold simultaneously with the frame of (2.2i), however, since in this case we are concerned only with materially-conserved quantities of up to second order then we do not need to consider differential consequences of (2.2i), only the equation itself. We choose to eliminate M_{zzz} between (2.2i) and (3.9) yielding

\begin{matrix}− M_{yz}\{Q_{x} + Q_{M}M_{x} + … + Q_{M_{yz}}M_{xyz}\} \\ + M_{xz}\{Q_{y} + Q_{M}M_{y} + … + Q_{M_{yz}}M_{yzz}\} \\ + \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + … + Q_{M_{yz}}M_{yzz}\} = 0. \end{matrix}  (3.10)

(Care is needed in choosing which term to eliminate since Q and its derivatives potentially depend on some of the terms in (2.2i); Q does not depend on M_{zzz} in this case.) Note the coefficient of Q_{M_{zz}} is identically zero owing to the cancellation of certain terms. Since (3.10) must hold for all solutions M(x, y, z) and Q does not depend on third-order derivatives of M, then the coefficients of products of powers of these third-order derivatives must each be identically zero: from the coefficients of M_{xxx}, M_{xxy}, M_{xxz}, M_{xyy}, M_{xyz}, M_{yyz}, M_{yyy}, M_{xzz} and M_{yzz}, respectively, we obtain

M_{yz}Q_{M_{xx}} = 0,  (3.11a)
M_{xz}Q_{M_{xx}} − M_{yz}Q_{M_{xy}} = 0,  (3.11b)
M_{x}Q_{M_{xx}} − yM_{yz}Q_{M_{xz}} = 0,  (3.11c)
M_{xz}Q_{M_{xy}} − M_{yz}Q_{M_{yy}} = 0,  (3.11d)
y(M_{xz}Q_{M_{xz}} − M_{yz}Q_{M_{yz}}) + M_{x}Q_{M_{xy}} = 0,  (3.11e)
yM_{xz}Q_{M_{yz}} + M_{x}Q_{M_{yy}} = 0,  (3.11f)
M_{xz}Q_{M_{yy}} = 0,  (3.11g)
M_{x}Q_{M_{xz}} = 0,  (3.11h)
M_{x}Q_{M_{yz}} = 0.  (3.11i)

In general M_{x}, M_{xz} and M_{yz} are nonzero, so from (3.11) we conclude that

Q_{M_{xx}} = Q_{M_{xy}} = Q_{M_{xz}} = Q_{M_{yy}} = Q_{M_{yz}} = 0.  (3.12)

Note that we have not shown that Q_{M_{zz}} = 0 (recall the cancellation of terms mentioned above); indeed it turns out that Q_{M_{zz}} ≠ 0.

Using (3.12) we find that (3.10) simplifies substantially, viz.

\begin{matrix}− M_{yz} \{Q_{x} + Q_{M}M_{x} + Q_{M_{x}}M_{xx} + Q_{M_{y}}M_{xy}\} \\ + M_{xz} \{Q_{y} + Q_{M}M_{x} + Q_{M_{x}}M_{xx} + Q_{M_{y}}M_{yy}\} \\ + \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + Q_{M_{x}}M_{xz} + Q_{M_{y}}M_{yz} + Q_{M_{z}}M_{zz}\} = 0. \end{matrix}  (3.13)

Since we have determined that Q depends on only one second-order derivative of M, M_{zz}, then the coefficients of products of powers of all other second-order derivatives mush each be identically zero. From the coefficients of M_{xx}M_{yz}, M_{yy}M_{xz}, M_{yz}M_{xy} and M_{xz}M_{xy} we see that

Q_{M_{x}} = Q_{M_{y}} = 0,  (3.14)

and then using this result, from the coefficients of M_{xz} and M_{yz} we find

Q_{y} + Q_{M}M_{y} = 0,  (3.15a)
Q_{x} + Q_{M}M_{x} = 0,  (3.15b)
\frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + Q_{M_{z}}M_{zz}\} = 0.  (3.15c)

From (3.14) we know that Q is independent of M_{x} and M_{y}, so that from (3.15a) and (3.15b) we find

Q_{x} = Q_{y} = Q_{M} = 0.  (3.16)

It remains to satisfy (3.15c). M_{x} is not in general zero and then since Q_{M} = 0 then we require the solution of Q_{z} +Q_{M_{z}}M_{zz} = 0; there are two cases to consider: (i), Q_{z} = Q_{M_{z}} = 0 and (ii), Q_{z}Q_{M_{z}} ≠ 0. In Case (i) the general solution is Q = Q(M_{zz}) and in Case (ii) it is Q = Q(M_{z} − zM_{zz}). Hence all second-order, materially-conserved quantities of system (2.2) are of the form

Q = Q(M_{zz}, M_{z} − zM_{zz}).  (3.17)

It is easy to check that all quantities of this form are indeed materially-conserved.

We now turn to the problem of determining all materially-conserved quantities of up to third order, i.e., we assume

Q = Q(x, y, z, M, M_{x}, M_{y}, M_{z}, M_{xx}, M_{xy}, … M_{zz}, M_{xxx}, M_{xxy}, …, M_{zzz}).  (3.18)

Substituting this into (3.7) we obtain

\begin{matrix}− M_{yz}\{Q_{x} + Q_{M}M_{x} … + Q_{M_{zzz}}M_{xzzz}\} \\ + M_{xz}\{Q_{y} + Q_{M}M_{y} … + Q_{M_{zzz}}M_{yzzz}\} \\ + \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} … + Q_{M_{zzz}}M_{zzzz}\} = 0 \end{matrix}  (3.19)

and this must hold simultaneously with the frame of (2.2i) (cf. (3.9), above). This time, since (3.19) contains derivatives of M of up to fourth order we must consider not only (2.2i) itself, but in addition its derivative with respect to each independent variable, x, y and z, viz.,

M_{xx}M_{zzz} + M_{x}M_{xzzz} + y(M_{xxz}M_{yzz} + … − M_{yz}M_{xxzz}) = 0,  (3.20a)
M_{xy}M_{zzz} + M_{x}M_{yzzz} + M_{xz}M_{yzz}
− M_{yz}M_{xzz} + y(M_{xyz}M_{yzz} + … − M_{yz}M_{xyzz}) = 0,  (3.20b)
M_{xz}M_{zzz} + M_{x}M_{zzzz} + y(M_{xzz}M_{yzz} + … − M_{yz}M_{xzzz}) = 0.  (3.20c)

respectively. The procedure is to use (2.2i), (3.20a), (3.20b) and (3.20c) to eliminate four terms from (3.19) and then continue in a similar way to the second order case above. Again we must proceed carefully as the derivatives of Q in (3.19) potentially depend upon all derivatives of M up to third order and we cannot eliminate implicit terms! We first choose to eliminate M_{xxzz}, M_{yyzz} and M_{zzzz} from (3.19) by using (3.20a), (3.20b) and (3.20c), respectively; we may do this without difficulty, Q being independent of fourth-order derivatives of M. We obtain

\begin{matrix}Q_{M_{xzz}}\{M_{xx}M_{zzz} + M_{x}M_{xzzz} + y(M_{xxz}M_{yzz} + M_{xz}M_{xyzz} − M_{xyz}M_{xzz})\} \\ + Q_{M_{yzz}}\{M_{yz}M_{xzz} … + y(M_{yyz}M_{xzz} + …)\} \\ − M_{yz}\{Q_{x} + Q_{M}M_{x} + Q_{M_{x}}M_{xx} + … + Q_{M_{zzz}}M_{xzzz}\} \\ + M_{xz}\{Q_{y} + Q_{M}M_{y} + Q_{M_{x}}M_{xy} + … + Q_{M_{zzz}}M_{yzzz}\} \\ + \frac{M_{x}}{y} \{\frac{Q_{M_{xxz}}}{M_{yz}} [M_{xx}M_{zzz} + … + y(M_{xxz}M_{yzz} + …)] \\ + \frac{Q_{M_{yyz}}}{M_{xz}} [M_{yz}M_{xzz} − … + y(M_{yyz}M_{xzz} + …)] \\ + \frac{Q_{M_{zzz}}}{M_{x}} [y(M_{yz}M_{xzzz} − M_{xz}M_{yzzz}) − M_{xz}M_{zzz}] \\ + Q_{z} + Q_{M}M_{z} + … + Q_{M_{yzz}}M_{yzzz}\} = 0. \end{matrix}  (3.21)

Note that owing to the cancellation of certain terms the coefficient of Q_{M_{zzz}} is zero. Since (3.21) must be satisfied for all solutions M(x, y, z) of (2.2i) then the coefficients of products of powers of each fourth-order derivative must be identically zero (in drawing conclusions from these identities we must take account of (2.2i), at least in principle, as we have considered only (3.20a), (3.20b) and (3.20c) so far). From the coefficients of M_{xxxx}, M_{yyyy}, M_{xxxy}, M_{xxxz}, M_{xxyy}, M_{xxyz}, M_{xyyy}, M_{xyyz}, M_{xzzz}, M_{yzzz}, M_{yyyz} and M_{xyzz} we obtain respectively

−M_{yz}Q_{M_{xxx}} = 0,  (3.22a)
−M_{xz}Q_{M_{yyy}} = 0,  (3.22b)
−M_{yz}Q_{M_{xxy}} + M_{xz}Q_{M_{xxx}} = 0,  (3.22c)

−M_{yz}Q_{M_{xxz}} = 0,  (3.22d)
−M_{yz}Q_{M_{xyy}} + M_{xz}Q_{M_{xxy}} = 0,  (3.22e)
−M_{yz}Q_{M_{xyz}} + M_{xz}Q_{M_{xxz}} + \frac{M_{x}}{y} Q_{M_{xxy}} = 0,  (3.22f)
−M_{yz}Q_{M_{yyy}} + M_{xz}Q_{M_{xyy}} = 0,  (3.22g)
−M_{yz}Q_{M_{yyz}} + M_{xz}Q_{M_{xyz}} + \frac{M_{x}}{y} Q_{M_{xyy}} = 0.  (3.22h)
− (1 − \frac{1}{y} ) M_{x}Q_{M_{xzz}} + \frac{M_{x}^{2}Q_{M_{xxz}}}{yM_{yz}} = 0.  (3.22i)
−M_{x}Q_{M_{yzz}} − \frac{M_{x}^{2}Q_{M_{yyz}}}{yM_{xz}} + \frac{M_{x}}{y} Q_{M_{yzz}} = 0.  (3.22j)
M_{xz}Q_{M_{yyz}} + \frac{M_{x}}{y} Q_{M_{yyy}} = 0.  (3.22k)

−yM_{xz}Q_{M_{xzz}} − M_{yz}Q_{M_{yzz}} + yM_{yz}Q_{M_{yzz}} + M_{xz}Q_{M_{xzz}}

+ (\frac{yM_{xz}Q_{M_{xxz}}}{M_{yz}} + \frac{yM_{yz}Q_{M_{yyz}}}{M_{xz}} + Q_{M_{xyz}}) \frac{M_{x}}{y} = 0.  (3.22l)

In general M_{x}, M_{xz} and M_{yz} are nonzero, so from (3.22) we conclude that

\begin{matrix}Q_{M_{xxx}} = Q_{M_{yyy}} = Q_{M_{xxy}} = Q_{M_{xxz}} \\ = Q_{M_{xyy}} = Q_{M_{xyz}} = Q_{M_{yyz}} = Q_{M_{xzz}} = Q_{M_{yzz}} = 0. \end{matrix}  (3.23)

Note that we have not shown that Q_{M_{zzz}} = 0 as no coefficient of a fourth-order derivative of M in (3.21) involves Q_{M_{zzz}} (recall the cancellation of terms mentioned above); indeed it turns out that in general Q_{M_{zzz}} ≠ 0. In the special case Q_{M_{zzz}} = 0 we reduce to the search for second-order materially-conserved quantities considered above. We proceed with the general case.

It turned out that in obtaining (3.23) the condition (2.2i) is not needed; since we are now sure that Q is independent of both M_{xzz} and M_{yzz} we can use (2.2i) to eliminate one of these third-order terms from (3.21). Choosing M_{xzz} and using (3.23) then (3.21) becomes

\begin{matrix}− M_{yz}\{Q_{x} + Q_{M}M_{x} + Q_{M_{x}}M_{xx} + Q_{M_{y}}M_{xy} + Q_{M_{z}}M_{xz} + Q_{M_{xx}}M_{xxz} \\ + Q_{M_{xy}}M_{xxy} + Q_{M_{xz}}M_{xzz} + Q_{M_{yy}}M_{xyy} + Q_{M_{yz}}M_{xyz} + Q_{M_{zz}}M_{xzz} \\ + \frac{Q_{M_{zz}}}{yM_{yz}} [M_{x}M_{zzz} + yM_{xz}M_{yzz}] \} \\ + M_{xz}\{Q_{y} + Q_{M}M_{y} + Q_{M_{x}}M_{xy} + Q_{M_{y}}M_{yy} + Q_{M_{z}}M_{yz} + Q_{M_{xx}}M_{xxy} \\ + Q_{M_{xy}}M_{xyy} + Q_{M_{xz}}M_{xyz} + Q_{M_{yy}}M_{yyy} + Q_{M_{yz}}M_{yzz} + Q_{M_{zz}}M_{yzz}\} \\ + \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + Q_{M_{x}}M_{xz} + Q_{M_{y}}M_{yz} + Q_{M_{z}}M_{zz} + Q_{M_{xx}}M_{xxz} \\ + Q_{M_{xy}}M_{xyz} + Q_{M_{xz}}M_{xzz} + Q_{M_{yy}}M_{yyz} + Q_{M_{yz}}M_{yzz} + Q_{M_{zz}}M_{zzz}\} \\ − \frac{1}{y} M_{xz}M_{zzz}Q_{M_{zzz}} = 0. \end{matrix} (3.24)

(We find that there are no terms involving Q_{M_{zz}} owing to cancellation of terms.) Since we have shown that Q is independent of all but one third-order derivative of M then the coefficients of different products of powers of third-order terms in (3.24), except those involving M_{zzz} must be identically zero, viz.,

−M_{yz}Q_{M_{xx}} = 0,  (3.25a)
−M_{yz}Q_{M_{xy}} + M_{xz}Q_{M_{xx}} = 0,  (3.25b)
−M_{yz}Q_{M_{xz}} + \frac{M_{x}}{y} Q_{M_{xx}} = 0,  (3.25c)
−M_{yz}Q_{M_{yy}} + M_{xz}Q_{M_{xy}} = 0,  (3.25d)
−M_{yz}Q_{M_{yz}} + M_{xz}Q_{M_{xz}} + \frac{M_{x}}{y} Q_{M_{yz}} = 0,  (3.25e)
\frac{M_{x}}{y} Q_{M_{yz}} = 0,  (3.25f)
M_{xz}Q_{M_{yy}} = 0,  (3.25g)
M_{xy}Q_{M_{yz}} + \frac{M_{x}}{y} Q_{M_{yy}} = 0,  (3.25h)
\frac{M_{x}}{y} Q_{M_{xz}} = 0.  (3.25i)

In general M_{x}, M_{xz} and M_{yz} are nonzero, so from (3.25) we conclude that

Q_{M_{xx}} = Q_{M_{xy}} = Q_{M_{xz}} = Q_{M_{yy}} = Q_{M_{yz}} = 0.  (3.26)

Now, since Q_{M_{xx}} = Q_{M_{yy}} = 0 then the coefficients of different products of powers of the corresponding second-order derivatives of M in (3.25) must be identically zero and so, from the coefficients of M_{xx} and M_{yy}, we obtain

−M_{yz}Q_{M_{x}} = 0,  (3.27a)
M_{xz}Q_{M_{y}} = 0,  (3.28b)

respectively, from which we conclude

Q_{M_{x}} = Q_{M_{y}} = 0.  (3.29)

Using results (3.26) and (3.29) then (3.24) simplifies considerably, viz.,

\begin{matrix}− M_{yz}(Q_{x} + Q_{M}M_{x}) \\ + M_{xz} \{Q_{y} + Q_{M}M_{y} − \frac{1}{y} M_{zzz}Q_{M_{zzz}}\} \\ − \frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z} + Q_{M_{z}}M_{zz}\} = 0. \end{matrix}  (3.30)

Finally, since Q_{M_{yz}} = Q_{M_{x}} = 0 then from the coefficients of M_{x}M_{yz} and M_{yz} we obtain

Q_{M} = Q_{x} = 0,  (3.31)

leaving

M_{xz} \{Q_{y} − \frac{1}{y} M_{zzz}Q_{M_{zzz}}\} + \frac{M_{x}}{y} \{Q_{z} + M_{zz}Q_{M_{z}}\} = 0.  (3.32)

Since in general M_{x}M_{xz} ≠ 0 then

Q_{y} − \frac{1}{y} M_{zzz}Q_{M_{zzz}} = 0,  (3.33a)
Q_{z} + M_{zz}Q_{M_{z}} = 0.  (3.33b)

The general solutions to these equations are respectively Q = Q_{a}(yM_{zzz}, c_{a}), Q = Q_{b}(M_{z} − zM_{zz}, M_{zz}, c_{b}), where c_{a} and c_{b} are constants of integration; therefore the general materially-conserved quantity of the system (2.2) is

Q = Q(yM_{zzz}, M_{z} − zM_{zz}, M_{zz}).  (3.34)

It is easy to check that all quantities of this form are indeed materially-conserved.

3.2. The Question of Further ���s

Using an essentially algorithmic method we have found three materially-conserved quantities of the system (2.2) or order 3 or less. Are there other such quantities, of higher order? Using the method, what happens when there are no mcqs corresponding to the ansatz being studied?

Recall our definition of a mcq, (3.7). In the context of (2.2) this becomes

\{−\frac{M_{yz}}{y} \frac{∂}{∂x} + \frac{M_{xz}}{y} \frac{∂}{∂y} + \frac{M_{x}}{y^{2}} \frac{∂}{∂z} \} Q = 0,  (3.35)

which may be written (informally) as

\frac{dx}{−M_{yz}/y} = \frac{dy}{M_{xz}/y} = \frac{dz}{M_{x}/y} = \frac{dQ}{0}.  (3.36)

The general solution of (3.36) includes exactly three “constants” of integration so that the general solution of (3.35) is necessarily of the form

Q = S(C_{1}, C_{2}, C_{3}),  (3.37)

where C_{1}, C_{2} and C_{3} are functions of y and derivatives of M, i.e., Q can at most be a function of three arguments. Hence there are no materially-conserved quantities of fourthorder or higher and (3.34) is indeed the most general one (we discuss this idea for a general fluid-dynamical system in Section 5).

It is worth looking briefly at how the method would show that there were no (further) mcqs corresponding to a particular ansatz, for example, for (2.2), that there are no fourthorder quantities (one will not always be as fortunate as for the system above!). Rather than consider ansätze of the form Q = Q(x, M, M_{x}, M_{xx}, M_{xxx}, M_{xxxx}), or considering an alternative system of governing equations, we can see what will happen by looking again at the computations of Section 3.1. When looking for mcqs of up to second order we noted that the coefficient of Q_{M_{zz}} was identically zero (see text between Equations (3.10) and (3.11)); also, when looking for mcqs of up to third order we noted that the coefficient of Q_{M_{zzz}} was identically zero (see text between Equations (3.21) and (3.22)). These fortuitous events lead to (3.15c), and to (3.33a) and (3.33b), the final systems of equations for mcqs, Q. Were this not the case then one would have obtained equations forcing the conclusion that Q_{M_{zz}} = Q_{M_{zzz}} = 0 (cf. Equations (3.11) and (3.22)), and our final equations would be, respectively,

\frac{M_{x}}{y} \{Q_{z} + Q_{M}M_{z}\} = 0,

and

Q_{y} = Q_{z} = 0,

from which we must conclude that there are no mcqs. In short, the fortuitous cancellation of certain terms from the computation leading to “zero coefficients” for corresponding derivatives of Q, as one proceeds, leads to the existence of mcqs — should such “zero coefficients” not appear then nor will mcqs.

4. Constructing Solutions

In this section we consider the construction of solutions of the governing equations, (2.2), from the functionally-related materially-conserved quantities, i.e., from (2.4). First we consider an example solution, in an ocean basin, which is of particular physical significance and illustrates the generality of solutions obtainable from (2.4). The example shows that solutions constructed by means of the mcqs satisfy all the boundary conditions which are usually applied to the problem. Previous analytical solutions are either unable to satisfy a full set of general boundary conditions (e.g., Salmon & Hollerbach, 1991; Hood, 1996; Hood & Williams, 1996) or impose a particular physical character on the solution (e.g., Robinson & Stommel, 1959; Welander, 1959; Needler 1967 and 1971; Hodnett, 1978); either way generality is lost. Secondly we consider how one might systematically determine the function relating the mcqs, M in (2.4), by generalising (2.2) in a physically-motivated way.

4.1. An Example — An Ocean Basin

We now describe how, using solutions determined from (2.4), one can determine the ocean dynamics and thermodynamics in the northern Atlantic Ocean (for a full description of the problem see (Hood & Williams, 1996, and references contained therein). Figure 4.1 shows a schematic west-east vertical section through the region of interest with the boundary conditions which are to be applied. It is not clear, either from a physical or mathematical

w_{E} = w(x, y) sea surface: z = B_{ρ} = 0  ρ = ρ_{S}(x, y) sea floor: ρ = u = v = w = B = 0 no flow

Figure 4.1. A schematic west-east vertical section through the northern Atlantic Ocean showing the boundary conditions usually applied to the model: at the seasurface, z = B_{ρ} = 0, the vertical velocity, w_{E} = w(x, y) and, optionally, the surface density distribution, ρ_{S} = ρ(x, y), are imposed; at depth, a level of no-motion may be assumed; the western boundary is passive, but at the eastern boundary an integral condition is applied (see text); in addition the the density profile may be prescribed at certain locations, (x, y).

standpoint, what the function M should be. However, previous work on special cases suggests we take the following approach.

First, following Killworth (1987), we change to a coordinate system in which we use density, ρ, as the “vertical” coordinate, rather than depth, z. Then, choosing B(x, y, ρ) = p + ρz, the Bernoulli function, as the dependent variable in place of M(x, y, z), the governing equations, (2.2), become

u = −\frac{B_{y}}{y} ,  v = \frac{B_{x}}{y} ,  z = B_{ρ},  w = uz_{x} + vz_{y},  (4.1i)
(B_{x}B_{yρρ} − B_{y}B_{xρρ})y − B_{x}B_{ρρ} = 0;  (4.1ii)

the potential vorticity is

q = \frac{y}{B_{ρρ}}  (4.2)

and so (2.4) becomes

B_{ρρ} = yB(B, ρ)  (4.3)

(a considerable simplification).

The simplest function B is linear giving

B_{ρρ} = f_{1}(ρ)B + f_{0}(ρ),  (4.4)

where both f_{1}(ρ) and f_{0}(ρ) are to be determined. In fact this simple case is of particular physical significance: Killworth (1987) considered the case f_{0} = 0 with the restricted boundary conditions ρ = ρ_{S}(y) (i.e., the surface density is a function of only latitude) with w = w_{E}(y) (the surface vertical velocity is also a function of latitude, only) and/or u(0, y, z) = 0 (eastern velocity at the eastern boundary is zero at all depths and latitudes); Salmon (1994) considered f_{1} = 0 and added time-dependence to the problem.

From (4.4) we expect to find two functions of integration, a(x, y) and b(x, y); we also have to fix f_{1}(ρ) and f_{0}(ρ). The boundary conditions for this problem are a subject of some debate; we therefore consider two cases. (2.2) are a hyperbolic system (Huang, 1988) so one cannot expect to apply conditions at all boundaries: the northern, southern and western boundaries are all passive. We apply conditions at the surface, bottom and an integral condition at the eastern boundary. (For more detailed discussion of boundary conditions for this problem see Killworth, 1987; Samelson & Vallis, 1997; Huang, 1988.)

Case i.: We impose the vertical velocity field at the surface (more accurately, below some surface layer, subject to different dominant dynamics, which we patch onto our solution domain), i.e., we prescribe w_{E}(x, y), where

w_{E}(x, y) = \frac{1}{y} (B_{x}B_{yρ} − B_{y}B_{xρ})_{B_{ρ}=0}.  (4.5a)

We also assume that flow below some depth, which we take to be ρ = 0, without loss of generality, is an order of magnitude weaker than in our solution domain; this is described by

B(x, y, 0) = 0  (4.5b)

(Killworth, 1987, Section 2).

Case ii.: Alternatively, we can prescribe both w and ρ at the surface, i.e.,

w_{E}(x, y) = \frac{1}{y} (B_{x}B_{yρ} − B_{y}B_{xρ})_{B_{ρ}=0},  (4.6a)
ρ_{S}(x, y) = ρ _{B_{ρ}=0}.  (4.6b)

In fact these conditions are not enough to fix a(x, y) and b(x, y). In each case we have imposed one algebraic constraint, so that a and b are no longer independent, together with one first-order partial differential equation with independent variables x and y. Consequently it remains to fix an arbitrary function of some known function of x and y resulting from the integration of the first-order pde. This can be done by placing an (integral) constraint on the flow or heat-flux through the eastern boundary. Finally, f_{1}(ρ) and f_{0}(ρ) may be set by prescribing the variation of density with depth at two points (x, y).

4.2. Systematically Finding M

Above we considered a particular case of M in (2.4) (i.e., B in (4.3)), motivated by other authors’ work. The question remain, how given a particular problem, would one find M systematically? An answer, it turns out, is to make the model which the governing equations (2.2) represents more realistic.

(2.2) represent an ideal system, i.e., a system in which diffusion, friction and other such effects are assumed negligible. This is done for simplicity and is only approximately true. What happens if we introduce some non-ideal effect, parametrised by κ, 0 < κ ≪ 1? The relation (4.3) will now be only approximately true. We suppose

B_{ρρ} = yB(ρ, B, κx, κy),  (4.7)

i.e., the non-ideal effect is balanced, mathematically, by a slow variation in the functional relation between the mcqs.

For illustration let us consider Fickian diffusion, so that, (4.1ii) becomes

(B_{x}B_{yρρ} − B_{y}B_{xρρ})y − B_{x}B_{ρρ} = κy^{2}B_{ρρρ},  (4.8)

and let us work with the special case studied by Salmon (1994), i.e.,

B_{ρρ} = yF(ρ).  (4.9)

Generalising (4.9) to include slow variation with x and y, and integrating we obtain

B(x, y, ρ) = yF (ρ, κx, κy) + ρb(x, y) + a(x, y),  F (ρ) = ∫ ^{ρ}∫ ^{ρ_{1}}\frac{1}{Q(ρ_{2})} dρ_{2} dρ_{1} (4.10)

(where a and b are functions of integration, to be determined), then substituting into (4.8) terms of O(1) cancel leaving

(ρb_{x} + a_{x})F_{Y ρρ} = y^{2}F_{ρρρ} + O(κ^{2}),  Y = κy.  (4.11)

This differential constraint does not fully determine F, or equivalently F, but perhaps points the way. This approach, for both Salmon’s special case, (4.9), and the general case is actively being considered by the author (Hood, 1997).

5. Discussion and Conclusions

In the preceding sections we have considered how one might use the conservation laws of the governing partial differential equations describing a system to construct exact analytical solutions; we have focused on the equations describing large-scale flow within the ocean (the geostrophic equations). For these equations we used a systematic method to recover three materially-conserved quantities which were necessarily functionally-related; this relation takes the form of a differential equation which is significantly simpler than the governing equations, the solution of which can be used as a quite general solution-ansatz for the governing pdes. We also proved that no more mcqs exist. Finally we considered an example, an ocean basin, which illustrated the generality of the solutions obtained and suggested a mechanism reflecting higher-order dynamics by which one might determine the function relating the quantities.

Generality of solutions It is important to emphasise that solutions obtained by means of relating materiallyconserved quantities are significantly more general than those obtained by other analytical methods, for these equations. In early work on (2.2) informed guesses where made for forms of solution sometimes based on similarity variables or partial separation of variables, e.g., Robinson & Stommel, 1959; Welander, 1959; Needler 1967 and 1971; Hodnett, 1978. These works mark significant analytical progress on a formidable system of equations. However, in each case there are limitations on boundary conditions which can be satisfied and the physical character of the solutions, perhaps motivated by the physics, is built into the solution method. An alternative, systematic (indeed algorithmic) method has been tried: Salmon & Hollerbach (1991) used Lie’s method (e.g., Bluman & Kumei, 1989; Olver, 1986; Stephani, 1990) to determine a set of point symmetries of (2.2), with Fickian diffusion added in (2.2c) (cf. (4.8)). This work was extended, for more general diffusion, using a more general method which obtained more symmetries by Hood (1996) and Hood & Williams (1996). Further, the time-dependent problem was studied (using Lie’s method) by Edwards (1996). These symmetries were used to construct classes of solutions which were general enough to carry out simple experiments to investigate particular physical processes, however, again, only a limited class of boundary conditions can be satisfied. In contrast, integration of (2.4) (equivalently (4.3)) yields the general solution of the governing equations, (2.2) (equivalently (4.1)), in the sense that there are two arbitrary functions of x and y and (at least) one of z (or ρ). As illustrated in Section 4.1 these solutions can satisfy the full set of boundary conditions that one would expect to impose.

Wider Application These ideas are certainly not restricted to the geostrophic equations. In general, given a system of pdes in n independent variables with n materially-conserved quantities known, one can construct a solution ansatz in exactly the same way. For example consider the 2D-Euler equations,

\frac{Du}{Dt} = g − \frac{1}{ρ}∇p.  (5.1)

Both vorticity, ω (where ωk = ω = ∇ × u), and the stream-function, ψ (where u = ∂ψ/∂y and v = −∂ψ/∂x, which automatically ensures mass conservation), are materially-conserved, whence

\frac{∂^{2}ψ}{∂x^{2}} + \frac{∂^{2}ψ}{∂y^{2}} = G(ψ),  (5.2)

for some function G. (5.2) alone determines the flow, given G. Certain functions have been considered (e.g., many workers have assumed a linear q-ψ relationship, but there is analytical, numerical and experimental evidence that this is not always possible, especially in the case of a “tripolar” vortex (van Heijst et al., 1991; Legras et al., 1988)); as in Section 4.2 one could add a small amount of some non-ideal effect to suggest what G might be.

It is likely that for any system of equations for which these ideas are fruitfully employed, the resulting analytical solutions will be significantly more general than those resulting from point-symmetries — stronger boundary conditions will be satisfied. Furthermore, one may be able to construct useful solutions even with fewer than n mcqs known — indeed, the most physically relevant case may involve fewer than n (e.g., Salmon, 1994).

In Section 3.2 we showed that the geostrophic equations, (2.2), have no more than three mcqs. This principle carries over to any fluid-dynamical system. An mcq must satisfy

\{u\frac{∂}{∂x} + v \frac{∂}{∂y} + w \frac{∂}{∂z} \} Q = 0  (5.3)

(for a steady system), which may be written (informally) as

\frac{dx}{u} = \frac{dy}{v} = \frac{dz}{w} = \frac{dQ}{0}.  (5.4)

In principle, one can always solve the governing equations for u(x, y, z), v(x, y, z) and w(x, y, z), and then integrate (5.4) for Q as a function of three “constants” of motion. This is true for any u(x, y, z), v(x, y, z) and w(x, y, z). It follows that a three-dimensional, steady system has three mcqs. This principle is readily extended to systems with more (or fewer) independent variables. However, to be of interest, a quantity must be materiallyconserved for a wide class of solutions of the governing equations (preferably all solutions), i.e., we want to be able to find mcqs without first solving the governing equations — indeed, we may well intend that knowledge of the mcqs is to help us solve them! If the equations describing a fluid-dynamical system can be put in potential form, cf. (2.2), then we have a much better chance of finding useful mcqs. In this case we can construct an equation such as (3.35). It may be possible to integrate this directly, though this is unlikely. Usually finding the mcqs by using the method described in Section 3.1 will be the most fruitful approach: given that the number of mcqs is known from the number of independent variables then one can consider candidates of increasing order, using a computer algebra system, until all these are found.

Relation of point-symmetry solutions and mcq solutions; Noether The question arises, what is the connection between solutions obtained by means of point symmetries (of the governing equations) and those constructed from materially-conserved quantities? Is one a subset of the other? Noether’s Theorem (e.g, Bluman & Kumei, 1989; Stephani, 1990) gives us some insight.

Solutions obtainable from material-conservation laws are a consequence of the existence of point-symmetries of the action integral of a system; in contrast, the usual pointsymmetry methods determine solutions which reflect symmetries of the governing equations themselves. Since the solutions constructed from mcqs are general in the sense stated above then one would usually expect solutions constructed from point symmetries to be special cases of the former — of course some of the point-symmetry solutions may well be singular.

Noether’s theorem suggests an alternative approach to the construction of conservation laws to that advocated in this paper: determine the point symmetries of the action integral and use these to construct the mcqs. This approach is also essentially algorithmic, however, one is limited to those systems for which a Lagrangian is known. This is a significant weakness from which the more direct approach of Section 3 does not suffer.

Backwards: Finally, one can turn the whole problem around. Given that one can determine the materially-conserved quantities of a fluid-dynamical system from the velocity field, an obvious question that arises is, given n mcqs of a system in n independent variables, can one determine the corresponding dynamics? Yes — but not uniquely. Let us return to the three mcqs given by (2.3). Any velocity field for which these three quantities are materiallyconserved must satisfy

uM_{xzz} + vM_{yzz} + wM_{zzz} = 0,  (5.5a)
u(M_{xz} − zM_{xzz}) + v(M_{yz} − zM_{yzz}) − wzM_{zzz} = 0,  (5.5b)
uyM_{xzzz} + v(M_{zzz} + yM_{yzzz}) + wyM_{zzzz} = 0,  (5.5c)

which is a 3 × 3 linear, algebraic system for u, v and w. Following some straightforward algebra we find a solution exists provided

M_{xz}M_{zzz}^{2} + yM_{zzzz}(M_{yz}M_{xzz} − M_{xz}M_{yzz}) + yM_{zzz}(M_{xz}M_{yzzz} − M_{yz}M_{xzzz}) = 0, (5.6)

and then the general solution, for u, v and w may be written,

M_{xz}u + M_{yz}v = 0,  (5.7a)
w = −\frac{M_{xzz}u + M_{yzz}v}{M_{zzz}} .  (5.7b)

An integrating factor for for (5.6) is easily found: multiplying through by M_{zzz}^{−2} and integrating we find

M_{x}M_{zzz} + y(M_{xz}M_{yzz} − M_{yz}M_{xzz}) = M_{0}(x, y)M_{zzz},  (5.8)

where M_{0}(x, y) is a function of integration. The dynamics corresponding to the mcqs (2.3) is therefore given by (5.7) and (5.8). (5.8) is in fact equivalent to (2.2i) (substitute M(x, y) + ∫^{x} M_{0}(x̂, y) dx̂ for M(x, y)); this equivalence reflects a point symmetry of (2.2i) (Salmon & Hollerbach, 1991). Nevertheless, a degree of freedom in the velocity field remains — the dynamics are not uniquely determined by the mcqs.

Acknowledgements

I would like to thank two anonymous referees for their comments on the first draft of this paper. This work was in part funded by grant no. R.D.F. ADMF5B from the University of Liverpool.

References

Batchelor, G.K., 1967, “An Introduction to Fluid Dynamics,” pub. C.U.P. Bluman G.W. and Kumei S., 1989, “Symmetries and Differential Equations,” Appl. Math. Sci., 81, Springer-Verlag, Berlin-Heidelberg-New York. Edwards N., 1996, “Unsteady similarity solutions and oscillating ocean gyres,” J. Mar. Res., 54, 793 – 826. Hodnett P.F., 1978, “On the advective model of the thermocline circulation,” J. Mar. Res., 36, 185 – 198. Hood, S., 1997, “Using materially-conserved quantities to determine the profile function in generalised two-layer models”, submitted to J. Mar. Res. Hood, S. and R.G. Williams, 1996, “On frontal and ventilated models of the main thermocline,” J. Mar. Res., 54, 211 – 238. Hood, S., 1996, “New similarity solutions of the thermocline equations — diffusion an arbitrary function of z,” J. Mar. Res., 54, 187 – 209. Huang, R.X., 1988, “On Boundary Value Problems of the Ideal-Fluid Thermocline,” J. Phys. Ocean., 18, 619 – 641. Ibragimov, N.H., 1994, “CRC Handbook of Lie Group Analysis of Differential Equations, Volume 1,” CRC Press Inc. Killworth, P.D., 1987, “A Continuously Stratified Nonlinear Ventilated Thermocline,” J. Phys. Ocean., 17, 1925 – 1943. Legras, B., P. Santangelo and R. Benzi, 1988, “High-resolution numerical experiments for forced twodimensional turbulence,” Europhys. Lett., 5, 37 – 42. Olver P.J., 1986, “Applications of Lie Groups to Differential Equations,” Graduate Texts Math, 107, Springer-Verlag, New York Needler, G.T., 1971, “Thermocline models with arbitrary barotropic flow,” Deep Sea Res., 18, 895– 903. Needler, G.T., 1967, “A model for thermohaline circulation in an ocean of finite depth,” J. Mar. Res., 25, 329 – 342. Pedlosky, J., 1987, “Geophysical Fluid Dynamics, 2nd. ed.,,” Springer-Verlag. Robinson, A.R. and H. Stommel, 1959, “The ocean thermocline and the associated thermohaline circulation,” Tellus, 11, 295 – 308. Samelson, R.M. and G.K. Vallis, 1997, “Large-scale circulation with small diapycnal diffusion: the twothermocline limit,” J. Mar. Res., 55, 223 – 275. Salmon, R., 1994, “Generalized two-layer models of ocean circulation,” J. Mar. Res., 52, 865 – 908. Salmon, R. and R. Hollerbach, 1991, “Similarity solutions of the thermocline equations,” J. Marine Research, 49, 249 – 280. Stephani, H., 1990, “Differential Equations, their Solution using Symmetries,” ed. M. MacCallum, Cambridge University Press, Cambridge van Heijst, G.J.F., R.C. Kloosterziel and C.W.M. Williams, 1991, “Laboratory experiments on the tripolar vortex in a rotating fluid,” J. Fluid. Mech., 225, 301 – 331. Welander, P., 1971, “Some exact solutions to the equations describing an ideal-fluid thermocline,” J. Mar. Res., 29, 60 – 68. Welander, P., 1959, “An Advective Model of the Ocean Thermocline,” Tellus, 11, 309 – 318.