# Geometric theory of non-regular separation of variables and the bi-Helmholtz
equation

Claudia Chanu
Dipartimento di Matematica, Università di Torino,
Torino, Italia
<EMAIL_ADDRESS>Basel Jayyusi
Department of Applied Mathematics, University of Waterloo,
Waterloo, Ontario, Canada
<EMAIL_ADDRESS>Raymond G McLenaghan
Department of Applied Mathematics, University of Waterloo,
Waterloo, Ontario, Canada
<EMAIL_ADDRESS>
(August 28, 2024)

###### Abstract

The geometric theory of additive separation of variables is applied to the
search for multiplicative separated solutions of the bi-Helmholtz equation. It
is shown that the equation does not admit regular separation in any coordinate
system in any pseudo-Riemannian space. The equation is studied in the four
coordinate systems in the Euclidean plane where the Helmholtz equation and
hence the bi-Helmholtz equation is separable. It is shown that the bi-Helmoltz
equation admits non-trivial non-regular separation in both Cartesian and polar
coordinates, while it possesses only trivial separability in parabolic and
elliptic-hyperbolic coordinates. The results are applied to the study of small
vibrations of a thin solid circular plate of uniform density which is governed
by the bi-Helmholtz equation.

##  1 Introduction

The solution of boundary value problems for the partial differential equations
of mathematical physics by the method of separation of variables is an
effective method that has been employed for almost two centuries [10, 9]. The
method, which assumes that the solution is a product of functions each one of
a single independent variable (product ansatz), reduces the partial
differential equation to a set of ordinary differential equations with a
corresponding set of separation constants. The application of the method to
the various equations of mathematical physics leads to certain ordinary
differential equations the solutions of which are studied in their own right
in special function theory. The special functions defined by an equation
arising from a given separable coordinate system define function spaces that
yield the solution of the given boundary problem under general boundary
values. The theory of this method for a general second order linear PDE with
variable coefficients under the product ansatz is well developed [1, 7, 8]
However, the theory for equations of higher order has only more recently been
studied [6] where general separation is considered and where the concept of
non-regular separation may be utilized. In this article a geometric
formulation of this theory [4] is applied to the bi-Helmholtz equation which
arises in the study of small vibrations of a thin, solid plate [11, 5].

The plan of the paper is as follows. In Section 2 the geometric theory of
separation of variables is reviewed and the concepts of regular and non-
regular separation defined. In Section 3 the theory of the previous section is
applied to the bi-Helmholtz equation. The main result of this section is the
proof of Theorem 1 which states that regular multiplicative separation for the
bi-Helmholtz equation $(\ref{eq_bilap})$ on any Riemannian or pseudo-
Riemannian $n$-dimensional manifold does not occur in any system of
coordinates. Conditions for non-regular separation to occur are given in
Proposition 2 and the concept of non-trivial non-regular separation is
introduced. Section 4 contains examples in $\mathbb{E}^{2}$ of coordinate
systems where the bi-Helmholtz equation admits trivial and non-trivial non-
regular separation. In Section 5 the results of the preceding section are
applied to the analysis of the small vibrations of a thin solid circular plate
which reproduce the classical results obtained by Rayleigh [11]. The
Conclusion is given in Section 6.

##  2 Geometric theory of separation of variables

One of the most known ansatz for solving a PDE depending on $n$ independent
variables $(q^{i})$ is the additive separation, that is the search of those
solutions $u(q^{i})$ which are written as sum of functions depending on a
single variable (in the given coordinate system $(q^{i})$):

$u=\sum_{i=1}^{n}S_{(i)}(q^{i})$

Moreover, we do not look for a single separated solution but we want to
determine a family (as big as possible) of separated solutions

$u=\sum_{i=1}^{n}S_{(i)}(q^{i},c_{\alpha}),$

where the real constant parameters $c_{\alpha}$ satisfy a suitable
completeness condition (see below). When the method of SoV works, it is
possible to split the PDE into $n$ ODEs of order $l$ which involve only one of
the functions $S_{(i)}$; these ODEs are known as _separated equations_.
Separability of a PDE strongly depends on the choice of the independent
variables and is destroyed by a general transformation of the $(q^{i})$.
However, we want to provide a geometrical interpretation of what we mean by
solving a PDE through the ansatz of separation of variables. We will
distinguish between two possible types of families of separated solutions: the
first type depends on the maximal number of parameters (we refer to this case
as ”free” [1] or ”regular” [6] separation of variables) and the second type
depending on less parameters. We refer to this case as ”non regular
separation” [6, 3].

We consider the separability of the $l$-th order PDE

$\mathcal{H}(q^{i},u,u_{i},\ldots,u_{i}^{(l)})=h\qquad(h\in\mathbb{R})$ (1)

in the coordinates $(q^{i})$ on the $n$-dimensional manifold $Q$. Here and in
the following we denote by $u=u(q^{i})$ the unknown function and set

$u_{i}=\frac{\partial u}{\partial q^{i}},\qquad
u^{(2)}_{i}=u_{ii}=\frac{\partial^{2}u}{(\partial q^{i})^{2}},\quad\ldots\quad
u_{i}^{(l)}=\frac{\partial^{l}u}{(\partial q^{i})^{l}}.$

For the sake of simplicity, we assume that the maximal order of derivatives
involved in $H$ is the same for each index $i=1,\ldots,n$. Let $Z$ be the
$(nl+1)$-dimensional space of the dependent variable and its separated
derivatives: coordinates on $Z$ are given by
$(u,u_{i},u_{i}^{(2)},\ldots,u_{i}^{(l)})$. We consider the trivial bundle
over $Q$, $M=Q\times Z$. In [1] free separation is defined as the existence of
an additively separated solution $u$ of (1), depending on $nl+1$ parameters
$(c_{A})$, satisfying the completeness condition

$\mathrm{rank}\left[\frac{\partial u}{\partial c_{A}}\bigg{|}\frac{\partial
u_{i}}{\partial c_{A}}\bigg{|}\ldots\bigg{|}\frac{\partial
u_{i}^{(l)}}{\partial c_{A}}\right]=nl+1.$

Free separation occurs if and only if the $n$ vector fields of the form

$D_{i}=\frac{\partial}{\partial q^{i}}+u_{i}\frac{\partial}{\partial
u}+u_{i}^{(2)}\frac{\partial}{\partial
u_{i}}+\ldots+u_{i}^{(l)}\frac{\partial}{\partial
u_{i}^{(l-1)}}+R_{i}\frac{\partial}{\partial u_{i}^{(l)}},$ (2)

where $R_{i}(q^{j},u,u_{j},\ldots,u_{j}^{(l)})$ are functions on $M$, are
commuting symmetries of (1) namely, the conditions

$\displaystyle D_{i}\mathcal{H}=0,$ (3) $\displaystyle[D_{i},D_{j}]=0$ (4)

are satisfied (see [1]). Eqs (3) determine the functions $R_{i}$, while (4)
are equivalent to

$D_{i}R_{j}=0\qquad(i\neq j).$ (5)

By expanding conditions (5), we get exactly the conditions given in [6] for
regular separation. Thus, free separation corresponds to regular separation.

######  Remark 1.

The geometrical interpretation of the free separation is summarized in the
following items (see [1]):

  * •

$\Delta=\hbox{\rm span}(D_{i})$ is an integrable distribution of rank $n$ on
$M$.

  * •

The foliation of $n$-dimensional integral manifolds of $\Delta$ is described
by a complete separated solution of $\mathcal{H}=h$:

$\left\\{\begin{array}[]{l}u=S=\sum_{i=1}^{n}S_{(i)}(q^{i},c_{\alpha})\\\
u_{i}=S_{(i)}^{\prime}\\\ u_{ii}=S_{(i)}^{\prime\prime}\\\
\cdots\end{array}\right.\quad\hbox{with }\partial_{i}S=S_{(i)}^{\prime},\
\partial^{2}_{ij}S=\delta_{ij}S_{(i)}^{\prime\prime},\ \ldots$

  * •

Free (regular) separated solutions depend on $nl+1$ parameters.

  * •

Completeness means that for any point $P_{0}$ with coordinates $(q_{0}^{i})$
there is a separated solution of $\mathcal{H}=h$ for each choice of the value
of $u$ and its $nl$ derivatives $u_{i},u_{i}^{2},\ldots,u_{i}^{l}$ at $P_{0}$.

######  Remark 2.

If $\mathcal{H}$ does not depend on $u$, then $u$ is defined up to an additive
constant and the relevant constants appearing in the separated solutions are
in fact $nl$ (the remaining one being the finial additive constant). In this
case we can eliminate the variable $u$ from the space $Z$ and (2) reduce to

$D_{i}=\frac{\partial}{\partial q^{i}}+u_{i}^{(2)}\frac{\partial}{\partial
u_{i}}+\ldots+u_{i}^{(l)}\frac{\partial}{\partial
u_{i}^{(l-1)}}+R_{i}\frac{\partial}{\partial u_{i}^{(l)}}.$ (6)

Moreover, even the separated equations do not depend on $S_{(i)}$. Thus, each
of $S_{(i)}$ is defined by its ODE up to an additive constant which can be
disregarded, since it only affect the inessential global additive constant of
$u$ . Hence, each $l$-th order ODE actually contribute with $l-1$ essential
integration constants. This means that the $nl$ constants which are involved
in the regular separated solutions splits into $n$ separation constants plus
$nl-n=n(l-1)$ integrating constants (see Remarks 4.8 and 4.9 in [1] for the
case of the Schrödinger equation).

We consider now the second type of separation, occurring when condition (5) is
not identically satisfied. This case was firstly considered by Kalnins and
Miller in [6] and called non-regular separation. The authors state that also
in this case separable solutions still may exist, but they will depend on less
than $nl+1$ parameters.

However, in [6] it is not specified under what conditions separable solutions
will exist and how to determine the number of the parameters involved.

In [3] one of us gives a geometric interpretation of the situation that
naturally leads to an effective definition of non-regular separation.

Let $N$ be a submanifold of $M=Q\times Z$ locally described by the $r$
equations

$f_{a}=0\qquad(a=1,\ldots,r).$

If the vectors $D_{i}$ commute on $N$, that is

$D_{i}R_{j}|_{N}=0\qquad(i\neq j,\ i,j=1,\ldots,n),$

and the vector fields $D_{i}$ are tangent to $N$, that is

$D_{i}f_{a}|_{N}=0\quad(i=1,\ldots,n,\ a=1,\ldots r),$

then we can restrict the distribution $\Delta$ generated by the $D_{i}$ to an
involutive distribution on $N$. Thus, on $N$ we have a complete separated
solution of $\mathcal{H}=h$ depending on $nl+1-r$ parameters $(c_{\alpha})$,
which can be considered as a constrained separated solution on $M$. Hence,
according to [3] we have the following

######  Definition 1.

The PDE $\mathcal{H}=h$ admits a non-regular or _c_ onstrained additive
separation on a submanifold $N$ of $M$ defined by the $r$ equations $f_{a}=0$,
if

  1. 1.

$u=\sum_{i}S^{(i)}(q^{i},c_{\alpha})$ is a solution of $\mathcal{H}=h$;

  2. 2.

$u$ depends on $nl+1-r$ parameters $(c_{\alpha})$ satisfying the completeness
conditions

$\hbox{rank}\left[\frac{\partial u}{\partial c_{\alpha}}\bigg{|}\frac{\partial
u_{i}}{\partial c_{\alpha}}\bigg{|}\ldots\bigg{|}\frac{\partial
u_{i}^{(l)}}{\partial c_{A}}\right]=nl+1-r;$

  3. 3.

$u$ and its derivatives satisfy $f_{\alpha}(q^{j},u,u_{i}...)=0$ for all
(admissible) values of the parameters $(c_{\alpha})$.

From the above discussion we get the following criterion for the constrained
separation

######  Proposition 1.

In a given coordinate system $(q^{i})$ equation $(\ref{eq})$ admits a
constrained separation on the submanifold $N$ defined by equations $f_{a}=0$
if and only if the vector fields $D_{i}$ $(\ref{Di})$ are symmetries of
$(\ref{eq})$, tangent to $N$ and commute on $N$, that is

$D_{i}\mathcal{H}=0,\quad D_{i}R_{j}|_{N}=0,\quad D_{i}f_{a}|_{N}=0\quad(i\neq
j=1,\ldots,n,\;a=1,\ldots,r).$

######  Remark 3.

If $D_{i}R_{j}$ are everywhere different from zero, then equation (1) has no
additive separable solutions. A possible choice for $N$ is the set of points
satisfying equations $D_{i}R_{j}=0$, but in many cases this set is not a well-
defined manifold or the vectors $D_{i}$ are not tangent to it.

######  Remark 4.

In analogy with Remark 1, the geometric interpretation of the constrained
separation on a submanifold $N$ is sketched in the following items:

  * •

$\Delta_{S}=\hbox{\rm span}(D_{i})$ is an integrable $n$-dimensional
distribution on $N$ whose integral manifolds are described by a constrained
complete separated solution of $\mathcal{H}=h$.

  * •

constrained separated solutions depend on $nl+1-r$ parameters.

  * •

completeness means that for any $P_{0}=(q^{i}_{0})$ there is a separated
solution of $\mathcal{H}=h$ for each choice of the value of $u$ and its $nl$
separated derivatives at $P_{0}$ satisfying $f_{\alpha}=0$ (i.e., the initial
condition belongs to $N$).

##  3 Application to the bi-Helmholtz equation

In order to apply the theory recalled in the previous section to the
multiplicative separation of the bi-Helmholtz equation, we should transform
the unknown function $\psi$ to its logarithm $u=\log\psi$ and then we should
write down the simplified equation for separated solutions. First of of all,
we need to express the bi-Helmholtz equation

$\Delta^{2}\psi=\lambda\psi,\qquad\lambda\in\mathbb{R}$ (7)

where $\Delta$ is the Laplace-Beltrami operator on a $n$-dimensional
Riemannian or pseudo-Riemannian manifold $(Q,\mathbf{g})$ with respect to a
coordinate system $(q^{i})$ in which the metric has contravariant components
$(g^{ij})$. For the moment we do not require any special property to the
coordinates (such as orthogonality, or specifying of the dimension $n$). It is
well known that

$\Delta\psi=g^{ij}\partial_{i}\partial_{j}\psi-\Gamma^{i}\partial_{i}\psi$ (8)

where $\partial_{i}$ is the partial derivative with respect to $q^{i}$ and

$\Gamma^{h}=g^{ij}\Gamma_{ij}^{h}=\frac{1}{2}\,g^{ij}g^{hk}(\partial_{i}g_{jk}+\partial_{j}g_{ik}-\partial_{k}g_{ij})$

are the contracted Christoffel symbols (see [1]).

If not explicitly written, Einstein summation notation is sytematically
employed: repeated upper and lower indices means summation.

By applying twice the Laplace operator we get

$\displaystyle\Delta^{2}\psi$ $\displaystyle=$ $\displaystyle
g^{ij}\partial_{ij}(g^{hk}\partial_{hk}\psi-\Gamma^{h}\partial_{h}\psi)-\Gamma^{i}\partial_{i}(g^{hk}\partial_{hk}\psi-\Gamma^{h}\partial_{h}\psi)=$
$\displaystyle=$ $\displaystyle
g^{ij}g^{hk}\partial_{ijhk}\psi+2(g^{ij}\partial_{j}g^{hk}-g^{hk}\Gamma^{i})\partial_{ihk}\psi+$
$\displaystyle(g^{ij}\partial_{ij}g^{hk}-2g^{jk}\partial_{j}\Gamma^{h}-\Gamma^{i}\partial_{i}g^{hk}+\Gamma^{h}\Gamma^{k})\partial_{hk}\psi+$
$\displaystyle(-g^{ij}\partial_{ij}\Gamma^{h}+\Gamma^{i}\partial_{i}\Gamma^{h})\partial_{h}\psi,$

where

$\partial_{i_{1}\ldots i_{p}}:=\partial_{i_{1}}\ldots\partial_{i_{p}}.$ (10)

In order to simplify notation we can write the bi-Laplacian as

$\Delta^{2}\psi=A^{ijkl}\partial_{ijkl}\psi+B^{ijk}\partial_{ijk}\psi+C^{ij}\partial_{ij}\psi+D^{i}\partial_{i}\psi,$
(11)

where

$\displaystyle A^{ijkl}=g^{(ij}g^{kl)}$ (12) $\displaystyle
B^{ijk}=2(g^{h(i}\partial_{h}g^{jk)}-g^{(ij}\Gamma^{k)})$ (13) $\displaystyle
C^{ij}=g^{kl}\partial_{kl}g^{ij}-2g^{k(i}\partial_{k}\Gamma^{j)}-\Gamma^{k}\partial_{k}g^{ij}+\Gamma^{i}\Gamma^{j}$
(14) $\displaystyle
D^{i}=-g^{jk}\partial_{jk}\Gamma^{i}+\Gamma^{j}\partial_{j}\Gamma^{i}=-\Delta\Gamma^{i}$
(15)

where $(\cdots)$ indicates symmetrization of the indices.

In order to pass from a multiplicative separated solution to an additive
separated solution we perform the change of unknown

$\psi=e^{u},\qquad u=\log\psi$

We denote the partial derivatives of $u$ by $u_{i}$, $u_{ij}$ etc. The link
between partial derivatives of $\psi$ and $u$ is given by

$\begin{array}[]{l}\partial_{i}\psi=e^{u}u_{i}\vphantom{\dfrac{1}{2}}\\\
\partial_{ij}\psi=e^{u}(u_{ij}+u_{i}u_{j})\vphantom{\dfrac{1}{2}}\\\
\partial_{ijk}\psi=e^{u}(u_{ijk}+3u_{(i}u_{jk)}+u_{i}u_{j}u_{k})\vphantom{\dfrac{1}{2}}\\\
\partial_{ijkl}\psi=e^{u}(u_{ijkl}+4u_{(i}u_{jkl)}+3u_{(ij}u_{kl)}+6u_{(i}u_{j}u_{kl)}+u_{i}u_{j}u_{k}u_{l})\vphantom{\dfrac{1}{2}}\end{array}$

Hence, the equation to which we want to apply the geometric theory of SoV is

$H(q^{i},u_{i},u_{ij},\ldots,u_{ijkl})=\lambda,\qquad\lambda\in\mathbb{R},$

where

$\begin{array}[]{l}H=A^{ijkl}u_{ijkl}+(B^{ijk}+4A^{ijkl}u_{l})u_{ijk}+(C^{ij}+3B^{ijk}u_{k})u_{ij}+\\\
\qquad
3A^{ijkl}u_{ij}u_{kl}+6A^{ijkl}u_{i}u_{j}u_{kl}+\vphantom{\dfrac{1}{2}}\\\
\qquad
A^{ijkl}u_{i}u_{j}u_{k}u_{l}+B^{ijk}u_{i}u_{j}u_{k}+C^{ij}u_{i}u_{j}+D^{i}u_{i},\end{array}$

If we are interested only in separated solutions

$u=\sum_{i}S_{(i)}(q^{i},c_{A}),\qquad A=1,\ldots,4n,$

where, since the PDE does not depend explicitly on $u$, the constants
$(c_{A})$ are at most $4n$, we replace $H$ by

$\begin{array}[]{l}H_{s}=(g^{ii})^{2}u_{i}^{(4)}+(4g^{ii}g^{ij}u_{j}+B^{iii})u_{i}^{(3)}+(g^{ii}g^{jj}+2(g^{ij})^{2})u_{i}^{(2)}u_{j}^{(2)}+\\\
\qquad\
(2(g^{ii}g^{hj}+2g^{ij}g^{ih})u_{j}u_{h}+(B^{iij}+B^{iji}+B^{jii})u_{j}+C^{ii})u_{i}^{(2)}+\vphantom{\dfrac{1}{2}}\\\
\qquad\
g^{ij}g^{hk}u_{i}u_{j}u_{h}u_{k}+B^{ihk}u_{i}u_{h}u_{k}+C^{ij}u_{i}u_{j}+D^{i}u_{i},\end{array}$
(16)

where $u_{i}^{(s)}$ is the $s$-th partial derivative of $u$ w.r.t. $q^{i}$.

The function (16) is a 4-th degree polynomial in the derivatives of $u$ to the
4-th order; by observing that the terms $B^{ijh},C^{ij},D^{i}$ contains first,
second, third order derivatives of the metric tensor respectively, we can say
that $H_{s}$ is homogeneous: the sum of the degree and the order of the
derivatives it is the same (four) in all terms.

For the further computation it is useful to write down the first partial
derivatives of $H_{s}$

$\displaystyle\frac{\partial H_{s}}{\partial u^{(4)}_{i}}$ $\displaystyle=$
$\displaystyle(g^{ii})^{2}$ (17) $\displaystyle\frac{\partial H_{s}}{\partial
u^{(3)}_{i}}$ $\displaystyle=$ $\displaystyle B^{iii}+4g^{ii}g^{ij}u_{j}$ (18)
$\displaystyle\frac{\partial H_{s}}{\partial u^{(2)}_{i}}$ $\displaystyle=$
$\displaystyle
2(g^{ii}g^{jj}+2(g^{ij})^{2})u_{j}^{(2)}+2(g^{ii}g^{hj}+2g^{ij}g^{ih})u_{j}u_{h}+$
$\displaystyle(B^{iij}+B^{iji}+B^{jii})u_{j}+C^{ii}$
$\displaystyle\frac{\partial H_{s}}{\partial u_{i}}$ $\displaystyle=$
$\displaystyle
4g^{jj}g^{ij}u_{j}^{(3)}+u_{j}^{(2)}(B^{jji}+B^{jij}+B^{ijj}+4(g^{jj}g^{hi}+2g^{ji}g^{jh})u_{h})+$
$\displaystyle
D^{i}+(C^{ij}+C^{ji})u_{j}+(B^{ijh}+B^{jhi}+B^{hij})u_{j}u_{h}+4g^{ij}g^{hk}u_{j}u_{h}u_{k}$

We compute the vector fields associated with the separation, which are of the
form (6):

$D_{i}=\partial_{i}+u_{i}^{(2)}\frac{\partial}{\partial
u_{i}}+u_{i}^{(3)}\frac{\partial}{\partial
u_{i}^{(2)}}+u_{i}^{(4)}\frac{\partial}{\partial
u_{i}^{(3)}}+R_{i}\frac{\partial}{\partial u_{i}^{(4)}},$

where the index $i$ is not summed and the functions
$R_{i}(q^{h},u_{h},...u_{h}^{(4)})$ are determined by the condition

$D_{i}(H_{s})=0.$

Hence, under the technical assumption that $g^{ii}\neq 0$, we have

$R_{i}=-\left(\frac{\partial H_{s}}{\partial
u_{i}^{(4)}}\right)^{-1}\left(\partial_{i}H_{s}+u_{i}^{(2)}\frac{\partial
H_{s}}{\partial u_{i}}+u_{i}^{(3)}\frac{\partial H_{s}}{\partial
u_{i}^{(2)}}+u_{i}^{(4)}\frac{\partial H_{s}}{\partial u_{i}^{(3)}}\right)=$
$=-\frac{1}{(g^{ii})^{2}}\left(\partial_{i}H_{s}+u_{i}^{(2)}\frac{\partial
H_{s}}{\partial u_{i}}+u_{i}^{(3)}\frac{\partial H_{s}}{\partial
u_{i}^{(2)}}+u_{i}^{(4)}\frac{\partial H_{s}}{\partial u_{i}^{(3)}}\right)$

The condition for free (or regular) separation is that

$D_{i}R_{j}=0\qquad\forall i\neq j$

that is

$\partial_{i}R_{j}+u_{i}^{(2)}\frac{\partial R_{j}}{\partial
u_{i}}+u_{i}^{(3)}\frac{\partial R_{j}}{\partial
u_{i}^{(2)}}+u_{i}^{(4)}\frac{\partial R_{i}}{\partial
u_{i}^{(3)}}+R_{i}\frac{\partial R_{j}}{\partial u_{i}^{(4)}}=0.$

We look for terms containing $u_{h}^{(3)}$ or $u_{h}^{(4)}$ and such that the
sum of the orders of the derivatives of $u$ is six i.e., the terms containing
$u^{(4)}_{h}u^{(2)}_{k}$ or $u^{(3)}_{h}u^{(3)}_{k}$ (for $k\neq h$) and with
the components of the metric tensor not derived. Since we have for $i\neq j$

$\frac{\partial R_{j}}{\partial
u_{i}^{(4)}}=-\frac{1}{(g^{jj})^{2}}\left(\frac{\partial^{2}H_{s}}{\partial
q^{j}\partial u_{i}^{(4)}}+u_{j}^{(2)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(4)}\partial u_{j}}+u_{j}^{(3)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(4)}\partial u_{j}^{(2)}}+u_{j}^{(4)}\frac{\partial^{2}H_{s}}{\partial
u_{j}^{(3)}\partial u_{i}^{(4)}}\right)$
$=-\frac{\partial_{j}(g^{ii})^{2}}{(g^{jj})^{2}},$

the addendum $R_{i}\frac{\partial R_{j}}{\partial u_{i}^{(4)}}$ does not
contain such a term, as well as $\partial_{i}R_{j}$. Moreover, we have

$\frac{\partial R_{j}}{\partial
u_{i}^{(3)}}=-\frac{1}{(g^{jj})^{2}}\left(\frac{\partial^{2}H_{s}}{\partial
q^{j}\partial u_{i}^{(3)}}+u_{j}^{(2)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(3)}\partial u_{j}}+u_{j}^{(3)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(3)}\partial u_{j}^{(2)}}+u_{j}^{(4)}\frac{\partial^{2}H_{s}}{\partial
u_{j}^{(3)}\partial u_{i}^{(3)}}\right)$
$=-\frac{1}{(g^{jj})^{2}}\left(\partial_{j}B^{iii}+4\partial_{j}(g^{ii}g^{ih})u_{h}+4u_{j}^{(2)}g^{ii}g^{ij}\right)$

Thus, in the addendum $u^{(4)}_{i}\frac{\partial R_{j}}{\partial
u_{i}^{(3)}},$ the function $D_{i}R_{j}$ contains the term

$-4\frac{g^{ii}g^{ij}}{(g^{jj})^{2}}u_{j}^{(2)}u_{i}^{(4)}$ (21)

For the derivative of $R_{j}$ w.r.t. $u_{i}^{(2)}$ we have

$\frac{\partial R_{j}}{\partial
u_{i}^{(2)}}=-\frac{1}{(g^{jj})^{2}}\left(\frac{\partial^{2}H_{s}}{\partial
q^{j}\partial u_{i}^{(2)}}+u_{j}^{(2)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(2)}\partial u_{j}}+u_{j}^{(3)}\frac{\partial^{2}H_{s}}{\partial
u_{i}^{(2)}\partial u_{j}^{(2)}}+u_{j}^{(4)}\frac{\partial^{2}H_{s}}{\partial
u_{j}^{(3)}\partial u_{i}^{(2)}}\right)=$
$=-\frac{1}{(g^{jj})^{2}}\left(\frac{\partial^{2}H_{s}}{\partial q^{j}\partial
u_{i}^{(2)}}+u_{j}^{(2)}\frac{\partial^{2}H_{s}}{\partial u_{i}^{(2)}\partial
u_{j}}+2(g^{ii}g^{jj}+2(g^{ij})^{2})u_{j}^{(3)}\right).$

Therefore, in $u^{(3)}_{i}\frac{\partial R_{j}}{\partial u_{i}^{(2)}}$, the
function $D_{i}R_{j}$ contains the term

$\frac{-2}{(g^{jj})^{2}}(g^{ii}g^{jj}+2(g^{ij})^{2})u_{j}^{(3)}u_{i}^{(3)}$
(22)

and, since $\frac{\partial H_{s}}{\partial u_{i}^{(2)}}$ does not depend on
$u_{h}^{(3)}$ and $u_{h}^{(4)}$, we disregard other addenda of $\frac{\partial
R_{j}}{\partial u_{i}^{(2)}}$. The last term to consider is

$\frac{\partial R_{j}}{\partial
u_{i}}=-\frac{1}{(g^{jj})^{2}}\left(\frac{\partial^{2}H_{s}}{\partial
q^{j}\partial u_{i}}+u_{j}^{(2)}\frac{\partial^{2}H_{s}}{\partial
u_{i}\partial u_{j}}+u_{j}^{(3)}\frac{\partial^{2}H_{s}}{\partial
u_{i}\partial u_{j}^{(2)}}+u_{j}^{(4)}\frac{\partial^{2}H_{s}}{\partial
u_{j}^{(3)}\partial u_{i}}\right)$
$=-\frac{1}{(g^{jj})^{2}}\left(4g^{jj}g^{ij}u_{j}^{(4)}+\textrm{terms with
derivatives of $u$ of order less than $4^{th}$}\right)$

Thus $D_{i}R_{j}$ contains the term

$-4\frac{g^{ij}}{g^{jj}}u_{j}^{(4)}u_{i}^{(2)}$ (23)

and no other terms in $u_{h}^{(4)}u_{k}^{(2)}$ or $u_{h}^{(3)}u_{k}^{(3)}$ are
present in $D_{i}R_{j}$ except for (21,22,23). It is easy to see that the
three terms cannot be all zero since $g^{ii}\neq 0,$ for all $i$. Hence, we
have proved that

######  Theorem 1.

Regular multiplicative separation for the bi-Helmholtz equation
$(\ref{eq_bilap})$ on any Riemannian or pseudo-Riemannian $n$-dimensional
manifold does not occur in any system of coordinates.

######  Remark 5.

This means that we cannot find a family of separated solutions that depends on
$4n$ parameters and such that the values of

$u_{i},\qquad u_{i}^{(2)},\qquad u_{i}^{(3)},\qquad u_{i}^{(4)}$ (24)

can be fixed in an arbitrary way. Going back to the multiplicatively separated
solutions $\prod_{i}\psi_{i}$ of (7), the impossibility that regular
separation occurs can be interpreted as follows: the one-to-one relationship
existing between (24) and

$\frac{\psi_{i}^{\prime}}{\psi_{i}}\qquad\frac{\psi_{i}^{\prime\prime}}{\psi_{i}}\qquad\frac{\psi_{i}^{\prime\prime\prime}}{\psi_{i}}\qquad\frac{\psi_{i}^{\prime\prime\prime\prime}}{\psi_{i}}$

(see [1] for the explicit link till order 2) implies that it is not possible
to assign these values in an arbitrary way at a point $q_{0}$ and finding a
separated solution satisfying these initial conditions.

However, smaller families of separated solutions may exist. In this case the
theory of non-regular separation (see [3]) may be usefully applied.

######  Proposition 2.

In any coordinate system $(q^{i})$ allowing regular separation for the
Helmholtz equation, there exists a submanifold $N$ of $M$ with dimension at
least $2n$ where bi-Helmholtz admits non-regular separation.

The statement follows form the fact that every solution of the Helmholtz
equation is also a solution of the bi-Helmholtz equation. Hence, the
$2n$-parameter families of the multiplicatively separated solutions of the
Helmholtz equation are multiplicatively separated solution of (7). Non-regular
separation separation is said to be non-trivial if there exists a bigger
family of separable solutions for the bi-Helmholtz equation than that for the
Helmholtz equation. Otherwise, the separation is said to be trivial. A
detailed analysis of some examples shows that non-trivial non-regular
separation is possible.

##  4 Examples

In this section we provide some explicit computations of the manifolds $N$
allowing non-regular separation of variables for the bi-Helmholtz equation.
For sake of simplicity, we restrict ourselves to the Euclidean plane.
Moreover, in order to be sure that there exists a $N$ where non-regular
separation occurs, we shall consider only coordinates that allow separation
for the Helmholtz equation. Several situations will be described. In two of
the four possible separable coordinate systems (Cartesian and polar) we show
that non-trivial non-regular separation of the bi-Helmholtz equation is
possible (sections 4.1 and 4.2). One the other hand, if we use parabolic or
elliptic-hyperbolic coordinates, the only possibility is trivial separation.

###  4.1 Cartesian coordinates on the plane

Let us examine the case of bi-Helmholtz equation on the plane in Cartesian
coordinates $(q^{1},q^{2})=(x,y)$. The only separability condition
$D_{1}R_{2}=0$ becomes

$-(2u^{(2)}_{2}u_{2}+u_{2}^{(3)})(2u^{(2)}_{1}u_{1}+u_{1}^{(3)})=0$ (25)

which is satisfied if one of the factors vanishes. However, this condition
alone does not define a submanifold such that $D_{1}$ and $D_{2}$ are tangent
to it: for example

$D_{1}(2u^{(2)}_{1}u_{1}+u_{1}^{(3)})=2(u_{1}^{(2)})^{2}+2u_{1}^{(3)}u_{1}+u_{1}^{(4)}$
(26)

Thus this condition must also be added in order to define a manifold where the
non-regular separation could occur. Let us call $N$ the submanifold defined by
the equations

$f_{1}=u_{1}^{(3)}+2u^{(2)}_{1}u_{1}=0,\qquad\qquad
f_{2}=u_{1}^{(4)}+2(u_{1}^{(2)})^{2}+2u_{1}^{(3)}u_{1}=0$ (27)

It is easy to check that $D_{i}f_{a}|_{N}=0$ ($i=1,2$ $a=1,2$) and that
$D_{i}R_{j}|_{N}=0$. Hence, on the 6-dimensional manifold $N$ non-regular
separation occurs in Cartesian coordinates. The four derivatives of $S_{2}(y)$
($u_{2},u_{2}^{(2)},u_{2}^{(3)},u_{2}^{(4)}$) can be arbitrarily assigned at
an initial point, while only two of the derivatives of $S_{1}(x)$
($u_{1},u_{1}^{(2)}$) are free and the remaining ones are determined by the
equations $f_{a}$ defining $N$.

Since the equations for $N$ only involve the dependence of $u$ on $q^{1}=x$,
they can be seen as separated equations for the function $S_{1}(q^{1})$: in
particular $f_{1}$ and $f_{2}$ (its differential consequence) imply that
$S_{1}$ has to satisfy

$\frac{d^{2}}{dx^{2}}S_{1}+(\frac{d}{dx}S_{1})^{2}=c_{1},\qquad
c_{1}\in\mathbb{R}$ (28)

It is interesting to remark that (28) means precisely that $S_{1}$ is a
separated solution of the Laplace equation. Indeed, by inserting the function
$\psi_{1}(x)=e^{S_{1}}$ in (28), we get

$\psi_{1}^{\prime\prime}=c_{1}\psi_{1}\qquad(c_{1}\in\mathbb{R}).$

However, the family of the separated solutions is bigger than the separated
solutions of the Helmholtz equation, since for $S_{2}$ we can get functions
which are not solutions of the Helmholtz equation but only of the bi-Helmholtz
one. Indeed, by (28) the separated equation for $H_{s}$ written in Cartesian
coordinates takes the form

$u_{2}^{(4)}+3(u_{2}^{(2)})^{2}+4u_{2}u_{2}^{(3)}+6(u_{2})^{2}u_{2}^{(2)}+(u_{2})^{4}+2c_{1}\left(u_{2}^{(2)}+(u_{2})^{2}\right)+c_{1}^{2}-\lambda=0,$
(29)

where $c_{1}=u_{1}^{(2)}+(u_{1})^{2}$ which is constant on the constraint
surface $N$. A simpler form for (29) is given in terms of
$\psi_{2}=\ln(S_{2})$: with respect to this new independent variable the ODE
becomes the linear 4-th order equation

$(\psi_{2})^{(4)}+2c_{1}\psi_{2}^{\prime\prime}+(c_{1}^{2}-\lambda)\psi_{2}=0.$

###  4.2 Polar coordinates

Let us examine the case of bi-Helmholtz equation on the plane in polar
coordinates $(q^{1},q^{2})=(r,\theta)$. The situation is similar, but is no
longer symmetric in the two variables because we pass from two ignorable
coordinates to one ignorable coordinate.

$(u_{1}^{(3)}r^{2}+2r(u_{1})^{2}+u_{1}^{(2)}r-u_{1}+2u_{1}^{(2)}u_{1}r^{2})(2u^{(2)}_{2}u_{2}+u_{2}^{(3)})=0$
(30)

As in the Cartesian case, if we consider the condition on $\theta$-depending
function, on the six-dimensional submanifold $N_{pol}$ defined by

$f_{1}=2u^{(2)}_{2}u_{2}+u_{2}^{(3)}=0,\qquad
f_{2}=u_{2}^{(4)}+2(u_{2}^{(2)})^{2}+2u_{2}^{(3)}u_{2}=0$ (31)

reduced separation occurs. Conditions (31) mean that the following separated
equation for $S_{2}(\theta)$ holds

$\frac{d^{2}}{d\theta^{2}}S_{2}+(\frac{d}{d\theta}S_{2})^{2}=c_{2},\qquad
c_{2}\in\mathbb{R}$ (32)

i.e., that $\psi_{2}(\theta)=e^{S_{2}}$ satisfies

$\psi_{2}^{\prime\prime}=c_{2}\psi_{2}$ (33)

As in the Cartesian case plugging (32) into the separated equation $H_{s}$
provides a separated equation which involves the variable $r$ only:

$\displaystyle
r^{4}u_{1}^{(4)}+\left(4r^{4}u_{1}+2r^{3}\right)u_{1}^{(3)}+3r^{4}(u_{1}^{(2)})^{2}+r^{2}\left(6r^{2}(u_{1})^{2}+6ru_{1}+2-1\right)u_{1}^{(2)}+$
$\displaystyle\\!+r^{4}(u_{1})^{4}+2r^{3}(u_{1})^{3}+(2c_{2}-1)(r^{2}(u_{1})^{2}-ru_{1})-\lambda
r^{4}+c_{2}^{2}-4c_{2}=0$ (34)

here $c_{2}=u_{2}^{(2)}+(u_{2})^{2}$ which is constant on the constraint
surface $N_{pol}$. The above equation becomes simpler passing to
$\psi_{1}=e^{S_{1}}$

$r^{4}(\psi_{1})^{(4)}+2r^{3}(\psi_{1})^{(3)}+r(2c_{2}-1)(r\psi_{1}^{\prime\prime}-\psi_{1}^{\prime})+(c_{2}^{2}+4c_{2}-\lambda
r^{4})\psi_{1}=0$ (35)

whose solutions will be analyzed in Section 5.

The condition for the $r$-depending part of (30) is equivalent to

$r^{2}(\frac{d^{2}}{dr^{2}}S_{1}+(\frac{d}{dr}S_{1})^{2})-r\frac{d}{dr}S_{1}=c_{1},\qquad
c_{1}\in\mathbb{R}$

However, this condition and its differential conditions defines a submanifold
which is tangent to the generator $D_{1}$ only for $\lambda=0$.

###  4.3 Parabolic and Elliptic-Hyperbolic coordinates

We consider the general case of Liouville coordinates on $\mathbb{R}^{2}$
which produces Cartesian, polar, parabolic and elliptic-hyperbolic coordinates
as special cases when we demand that the Gaussian curvature vanishes. This
coordinate system is characterized by the following metric

$\displaystyle g=(f(u)+g(v))(du^{2}+dv^{2}),$ (36)

for some arbitrary smooth functions $f$ and $g$. We first consider the
Helmholtz equation

$\displaystyle\Delta\psi=\gamma\psi,$ (37)

where we recall that $\Delta$ denotes the Laplace-Beltrami operator (8). In
Liouville coordinates it takes the form

$\displaystyle\Delta\psi=\frac{1}{f+g}\left(\partial_{u}^{2}\psi+\partial_{v}^{2}\psi\right)$
(38)

If we assume $\psi(u,v)=U(u)V(v)$ the Helmholtz equation may be written as

$\displaystyle\left(\frac{U^{\prime\prime}}{U}-\gamma
f\right)+\left(\frac{V^{\prime\prime}}{V}-\gamma g\right)=0$ (39)

We observe that the equation separates for any smooth functions $f$ and $g$.

$\displaystyle U^{\prime\prime}-(c+\gamma f)U=0\hskip
56.9055ptV^{\prime\prime}+(c-\gamma g)V=0,$ (40)

where $c$ denotes the separation constant.

We now turn our attention to the bi-Helmoltz equation. The Gaussian curvature
in this coordinate system has the following form

$\displaystyle
K=-\frac{1}{2(f+g)^{2}}\left(f^{\prime\prime}+g^{\prime\prime}-\frac{f^{\prime
2}+g^{\prime 2}}{f+g}\right)$ (41)

From here on we will set $K=0$ in order to consider the parabolic and
elliptic-hyperbolic cases in which we also have that $f^{\prime}\neq
0,g^{\prime}\neq 0$. We use the Laplace-Beltrami (38) to write the bi-Laplace
operator as

$\displaystyle\Delta^{2}\psi=\frac{\Delta_{0}^{2}\psi}{(f+g)^{2}}-\frac{2}{(f+g)^{3}}(f^{\prime}\partial_{u}\Delta_{0}\psi+g^{\prime}\partial_{v}\Delta_{0}\psi)+\frac{f^{\prime\prime}+g^{\prime\prime}}{(f+g)^{3}}\Delta_{0}\psi$
(42)

From here we can write the bi-Helmholtz equation $\Delta^{2}\psi=\lambda\psi$,
using (42). The product ansatz $\psi(u,v)=U(u)V(v)$ implies that (7) takes the
form

$\displaystyle
f\frac{U^{(4)}}{U}+g\frac{V^{(4)}}{V}-2f^{\prime}\frac{U^{(3)}}{U}-2g^{\prime}\frac{V^{(3)}}{V}+f^{\prime\prime}\frac{U^{\prime\prime}}{U}+g^{\prime\prime}\frac{V^{\prime\prime}}{V}$
$\displaystyle+$
$\displaystyle\frac{V^{\prime\prime}}{V}\left(f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}\right)+\frac{U^{\prime\prime}}{U}\left(g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}\right)$
$\displaystyle-$ $\displaystyle\lambda(f+g)^{3}=0$ (43)

We observe that the above equation is not separable and hence does not admit
regular multiplicative separation. This result is consistent with the
conclusion of Theorem 1. However, non-regular separation may be possible. By
taking derivatives of (4.3) we can obtain necessary conditions for separation
(constraint equations). We find a first necessary condition by applying
$\partial^{2}_{u,v}$ which yields

$\displaystyle
f^{\prime}\left(\frac{V^{(4)}}{V}\right)^{\prime}+g^{\prime}\left(\frac{U^{(4)}}{U}\right)^{\prime}+\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\left(f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}\right)^{\prime}$
$\displaystyle+\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}\left(g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}\right)^{\prime}-6\lambda
f^{\prime}g^{\prime}(f+g)=0$ (44)

We again note that this condition is not separable. We proceed by dividing by
$f^{\prime}g^{\prime}$ and applying $\partial^{2}_{u,v}$ to get the separable
condition

$\displaystyle\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}\left(\frac{\left(f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}\right)^{\prime}}{f^{\prime}}\right)^{\prime}$
$\displaystyle+\left(\frac{1}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}\right)^{\prime}\left(\frac{\left(g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}\right)^{\prime}}{g^{\prime}}\right)^{\prime}=0.$
(45)

From here we obtain several cases: (i) neither
$\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}$
nor
$\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}$
vanish, from which we obtain $\lambda=0$ (ii) only one vanishes which also
implies $\lambda=0$ (iii) both vanish which gives us a solution to the
Helmholtz equation in Liouville coordinates. Hence, we do not get any non-
trivial separated solutions. The proof in case (i) is given below and in
Appendix A while that for cases (ii) and (iii) is given in Appendix B.

Case(i): Dividing (4.3) by
$\left(\frac{1}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}\right)^{\prime}$
$\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}$
we obtain a separable equation. Separating and integrating (see appendix A for
this calculation) we find the following equations

$\displaystyle(2g-C)\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}+g^{\prime\prime}-C_{1}g-C_{2}$
$\displaystyle=0,$ (46)
$\displaystyle(2f+C)\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}+f^{\prime\prime}-D_{1}f-D_{2}$
$\displaystyle=0,$ (47)

where $C,C_{1},C_{2},D_{1}.D_{2}\in\mathbb{R}$. We can use these conditions in
(4.3) to obtain

$\displaystyle\frac{1}{g^{\prime}}\left(\frac{V^{(4)}}{V}\right)^{\prime}+\frac{1}{f^{\prime}}\left(\frac{U^{(4)}}{U}\right)^{\prime}+\frac{D_{1}}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}+\frac{C_{1}}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}-6\lambda(f+g)=0.$
(48)

Notice that this condition is now separable. And using this condition we can
integrate back to get a simplified form of equation (4.3)

$\displaystyle
f\frac{U^{(4)}}{U}+g\frac{V^{(4)}}{V}-2f^{\prime}\frac{U^{(3)}}{U}-2g^{\prime}\frac{V^{(3)}}{V}+f^{\prime\prime}\frac{U^{\prime\prime}}{U}+g^{\prime\prime}\frac{V^{\prime\prime}}{V}-\lambda
f^{3}-\lambda g^{3}+$
$\displaystyle+f\frac{V^{(4)}}{V}+g\frac{U^{(4)}}{U}+(C_{1}g+C_{2})\frac{U^{\prime\prime}}{U}+(D_{1}f+D_{2})\frac{V^{\prime\prime}}{V}-3\lambda
f^{2}g-3\lambda fg^{2}=0.$ (49)

To make sure we account for all the constraints properly we also need to use
the derivatives of equations (46) and (47). There are other constraints coming
from the additional requirement $K=0$, for the interest of brevity we will
include these calculations Appendix A.

After using all the constraint equations we end up with the following

$\displaystyle\frac{(C_{2}-D_{2}-(C_{1}+D_{1})f)}{C_{1}+D_{1}}\left(\alpha
f+\beta+3\lambda
f^{2}-\frac{(k+C_{1})(D_{1}-k)f+C_{1}(D_{2}-D)+D(D_{1}-k)}{2f+C}\right.$
$\displaystyle+2\frac{f^{\prime
2}}{(2f+C)^{2}}(D_{1}-k)-\left.\frac{(D_{1}f+D_{2})^{2}-f^{\prime\prime
2}}{(2f+C)^{2}}\right)-2\frac{(D_{1}-k)}{2f+C}f^{\prime 2}$
$\displaystyle+\left((k-C_{1})f+C_{2}+D\right)\left(\frac{(D_{1}-k)f+D_{2}-D}{2f+C}\right)+2\lambda
f^{3}+\alpha f^{2}+(\beta+\gamma)f=\delta$ (50)

Where from the condition $K=0$ we have $f^{\prime 2}=kf^{2}+2Df-k$ and
$g^{\prime 2}=-kg^{2}+2Dg+k$. Thus the above equation simplifies to a
polynomial in $f$ after multiplication by $(2f+C)^{2}$. Furthermore, since
$f^{\prime}\neq 0$, the set $\\{1,f,f^{2},\ldots,f^{n}\\}$ is linearly
independent. This implies that the coefficients of the different powers of $f$
must all vanish. The coefficient of the highest power of $f$ the is $8\lambda$
the vanishing of which implies that $\lambda=0$. In the case that
$C_{1}+D_{1}=0$ (A) is a polynomial in $f$ with highest order term $3\lambda
f^{2}$ which also implies $\lambda=0$.

##  5  The circular vibrating plate

As mentioned in Section 3 one can pass between multiplicative and additive
separation by the change of variable $\psi=e^{u}$. In this section we use the
framework of multiplicative separation.

We apply the results for polar coordinates obtained in Section 4 to study
vibrations of a thin, solid plate of uniform material, constant thickness, and
areal density $\rho$ [11, 5].

The equation of motion for the plate under the assumption of small
oscillations is given by [5]

$\displaystyle\rho\ddot{\psi}+c\Delta^{2}\psi=0,$ (51)

where $\psi$ is the lateral elevation at any point on the region $\Omega$ in
the interior of the plate and $c$ is some positive constant depending on the
material. If the plate is clamped at the edge he boundary conditions on $\psi$
are $\psi\lvert_{\partial\Omega}=\psi_{n}\lvert_{\partial\Omega}=0$, where the
subscript $n$ denotes differentiation with respect to the normal to the
boundary. For the case of a circular plate $\Omega$ is a disk of some radius
$a$. We obtain the solution of (51) by utilizing a separation of variables
approach rather that the Fourier series method employed by Rayleigh [11].

Consider solutions of (51) of the form $\psi=w(r,\theta)T(t)$ which yield the
separated equations

$\displaystyle\frac{\Delta^{2}w}{w}=-\frac{\rho}{c}\frac{\ddot{T}}{T}=k^{4}.$
(52)

The separation constant is assumed to be positive (hence written as $k^{4}$
for convenience), since we expect oscillatory behaviour in time. Notice that
this separation is regular since there are no additional constraints on the
separated functions. The time part of (52) has the form

$\displaystyle\ddot{T}+\omega^{2}T=0,$ (53)

where $\omega^{2}=\frac{ck^{4}}{\rho}$. The general solution is given by

$\displaystyle T(t)=G\cos(\omega t)+H\sin(\omega t),$ (54)

where $G$ and $H$ are arbitrary constants.

The spatial part (52) satisfies the bi-Helmholtz equation,

$\displaystyle\Delta^{2}w=k^{4}w$ (55)

Written out explicitly in polar coordinates this equation takes the form

$\displaystyle
w_{rrrr}+\frac{1}{r^{4}}w_{\theta\theta\theta\theta}+\frac{2}{r}w_{rrr}-\frac{1}{r^{2}}w_{rr}+\frac{1}{r^{3}}w_{r}+\frac{4}{r^{4}}w_{\theta\theta}-\frac{2}{r^{3}}w_{r\theta\theta}+\frac{2}{r^{2}}w_{rr\theta\theta}=k^{4}w,$
(56)

where the subscript $r$ denotes $\frac{d}{dr}$. We now consider product
solutions of the form

$\displaystyle w(r,\theta)=R(r)\Theta(\theta),$ (57)

which yields the separated equations

$\displaystyle\Theta^{\prime\prime}=-l^{2}\Theta,$ (58) $\displaystyle
D^{4}R+\frac{2}{r}D^{3}R-\frac{1}{r^{2}}(1+2l^{2})D^{2}R+\frac{1}{r^{3}}(1+2l^{2})DR-\frac{l^{2}}{r^{4}}(4-l^{2})R=k^{4}R,$
(59)

where $l$ is a separation constant and $D=\frac{d}{dr}$.

The solution of the differential equation (58) satisfied by the angular
function $\Theta$ has the form

$\displaystyle\Theta_{l}(\theta)=E_{l}\cos(l\theta)+F_{l}\sin(l\theta),$ (60)

where $E_{l}$ and $F_{l}$ are arbitrary constants. In order that the function
$\Theta$ be single valued on $\Omega$, the constant $l$ must be a positive
integer

$l=n\quad(n=1,2,\ldots).$

We now turn our attention to (59). It may be shown that the left-hand-side may
be factored in two ways to yield

$\displaystyle\left(D^{2}+\frac{1}{r}D-\frac{n^{2}}{r^{2}}+k^{2}\right)\left(D^{2}+\frac{1}{r}D-\frac{n^{2}}{r^{2}}-k^{2}\right)R$
$\displaystyle=0,$ (61)
$\displaystyle\left(D^{2}+\frac{1}{r}D-\frac{n^{2}}{r^{2}}-k^{2}\right)\left(D^{2}+\frac{1}{r}D-\frac{n^{2}}{r^{2}}+k^{2}\right)R$
$\displaystyle=0.$ (62)

This result shows that the differential operators which appear inside the
brackets (the Bessel and modified Bessel operators) commute. Since we know
bases of the solution spaces to the Bessel and modified Bessel equations, the
general solution of the radial equation may be written as

$\displaystyle R_{n}(r)$
$\displaystyle=A_{n}J_{n}(kr)+B_{n}Y_{n}(kr)+C_{n}I_{n}(kr)+D_{n}K_{n}(kr),$
(63)

where $J_{n},Y_{n},I_{n},K_{n}$ are respectively the Bessel and modified
Bessel functions of the first and second kind and $A_{n},B_{n},C_{n},D_{n}$
are arbitrary constants.

The Bessel functions of the second kind are inadmissible since they are
singular at the origin so we take $B_{n}=D_{n}=0$.

The boundary conditions tell us that $R_{n}(a)=R_{n}^{\prime}(a)=0$.

$\displaystyle R_{n}(a)$ $\displaystyle=A_{n}J_{n}(ka)+C_{n}I_{n}(ka)=0$ (64)
$\displaystyle R_{n}^{\prime}(a)$
$\displaystyle=kA_{n}J_{n}^{\prime}(ka)+kC_{n}I_{n}^{\prime}(ka)=0$ (65)

To solve (64) we either need $ka=j_{n,m}$ (the $m$th zero of $J_{n}$) and
$C_{n}=0$ in which case we retrieve the vibrating membrane solution, or
$A_{n}=-\frac{I_{n}(ka)}{J_{n}(ka)}C_{n}$ (we can absorb $C_{n}$ into the
definition of $E_{n},F_{n}$ to simplify). Since we are interested in solutions
other than the ones for the vibrating membrane we take the latter in which
case (65) becomes

$\displaystyle
I_{n}^{\prime}(ka)-\frac{I_{n}(ka)}{J_{n}(ka)}J_{n}^{\prime}(ka)$
$\displaystyle=0$ (66)
$\displaystyle\frac{I_{n}^{\prime}(ka)}{I_{n}(ka)}-\frac{J_{n}^{\prime}(ka)}{J_{n}(ka)}$
$\displaystyle=0$ (67)

We determine the roots of the above equation numerically to get a condition
$ka=l_{n,m},$ where $l_{n,m}$ is the $m$th root of equation (67) for some
fixed $m$

$\displaystyle
R_{n,m}(r)=C_{n}\left(I_{n}\left(\frac{l_{n,m}r}{a}\right)-\frac{I_{n}(l_{n,m})}{J_{n}(l_{n,m})}J_{n}\left(\frac{l_{n,m}r}{a}\right)\right)$
(68)

Thus the general solution for $u$ has the form

$\displaystyle\psi=\sum\limits_{m,n}\left(E_{n}\cos(n\theta)+F_{n}\sin(n\theta)\right)\left(G_{n,m}\cos(\omega_{n,m}t)+H_{n,m}\sin(\omega_{n,m}t)\right)R_{n,m}$
(69)

Where $\omega_{n,m}=c^{2}\frac{l_{n,m}^{2}}{a^{2}}$. See Rayleigh [11] for
further study of this solution.

##  6 Conclusion

The main idea of the paper is to apply the technique of the regular and non-
regular separation to the search of the multiplicative separated solutions of
the bi-Helmholtz equation, which is a classical 4-th order PDE equation of
Mathematical Physics which was solved in particular cases by means of this
technique ([11]). Nevertheless, deeper studies about separability usually deal
with first and second order PDE’s (as well as the Hamilton-Jacobi equation and
the stationary Schrödinger equation).

The choice of the bi-Helmholtz equation relies on two aspects: the existence
of a physical application (it is not a simple toy model) and the fact that
this kind of equation seemed to provide a good example for a deeper
understanding of the non-regular separation method.

Non-regular separation appears naturally in the study of separability of the
bi-Helmholtz equation, since we prove (Section 3) that regular separation
never occurs for this equation, but on the other hand we already know the
existence of a family of separated solutions (the trivial ones, arising from
the separability of the standard Helmholtz equation). In the examples on the
Euclidean plane, we show that, in two of the four separable coordinate systems
(Cartesian and polar coordinates), non-trivial separated solutions can be
determined, while in the remaining coordinate systems (parabolic and elliptic-
hyperbolic coordinates) the only possibles separated solutions are proved to
be the trivial ones.

The study of the separated solutions is done from the simplest geometrical
view point: the search for a submanifold, as large as possible, where the
separability conditions are satisfied, without assuming any particular
structure for the separated solution (such as as side conditions or
generalizations of Stäckel matrices as is done in [7]).

This paper represents the first step in a program to analyse interesting
examples of non-regular separation As further extensions of the present paper,
one could consider examples in higher dimensions or in different Riemannian
manifolds, to get more physical applications. Indeed, a wider collection of
examples is the natural starting point of a study to gain a deeper insight
into the geometric conditions (for instance in terms of symmetry operators) of
the non-regular separation of the bi-Helmholtz equation, as has been done for
the fixed energy R-separation for the Schrödinger equation [4, 7].

## Acknowledgements

R. G. McLenaghan wishes to acknowledge financial support from the Natural
Sciences and Engineering Research Council of Canada in the form of and a
Discovery Grant. He would also like to thank the Dipartimento di Matematica of
University of Torino (Italy) for hospitality during several visits when part
of this research was undertaken. C. M. Chanu acknowledges the contribution of
the local research project Metodi Geometrici in Fisica Matematica e
Applicazioni (2019) of Dipartimento di Matematica of University of Torino.

## Appendix A

Dividing (4.3) by
$\left(\frac{U^{(3)}}{f^{\prime}U}-\frac{U^{\prime\prime}U^{\prime}}{g^{\prime}U^{2}}\right)^{\prime}\left(\frac{V^{(3)}}{g^{\prime}V}-\frac{V^{\prime\prime}V^{\prime}}{g^{\prime}V^{2}}\right)^{\prime}$,
we see that we can separate as follows

$\displaystyle\frac{\left(\frac{\left(g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}\right)^{\prime}}{g^{\prime}}\right)^{\prime}}{\left(\frac{V^{(3)}}{g^{\prime}V}-\frac{V^{\prime\prime}V^{\prime}}{g^{\prime}V^{2}}\right)^{\prime}}=-\frac{\left(\frac{\left(f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}\right)^{\prime}}{f^{\prime}}\right)^{\prime}}{\left(\frac{U^{(3)}}{f^{\prime}U}-\frac{U^{\prime\prime}U^{\prime}}{f^{\prime}U^{2}}\right)^{\prime}}=C$
(70)

Where $C\in\mathbb{R}$. Separating and integrating

$\displaystyle\frac{\left(g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}\right)^{\prime}}{g^{\prime}}$
$\displaystyle=\frac{C}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}+C_{1}$
(71)
$\displaystyle\frac{\left(f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}\right)^{\prime}}{f^{\prime}}$
$\displaystyle=-\frac{C}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}+D_{1}$
(72)

The case where one of
$\left(\frac{U^{(3)}}{f^{\prime}U}-\frac{U^{\prime\prime}U^{\prime}}{g^{\prime}U^{2}}\right)^{\prime},\,\left(\frac{V^{(3)}}{g^{\prime}V}-\frac{V^{\prime\prime}V^{\prime}}{g^{\prime}V^{2}}\right)^{\prime}$
vanish corresponds to setting $C=0$ in one of the above equations. Multiplying
through by $f^{\prime}$ and $g^{\prime}$ respectively and integrating once
again we obtain

$\displaystyle
g^{\prime\prime}+2g\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}$
$\displaystyle=C\frac{V^{\prime\prime}}{V}+C_{1}g+C_{2}$ (73) $\displaystyle
f^{\prime\prime}+2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}$
$\displaystyle=-C\frac{U^{\prime\prime}}{U}+D_{1}f+D_{2}$ (74)

Where $C_{1},C_{2},D_{1}.D_{2}\in\mathbb{R}$. Or more compactly

$\displaystyle(2g-C)\frac{V^{\prime\prime}}{V}-2g^{\prime}\frac{V^{\prime}}{V}+g^{\prime\prime}-C_{1}g-C_{2}$
$\displaystyle=0$ (75)
$\displaystyle(2f+C)\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}+f^{\prime\prime}-D_{1}f-D_{2}$
$\displaystyle=0$ (76)

We can now separate (48)

$\displaystyle\frac{U^{(4)}}{U}+C_{1}\frac{U^{\prime\prime}}{U}-3\lambda
f^{2}$ $\displaystyle=\alpha f+\beta$ (77)
$\displaystyle\frac{V^{(4)}}{V}+D_{1}\frac{V^{\prime\prime}}{V}-3\lambda
g^{2}$ $\displaystyle=-\alpha g+\gamma$ (78)

The condition $K=0$ separates into

$\displaystyle\frac{f^{(3)}}{f^{\prime}}=-\frac{g^{(3)}}{g^{\prime}}=k$ (79)

Where $k\in\mathbb{R}$ is a separation constant. Integrating once

$\displaystyle f^{\prime\prime}=kf+D\hskip 28.45274ptg^{\prime\prime}=-kg+D$
(80)

These are extra conditions that have to be taken into account. We also need to
account for the derivatives of equation (75) and (76)

$\displaystyle\frac{U^{\prime\prime}}{U}$
$\displaystyle=\frac{2f^{\prime}}{2f+C}\frac{U^{\prime}}{U}+\frac{D_{1}f+D_{2}-f^{\prime\prime}}{2f+C}$
(81) $\displaystyle\frac{U^{(3)}}{U}$
$\displaystyle=\frac{D_{1}f+D_{2}+f^{\prime\prime}}{2f+C}\frac{U^{\prime}}{U}+\frac{D_{1}-k}{2f+C}f^{\prime}$
(82) $\displaystyle\frac{U^{(4)}}{U}$
$\displaystyle=\frac{D_{1}}{2f+C}f^{\prime\prime}+\frac{D_{1}f+D_{2}+f^{\prime\prime}}{2f+C}\frac{U^{\prime\prime}}{U}-\frac{2f^{\prime}}{2f+C}\frac{U^{(3)}}{U}$
(83)
$\displaystyle=\frac{(D_{1}-k)f^{\prime\prime}}{2f+C}+\frac{(D_{1}f+D_{2})^{2}-f^{\prime\prime
2}}{(2f+C)^{2}}-\frac{2f^{\prime 2}}{(2f+C)^{2}}(D_{1}-k)$ (84)
$\displaystyle\frac{V^{\prime\prime}}{V}$
$\displaystyle=\frac{2g^{\prime}}{2g-C}\frac{V^{\prime}}{V}+\frac{C_{1}g+C_{2}+g^{\prime\prime}}{2g-C}$
(85) $\displaystyle\frac{V^{(3)}}{V}$
$\displaystyle=\frac{C_{1}g+C_{2}+g^{\prime\prime}}{2g-C}\frac{V^{\prime}}{V}+\frac{C_{1}+k}{2g-C}g^{\prime}$
(86) $\displaystyle\frac{V^{(4)}}{V}$
$\displaystyle=\frac{C_{1}}{2g-C}g^{\prime\prime}+\frac{C_{1}g+C_{2}-g^{\prime\prime}}{2g-C}\frac{V^{\prime\prime}}{V}-\frac{2g^{\prime}}{2g-C}\frac{V^{(3)}}{V}$
(87)
$\displaystyle=\frac{(C_{1}+k)g^{\prime\prime}}{2g+C}+\frac{(C_{1}g+C_{2})^{2}-g^{\prime\prime
2}}{(2g-C)^{2}}-\frac{2g^{\prime 2}}{(2g-C)^{2}}(C_{1}+k)$ (88)

Using (81) to eliminate higher derivatives in (77)

$\displaystyle(C_{1}+D_{1})\frac{2f^{\prime}}{2f+C}\frac{U^{\prime}}{U}-2\frac{f^{\prime
2}}{(2f+C)^{2}}(D_{1}-k)+(D_{1}-k)\frac{f^{\prime\prime}}{2f+C}+\frac{(D_{1}f+D_{2})^{2}-f^{\prime\prime
2}}{(2f+C)^{2}}+$ $\displaystyle+$ $\displaystyle
C_{1}\frac{D_{1}f+D_{2}-f^{\prime\prime}}{2f+C}-3\lambda f^{2}=\alpha f+\beta$
(89)
$\displaystyle(C_{1}+D_{1})\frac{2g^{\prime}}{2g-C}\frac{V^{\prime}}{V}-2\frac{g^{\prime
2}}{(2g-C)^{2}}(C_{1}+k)+(C_{1}+k)\frac{g^{\prime\prime}}{2g-C}+\frac{(C_{1}g+C_{2})^{2}-g^{\prime\prime
2}}{(2g-C)^{2}}+$ $\displaystyle+$ $\displaystyle
D_{1}\frac{C_{1}g+C_{2}-g^{\prime\prime}}{2g-C}-3\lambda g^{2}=-\alpha
g+\gamma$ (90)

Substituting (77) into (4.3) we have the following

$\displaystyle-$ $\displaystyle
2f^{\prime}\frac{U^{(3)}}{U}-2g^{\prime}\frac{V^{(3)}}{V}+\left((k-C_{1})f+C_{2}+D\right)\frac{U^{\prime\prime}}{U}+\left(D_{2}+D-(k+D_{1})g\right)\frac{V^{\prime\prime}}{V}$
$\displaystyle+$ $\displaystyle
2\lambda(f^{3}+g^{3})+\alpha(f^{2}-g^{2})+(\beta+\gamma)(f+g)=0$ (91)

Separating this equation we have

$\displaystyle-2f^{\prime}\frac{U^{(3)}}{U}+\left((k-C_{1})f+C_{2}+D\right)\frac{U^{\prime\prime}}{U}+2\lambda
f^{3}+\alpha f^{2}+(\beta+\gamma)f$ $\displaystyle=\delta$ (92)
$\displaystyle-2g^{\prime}\frac{V^{(3)}}{V}+\left(-(k+D_{1})g+D_{2}+D\right)\frac{V^{\prime\prime}}{V}+2\lambda
g^{3}-\alpha g^{2}+(\beta+\gamma)g$ $\displaystyle=-\delta$ (93)

For some $\delta\in\mathbb{R}$. Eliminating the derivatives from
$(\ref{E100})$ we have the following

$\displaystyle
2f^{\prime}\frac{C_{2}-D_{2}-(C_{1}+D_{1})f}{2f+C}\frac{U^{\prime}}{U}-2\frac{(D_{1}-k)}{2f+C}f^{\prime
2}$ (94)
$\displaystyle+\left((k-C_{1})f+C_{2}+D\right)\left(\frac{(D_{1}-k)f+D_{2}-D}{2f+C}\right)$
$\displaystyle+$ $\displaystyle 2\lambda f^{3}+\alpha
f^{2}+(\beta+\gamma)f=\delta$ (95)

Isolating for $\frac{2f^{\prime}}{2f+C}\frac{U^{\prime}}{U}$ from
$(\ref{FODE1})$ assuming $C_{1}+D_{1}\neq 0$ and using this to eliminate
derivatives in (94)

$\displaystyle\frac{(C_{2}-D_{2}-(C_{1}+D_{1})f)}{C_{1}+D_{1}}\left(\alpha
f+\beta+3\lambda
f^{2}-\frac{(k+C_{1})(D_{1}-k)f+C_{1}(D_{2}-D)+D(D_{1}-k)}{2f+C}\right.$
$\displaystyle+2\frac{f^{\prime
2}}{(2f+C)^{2}}(D_{1}-k)-\left.\frac{(D_{1}f+D_{2})^{2}-f^{\prime\prime
2}}{(2f+C)^{2}}\right)-2\frac{(D_{1}-k)}{2f+C}f^{\prime 2}$
$\displaystyle+\left((k-C_{1})f+C_{2}+D\right)\left(\frac{(D_{1}-k)f+D_{2}-D}{2f+C}\right)+2\lambda
f^{3}+\alpha f^{2}+(\beta+\gamma)f=\delta$ (96)

When $C_{1}+D_{1}=0$, (A) implies $\lambda=0$.

## Appendix B

Case (ii): In this case we have

$\displaystyle\left(\frac{1}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}\right)^{\prime}=0,\quad\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}\neq
0.$ (97)

Integrating we get

$\displaystyle\frac{U^{\prime\prime}}{U}=\alpha_{5}f+\alpha_{6},$ (98)

where $\alpha_{5},\alpha_{6}$ are constants. Equation (4.3) implies that

$\displaystyle
2f\frac{U^{\prime\prime}}{U}-2f^{\prime}\frac{U^{\prime}}{U}+f^{\prime\prime}=\alpha_{3}f+\alpha_{4}$
(99)

$\alpha_{3},\alpha_{4}$ are constants. Substitution for
$\frac{U^{\prime\prime}}{U}$ from (98) yields

$\displaystyle
2f^{\prime}U^{\prime}=(f^{\prime\prime}+2\alpha_{5}f^{2}+(2\alpha_{6}-\alpha_{3})f-\alpha_{4})U.$
(100)

Differentiation of the above equation followed substitution for
$U^{\prime\prime}$ from (98) and $2f^{\prime}U^{\prime}$ from (100) yields
after simplification

$\displaystyle(2\alpha_{5}f^{2}+(2\alpha_{6}-\alpha_{3})f-\alpha_{4}^{2})^{2}-f^{\prime\prime
2}+2f^{\prime}f^{(3)}+4\alpha_{5}ff^{\prime 2}-\alpha_{3}f^{\prime 2}=0$ (101)

Using the relations between $f$ and its derivatives

$\displaystyle f^{\prime\prime}=kf+D,\quad f^{\prime 2}=kf^{2}+2Df+\Lambda$
(102)

Equation (101) becomes a polynomial in $f$, the coefficient of the highest
power of $f$ is $4\alpha_{5}^{2}$ which implies that $\alpha_{5}=0$. In view
of the above equation (4.3) separates, the compatibility of the separated
equation for $U$ with (98) and (100) gives us that $\lambda=\alpha_{5}^{2}$,
thus we conclude that $\lambda=0$ in this case as well.

Case (iii): both of the following conditions hold:

$\displaystyle\left(\frac{1}{f^{\prime}}\left(\frac{U^{\prime\prime}}{U}\right)^{\prime}\right)^{\prime}=0.$
(103)
$\displaystyle\left(\frac{1}{g^{\prime}}\left(\frac{V^{\prime\prime}}{V}\right)^{\prime}\right)^{\prime}=0.$
(104)

The solutions of (103) and (104) are given by

$\displaystyle\frac{U^{\prime\prime}}{U}=\alpha_{5}f+\alpha_{6}$ (105)
$\displaystyle\frac{V^{\prime\prime}}{V}=(\beta_{5}f+\beta_{6})$ (106)

Computing derivatives of the above equations:

$\displaystyle
U^{(3)}=(\alpha_{5}f+\alpha_{6})U^{\prime}+\alpha_{5}f^{\prime}U$ (107)
$\displaystyle
U^{(4)}=2\alpha_{5}f^{\prime}U^{\prime}+(\alpha_{5}f^{\prime\prime}+(\alpha_{5}f+\alpha_{6})^{2})U$
(108) $\displaystyle
V^{(3)}=(\beta_{5}g+\beta_{6})V^{\prime}+\beta_{5}g^{\prime}V$ (109)
$\displaystyle
V^{(4)}=2\beta_{5}g^{\prime}V^{\prime}+(\beta_{5}g^{\prime\prime}+(\beta_{5}g+\beta_{6})^{2})V$
(110)

With the use of the above derivatives the integrability condition (4.3)
separates to yield the following equations:

$\displaystyle
2(\beta_{5}-\alpha_{5})f^{\prime}U^{\prime}=((\alpha_{5}+\beta_{5})f^{\prime\prime}+\alpha_{5}(\alpha_{5}+2\beta_{5})f^{2}$
(111) $\displaystyle+(2\alpha_{6}(\alpha_{5}+\beta_{5})-\alpha)f-3\lambda
f^{2}+\alpha_{6}^{2}-\alpha_{7})U,$ (112) $\displaystyle
2(\alpha_{5}-\beta_{5})g^{\prime}V^{\prime}=((\alpha_{5}+\beta_{5})g^{\prime\prime}+\beta_{5}(\alpha_{5}+2\beta_{5})g^{2}$
(113) $\displaystyle+(2\beta_{6}(\alpha_{5}+\beta_{5})+\alpha)g-3\lambda
g^{2}+\alpha_{6}^{2}-\alpha_{7})V,$ (114)

where $\alpha$ is the separation constant. If $\beta_{5}=\alpha_{5}$, (111)
and (113) imply that

$\displaystyle\lambda=\alpha_{5}^{2},\quad\beta_{5}=\alpha_{5},\quad\beta_{6}=-\alpha_{6}$
(115)

We conclude that (115) implies that $\phi(u,v)=U(u)V(v)$ defines a separable
solution of the Helmholtz equation (See (40).). If $\beta_{5}\neq\alpha_{5}$,
one differentiates (111) and (113) and eliminates all derivatives of $U$ and
$V$. One obtains polynomial equations in $f$ and $g$ which imply that
$\beta_{5}^{2}=\alpha_{5}^{2}$. The case $\beta_{5}=\alpha_{5}$ has already
been considered. The case $\beta_{5}=-\alpha_{5}$, yields
$3\lambda=-\alpha_{5}^{2}$, which is un-physical. This completes the proof of
Case (iii).

## References

  * [1] S. Benenti, C. Chanu, and G. Rastelli ”Remarks on the connection between the additive separation of the Hamilton–Jacobi equation and the multiplicative separation of the Schr ödinger equation”. I. The completeness and Robert-son conditions”. In: Journal of Mathematical Physics 43.11 (2002), pp. 5183–5222.doi:10.1063/1.1506180.
  * [2] M. Chanachowicz and C. M. Chanu and R. G. McLenaghan ”$R$-separation of variables for the conformally invariant Laplace-Beltrami equation” Journal of Geometry and Physics 59 (2009) 876-884.
  * [3] C. M. Chanu “Geometry of non-regular separation”. In: Symmetries and Overdetermined Systems of Partial Differential Equations. Ed. by M. East-wood and W. Miller Jr. The IMA Volumes in Mathematics and its Applications vols. 144. Springer, 2008.
  * [4] C. Chanu, L. and G. Rastelli “Fixed Energy $R$-separation for the Scrödinger equation”. In: International Journal of Geometric Methods in Modern Physics 3(3) (2006), 489-508.
  * [5] G. F. D. Duff and D. Naylor. Differential Equations of Applied Mathematics.John Wiley & Sons, 1966.
  * [6] E. G Kalnins and W. Miller Jr. Intrinsic characterization of variable separation for the partial differential equations of mechanics. In: Proceeding of IUTAM-ISIMM Symposium on Modern Developments in Analytical Mechanics. 172. Atti Accad. Sci. Torino, 1983.
  * [7] E. G. Kalnins, J. M. Kress and W. Miller Jr. Separation of Variables and Superintegrability. IOP, 2018.isbn: 978-7503-1314-8.doi:10.1088/978-0-7503-1314-8.
  * [8] W. Miller Jr. Symmetries and Separation of Variables, Addison-Wesley, 1977.
  * [9] P. Moon and D. E. Spencer Field Theory Handbook Springer-Verlag, Berlin, 1961.
  * [10] P. M. Morse and H. Feshbach Methods of Theoretical Physics, Vols. 1 and 2, McGraw-Hill, 1953 (cit.on pp. 1).
  * [11] J. W. Strutt, Baron Raleigh The Theory of Sound. second ed. Vol. 1. Dover, 1945.


FPN short=FPN, long=feature pyramid network PAN short=PAN, long=path
aggregation network NMS short=NMS, long=non-maximum suppression GSI short=GSI,
long=Gaussian-smoothed interpolation FPs short=FPs, long=false positives FNs
short=FNs, long=false negatives IDs short=IDs, long=ID switches

# PieTrack: An MOT solution based on synthetic data training and self-
supervised domain adaptation

Yirui Wang  Shenghua He  Youbao Tang  Jingyu Chen  Honghao Zhou  Sanliang Hong
Junjie Liang  Yanxin Huang  Ning Zhang  Ruei-Sung Lin  Mei Han
PAII Inc. & Ping An Bank Co  Ltd.

###### Abstract

In order to cope with the increasing demand for labeling data and privacy
issues with human detection, synthetic data has been used as a substitute and
showing promising results in human detection and tracking tasks. We
participate in the 7th Workshop on Benchmarking Multi-Target Tracking (BMTT),
themed on ”How Far Can Synthetic Data Take us”? Our solution, PieTrack, is
developed based on synthetic data without using any pre-trained weights. We
propose a self-supervised domain adaptation method that enables mitigating the
domain shift issue between the synthetic (, _MOTSynth_) and real data (,
_MOT17_) without involving extra human labels. By leveraging the proposed
multi-scale ensemble inference, we achieved a final HOTA score of 58.7 on the
MOT17 testing set, ranked third place in the challenge.

##  1 Introduction

Multiple object tracking (MOT) is an active research area with profound
applications. With more advanced and bigger deep learning-based models, the
need of labeled real data has also been increased greatly. This is especially
the case in human-centered applications such as detection, tracking, action
recognition and spatial-temporal localization, etc. Although labeling people
in clear background with large-scale bounding-box is not difficult, precising
labelling of people in small-scale with occlusion is not easy and could be
error-prone. Meanwhile, there has been an increasing awareness of the privacy
protection and avoid exposing real human faces for deep learning applications.
This is such an important and ethical precaution but difficult in practice due
to the face obfuscation (blurring, mosaicing, etc.) still far from practical
in most human detection-related applications [27].

On the other hand, the synthetic data has served as a good alternative in
certain applications such as autonomous-driven vehicles [14], with virtual
driving datasets such as virtual KITTI [10], GTA5[21], VIPER [20]. Part of the
success with the virtual driving datasets is due to the rigid body of cars
which does not have flexible limbs like human being.

Last year (2021), thanks to the advanced game engine with physics simulation
as well as human body simulation, MOTSynth, a large-scale synthetic data has
been released for the multiple object tracking [8]. This dataset benefits the
community in multiple areas of research and applications including pedestrian
detection, tracking, re-identification, and instance segmentation. Such a
synthetic data drew our attention because of its advantages in both labeling
effort and privacy preservation. In particular, we are delighted to be part of
the 7th BMTT Challenges “MOTSynth-MOT-CVPR22”, and our participating approach
has the following highlights:

  * –

A two-stage detection and association strategy is adopted, where our detector
is trained using synthetic data from scratch without pre-trained weights.

  * –

Self-supervised domain adaptation approach is used to fill the gap between the
synthetic trained detector and the real MOT data, by using only unlabeled
MOT17 training split.

  * –

Multi-scale ensemble inference by taking different resolutions as input is
used. This strategy helps to further boost the performance.

##  2 Related Work

Tracking-by-Detection. Powerful detection capability thanks to the deep-
learning models, enables better tracking performance in the MOT domain. One-
stage approaches [1, 29, 18], tackle the detection and the tracking tasks
together in an end-to-end fashion and achieves high-speed process. However,
this group of approaches suffers a lower performance than the two-stage
because of its compromise between detection and tracking branches of the deep
learning model. Recent developments of YOLO-based approaches [19, 3, 11]
provides high-speed, high-accuracy detectors. These advances facilitate the
two-stage approaches by striking a balance between the detection speed and the
data association [28].

Synthetic Training Data. Autonomous driving applications have adopted
synthetic data in its model training since early time. Virtual KITTI [10]
cloned the real data to virtual proxy for detection validation. Omni-MOT [23]
provided a dataset, focuses on car-based detection and tracking. Parallel
Domain is a commercial content provider for MOT of cars [24]. In the people-
centered MOT domain, JTA (Joint Track Auto) [9] proposed a dataset on CG human
poses for estimation and tracking. MOTSynth [8] is the most comprehensive
dataset for pedestrian-based MOT and have drawn a lot attention since its
release. We are also drawn to this dataset and participated the Challenge
using the MOTSynth. Next section presents our methodology for the challenge
with detailed explanation and illustration.

##  3 Method

In this work, we adopt the tracking by detection paradigm and propose an
iterative approach to alleviate the domain gap between the synthetic and the
real data. To follow the competition rule that no human labels for the MOT17
dataset can be used for fine-tuning, we propose a self-supervised pseudo-label
mining and refining procedure on the MOT17 training set to shrink the domain
gap. And we also adopt an multi-scale ensemble strategy to improve the
robustness during the inference. The overview of the training framework is
shown in Fig. 1.

###  3.1 Multi-Object Detection and Tracking

Object Detection Network. The proposed PieTrack employs YOLOX [12] detection
network to detect object of interests (, pedestrian) in each frame. The
network is built upon YOLOv3 [19] architecture with modified CSPNet [25]
backbone and PAN [15] head to enhance gradient propagation and feature
hierarchy. YOLOX decouples the classification and regression tasks into
separate branches as they are prone to conflict with each other [22, 26]. The
detector adopts an anchor-free design to avoid domain-specific anchor
configuration and reduce computational cost during inference [12]. To
alleviate the imbalance between positive and negative samples, YOLOX consider
each $3\times 3$ area around the true object centers as positive samples.
Moreover, a simplified OTA [11] strategy, named SimOTA, is proposed to improve
the training efficiency without sacrificing the performance. The network can
be customized to different sizes and trained in an end-to-end manner. In our
work, we adopt YOLOX-X as our detection network and we follow [12] to train
the detector from scratch without using any pretrained weights.

Data Association and Tracking. We adopt a two-stage data association method
proposed in [28] with modifications to partially mitigate the domain shift
issue. In the first stage, the association is performed between the high
confidence detection boxes and all tracklets. In the second stage, the
association is further performed between the unmatches less confident
detection boxes and the remain tracklets. Different from the original strategy
in [28] that excludes lost tracklets from the second association, we propose
to consider both tracked and lost tracklets for sensitivity consideration.
Because of domain shift issue, we observed significantly increased false
negatives. By considering the lost tracklets in both stage, we empirically
find it effectively increases the tracking sensitivity by at least $1.1\%$. In
addition, we employ NSA Kalman algorithm [7] as the dynamic motion model and
use GSI to remediate missing tracklets [6].

###  3.2 Self-supervised Domain Adaptation

Motivations. Domain shift issue is the most challenging part in this task. It
is noticeable that patterns, such as pedestrian’s appearance, lighting
conditions, distractors, , are significantly different between the source
domain ( MOTSynth) and the target domain ( MOT17), casting a shadow on
achieving high detection/tracking performance on the testing scenario.
Therefore, in our work, we focus on leveraging self-supervised domain
adaptation method to boost the detection performance by using only unlabeled
MOT17 training split.

Figure 1: Overview of the iterative domain adaptation pipeline. Firstly, the
detector is trained on pure MOTSynth dataset. Then we perform the iterative
domain adaptation by repeatedly generating and training with the pseudo
bounding-box labels of the MOT17 training set.

Domain Adaptation via Iterative Refinement Domain adaptation is crucial for
acquiring good generalization ability on unseen data, especially when there is
a domain gap between the training and testing scenario. Recently, training
from synthetic data draws many attentions, as it is an efficient way to
collect large-scale dataset. In our work, we propose to leverage a simple yet
effective self-supervised iterative refinement approach to address the domain
shift issue. We firstly train the detector for $80$ epochs on the whole
MOTSynth dataset (with a subsampling ratio of 10), and we validate the
performance on the MOT17 validation half to determine the best model. Then, we
start the iterative domain adaptation by firstly using a pre-defined detection
confidence threshold ($T=0.5$ for the initial iteration, and $T=0.1$ for the
rest) to generate pseudo bounding-box labels for the MOT17 training set, and
we fine-tune the detector with the pseudo-labels for another $N$ ($N=40$ for
the initial iteration, and $N=20$ for the rest) epochs. We repeat this process
until there is no improvement when validating the model with MOT17 training
set true labels. We find this process significantly reduce false negatives and
noticeably improve the overall tracking performance, as demonstrated in the
Sec. 4.

###  3.3 Multi-scale Ensemble Inference

To improve the robustness of detection results during inference, similar to
[4, 5, 13], we employ an ensemble strategy to further boost the model
performance. We take three different resolutions as inputs (, $640\times
1152$, $800\times 1440$, and $1280\times 2304$), and augment each input by
horizontal flipping. This creates three pairs of inputs for the detection
network. We design a strategy to make the input resolution compatible with the
predicted box size. To be specific, we only allow the lowest resolution input
to predict boxes that larger than $32\times 32$, and restrict the largest
input to produce boxes that smaller than $96\times 96$. All box predictions
from the medium resolution inputs are used. And we collect all the legal
predictions together and perform NMS at once. The effectiveness of this multi-
scale ensemble strategy is demonstrated in Sec. 4.

##  4 Experiments

###  4.1 Dataset and Evaluation Metrics

MOTSynth [8] is a large-scale synthetic dataset for object detection and
tracking. It consists of $768$ sequences and results in $1,382$k frames. We
employ MOTSynth as the only source of labeled data for the development. For
parameter tuning and methodology validation, we use MOT17 [17] training split
as our validation set. We primarily use HOTA [16] and MOTA [2] to assess the
overall tracking performance, and also take FPs, FNs, and IDs as references.

###  4.2 Ablation Study

We conduct ablation study to reveal the effectiveness of the proposed
iterative domain adaptation and multi-scale ensemble inference. As
demonstrated in Table 1, the baseline method that trained on the MOTSynth
labeled dataset produces large amount of false negatives when testing on the
MOT17 training set. We hypotheses that this poor performance stems from the
domain shift issue as described in Sec. 3.2. Therefore, we further conduct
iterative domain adaptation to improve the sensitivity of the model. As can be
seen, we repeat this process twice and obtain consistent and significant
improvements over the baseline method. By introducing the multi-scale ensemble
inference, we further boost the HOTA score to $57.68$.

Table 1: Ablation study on the iterative domain adaptation and the multi-scale ensemble learning on MOT17 training set. Method |  HOTA |  MOTA |  FP |  FN |  IDs
---|---|---|---|---|---
Baseline |  52.87 |  60.50 |  10553 |  33258 |  514
Iter.1 |  55.66 |  63.9 |  11426 |  28643 |  517
Iter.2 |  56.82 |  63.8 |  11574 |  28569 |  456
Ensemble |  57.68 |  64.30 |  15528 |  23985 |  534

###  4.3 MOTSynth Challenge Performance

We use the best model obtained by the iterative domain adaptation to perform
inference the MOT17 testing split. We submit the results to the MOTChallenge
platform, which yields a HOTA score of 58.7.

##  5 Conclusion

We present our approach for 7th BMTT’s Challenge on MOTSynth-MOT-CVPR22. We
adopt tracking-by-detection paradigm with two-stage approach, where detector
is trained from scratch using MOTSynth data. Domain adaptation uses self-
supervised fashion where only the unlabeled MOT17 training data is used. We
also further boosts our performance with multiple-scale ensemble inference. We
achieved a HOTA score of 58.7 as the final ranking result showing on the
challenge website.

## References

  * [1] Philipp Bergmann, Tim Meinhardt, and Laura Leal-Taixé.  Tracking without bells and whistles.  In The IEEE International Conference on Computer Vision (ICCV), October 2019.
  * [2] Keni Bernardin and Rainer Stiefelhagen.  Evaluating multiple object tracking performance: the clear mot metrics.  EURASIP Journal on Image and Video Processing, 2008:1–10, 2008\.
  * [3] Alexey Bochkovskiy, Chien-Yao Wang, and Hong-Yuan Mark Liao.  Yolov4: Optimal speed and accuracy of object detection.  arXiv preprint arXiv:2004.10934, 2020.
  * [4] Ángela Casado-García and Jónathan Heras.  Ensemble methods for object detection.  In ECAI 2020, pages 2688–2695. IOS Press, 2020.
  * [5] Chi-Tung Cheng, Yirui Wang, Huan-Wu Chen, Po-Meng Hsiao, Chun-Nan Yeh, Chi-Hsun Hsieh, Shun Miao, Jing Xiao, Chien-Hung Liao, and Le Lu.  A scalable physician-level deep learning algorithm detects universal trauma on pelvic radiographs.  Nature communications, 12(1):1–10, 2021.
  * [6] Yunhao Du, Yang Song, Bo Yang, and Yanyun Zhao.  Strongsort: Make deepsort great again.  arXiv preprint arXiv:2202.13514, 2022.
  * [7] Yunhao Du, Junfeng Wan, Yanyun Zhao, Binyu Zhang, Zhihang Tong, and Junhao Dong.  Giaotracker: A comprehensive framework for mcmot with global information and optimizing strategies in visdrone 2021.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 2809–2819, 2021.
  * [8] Matteo Fabbri, Guillem Brasó, Gianluca Maugeri, Orcun Cetintas, Riccardo Gasparini, Aljoša Ošep, Simone Calderara, Laura Leal-Taixé, and Rita Cucchiara.  Motsynth: How can synthetic data help pedestrian detection and tracking?  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 10849–10859, 2021.
  * [9] Matteo Fabbri, Fabio Lanzi, Simone Calderara, Andrea Palazzi, Roberto Vezzani, and Rita Cucchiara.  Learning to detect and track visible and occluded body joints in a virtual world.  In Proceedings of the European conference on computer vision (ECCV), pages 430–446, 2018.
  * [10] Adrien Gaidon, Qiao Wang, Yohann Cabon, and Eleonora Vig.  Virtual worlds as proxy for multi-object tracking analysis.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 4340–4349, 2016.
  * [11] Zheng Ge, Songtao Liu, Zeming Li, Osamu Yoshie, and Jian Sun.  Ota: Optimal transport assignment for object detection.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 303–312, 2021.
  * [12] Zheng Ge, Songtao Liu, Feng Wang, Zeming Li, and Jian Sun.  Yolox: Exceeding yolo series in 2021.  arXiv preprint arXiv:2107.08430, 2021.
  * [13] Chen-I Hsieh, Kang Zheng, Chihung Lin, Ling Mei, Le Lu, Weijian Li, Fang-Ping Chen, Yirui Wang, Xiaoyun Zhou, Fakai Wang, et al.  Automated bone mineral density prediction and fracture risk assessment using plain radiographs via deep learning.  Nature communications, 12(1):1–9, 2021.
  * [14] Hou-Ning Hu, Qi-Zhi Cai, Dequan Wang, Ji Lin, Min Sun, Philipp Krähenbühl, Trevor Darrell, and Fisher Yu.  Joint monocular 3d vehicle detection and tracking.  2019\.
  * [15] Shu Liu, Lu Qi, Haifang Qin, Jianping Shi, and Jiaya Jia.  Path aggregation network for instance segmentation.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 8759–8768, 2018.
  * [16] Jonathon Luiten, Aljosa Osep, Patrick Dendorfer, Philip Torr, Andreas Geiger, Laura Leal-Taixé, and Bastian Leibe.  Hota: A higher order metric for evaluating multi-object tracking.  International journal of computer vision, 129(2):548–578, 2021\.
  * [17] Anton Milan, Laura Leal-Taixé, Ian Reid, Stefan Roth, and Konrad Schindler.  Mot16: A benchmark for multi-object tracking.  arXiv preprint arXiv:1603.00831, 2016.
  * [18] Jinlong Peng, Changan Wang, Fangbin Wan, Yang Wu, Yabiao Wang, Ying Tai, Chengjie Wang, Jilin Li, Feiyue Huang, and Yanwei Fu.  Chained-tracker: Chaining paired attentive regression results for end-to-end joint multiple-object detection and tracking.  In Proceedings of the European Conference on Computer Vision, 2020\.
  * [19] Joseph Redmon and Ali Farhadi.  Yolov3: An incremental improvement.  arXiv preprint arXiv:1804.02767, 2018.
  * [20] Stephan R Richter, Zeeshan Hayder, and Vladlen Koltun.  Playing for benchmarks.  In Proceedings of the IEEE International Conference on Computer Vision, pages 2213–2222, 2017.
  * [21] Stephan R Richter, Vibhav Vineet, Stefan Roth, and Vladlen Koltun.  Playing for data: Ground truth from computer games.  In European conference on computer vision, pages 102–118. Springer, 2016.
  * [22] Guanglu Song, Yu Liu, and Xiaogang Wang.  Revisiting the sibling head in object detector.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 11563–11572, 2020.
  * [23] ShiJie Sun, Naveed Akhtar, XiangYu Song, HuanSheng Song, Ajmal Mian, and Mubarak Shah.  Simultaneous detection and tracking with motion modelling for multiple object tracking.  In European Conference on Computer Vision, pages 626–643. Springer, 2020.
  * [24] Pavel Tokmakov, Jie Li, Wolfram Burgard, and Adrien Gaidon.  Learning to track with object permanence.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 10860–10869, 2021.
  * [25] Chien-Yao Wang, Hong-Yuan Mark Liao, Yueh-Hua Wu, Ping-Yang Chen, Jun-Wei Hsieh, and I-Hau Yeh.  Cspnet: A new backbone that can enhance learning capability of cnn.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition workshops, pages 390–391, 2020.
  * [26] Yue Wu, Yinpeng Chen, Lu Yuan, Zicheng Liu, Lijuan Wang, Hongzhi Li, and Yun Fu.  Rethinking classification and localization for object detection.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 10186–10195, 2020.
  * [27] Kaiyu Yang, Jacqueline Yau, Li Fei-Fei, Jia Deng, and Olga Russakovsky.  A study of face obfuscation in imagenet.  arXiv preprint arXiv:2103.06191, 2021.
  * [28] Yifu Zhang, Peize Sun, Yi Jiang, Dongdong Yu, Zehuan Yuan, Ping Luo, Wenyu Liu, and Xinggang Wang.  Bytetrack: Multi-object tracking by associating every detection box.  arXiv preprint arXiv:2110.06864, 2021.
  * [29] Xingyi Zhou, Vladlen Koltun, and Philipp Krähenbühl.  Tracking objects as points.  ECCV, 2020.