Difference between revisions of "Meylan 2002b"

From WikiWaves
Jump to navigationJump to search
Line 10: Line 10:
  
 
[[Category:Reference]]
 
[[Category:Reference]]
 +
 +
 +
 +
= Spectral Solution of Time Dependent Shallow Water Hydroelasticity =
 +
\author{M\ls I\ls C\ls H\ls A\ls E\ls L\ns H.\ns M\ls E\ls Y\ls L\ls A\ls N<math>
 +
^{1}</math>}
 +
?? and in revised form ??
 +
 +
 +
 +
The spectral theory of a thin plate floating on shallow water is derived and
 +
used to solve the time-dependent motion. This theory is based on an energy
 +
inner product in which the evolution operator becomes unitary. Two solution
 +
methods are presented. In the first, the solution is expanded in the
 +
eigenfunctions of a self-adjoint operator, which are the incoming wave
 +
solutions for a single frequency. In the second, the scattering theory of
 +
Lax-Phillips is used. The Lax-Phillips scattering solution is suitable for
 +
calculating only the free motion of the plate. However, it determines the
 +
modes of vibration of the plate-water system. These modes, both oscillate
 +
and decay, are found by a complex search algorithm based contour
 +
integration. As well as an application to modelling floating runways, the
 +
spectral-theory for a floating thin plate on shallow water is a solvable
 +
model for more complicated hydroelastic systems.
 +
 +
 +
==Introduction==
 +
 +
Hydroelasticity is the study of immersed or floating elastic bodies in a
 +
fluid. It has a wide range of applications including very large floating
 +
structures ([[Kashiwagihydro98]]), ships ([[Bishop]]), breakwaters (
 +
[[Stoker]]) and sea ice ([[Squire_Review]]). One of the best studied
 +
hydroelastic models is the linear floating thin plate ([[Newman_deform]],
 +
[[ohmatsuVLFS]], [[Kagemoto97]], and [[Kashiwagihydro98]]) because
 +
it models many physical systems, such as a floating runway or an ice floe (
 +
[[jgrfloecirc]]).
 +
 +
The time-dependence in linear hydroelastic problems is usually removed by
 +
solving for a single frequency which we will refer to as the \emph{single
 +
frequency solution}. The solution is normally found by expanding the elastic
 +
body motion in the free modes of vibration and solving the fluid equations
 +
using a Green function ([[Bishop]]). This is analogous to solving for a
 +
rigid floating body using the six rigid modes ([[Sarp_Isa]]). While
 +
alternative methods have been developed ([[Kashiwagibspline]], \cite
 +
{ohmatsuVLFS}, and [[Kagemoto97]]), these are based on exploiting some
 +
property such as a simple geometry or high relative stiffness.
 +
 +
In contrast to the single frequency solutions, solving time-dependent linear
 +
hydroelastic systems remains a major challenge. [[Kashiwagitime]] and
 +
[[Endotime]] have developed a time-stepping procedure; however, this
 +
method results in error growth in time. Since the problem is linear it is
 +
solvable by a spectral method which eliminates the long-time growth of
 +
errors. Furthermore, such a method provides information about the behaviour
 +
of the solution, such as the decay constant of the motion. However, the
 +
spectral theory for linear hydroelasticity has not been developed. For this
 +
reason, spectral type solutions such as [[Ohmatsutimedep]], based on a
 +
Fourier expansion of the single frequency solution, have only solved the
 +
problem in restricted circumstances.
 +
 +
A floating thin plate on shallow water is considered in this paper. This
 +
problem has been chosen for the following reasons: while the single
 +
frequency solution is straightforward ([[Stoker]]), the time-dependent
 +
solution has never been calculated; recently [[OhkusuNamba]], \cite
 +
{OhkusuISOPE}, and [[Ertekinshallow1999]] used it to model a floating
 +
runway; the spectral-theory developed here is a solvable model for more
 +
complex hydroelastic systems.
 +
 +
The spectral theory for a thin plate on shallow water is based on an inner
 +
product which gives the energy of the plate-water system. With respect to
 +
this inner product the evolution operator becomes unitary. Two different
 +
solution methods are derived from this spectral theory. The first method is
 +
based on an expansion of the solution in eigenfunctions of a self-adjoint
 +
operator. These are the single frequency solutions. The second solution
 +
method is based on the scattering theory of Lax-Phillips ([[laxphilips]]
 +
). It provides the solution in terms of a countable number of modes which
 +
have both an oscillation and a decay. These modes are important as they can
 +
be used to characterise the response of the system. With the exception of
 +
[[Hazard]], they have not been investigated for hydroelastic systems.
 +
 +
==Formulation: A Thin Plate on Shallow Water==
 +
 +
Figure \ref{shallow} shows a schematic diagram of the problem. The plate is
 +
infinite in the <math>y</math> direction, so that only the <math>x</math> and <math>z</math> directions are
 +
considered. The <math>x</math> direction is horizontal, the positive <math>z</math> axis points
 +
vertically up, and the plate covers the region <math>-b\leqslant x\leqslant b.</math>
 +
The water is of uniform depth <math>h</math> which is small<math>\ </math>enough that the water
 +
may be approximated as shallow ([[Stoker]]). The amplitudes are assumed
 +
small enough that the linear theory is appropriate, and the plate is
 +
sufficiently thin that the shallow draft approximation may be made (\cite
 +
{OhkusuISOPE}). The solution could be extended to waves incident at an angle
 +
on a infinite two dimensional plate, as described in [[OhkusuISOPE]], but
 +
to keep the treatment straightforward this is not done here.
 +
 +
The mathematical description of the problem follows from [[Stoker]]. The
 +
kinematic condition is
 +
<center><math>
 +
\partial _{t}\zeta =-h\partial _{x}^{2}\phi ,  (1)
 +
</math></center>
 +
where <math>\phi </math> is the velocity potential of the water (averaged over the
 +
depth) and <math>\zeta </math> is the displacement of the water surface or the plate
 +
(from the shallow draft approximation). The equation derived by equating the
 +
pressure at the free surface is
 +
<center><math>
 +
-\rho g\zeta -\rho \partial _{t}\phi =\left\{
 +
\begin{matrix}{c}
 +
0,\;\;x\notin (-b,b), \\
 +
D\partial _{x}^{4}\zeta +\rho ^{\prime }d\partial _{t}^{2}\zeta ,\;\;x\in
 +
(-b,b),
 +
\end{matrix}
 +
\right.  (2)
 +
</math></center>
 +
where <math>D</math> is the bending rigidity of the plate per unit length, <math>\rho </math> is
 +
the density of water, <math>\rho ^{\prime }</math> is the average density of the plate,
 +
<math>g</math> is the acceleration due to gravity, and <math>d</math> is the thickness of the plate
 +
<math>.</math> At the ends of the plate the free edge boundary conditions
 +
<center><math>
 +
\lim_{x\downarrow -b}\partial _{x}^{2}\zeta =\lim_{x\uparrow b}\partial
 +
_{x}^{2}\zeta =\lim_{x\downarrow -b}\partial _{x}^{3}\zeta =\lim_{x\uparrow
 +
b}\partial _{x}^{3}\zeta =0  (3)
 +
</math></center>
 +
are applied since these are common in offshore engineering applications (
 +
[[OhkusuISOPE]]). However the theory which will be developed applies
 +
equally to any of the energy-conserving edge conditions such as clamped or
 +
pinned and there is no need for the boundary conditions to be symmetric.
 +
Equation (3) gives the following implied boundary
 +
conditions for <math>\phi </math>
 +
<center><math>
 +
\lim_{x\downarrow -b}\partial _{x}^{4}\phi =\lim_{x\uparrow b}\partial
 +
_{x}^{4}\phi =\lim_{x\downarrow -b}\partial _{x}^{5}\phi =\lim_{x\uparrow
 +
b}\partial _{x}^{5}\phi =0  (4)
 +
</math></center>
 +
which will be used subsequently.
 +
 +
Non-dimensional variables are now introduced. The space variables are
 +
non-dimensionalised using the water depth <math>h,</math> and the time variables are
 +
non-dimensionalised using <math>\sqrt{h/g}</math>. The non-dimensional variables are
 +
<center><math>
 +
\bar{x}=\frac{x}{h},\;\bar{t}=t\sqrt{\frac{g}{h}},\;\bar{\zeta}=\frac{\zeta
 +
}{h},\;\mathrm{and}\;\bar{\phi}=\frac{\phi }{h^{2}\sqrt{g/h}}.
 +
</math></center>
 +
In these new variables, (1) and (2) become
 +
<center><math>
 +
\partial _{\bar{t}}\bar{\zeta}=-\partial _{\bar{x}}^{2}\bar{\phi}
 +
(5)
 +
</math></center>
 +
and
 +
<center><math>
 +
-\bar{\zeta}-\partial _{\bar{t}}\bar{\phi}=\left\{
 +
\begin{matrix}{c}
 +
0,\;\;\bar{x}\notin (-\bar{b},\bar{b}), \\
 +
\beta \partial _{\bar{x}}^{4}\bar{\zeta}+\gamma \partial _{\bar{t}}^{2}\bar{
 +
\zeta},\;\;\bar{x}\in (-\bar{b},\bar{b}),
 +
\end{matrix}
 +
\right.  (6)
 +
</math></center>
 +
where <math>\beta </math> and <math>\gamma </math> are
 +
<center><math>
 +
\beta =\frac{D}{\rho gh^{4}}\;\;\mathrm{and\ \ }\gamma =\frac{\rho ^{\prime
 +
}d}{\rho h}.
 +
</math></center>
 +
For clarity the overbar is dropped from now on.
 +
 +
The main change in extending the formulation to water of finite depth is
 +
that the velocity potential will be governed by Laplace's equation. This
 +
makes the solution of the problem much more computationally demanding since
 +
Laplace's equation must be solved by a numerical method, for example the
 +
boundary element method. Furthermore, the extension of the spectral theory,
 +
which will be developed here for shallow water, to water of finite depth is
 +
non-trivial and remains a subject for further work.
 +
 +
===Neglecting the inertia term.===
 +
 +
It can be assumed that <math>\left| \gamma \partial _{t}^{2}\zeta \right| \ll
 +
\left| \zeta \right| </math> for the following reasons ([[OhkusuISOPE]]). If we
 +
consider a mode of the displacement <math>\zeta =ae^{i\omega t}</math> (where <math>a</math> is
 +
the amplitude) then <math>\partial _{t}^{2}\zeta =-\omega ^{2}ae^{i\omega t}.</math>
 +
For each frequency, <math>\omega ,</math> there is a corresponding wavelength <math>\lambda
 +
=2\pi /\omega .</math> In the non-dimensional variables the wave speed and water
 +
depth are both unity. Since the water is shallow the wavelength <math>\lambda \gg
 +
1<math> and thus </math>\omega \ll 1.</math> It follows that any shallow water mode must
 +
satisfy <math>\left| \partial _{t}^{2}\zeta \right| \ll \left| \zeta \right| .</math>
 +
Also, <math>\gamma \ll 1</math> since <math>\rho ^{\prime }<\rho </math> (otherwise the plate
 +
would sink) and <math>d\ll h</math> (otherwise the submergence of the plate would not
 +
be negligible). Therefore,
 +
<center><math>
 +
\left| \gamma \partial _{t}^{2}\zeta \right| \ll \left| \zeta \right| ,
 +
</math></center>
 +
and we assume in what follows that the inertia, <math>\gamma \partial
 +
_{t}^{2}\zeta ,</math> is zero.
 +
 +
It should be noted that the inclusion of the inertia term in the spectral
 +
theory which will be developed is difficult because it introduces a time
 +
dependence in the energy inner product.
 +
 +
==The Energy Inner Product==
 +
 +
While equations (5) and (6) are not
 +
complicated they cannot be solved in a simple manner. It is not possible to
 +
Fourier transform in space because of the spatial discontinuity of the
 +
differential equations. The Weiner-Hopf technique cannot be used because the
 +
discontinuities divide the space into three regions. A Laplace
 +
transformation in time can be applied but this leads to non-trivial
 +
equations involving a spatially discontinuous differential equation subject
 +
to arbitrary initial conditions. However, straightforward solutions can be
 +
derived using spectral theory.
 +
 +
The spectral-theory solution of equations (5) and (\ref
 +
{displacement2}) is based on the spectral theory for a unitary operator
 +
(essentially, an operator is unitary if the adjoint is also the inverse). We
 +
therefore require an inner product in which the evolution operator is
 +
unitary. This inner product, since the system is conservative, is derived
 +
from the energy. The potential and displacement both contribute to this
 +
energy and we combine them in a two component vector, <math>U\left( x,t\right) </math>,
 +
given by
 +
<center><math>
 +
U\left( x,t\right) =\left(
 +
\begin{matrix}{c}
 +
\phi (x,t) \\
 +
i\zeta (x,t)
 +
\end{matrix}
 +
\right) .  (7)
 +
</math></center>
 +
The energy consists of the kinetic energy of the water (<math>\propto \left| \phi
 +
_{t}^{2}\right| <math>), the potential energy of the water (</math>\propto \left| \phi
 +
^{2}\right| </math>), and the energy of the plate. The energy inner product for
 +
the two vectors
 +
<center><math>
 +
U_{1}=\left(
 +
\begin{matrix}{c}
 +
\phi _{1} \\
 +
i\zeta _{1}
 +
\end{matrix}
 +
\right) \;\;\mathrm{and\ \ }U_{2}=\left(
 +
\begin{matrix}{c}
 +
\phi _{2} \\
 +
i\zeta _{2}
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
is given by
 +
<center><math>
 +
\left\langle U_{1},U_{2}\right\rangle _{\mathcal{H}}=\left\langle \partial
 +
_{x}\phi _{1},\partial _{x}\phi _{2}\right\rangle +\left\langle \left(
 +
1+\beta \left( H\left( x-b\right) -H\left( x+b\right) \right) \partial
 +
_{x}^{4}\right) i\zeta _{1},i\zeta _{2}\right\rangle ,
 +
(8)
 +
</math></center>
 +
where <math>H</math> is the Heaviside function. The subscript <math>\mathcal{H}</math> is used to
 +
denote the special inner product and the angle brackets without the <math>
 +
\mathcal{H}</math> denote the standard inner product, i.e.
 +
<center><math>
 +
\left\langle f\left( x\right) ,g\left( x\right) \right\rangle =\int_{-\infty
 +
}^{\infty }f\left( x\right) g^{\ast }\left( x\right) dx.
 +
</math></center>
 +
We now write (5) and (6) as
 +
<center><math>\begin{matrix}
 +
\frac{1}{i}\partial _{t}U &=&\mathcal{P}U  (9) \\
 +
U\left( x,t\right) _{t=0} &=&U_{0}\left( x\right) =\left(
 +
\begin{matrix}{c}
 +
\phi _{0}(x) \\
 +
i\zeta _{0}(x)
 +
\end{matrix}
 +
\right)  \nonumber
 +
\end{matrix}</math></center>
 +
where the operator <math>\mathcal{P}</math> is
 +
<center><math>
 +
\mathcal{P=}\left(
 +
\begin{matrix}{cc}
 +
0 & 1+\beta \left( H\left( x-b\right) -H\left( x+b\right) \right) \partial
 +
_{x}^{4} \\
 +
-\partial _{x}^{2} & 0
 +
\end{matrix}
 +
\right) .
 +
</math></center>
 +
<math>\mathcal{P}</math> is self-adjoint with respect to the inner product (\ref
 +
{energyinnerprod}) since <math>\mathcal{P}</math> satisfies
 +
<center><math>
 +
\left\langle \mathcal{P}U_{1},U_{2}\right\rangle _{\mathcal{H}}=\left\langle
 +
U_{1},\mathcal{P}U_{2}\right\rangle _{\mathcal{H}}
 +
</math></center>
 +
from integration by parts and the boundary conditions at the end of the
 +
plate (3). We can express the solution to (\ref
 +
{selfadjoint2}) as
 +
<center><math>
 +
U\left( x,t\right) =e^{i\mathcal{P}t}U_{0}  (10)
 +
</math></center>
 +
where <math>e^{i\mathcal{P}t}</math> is a unitary evolution operator.
 +
 +
==The self-adjoint solution method(11)==
 +
 +
In this section, a solution for the time dependent motion of the plate-water
 +
system is developed using the theory of self-adjoint operators. To evaluate
 +
equation (10) we require a method to calculate the evolution
 +
operator <math>e^{i\mathcal{P}t}</math>. This can be accomplished by using the
 +
eigenfunctions of the operator <math>\mathcal{P},</math> which are the single frequency
 +
solutions.
 +
 +
===Finding the eigenfunctions(12)===
 +
 +
Since <math>\mathcal{P}</math> is self-adjoint, the eigenvalues, <math>\lambda ,</math> must be
 +
real and therefore<math>\ </math>the eigenfunctions of <math>\mathcal{P}</math> are oscillatory
 +
exponentials outside the region of water covered by the plate. Furthermore,
 +
since the plate is finite, the spectrum (set of eigenvalues) is the entire
 +
real numbers. As is expected for two-component systems, there are two
 +
eigenfunctions associated with each eigenvalue <math>\lambda </math>. We choose
 +
incoming waves from the left (<math>\Phi ^{>})</math> and the right (<math>\Phi ^{<})</math> of
 +
unit amplitude as a basis for the eigenspace since they are the standard
 +
single frequency solutions. They have the following asymptotics,
 +
<center><math>
 +
\lim_{x\rightarrow -\infty }\Phi ^{>}=\left(
 +
\begin{matrix}{c}
 +
e^{i\lambda x} \\
 +
\lambda e^{i\lambda x}
 +
\end{matrix}
 +
\right) +S_{12}\left(
 +
\begin{matrix}{c}
 +
e^{-i\lambda x} \\
 +
\lambda e^{-i\lambda x}
 +
\end{matrix}
 +
\right) \;\;\mathrm{and\ \ }\lim_{x\rightarrow \infty }\Phi
 +
^{>}=S_{11}\left(
 +
\begin{matrix}{c}
 +
e^{i\lambda x} \\
 +
\lambda e^{i\lambda x}
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
and
 +
<center><math>
 +
\lim_{t\rightarrow -\infty }\Phi ^{<}=S_{22}\left(
 +
\begin{matrix}{c}
 +
e^{-i\lambda x} \\
 +
\lambda e^{-i\lambda x}
 +
\end{matrix}
 +
\right) \;\;\mathrm{and\ \ }\lim_{x\rightarrow \infty }\Phi ^{>}=\left(
 +
\begin{matrix}{c}
 +
e^{-i\lambda x} \\
 +
\lambda e^{-i\lambda x}
 +
\end{matrix}
 +
\right) +S_{21}\left(
 +
\begin{matrix}{c}
 +
e^{i\lambda x} \\
 +
\lambda e^{i\lambda x}
 +
\end{matrix}
 +
\right) ,
 +
</math></center>
 +
where <math>S_{11}</math>, <math>S_{12,}</math> <math>S_{21},</math> and <math>S_{22}</math> are the reflection and
 +
transmission coefficients (which must be determined). These eigenfunctions,
 +
which are analogous to the Jost solutions of Schrodinger's equation (\cite
 +
{Chadan89}), will be used to calculated the time-dependent solution.
 +
 +
We find the eigenfunction <math>\Phi ^{>}(\lambda ,x)</math> by solving (\ref
 +
{kinematic2}) and (6) in each region. The two components, <math>
 +
\phi ^{>}\left( \lambda ,x\right) <math> and </math>i\zeta ^{>}\left( \lambda ,x\right)
 +
,</math> are given by
 +
<center><math>
 +
\phi ^{>}\left( \lambda ,x\right) =\left\{
 +
\begin{matrix}{c}
 +
e^{-i\lambda x}+S_{11}\left( \lambda \right) e^{i\lambda x},\;\;x<-b, \\
 +
\sum\limits_{j=1}^{6}\alpha _{j}e^{\mu _{j}\left( \lambda \right)
 +
x},\;\;-b<x<b, \\
 +
S_{12}\left( \lambda \right) e^{-i\lambda x},\;\;x>b,
 +
\end{matrix}
 +
\right.  (13)
 +
</math></center>
 +
and
 +
<center><math>
 +
i\zeta ^{>}\left( \lambda ,x\right) =\left\{
 +
\begin{matrix}{c}
 +
\lambda e^{-i\lambda x}+\lambda S_{11}\left( \lambda \right) e^{i\lambda
 +
x},\;\;x<-b, \\
 +
\frac{-1}{\lambda }\sum\limits_{j=1}^{6}\mu _{j}\left( \lambda \right)
 +
^{2}\alpha _{j}e^{\mu _{j}\left( \lambda \right) x},\;\;-b<x<b, \\
 +
\lambda S_{12}\left( \lambda \right) e^{-i\lambda x},\;\;x>b,
 +
\end{matrix}
 +
\right.
 +
</math></center>
 +
where the coefficients <math>\mu _{j}\left( \lambda \right) </math> are the six roots
 +
of the equation
 +
<center><math>
 +
\beta \mu ^{6}+\mu ^{2}+\lambda ^{2}=0.  (14)
 +
</math></center>
 +
The values of <math>S_{11}\left( \lambda \right) ,</math> <math>S_{12}\left( \lambda \right)
 +
,<math> and </math>\alpha _{j}</math> are found from the boundary conditions (4)
 +
and the continuity of <math>\phi </math> and <math>\partial _{x}\phi </math> at <math>x=\pm b.</math>
 +
Therefore, to find the eigenfunction <math>\Phi ^{>}\left( \lambda ,x\right) ,</math>
 +
we solve the 8 by 8 linear system
 +
<center><math>
 +
\mathbf{M}\vec{a}\mathbf{=}\vec{b},  (15)
 +
</math></center>
 +
where <math>\mathbf{M}</math> is the matrix
 +
<center><math>
 +
\mathbf{M}=\left(
 +
\begin{matrix}{cccccccc}
 +
\mu _{1}^{4}e^{-\mu _{1}b} & \mu _{2}^{4}e^{-\mu _{2}b} & \mu
 +
_{3}^{4}e^{-\mu _{3}b} & \mu _{4}^{4}e^{-\mu _{4}b} & \mu _{5}^{4}e^{-\mu
 +
_{5}b} & \mu _{6}^{4}e^{-\mu _{6}b} & 0 & 0 \\
 +
\mu _{1}^{5}e^{-\mu _{1}b} & \mu _{2}^{5}e^{-\mu _{2}b} & \mu
 +
_{3}^{5}e^{-\mu _{3}b} & \mu _{4}^{5}e^{-\mu _{4}b} & \mu _{5}^{5}e^{-\mu
 +
_{5}b} & \mu _{6}^{5}e^{-\mu _{6}b} & 0 & 0 \\
 +
\mu _{1}^{4}e^{\mu _{1}b} & \mu _{2}^{4}e^{\mu _{2}b} & \mu _{3}^{4}e^{\mu
 +
_{3}b} & \mu _{4}^{4}e^{\mu _{4}b} & \mu _{5}^{4}e^{\mu _{5}b} & \mu
 +
_{6}^{4}e^{\mu _{6}b} & 0 & 0 \\
 +
\mu _{1}^{5}e^{\mu _{1}b} & \mu _{2}^{5}e^{\mu _{2}b} & \mu _{3}^{5}e^{\mu
 +
_{3}b} & \mu _{4}^{5}e^{\mu _{4}b} & \mu _{5}^{5}e^{\mu _{5}b} & \mu
 +
_{6}^{5}e^{\mu _{6}b} & 0 & 0 \\
 +
e^{-\mu _{1}b} & e^{-\mu _{2}b} & e^{-\mu _{3}b} & e^{-\mu _{4}b} & e^{-\mu
 +
_{5}b} & e^{-\mu _{6}b} & -e^{-i\lambda b} & 0 \\
 +
\mu _{1}e^{-\mu _{1}b} & \mu _{2}e^{-\mu _{2}b} & \mu _{3}e^{-\mu _{3}b} &
 +
\mu _{4}e^{-\mu _{4}b} & \mu _{5}e^{-\mu _{5}b} & \mu _{6}e^{-\mu _{6}b} &
 +
-i\lambda e^{-i\lambda b} & 0 \\
 +
e^{\mu _{1}b} & e^{\mu _{2}b} & e^{\mu _{3}b} & e^{\mu _{4}b} & e^{\mu _{5}b}
 +
& e^{\mu _{6}b} & 0 & -e^{-i\lambda b} \\
 +
\mu _{1}e^{\mu _{1}b} & \mu _{2}e^{\mu _{2}b} & \mu _{3}e^{\mu _{3}b} & \mu
 +
_{4}e^{\mu _{4}b} & \mu _{5}e^{\mu _{5}b} & \mu _{6}e^{\mu _{6}b} & 0 &
 +
i\lambda e^{-i\lambda b}
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
and <math>\vec{a}</math> and <math>\vec{b}</math> are given by
 +
<center><math>
 +
\vec{a}=\left(
 +
\begin{matrix}{c}
 +
\alpha _{1} \\
 +
\alpha _{2} \\
 +
\alpha _{3} \\
 +
\alpha _{4} \\
 +
\alpha _{5} \\
 +
\alpha _{6} \\
 +
S_{11} \\
 +
S_{12}
 +
\end{matrix}
 +
\right) ,\;\;\;\;\;\vec{b}=\left(
 +
\begin{matrix}{c}
 +
0 \\
 +
0 \\
 +
0 \\
 +
0 \\
 +
e^{i\lambda b} \\
 +
-i\lambda e^{i\lambda b} \\
 +
0 \\
 +
0
 +
\end{matrix}
 +
\right) .
 +
</math></center>
 +
Note that the coefficients <math>S_{11}</math> and <math>S_{12}</math> are contained in <math>\vec{a}.</math>
 +
 +
The eigenfunctions for the wave propagating from the right <math>\Phi ^{<}\left(
 +
\lambda ,x\right) <math> are found similarly. Since </math>S_{11}</math> represents the
 +
amplitude of the reflected wave and <math>S_{12}</math> represents the amplitude of the
 +
transmitted wave, conservation of energy requires that <math>\left| S_{11}\right|
 +
^{2}+\left| S_{12}\right| ^{2}=1.</math> Similarly, since the boundary conditions
 +
are symmetric <math>S_{22}\left( \lambda \right) =S_{11}\left( \lambda \right) </math>
 +
and <math>S_{12}\left( \lambda \right) =S_{21}\left( \lambda \right) .</math>
 +
 +
===Solution with the eigenfunctions===
 +
 +
Equation (10) can be solved by a generalised Fourier transform
 +
based on the eigenfunctions of the operator <math>\mathcal{P}</math>. The
 +
eigenfunctions are orthogonal since <math>\mathcal{P}</math> is self-adjoint, but they
 +
must be normalised. This is accomplished by using the following identity
 +
<center><math>
 +
\int_{0}^{\infty }e^{i\left( \lambda _{1}-\lambda _{2}\right) t}dt=\pi
 +
\delta \left( \lambda _{1}-\lambda _{2}\right)
 +
</math></center>
 +
where <math>\delta </math> is the Dirac delta function. Therefore
 +
<center><math>\begin{matrix}
 +
\left\langle \Phi ^{>}\left( x,\lambda _{1}\right) ,\Phi ^{>}\left(
 +
x,\lambda _{2}\right) \right\rangle _{\mathcal{H}} &=&\pi \delta \left(
 +
\lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}\left( 1+\left|
 +
S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}\right) \\
 +
&&+\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}\lambda
 +
_{2}\left( 1+\left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}\right) \\
 +
&=&4\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2},
 +
\nonumber
 +
\end{matrix}</math></center>
 +
using the condition <math>\left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}=1.</math>
 +
Similarly
 +
<center><math>
 +
\left\langle \Phi ^{<}\left( x,\lambda _{1}\right) ,\Phi ^{<}\left(
 +
x,\lambda _{2}\right) \right\rangle _{\mathcal{H}}=4\pi \delta \left(
 +
\lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}
 +
</math></center>
 +
and
 +
<center><math>\begin{matrix}
 +
\left\langle \Phi ^{>}\left( x,\lambda _{1}\right) ,\Phi ^{<}\left(
 +
x,\lambda _{2}\right) \right\rangle _{\mathcal{H}} &=&2\pi \delta \left(
 +
\lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}\left( S_{11}S_{21}^{\ast
 +
}+S_{12}S_{22}^{\ast }\right) \\
 +
&=&0.  \nonumber
 +
\end{matrix}</math></center>
 +
 +
The generalised Fourier transform which solves the evolution equation (\ref
 +
{unitary}) is
 +
<center><math>\begin{matrix}
 +
U\left( x,t\right) &=&\int_{-\infty }^{\infty }\left\langle U_{0}\left(
 +
x\right) ,\frac{\Phi ^{>}\left( x,\lambda \right) }{4\pi \lambda ^{2}}
 +
\right\rangle _{\mathcal{H}}\Phi ^{>}\left( x,\lambda \right) e^{i\lambda
 +
t}d\lambda  (16) \\
 +
&&+\int_{-\infty }^{\infty }\left\langle U_{0}\left( x\right) ,\frac{\Phi
 +
^{<}\left( x,\lambda \right) }{4\pi \lambda ^{2}}\right\rangle _{\mathcal{H}
 +
}\Phi ^{<}\left( x,\lambda \right) e^{i\lambda t}d\lambda .  \nonumber
 +
\end{matrix}</math></center>
 +
Equation (16) is the cornerstone of the approach. The
 +
integral in equation (16) can be calculated by the fast
 +
Fourier transform while the inner product can calculated by the fast Fourier
 +
transform if the initial condition <math>U_{0}</math> is zero underneath the plate (<math>
 +
-b<x<b).</math>
 +
 +
===Numerical Calculations===
 +
 +
The intention of this paper is to develop the solution methods rather than
 +
describe the physics of the motion and therefore only a few solutions are
 +
presented. From [[OhkusuISOPE]] we assume the plate stiffness is <math>\beta
 +
=2\times 10^{4}<math> and the plate length is </math>b=50</math> throughout. These values are
 +
typical for a floating runway. Figures \ref{incomingspecplot1} and \ref
 +
{incomingspecplot2} show the displacement and potential, respectively, for a
 +
pulse travelling to the right at the times <math>t=0,</math> 30, 60, 90, 120, 150, 180,
 +
210, and 240. The incoming wave pulse was chosen to be a Gaussian in
 +
potential centered at <math>x=-125</math> and sufficiently sharp to be negligible under
 +
the plate,
 +
<center><math>
 +
U_{0}\left( x\right) =\left(
 +
\begin{matrix}{c}
 +
\phi \left( x\right) \\
 +
i\phi ^{\prime }\left( x\right)
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
where
 +
<center><math>
 +
\phi \left( x\right) =\left\{
 +
\begin{matrix}{c}
 +
e^{-\tfrac{(x+125)^{2}}{350}},\;\;x<-50, \\
 +
0,\;\;x>-50.
 +
\end{matrix}
 +
\right.
 +
</math></center>
 +
At <math>t=0</math> the plate is initially at rest and the wave is to the left of the
 +
plate propagating towards it. From <math>t=30</math> the wave has reached the plate and
 +
the plate begins to undergo a complex bending motion in response to the
 +
incoming wave. The response of the plate in turn induces waves in the
 +
surrounding water which propagate away from the plate to the left and right.
 +
The final picture, <math>t=240</math>, shows the plate at rest with waves now
 +
propagating away from it. The majority of the wave energy has passed under
 +
the plate and continues to propagate to the right. However, the shape of the
 +
outgoing wave profile is markedly different from the incoming wave profile.
 +
Also, there is a significant reflected wave propagating away from the plate
 +
to the left.
 +
 +
Figures \ref{spectral1} and \ref{spectral2} show the evolution of the plate
 +
from an initial displacement in the absence of wave forcing for the times <math>
 +
t=0,</math> 20, 40, 60, 80, 100, 120, 140, and 160. Only the plate displacement is
 +
initially non-zero so that
 +
<center><math>
 +
U_{0}\left( x\right) =\left(
 +
\begin{matrix}{c}
 +
0 \\
 +
i\zeta \left( x\right)
 +
\end{matrix}
 +
\right) .
 +
</math></center>
 +
Figure \ref{spectral1} shows the motion for the symmetric initial plate
 +
displacement
 +
<center><math>
 +
\zeta \left( x\right) =e^{-\tfrac{x^{2}}{350}}.
 +
</math></center>
 +
As the plate evolves the plate vibrates, straightens, and the amplitude
 +
decays. The decay is due to the radiation of energy by the waves generated
 +
in the surrounding water. A complex wave train is produced by the plate
 +
motion and can be seen propagating away from the plate. Figure \ref
 +
{spectral2} shows the motion for the non-symmetric initial plate
 +
displacement
 +
<center><math>
 +
\zeta \left( x\right) =e^{-\tfrac{\left( x-50\right) ^{2}}{350}}.
 +
</math></center>
 +
Again as the plate evolves it straightens, vibrates, and decays and induces
 +
waves in the surrounding water.
 +
 +
==The Lax-Phillips Scattering Solution Method.==
 +
 +
In this section, a solution to the time-dependent motion of the plate-water
 +
system is developed using the Lax-Phillips scattering theory (\cite
 +
{laxphilips}). This solution method will only solve for an initial condition
 +
which is zero outside the plate, i.e. <math>U_{0}\left( x\right) =0</math> if <math>\left|
 +
x\right| >b</math>. However, it calculates the solution by an expansion in a
 +
countable number of modes.
 +
 +
===Lax-Phillips Scattering(17)===
 +
 +
The Lax-Phillips scattering theory will be briefly outlined here for our
 +
specific problem. The Hilbert space <math>\mathcal{H}</math>\ is decomposed into three
 +
subspaces. The incoming space, denoted by <math>D_{-},</math> consists of all waves
 +
travelling towards the plate, either from the left or the right, as
 +
appropriate. The outgoing subspace, denoted by <math>D_{+},</math> consists of all
 +
waves travelling away from the plate, again either to the left or right, as
 +
appropriate. What remains is the scattering space, denoted by <math>K,\ </math>
 +
consisting of the potential and displacement under the plate.
 +
 +
To apply the Lax-Phillips scattering the following conditions are required: <math>
 +
D_{-}<math> and </math>D_{+}</math> must be orthogonal; the incoming subspace must span the
 +
entire space under temporal evolution. For our system, the first condition
 +
follows from the inner product and the second condition follows from the
 +
simple structure of the eigenfunctions of the operator <math>\mathcal{P}.</math> From
 +
the Lax-Phillips scattering theory, since these conditions hold, the
 +
equation of motion for the plate in the absence of incoming waves can be
 +
written
 +
<center><math>
 +
\frac{1}{i}\partial _{t}U=\mathcal{B}U,  (18)
 +
</math></center>
 +
where <math>\mathcal{B}</math> is a non-self-adjoint operator. <math>\mathcal{B}</math> is related
 +
to <math>\mathcal{P}</math> by
 +
<center><math>
 +
e^{i\mathcal{B}t}=P_{K}\left. e^{i\mathcal{P}t}\right| _{K}
 +
</math></center>
 +
where <math>P_{K}</math> is the projection onto the subspace <math>K</math> and <math>\left. {}\right|
 +
_{K}<math> denotes a restriction to </math>K.<math> Therefore </math>e^{i\mathcal{B}t}</math> is the
 +
restricted to <math>K</math> of the evolution of an initial condition which is zero
 +
outside <math>K.</math> It should be noted that the equality in equation (\ref
 +
{nonselffirst}) is in general only true asymptotically. However the
 +
numerical results show for our case we have equality for all times.
 +
 +
From the Lax-Phillips scattering theory, equation (18) can
 +
be solved by finding the eigenvalues (sometimes referred to as scattering
 +
frequencies or resonances) and eigenfunctions of <math>\mathcal{B}</math>. The
 +
eigenvalues of <math>\mathcal{B}</math> occur at the singularities of the analytic
 +
extension to <math>\mathbb{C}</math> of the scattering matrix, <math>\mathbf{S}(\lambda ).</math>
 +
This is given by
 +
<center><math>
 +
\mathbf{S}(\lambda )=\left(
 +
\begin{matrix}{cc}
 +
S_{11}\left( \lambda \right) & S_{12}\left( \lambda \right) \\
 +
S_{21}\left( \lambda \right) & S_{22}\left( \lambda \right)
 +
\end{matrix}
 +
\right)  (19)
 +
</math></center>
 +
where <math>S_{11}\left( \lambda \right) </math>, <math>S_{12}\left( \lambda \right) ,</math> <math>
 +
S_{21}\left( \lambda \right) ,<math> and </math>S_{22}\left( \lambda \right) </math> are the
 +
scattered wave coefficients found from the single frequency solutions in
 +
section 12. As a consequence of the Lax-Phillips scattering
 +
theory the scattering matrix is unitary for real <math>\lambda </math> and the
 +
singularities must all lie in the upper complex plane (<math>{Im}\left(
 +
\lambda \right) >0).</math> Once the singularities have been found, the
 +
eigenfunctions can be calculated. They are not orthogonal, since <math>\mathcal{B}
 +
</math> is a non-self-adjoint operator, but a biorthogonal system can be formed
 +
using the eigenfunctions of the adjoint operator, <math>\mathcal{B}^{\ast }.</math>
 +
 +
The eigenfunctions of <math>\mathcal{B}</math> are the modes of vibration for the
 +
plate-water system. These modes have a decay as well as an oscillation due
 +
to the radiation of energy into the surrounding water. The frequency of the
 +
oscillation is determined by the real part of the eigenvalue and the rate of
 +
decay is determined by the imaginary part of the eigenvalue.
 +
 +
While the eigenvalues of <math>\mathcal{B}</math> occur precisely at the singularities
 +
of the solution found by a Laplace transformation in time the Lax-Phillips
 +
scattering theory solution has three major advantages over the Laplace
 +
transform solution: the eigenvalues (singularities) can be found using the
 +
scattering matrix; the difficult equations in the Laplace space involving
 +
the initial condition do not need to be solved; the contribution of the
 +
singularity (the residue) can be found directly from the inner product of
 +
the initial condition with the corresponding eigenfunction of the adjoint
 +
operator, <math>\mathcal{B}^{\ast }.</math>
 +
 +
===Finding the Singularities of the Scattering Matrix===
 +
 +
While the analytic extension of the scattering matrix is straightforward,
 +
the linear system (15) is solved for complex <math>\lambda ,</math>
 +
finding the singularities of the scattering matrix is non-trivial<math>.</math> The
 +
difficulty lies in the fact that we must search the complex plane for the
 +
singularities with no ''a priori''  knowledge about their location. We use
 +
a complex search algorithm based contour integration. The determinant of the
 +
scattering matrix is integrated around the contour of a region of the
 +
complex plane. If the value of this integral is zero, then the region is
 +
assumed to contain no singularities and the search is terminated (the
 +
possibility that the contribution of two or more singularities might cancel
 +
can be treated by considering further integrals, such as the variation of
 +
the argument around the contour). If the value of the integral is not zero,
 +
then the region must contain singularities and it is then divided into
 +
subregions and the search is repeated. Once the singularities have been
 +
located sufficiently well they are used as seeds for Newton's method and
 +
found to high accuracy.
 +
 +
If the eigenvalues have to be found for different parameter values then a
 +
homotopy, or continuation, method can be used, which avoids the slow complex
 +
search method. This method uses the known locations of the eigenvalues for
 +
one parameter value to determine the eigenvalues for a new parameter value
 +
by taking sufficiently small steps that Newton's method can be used with the
 +
previous solution as a seed. Unfortunately, a homotopy method requires the
 +
solution of the eigenvalues for at least one parameter value as an initial
 +
seed and this must be accomplished by a complex search algorithm.
 +
 +
The position of the eigenvalues for <math>\beta =2\times 10^{4}</math> and <math>b=50</math> are
 +
shown in Figure \ref{spectrum}. They are denoted by <math>\lambda _{n},</math> where <math>
 +
n\in \mathbb{Z}<math>, and ordered by increasing real part, with </math>n=0</math>
 +
corresponding the eigenvalue with smallest absolute real part. From the
 +
picture and on physical grounds, it seems likely that there exist
 +
asymptotics for the eigenvalues, however this theory is not developed here.
 +
 +
===Eigenfunctions===
 +
 +
The eigenfunctions of <math>\mathcal{B}</math> associated with the eigenvalue <math>\lambda
 +
_{n}<math> are denoted by </math>\Phi ^{+}(\lambda _{n},x),<math> and those of </math>\mathcal{B}
 +
^{\ast }<math> (the adjoint of </math>\mathcal{B})<math> associated with the eigenvalue </math>
 +
\lambda _{n}^{\ast }<math> are denoted by </math>\hat{\Phi}^{+}\left( \lambda
 +
_{n}^{\ast },x\right) </math>. That is,
 +
<center><math>
 +
\mathcal{B}\Phi ^{+}\left( \lambda _{n},x\right) =\lambda _{n}\Phi
 +
^{+}\left( \lambda _{n},x\right)
 +
</math></center>
 +
and
 +
<center><math>
 +
\mathcal{B}^{\ast }\hat{\Phi}^{+}\left( \lambda _{n}^{\ast },x\right)
 +
=\lambda _{n}^{\ast }\hat{\Phi}^{+}\left( \lambda _{n}^{\ast },x\right) .
 +
</math></center>
 +
The eigenfunction <math>\Phi ^{+}\left( \lambda _{n},x\right) </math> can be written
 +
<center><math>
 +
\Phi ^{+}\left( \lambda _{n},x\right) =\left(
 +
\begin{matrix}{c}
 +
\phi ^{+}\left( \lambda _{n},x\right) \\
 +
i\zeta ^{+}\left( \lambda _{n},x\right)
 +
\end{matrix}
 +
\right) =\left(
 +
\begin{matrix}{c}
 +
\sum\limits_{j=1}^{6}\alpha _{j}e^{\mu _{j}\left( \lambda _{n}\right) x} \\
 +
\sum\limits_{j=1}^{6}-\frac{\alpha _{j}\mu _{j}\left( \lambda _{n}\right)
 +
^{2}}{\lambda _{n}}e^{\mu _{j}\left( \lambda _{n}\right) x}
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
where <math>\mu _{j}\left( \lambda \right) </math> are the six roots of equation (\ref
 +
{cubicinlambda}) and the coefficients <math>\alpha _{j}</math> are found from the
 +
boundary conditions at the end of the plate (4) and the
 +
following condition. Since the scattering matrix is singular, there are no
 +
incoming wave from either direction. We use the condition that there is no
 +
incoming wave at <math>x=-b</math> and that the outgoing wave is of unit amplitude,
 +
i.e.
 +
<center><math>
 +
\phi ^{+}\left( \lambda _{n},-b\right) =e^{i\lambda _{n}b},\;=\ = \left.
 +
\partial _{x}\phi ^{+}\left( \lambda _{n},x\right) \right| _{x=-b}=i\lambda
 +
_{n}e^{i\lambda _{n}b}.
 +
</math></center>
 +
We do not use the condition that there is no outgoing wave at <math>x=b</math> because
 +
the system will become over determined. Therefore, the coefficients <math>\alpha
 +
_{j}</math> satisfy the linear equation
 +
<center><math>
 +
\mathbf{M}\vec{a}\mathbf{=}\vec{b},
 +
</math></center>
 +
where
 +
<center><math>
 +
\mathbf{M}=\left(
 +
\begin{matrix}{cccccc}
 +
\mu _{1}^{4}e^{-\mu _{1}b} & \mu _{2}^{4}e^{-\mu _{2}b} & \mu
 +
_{3}^{4}e^{-\mu _{3}b} & \mu _{4}^{4}e^{-\mu _{4}b} & \mu _{5}^{4}e^{-\mu
 +
_{5}b} & \mu _{6}^{4}e^{-\mu _{6}b} \\
 +
\mu _{1}^{5}e^{-\mu _{1}b} & \mu _{2}^{5}e^{-\mu _{2}b} & \mu
 +
_{3}^{5}e^{-\mu _{3}b} & \mu _{4}^{5}e^{-\mu _{4}b} & \mu _{5}^{5}e^{-\mu
 +
_{5}b} & \mu _{6}^{5}e^{-\mu _{6}b} \\
 +
\mu _{1}^{4}e^{\mu _{1}b} & \mu _{2}^{4}e^{\mu _{2}b} & \mu _{3}^{4}e^{\mu
 +
_{3}b} & \mu _{4}^{4}e^{\mu _{4}b} & \mu _{5}^{4}e^{\mu _{5}b} & \mu
 +
_{6}^{4}e^{\mu _{6}b} \\
 +
\mu _{1}^{5}e^{\mu _{1}b} & \mu _{2}^{5}e^{\mu _{2}b} & \mu _{3}^{5}e^{\mu
 +
_{3}b} & \mu _{4}^{5}e^{\mu _{4}b} & \mu _{5}^{5}e^{\mu _{5}b} & \mu
 +
_{6}^{5}e^{\mu _{6}b} \\
 +
e^{-\mu _{1}b} & e^{-\mu _{2}b} & e^{-\mu _{3}b} & e^{-\mu _{4}b} & e^{-\mu
 +
_{5}b} & e^{-\mu _{6}b} \\
 +
\mu _{1}e^{-\mu _{1}b} & \mu _{2}e^{-\mu _{2}b} & \mu _{3}e^{-\mu _{3}b} &
 +
\mu _{4}e^{-\mu _{4}b} & \mu _{5}e^{-\mu _{5}b} & \mu _{6}e^{-\mu _{6}b}
 +
\end{matrix}
 +
\right)
 +
</math></center>
 +
and <math>\vec{a}</math> and <math>\vec{b}</math> are given by
 +
<center><math>
 +
\vec{a}=\left(
 +
\begin{matrix}{c}
 +
\alpha _{1} \\
 +
\alpha _{2} \\
 +
\alpha _{3} \\
 +
\alpha _{4} \\
 +
\alpha _{5} \\
 +
\alpha _{6}
 +
\end{matrix}
 +
\right) ,\;\;\;\;\;\vec{b}=\left(
 +
\begin{matrix}{c}
 +
0 \\
 +
0 \\
 +
0 \\
 +
0 \\
 +
e^{i\lambda b} \\
 +
i\lambda e^{i\lambda b}
 +
\end{matrix}
 +
\right) .
 +
</math></center>
 +
The eigenfunctions for the adjoint operator are found similarly.
 +
 +
Figure \ref{eigfunctions} shows the real and imaginary parts of the
 +
eigenfunctions of <math>\mathcal{B}</math> for <math>n=1,</math> <math>3</math>, 5, and <math>7,</math> again with <math>
 +
\beta =2\times 10^{4}<math> and </math>b=50.</math> While the eigenfunctions do not have a
 +
simple shape, increasing oscillation is apparent as <math>n</math> increases.
 +
 +
===Inner products===
 +
 +
A biorthogonal system with respect to the energy inner product (\ref
 +
{energyinnerprod}) is formed by the eigenfunctions of <math>\mathcal{B}</math>, <math>\Phi
 +
^{+}\left( \lambda _{n},x\right) ,<math> and the eigenfunctions of </math>\mathcal{B}
 +
^{\ast },<math> </math>\hat{\Phi}^{+}\left( \lambda _{n},x\right) </math>. To normalise the
 +
biorthogonal system, the inner product of <math>\Phi ^{+}\left( \lambda
 +
_{n},x\right) <math> and </math>\hat{\Phi}^{+}\left( \lambda _{n},x\right) </math> has to be
 +
determined. From the definition of the energy inner product (\ref
 +
{energyinnerprod})
 +
<center><math>\begin{matrix}
 +
\left\langle \Phi \left( \lambda _{m},x\right) ,\hat{\Phi}\left( \lambda
 +
_{n},x\right) \right\rangle _{\mathcal{H}} &=&\int_{-b}^{b}\partial _{x}\phi
 +
^{+}\left( \lambda _{m},x\right) \left( \partial _{x}\hat{\phi}^{+}\left(
 +
\lambda _{n}^{\ast },x\right) \right) ^{\ast }dx  (20) \\
 +
&&+\int_{-b}^{b}(1+P)i\zeta ^{+}\left( \lambda _{m}^{\ast },x\right) \left( i
 +
\hat{\zeta}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx.
 +
\nonumber
 +
\end{matrix}</math></center>
 +
The two integrals in (20) are considered separately. The
 +
first is
 +
<center><math>\begin{matrix}
 +
&&\int_{-b}^{b}\partial _{x}\phi ^{+}\left( \lambda _{m},x\right) \left(
 +
\partial _{x}\hat{\phi}^{+}\left( \lambda _{n}^{\ast },x\right) \right)
 +
^{\ast }dx \\
 +
&=&\int_{-b}^{b}\left( \sum_{j=1}^{6}\mu _{j}\left( \lambda _{m}\right)
 +
\alpha _{j}e^{\mu _{j}\left( \lambda _{m}\right) x}\right) \left(
 +
\sum_{k=1}^{6}\mu _{k}\left( \lambda _{m}\right) \alpha _{k}e^{\mu
 +
_{k}\left( \lambda _{n}\right) x}\right) dx \\
 +
&=&\sum_{j=1}^{6}\sum_{k=1}^{6}\int_{-b}^{b}-\mu _{j}\left( \lambda
 +
_{m}\right) ^{2}\alpha _{j}e^{\mu _{j}\left( \lambda _{m}\right) x}\alpha
 +
_{k}e^{\mu _{k}\left( \lambda _{n}\right) x}dx \\
 +
&=&\sum_{j=1}^{6}\sum_{k=1}^{6}-\mu _{j}\left( \lambda _{m}\right)
 +
^{2}\alpha _{j}\alpha _{k}\left( \frac{e^{\left( \mu _{j}\left( \lambda
 +
_{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu
 +
_{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}
 +
}{\mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) }
 +
\right)  \nonumber
 +
\end{matrix}</math></center>
 +
and the second is
 +
<center><math>\begin{matrix}
 +
&&\int_{-b}^{b}(1+P)i\zeta ^{+}\left( \lambda _{m}^{\ast },x\right) \left( i
 +
\hat{\zeta}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx \\
 +
&=&\int_{-b}^{b}\left( \sum_{j=1}^{6}-\frac{\alpha _{j}}{\lambda _{m}}
 +
\left( \mu _{j}\left( \lambda _{m}\right) ^{2}+\beta \mu _{j}\left( \lambda
 +
_{m}\right) ^{6}\right) e^{\mu _{j}\left( \lambda _{m}\right) x}\right)
 +
\left( \sum_{k=1}^{6}-\frac{\alpha _{k}}{\lambda _{n}}\mu _{k}\left(
 +
\lambda _{n}\right) ^{2}e^{\mu _{k}\left( \lambda _{n}\right) x}\right) dx \\
 +
&=&\sum_{j=1}^{6}\sum_{k=1}^{6}\int_{-b}^{b}\frac{\alpha _{j}\alpha _{k}}{
 +
\lambda _{m}\lambda _{n}}\left( \mu _{j}\left( \lambda _{m}\right)
 +
^{2}+\beta \mu _{j}\left( \lambda _{m}\right) ^{6}\right) \mu _{k}\left(
 +
\lambda _{n}\right) ^{2}e^{\mu _{j}\left( \lambda _{m}\right) x}e^{\mu
 +
_{k}\left( \lambda _{n}\right) x}dx \\
 +
&=&\sum_{j=1}^{6}\sum_{k=1}^{6}\frac{\alpha _{j}\alpha _{k}}{\lambda
 +
_{m}\lambda _{n}}\left( \mu _{j}\left( \lambda _{m}\right) ^{2}+\beta \mu
 +
_{j}\left( \lambda _{m}\right) ^{6}\right) \mu _{k}\left( \lambda
 +
_{n}\right) ^{2}\left( \frac{e^{\left( \mu _{j}\left( \lambda _{m}\right)
 +
+\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu _{j}\left(
 +
\lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}}{\mu
 +
_{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) }\right)
 +
\\
 +
&=&\sum_{j=1}^{6}\sum_{k=1}^{6}\frac{\alpha _{j}\alpha _{k}}{\lambda
 +
_{m}\lambda _{n}}\left( -\lambda _{m}^{2}\right) \mu _{k}\left( \lambda
 +
_{n}\right) ^{2}\left( \frac{e^{\left( \mu _{j}\left( \lambda _{m}\right)
 +
+\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu _{j}\left(
 +
\lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}}{\mu
 +
_{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) }\right) .
 +
\nonumber
 +
\end{matrix}</math></center>
 +
Therefore the calculation of the inner product in equation (\ref
 +
{innerprodall}) does not require numerical integration.
 +
 +
===Solution===
 +
 +
By solving (18) using the eigenfunctions of <math>\mathcal{B}</math>
 +
and <math>\mathcal{B}^{\ast }</math> we find the evolution of the plate, from an
 +
initial displacement <math>U_{0}(x),</math> is
 +
<center><math>
 +
U\left( x,t\right) =\sum_{n=-\infty }^{\infty }e^{i\lambda _{n}t}\frac{
 +
\left\langle U_{0}\left( x\right) ,\hat{\Phi}\left( \lambda _{n},x\right)
 +
\right\rangle _{\mathcal{H}}}{\left\langle \Phi \left( \lambda _{n},x\right)
 +
,\hat{\Phi}\left( \lambda _{n},x\right) \right\rangle _{\mathcal{H}}}\Phi
 +
\left( \lambda _{n},x\right) .  (21)
 +
</math></center>
 +
The inner product of <math>U_{0}</math> with the eigenfunction <math>\hat{\Phi}\left(
 +
\lambda _{n},x\right) </math> is the only quantity left to compute in (\ref
 +
{evolutionB}). This inner product is written
 +
<center><math>
 +
\left\langle U_{0}\left( x\right) ,\hat{\Phi}\left( \lambda _{n},x\right)
 +
\right\rangle _{\mathcal{H}}=\sum\limits_{j=1}^{6}\left(
 +
\begin{matrix}{c}
 +
\alpha _{j}\dint_{-b}^{b}-\mu _{j}\left( \lambda _{n}\right) ^{2}e^{\mu
 +
_{j}\left( \lambda _{n}\right) x}\phi _{0}\left( x\right) dx+\left. \alpha
 +
_{j}\mu _{j}\left( \lambda _{n}\right) e^{\mu _{j}\left( \lambda _{n}\right)
 +
x}\phi _{0}\left( x\right) \right| _{-b}^{b} \\
 +
\alpha _{j}\lambda _{n}\dint_{-b}^{b}e^{\mu _{j}\left( \lambda _{n}\right)
 +
x}i\zeta _{0}\left( x\right) dx
 +
\end{matrix}
 +
\right)  (22)
 +
</math></center>
 +
and the integrals must be evaluated by numerical integration. The solutions
 +
calculated using the Lax-Phillips scattering theory are identical to those
 +
found using the self-adjoint operator method and for this reason no further
 +
figures are shown.\pagebreak
 +
 +
==Conclusions==
 +
 +
The spectral theory of a linear thin plate floating on shallow water has
 +
been derived. Two spectral-theory solutions have been presented which
 +
determine the time-dependent motion of the thin plate. The first method was
 +
based on self-adjoint operator theory and the second method was based on the
 +
Lax-Phillips scattering. The self-adjoint method solved both the wave
 +
forcing and the free plate problem while the Lax-Phillips method only solved
 +
for a free plate. The eigenfunctions for the self-adjoint method are
 +
orthogonal and the eigenvalues are continuous and consist of all <math>\mathbb{R}
 +
, </math> which makes the calculation of the eigenvalues trivial. The Lax-Phillips
 +
method has discrete eigenvalues which must be calculated numerically and the
 +
system of eigenfunctions is biorthogonal. The advantage of the Lax-Phillips
 +
method is that the modes of vibration of the plate-water system and their
 +
frequency and rate of decay are found. While the relative speeds of the two
 +
methods depends of the exact way in which they are implemented, the
 +
Lax-Phillips method should be considerably faster if the eigenvalues have
 +
been determined.
 +
 +
The development of a spectral theory for more complicated hydroelastic
 +
problems remains a major challenge. While this theory must be more
 +
complicated than that presented here, many features can be expected to
 +
remain. For example, [[ohmatsuVLFS]] has shown that the single frequency
 +
solutions can be used to solve certain time-dependent problems and \cite
 +
{Hazard} have shown that modes, in which the solution can be expanded, exist
 +
for other hydroelastic systems.
 +
 +
{\Large Acknowledgments}
 +
 +
\begin{acknowledgment}
 +
I would like to thank the anonymous reviewers, Dr. Kathy Ruggerio, and Prof.
 +
James Sneyd for their very helpful comments. Also, Prof. Boris Pavlov for
 +
explaining the Lax-Phillips scattering. \pagebreak
 +
\end{acknowledgment}
 +
 +
\bibliographystyle{jfm}
 +
\bibliography{mike,others}
 +
\pagebreak
 +
 +
\begin{center}
 +
{\huge Figure Captions}
 +
\end{center}
 +
 +
\textsc{Figure} 1. Schematic diagram of a thin plate floating on shallow
 +
water and the coordinates and dimensions of the problem.
 +
 +
\textsc{Figure} 2. The evolution of the displacement due to a pulse
 +
travelling to the right for the times shown. The plate occupies the region <math>
 +
-50\leq x\leq 50<math> and is shown by the bold line. </math>\beta =2\times
 +
10^{4},b=50. </math>
 +
 +
\textsc{Figure} 3. The evolution of the potential due to a pulse travelling
 +
to the right for the times shown. The plate occupies the region <math>-50\leq
 +
x\leq 50<math> and is shown by the bold line. </math>\beta =2\times 10^{4},b=50.</math>
 +
 +
\textsc{Figure} 4. The evolution of the displacement for a plate released at
 +
<math>t=0</math> for the times shown. The plate occupies the region <math>-50\leq x\leq 50</math>
 +
and is shown by the bold line. <math>\beta =2\times 10^{4},b=50.</math>
 +
 +
\textsc{Figure} 5. The evolution of the displacement for a plate released at
 +
<math>t=0</math> for the times shown. The plate occupies the region <math>-50\leq x\leq 50</math>
 +
and is shown by the bold line. <math>\beta =2\times 10^{4},b=50.</math>
 +
 +
\textsc{Figure} 6. The location of the first 19 eigenvalues <math>\lambda _{n}</math>
 +
of <math>\mathcal{B}</math> for <math>\beta =2\times 10^{4},\;b=50.</math>
 +
 +
\textsc{Figure} 7. The real (solid) and imaginary (dashed) parts of the
 +
resonance eigenfunctions for <math>n=1,</math> <math>3</math>, 5, and <math>7</math> as shown. <math>\beta
 +
=2\times 10^{4},\;b=50.</math>
 +
 +
 +
\pagebreak
 +
 +
\FRAME{ftFU}{5.5434in}{1.6544in}{0pt}{\Qcb{}}{\Qlb{shallow}}{shallow.eps}{
 +
\special{language "Scientific Word";type "GRAPHIC";display "ICON";valid_file
 +
"F";width 5.5434in;height 1.6544in;depth 0pt;original-width
 +
6.9297in;original-height 1.6544in;cropleft "0";croptop "1";cropright
 +
"1";cropbottom "0";filename 'shallow.eps';file-properties "XNPEU";}}
 +
 +
Some nonsense\pagebreak
 +
 +
\FRAME{ftF}{4.6164in}{3.6296in}{0pt}{}{\Qlb{incomingspecplot1}}{
 +
incomingspecplot1.eps}{\special{language "Scientific Word";type
 +
"GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width
 +
4.6164in;height 3.6296in;depth 0pt;original-width 6.8632in;original-height
 +
5.3869in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename
 +
'incomingspecplot1.eps';file-properties "XNPEU";}}
 +
 +
\FRAME{ftF}{4.5602in}{3.6288in}{0pt}{}{\Qlb{incomingspecplot2}}{
 +
incomingspecplot2.eps}{\special{language "Scientific Word";type
 +
"GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width
 +
4.5602in;height 3.6288in;depth 0pt;original-width 6.7793in;original-height
 +
5.386in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename
 +
'incomingspecplot2.eps';file-properties "XNPEU";}}
 +
 +
\FRAME{fthFU}{4.6155in}{3.6357in}{0pt}{\Qcb{{}}}{\Qlb{spectral1}}{
 +
spectral1.eps}{\special{language "Scientific Word";type
 +
"GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width
 +
4.6155in;height 3.6357in;depth 0pt;original-width 6.8493in;original-height
 +
5.3852in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename
 +
'spectral1.eps';file-properties "XNPEU";}}
 +
 +
\FRAME{fthFU}{4.5887in}{3.608in}{0pt}{\Qcb{{}}}{\Qlb{spectral2}}{
 +
spectral2.eps}{\special{language "Scientific Word";type
 +
"GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width
 +
4.5887in;height 3.608in;depth 0pt;original-width 6.8493in;original-height
 +
5.3852in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename
 +
'spectral2.eps';file-properties "XNPEU";}}
 +
 +
\FRAME{ftFU}{4.5982in}{3.7853in}{0pt}{\Qcb{{}}}{\Qlb{spectrum}}{spectrum.eps
 +
}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio
 +
TRUE;display "ICON";valid_file "F";width 4.5982in;height 3.7853in;depth
 +
0pt;original-width 6.8165in;original-height 5.5486in;cropleft "0";croptop
 +
"1";cropright "1";cropbottom "0";filename 'spectrum.eps';file-properties
 +
"XNPEU";}}
 +
 +
\FRAME{ftbpFU}{4.6164in}{3.7377in}{0pt}{\Qcb{{}}}{\Qlb{eigfunctions}}{
 +
eigfunctions.eps}{\special{language "Scientific Word";type
 +
"GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width
 +
4.6164in;height 3.7377in;depth 0pt;original-width 6.845in;original-height
 +
5.5365in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename
 +
'eigfunctions.eps';file-properties "XNPEU";}}

Revision as of 07:48, 25 July 2006

Michael H Meylan, Spectral Solution of Time Dependent Shallow Water Hydroelasticity, J. of Fluid Mechanics, 454, pp 387-402, 2002.

Time-Dependent Linear Water Wave problem of a Floating Elastic Plate on Shallow Depth. The solution was found using a Generalised Eigenfunction Expansion and as a sum over Scattering Frequencies.


Spectral Solution of Time Dependent Shallow Water Hydroelasticity

\author{M\ls I\ls C\ls H\ls A\ls E\ls L\ns H.\ns M\ls E\ls Y\ls L\ls A\ls N[math]\displaystyle{ ^{1} }[/math]} ?? and in revised form ??


The spectral theory of a thin plate floating on shallow water is derived and used to solve the time-dependent motion. This theory is based on an energy inner product in which the evolution operator becomes unitary. Two solution methods are presented. In the first, the solution is expanded in the eigenfunctions of a self-adjoint operator, which are the incoming wave solutions for a single frequency. In the second, the scattering theory of Lax-Phillips is used. The Lax-Phillips scattering solution is suitable for calculating only the free motion of the plate. However, it determines the modes of vibration of the plate-water system. These modes, both oscillate and decay, are found by a complex search algorithm based contour integration. As well as an application to modelling floating runways, the spectral-theory for a floating thin plate on shallow water is a solvable model for more complicated hydroelastic systems.


Introduction

Hydroelasticity is the study of immersed or floating elastic bodies in a fluid. It has a wide range of applications including very large floating structures (Kashiwagihydro98), ships (Bishop), breakwaters ( Stoker) and sea ice (Squire_Review). One of the best studied hydroelastic models is the linear floating thin plate (Newman_deform, ohmatsuVLFS, Kagemoto97, and Kashiwagihydro98) because it models many physical systems, such as a floating runway or an ice floe ( jgrfloecirc).

The time-dependence in linear hydroelastic problems is usually removed by solving for a single frequency which we will refer to as the \emph{single frequency solution}. The solution is normally found by expanding the elastic body motion in the free modes of vibration and solving the fluid equations using a Green function (Bishop). This is analogous to solving for a rigid floating body using the six rigid modes (Sarp_Isa). While alternative methods have been developed (Kashiwagibspline, \cite {ohmatsuVLFS}, and Kagemoto97), these are based on exploiting some property such as a simple geometry or high relative stiffness.

In contrast to the single frequency solutions, solving time-dependent linear hydroelastic systems remains a major challenge. Kashiwagitime and Endotime have developed a time-stepping procedure; however, this method results in error growth in time. Since the problem is linear it is solvable by a spectral method which eliminates the long-time growth of errors. Furthermore, such a method provides information about the behaviour of the solution, such as the decay constant of the motion. However, the spectral theory for linear hydroelasticity has not been developed. For this reason, spectral type solutions such as Ohmatsutimedep, based on a Fourier expansion of the single frequency solution, have only solved the problem in restricted circumstances.

A floating thin plate on shallow water is considered in this paper. This problem has been chosen for the following reasons: while the single frequency solution is straightforward (Stoker), the time-dependent solution has never been calculated; recently OhkusuNamba, \cite {OhkusuISOPE}, and Ertekinshallow1999 used it to model a floating runway; the spectral-theory developed here is a solvable model for more complex hydroelastic systems.

The spectral theory for a thin plate on shallow water is based on an inner product which gives the energy of the plate-water system. With respect to this inner product the evolution operator becomes unitary. Two different solution methods are derived from this spectral theory. The first method is based on an expansion of the solution in eigenfunctions of a self-adjoint operator. These are the single frequency solutions. The second solution method is based on the scattering theory of Lax-Phillips (laxphilips ). It provides the solution in terms of a countable number of modes which have both an oscillation and a decay. These modes are important as they can be used to characterise the response of the system. With the exception of Hazard, they have not been investigated for hydroelastic systems.

Formulation: A Thin Plate on Shallow Water

Figure \ref{shallow} shows a schematic diagram of the problem. The plate is infinite in the [math]\displaystyle{ y }[/math] direction, so that only the [math]\displaystyle{ x }[/math] and [math]\displaystyle{ z }[/math] directions are considered. The [math]\displaystyle{ x }[/math] direction is horizontal, the positive [math]\displaystyle{ z }[/math] axis points vertically up, and the plate covers the region [math]\displaystyle{ -b\leqslant x\leqslant b. }[/math] The water is of uniform depth [math]\displaystyle{ h }[/math] which is small[math]\displaystyle{ \ }[/math]enough that the water may be approximated as shallow (Stoker). The amplitudes are assumed small enough that the linear theory is appropriate, and the plate is sufficiently thin that the shallow draft approximation may be made (\cite {OhkusuISOPE}). The solution could be extended to waves incident at an angle on a infinite two dimensional plate, as described in OhkusuISOPE, but to keep the treatment straightforward this is not done here.

The mathematical description of the problem follows from Stoker. The kinematic condition is

[math]\displaystyle{ \partial _{t}\zeta =-h\partial _{x}^{2}\phi , (1) }[/math]

where [math]\displaystyle{ \phi }[/math] is the velocity potential of the water (averaged over the depth) and [math]\displaystyle{ \zeta }[/math] is the displacement of the water surface or the plate (from the shallow draft approximation). The equation derived by equating the pressure at the free surface is

[math]\displaystyle{ -\rho g\zeta -\rho \partial _{t}\phi =\left\{ \begin{matrix}{c} 0,\;\;x\notin (-b,b), \\ D\partial _{x}^{4}\zeta +\rho ^{\prime }d\partial _{t}^{2}\zeta ,\;\;x\in (-b,b), \end{matrix} \right. (2) }[/math]

where [math]\displaystyle{ D }[/math] is the bending rigidity of the plate per unit length, [math]\displaystyle{ \rho }[/math] is the density of water, [math]\displaystyle{ \rho ^{\prime } }[/math] is the average density of the plate, [math]\displaystyle{ g }[/math] is the acceleration due to gravity, and [math]\displaystyle{ d }[/math] is the thickness of the plate [math]\displaystyle{ . }[/math] At the ends of the plate the free edge boundary conditions

[math]\displaystyle{ \lim_{x\downarrow -b}\partial _{x}^{2}\zeta =\lim_{x\uparrow b}\partial _{x}^{2}\zeta =\lim_{x\downarrow -b}\partial _{x}^{3}\zeta =\lim_{x\uparrow b}\partial _{x}^{3}\zeta =0 (3) }[/math]

are applied since these are common in offshore engineering applications ( OhkusuISOPE). However the theory which will be developed applies equally to any of the energy-conserving edge conditions such as clamped or pinned and there is no need for the boundary conditions to be symmetric. Equation (3) gives the following implied boundary conditions for [math]\displaystyle{ \phi }[/math]

[math]\displaystyle{ \lim_{x\downarrow -b}\partial _{x}^{4}\phi =\lim_{x\uparrow b}\partial _{x}^{4}\phi =\lim_{x\downarrow -b}\partial _{x}^{5}\phi =\lim_{x\uparrow b}\partial _{x}^{5}\phi =0 (4) }[/math]

which will be used subsequently.

Non-dimensional variables are now introduced. The space variables are non-dimensionalised using the water depth [math]\displaystyle{ h, }[/math] and the time variables are non-dimensionalised using [math]\displaystyle{ \sqrt{h/g} }[/math]. The non-dimensional variables are

[math]\displaystyle{ \bar{x}=\frac{x}{h},\;\bar{t}=t\sqrt{\frac{g}{h}},\;\bar{\zeta}=\frac{\zeta }{h},\;\mathrm{and}\;\bar{\phi}=\frac{\phi }{h^{2}\sqrt{g/h}}. }[/math]

In these new variables, (1) and (2) become

[math]\displaystyle{ \partial _{\bar{t}}\bar{\zeta}=-\partial _{\bar{x}}^{2}\bar{\phi} (5) }[/math]

and

[math]\displaystyle{ -\bar{\zeta}-\partial _{\bar{t}}\bar{\phi}=\left\{ \begin{matrix}{c} 0,\;\;\bar{x}\notin (-\bar{b},\bar{b}), \\ \beta \partial _{\bar{x}}^{4}\bar{\zeta}+\gamma \partial _{\bar{t}}^{2}\bar{ \zeta},\;\;\bar{x}\in (-\bar{b},\bar{b}), \end{matrix} \right. (6) }[/math]

where [math]\displaystyle{ \beta }[/math] and [math]\displaystyle{ \gamma }[/math] are

[math]\displaystyle{ \beta =\frac{D}{\rho gh^{4}}\;\;\mathrm{and\ \ }\gamma =\frac{\rho ^{\prime }d}{\rho h}. }[/math]

For clarity the overbar is dropped from now on.

The main change in extending the formulation to water of finite depth is that the velocity potential will be governed by Laplace's equation. This makes the solution of the problem much more computationally demanding since Laplace's equation must be solved by a numerical method, for example the boundary element method. Furthermore, the extension of the spectral theory, which will be developed here for shallow water, to water of finite depth is non-trivial and remains a subject for further work.

Neglecting the inertia term.

It can be assumed that [math]\displaystyle{ \left| \gamma \partial _{t}^{2}\zeta \right| \ll \left| \zeta \right| }[/math] for the following reasons (OhkusuISOPE). If we consider a mode of the displacement [math]\displaystyle{ \zeta =ae^{i\omega t} }[/math] (where [math]\displaystyle{ a }[/math] is the amplitude) then [math]\displaystyle{ \partial _{t}^{2}\zeta =-\omega ^{2}ae^{i\omega t}. }[/math] For each frequency, [math]\displaystyle{ \omega , }[/math] there is a corresponding wavelength [math]\displaystyle{ \lambda =2\pi /\omega . }[/math] In the non-dimensional variables the wave speed and water depth are both unity. Since the water is shallow the wavelength [math]\displaystyle{ \lambda \gg 1\lt math\gt and thus }[/math]\omega \ll 1.</math> It follows that any shallow water mode must satisfy [math]\displaystyle{ \left| \partial _{t}^{2}\zeta \right| \ll \left| \zeta \right| . }[/math] Also, [math]\displaystyle{ \gamma \ll 1 }[/math] since [math]\displaystyle{ \rho ^{\prime }\lt \rho }[/math] (otherwise the plate would sink) and [math]\displaystyle{ d\ll h }[/math] (otherwise the submergence of the plate would not be negligible). Therefore,

[math]\displaystyle{ \left| \gamma \partial _{t}^{2}\zeta \right| \ll \left| \zeta \right| , }[/math]

and we assume in what follows that the inertia, [math]\displaystyle{ \gamma \partial _{t}^{2}\zeta , }[/math] is zero.

It should be noted that the inclusion of the inertia term in the spectral theory which will be developed is difficult because it introduces a time dependence in the energy inner product.

The Energy Inner Product

While equations (5) and (6) are not complicated they cannot be solved in a simple manner. It is not possible to Fourier transform in space because of the spatial discontinuity of the differential equations. The Weiner-Hopf technique cannot be used because the discontinuities divide the space into three regions. A Laplace transformation in time can be applied but this leads to non-trivial equations involving a spatially discontinuous differential equation subject to arbitrary initial conditions. However, straightforward solutions can be derived using spectral theory.

The spectral-theory solution of equations (5) and (\ref {displacement2}) is based on the spectral theory for a unitary operator (essentially, an operator is unitary if the adjoint is also the inverse). We therefore require an inner product in which the evolution operator is unitary. This inner product, since the system is conservative, is derived from the energy. The potential and displacement both contribute to this energy and we combine them in a two component vector, [math]\displaystyle{ U\left( x,t\right) }[/math], given by

[math]\displaystyle{ U\left( x,t\right) =\left( \begin{matrix}{c} \phi (x,t) \\ i\zeta (x,t) \end{matrix} \right) . (7) }[/math]

The energy consists of the kinetic energy of the water ([math]\displaystyle{ \propto \left| \phi _{t}^{2}\right| \lt math\gt ), the potential energy of the water ( }[/math]\propto \left| \phi ^{2}\right| </math>), and the energy of the plate. The energy inner product for the two vectors

[math]\displaystyle{ U_{1}=\left( \begin{matrix}{c} \phi _{1} \\ i\zeta _{1} \end{matrix} \right) \;\;\mathrm{and\ \ }U_{2}=\left( \begin{matrix}{c} \phi _{2} \\ i\zeta _{2} \end{matrix} \right) }[/math]

is given by

[math]\displaystyle{ \left\langle U_{1},U_{2}\right\rangle _{\mathcal{H}}=\left\langle \partial _{x}\phi _{1},\partial _{x}\phi _{2}\right\rangle +\left\langle \left( 1+\beta \left( H\left( x-b\right) -H\left( x+b\right) \right) \partial _{x}^{4}\right) i\zeta _{1},i\zeta _{2}\right\rangle , (8) }[/math]

where [math]\displaystyle{ H }[/math] is the Heaviside function. The subscript [math]\displaystyle{ \mathcal{H} }[/math] is used to denote the special inner product and the angle brackets without the [math]\displaystyle{ \mathcal{H} }[/math] denote the standard inner product, i.e.

[math]\displaystyle{ \left\langle f\left( x\right) ,g\left( x\right) \right\rangle =\int_{-\infty }^{\infty }f\left( x\right) g^{\ast }\left( x\right) dx. }[/math]

We now write (5) and (6) as

[math]\displaystyle{ \begin{matrix} \frac{1}{i}\partial _{t}U &=&\mathcal{P}U (9) \\ U\left( x,t\right) _{t=0} &=&U_{0}\left( x\right) =\left( \begin{matrix}{c} \phi _{0}(x) \\ i\zeta _{0}(x) \end{matrix} \right) \nonumber \end{matrix} }[/math]

where the operator [math]\displaystyle{ \mathcal{P} }[/math] is

[math]\displaystyle{ \mathcal{P=}\left( \begin{matrix}{cc} 0 & 1+\beta \left( H\left( x-b\right) -H\left( x+b\right) \right) \partial _{x}^{4} \\ -\partial _{x}^{2} & 0 \end{matrix} \right) . }[/math]

[math]\displaystyle{ \mathcal{P} }[/math] is self-adjoint with respect to the inner product (\ref {energyinnerprod}) since [math]\displaystyle{ \mathcal{P} }[/math] satisfies

[math]\displaystyle{ \left\langle \mathcal{P}U_{1},U_{2}\right\rangle _{\mathcal{H}}=\left\langle U_{1},\mathcal{P}U_{2}\right\rangle _{\mathcal{H}} }[/math]

from integration by parts and the boundary conditions at the end of the plate (3). We can express the solution to (\ref {selfadjoint2}) as

[math]\displaystyle{ U\left( x,t\right) =e^{i\mathcal{P}t}U_{0} (10) }[/math]

where [math]\displaystyle{ e^{i\mathcal{P}t} }[/math] is a unitary evolution operator.

The self-adjoint solution method(11)

In this section, a solution for the time dependent motion of the plate-water system is developed using the theory of self-adjoint operators. To evaluate equation (10) we require a method to calculate the evolution operator [math]\displaystyle{ e^{i\mathcal{P}t} }[/math]. This can be accomplished by using the eigenfunctions of the operator [math]\displaystyle{ \mathcal{P}, }[/math] which are the single frequency solutions.

Finding the eigenfunctions(12)

Since [math]\displaystyle{ \mathcal{P} }[/math] is self-adjoint, the eigenvalues, [math]\displaystyle{ \lambda , }[/math] must be real and therefore[math]\displaystyle{ \ }[/math]the eigenfunctions of [math]\displaystyle{ \mathcal{P} }[/math] are oscillatory exponentials outside the region of water covered by the plate. Furthermore, since the plate is finite, the spectrum (set of eigenvalues) is the entire real numbers. As is expected for two-component systems, there are two eigenfunctions associated with each eigenvalue [math]\displaystyle{ \lambda }[/math]. We choose incoming waves from the left ([math]\displaystyle{ \Phi ^{\gt }) }[/math] and the right ([math]\displaystyle{ \Phi ^{\lt }) }[/math] of unit amplitude as a basis for the eigenspace since they are the standard single frequency solutions. They have the following asymptotics,

[math]\displaystyle{ \lim_{x\rightarrow -\infty }\Phi ^{\gt }=\left( \begin{matrix}{c} e^{i\lambda x} \\ \lambda e^{i\lambda x} \end{matrix} \right) +S_{12}\left( \begin{matrix}{c} e^{-i\lambda x} \\ \lambda e^{-i\lambda x} \end{matrix} \right) \;\;\mathrm{and\ \ }\lim_{x\rightarrow \infty }\Phi ^{\gt }=S_{11}\left( \begin{matrix}{c} e^{i\lambda x} \\ \lambda e^{i\lambda x} \end{matrix} \right) }[/math]

and

[math]\displaystyle{ \lim_{t\rightarrow -\infty }\Phi ^{\lt }=S_{22}\left( \begin{matrix}{c} e^{-i\lambda x} \\ \lambda e^{-i\lambda x} \end{matrix} \right) \;\;\mathrm{and\ \ }\lim_{x\rightarrow \infty }\Phi ^{\gt }=\left( \begin{matrix}{c} e^{-i\lambda x} \\ \lambda e^{-i\lambda x} \end{matrix} \right) +S_{21}\left( \begin{matrix}{c} e^{i\lambda x} \\ \lambda e^{i\lambda x} \end{matrix} \right) , }[/math]

where [math]\displaystyle{ S_{11} }[/math], [math]\displaystyle{ S_{12,} }[/math] [math]\displaystyle{ S_{21}, }[/math] and [math]\displaystyle{ S_{22} }[/math] are the reflection and transmission coefficients (which must be determined). These eigenfunctions, which are analogous to the Jost solutions of Schrodinger's equation (\cite {Chadan89}), will be used to calculated the time-dependent solution.

We find the eigenfunction [math]\displaystyle{ \Phi ^{\gt }(\lambda ,x) }[/math] by solving (\ref {kinematic2}) and (6) in each region. The two components, [math]\displaystyle{ \phi ^{\gt }\left( \lambda ,x\right) \lt math\gt and }[/math]i\zeta ^{>}\left( \lambda ,x\right) ,</math> are given by

[math]\displaystyle{ \phi ^{\gt }\left( \lambda ,x\right) =\left\{ \begin{matrix}{c} e^{-i\lambda x}+S_{11}\left( \lambda \right) e^{i\lambda x},\;\;x\lt -b, \\ \sum\limits_{j=1}^{6}\alpha _{j}e^{\mu _{j}\left( \lambda \right) x},\;\;-b\lt x\lt b, \\ S_{12}\left( \lambda \right) e^{-i\lambda x},\;\;x\gt b, \end{matrix} \right. (13) }[/math]

and

[math]\displaystyle{ i\zeta ^{\gt }\left( \lambda ,x\right) =\left\{ \begin{matrix}{c} \lambda e^{-i\lambda x}+\lambda S_{11}\left( \lambda \right) e^{i\lambda x},\;\;x\lt -b, \\ \frac{-1}{\lambda }\sum\limits_{j=1}^{6}\mu _{j}\left( \lambda \right) ^{2}\alpha _{j}e^{\mu _{j}\left( \lambda \right) x},\;\;-b\lt x\lt b, \\ \lambda S_{12}\left( \lambda \right) e^{-i\lambda x},\;\;x\gt b, \end{matrix} \right. }[/math]

where the coefficients [math]\displaystyle{ \mu _{j}\left( \lambda \right) }[/math] are the six roots of the equation

[math]\displaystyle{ \beta \mu ^{6}+\mu ^{2}+\lambda ^{2}=0. (14) }[/math]

The values of [math]\displaystyle{ S_{11}\left( \lambda \right) , }[/math] [math]\displaystyle{ S_{12}\left( \lambda \right) ,\lt math\gt and }[/math]\alpha _{j}</math> are found from the boundary conditions (4) and the continuity of [math]\displaystyle{ \phi }[/math] and [math]\displaystyle{ \partial _{x}\phi }[/math] at [math]\displaystyle{ x=\pm b. }[/math] Therefore, to find the eigenfunction [math]\displaystyle{ \Phi ^{\gt }\left( \lambda ,x\right) , }[/math] we solve the 8 by 8 linear system

[math]\displaystyle{ \mathbf{M}\vec{a}\mathbf{=}\vec{b}, (15) }[/math]

where [math]\displaystyle{ \mathbf{M} }[/math] is the matrix

[math]\displaystyle{ \mathbf{M}=\left( \begin{matrix}{cccccccc} \mu _{1}^{4}e^{-\mu _{1}b} & \mu _{2}^{4}e^{-\mu _{2}b} & \mu _{3}^{4}e^{-\mu _{3}b} & \mu _{4}^{4}e^{-\mu _{4}b} & \mu _{5}^{4}e^{-\mu _{5}b} & \mu _{6}^{4}e^{-\mu _{6}b} & 0 & 0 \\ \mu _{1}^{5}e^{-\mu _{1}b} & \mu _{2}^{5}e^{-\mu _{2}b} & \mu _{3}^{5}e^{-\mu _{3}b} & \mu _{4}^{5}e^{-\mu _{4}b} & \mu _{5}^{5}e^{-\mu _{5}b} & \mu _{6}^{5}e^{-\mu _{6}b} & 0 & 0 \\ \mu _{1}^{4}e^{\mu _{1}b} & \mu _{2}^{4}e^{\mu _{2}b} & \mu _{3}^{4}e^{\mu _{3}b} & \mu _{4}^{4}e^{\mu _{4}b} & \mu _{5}^{4}e^{\mu _{5}b} & \mu _{6}^{4}e^{\mu _{6}b} & 0 & 0 \\ \mu _{1}^{5}e^{\mu _{1}b} & \mu _{2}^{5}e^{\mu _{2}b} & \mu _{3}^{5}e^{\mu _{3}b} & \mu _{4}^{5}e^{\mu _{4}b} & \mu _{5}^{5}e^{\mu _{5}b} & \mu _{6}^{5}e^{\mu _{6}b} & 0 & 0 \\ e^{-\mu _{1}b} & e^{-\mu _{2}b} & e^{-\mu _{3}b} & e^{-\mu _{4}b} & e^{-\mu _{5}b} & e^{-\mu _{6}b} & -e^{-i\lambda b} & 0 \\ \mu _{1}e^{-\mu _{1}b} & \mu _{2}e^{-\mu _{2}b} & \mu _{3}e^{-\mu _{3}b} & \mu _{4}e^{-\mu _{4}b} & \mu _{5}e^{-\mu _{5}b} & \mu _{6}e^{-\mu _{6}b} & -i\lambda e^{-i\lambda b} & 0 \\ e^{\mu _{1}b} & e^{\mu _{2}b} & e^{\mu _{3}b} & e^{\mu _{4}b} & e^{\mu _{5}b} & e^{\mu _{6}b} & 0 & -e^{-i\lambda b} \\ \mu _{1}e^{\mu _{1}b} & \mu _{2}e^{\mu _{2}b} & \mu _{3}e^{\mu _{3}b} & \mu _{4}e^{\mu _{4}b} & \mu _{5}e^{\mu _{5}b} & \mu _{6}e^{\mu _{6}b} & 0 & i\lambda e^{-i\lambda b} \end{matrix} \right) }[/math]

and [math]\displaystyle{ \vec{a} }[/math] and [math]\displaystyle{ \vec{b} }[/math] are given by

[math]\displaystyle{ \vec{a}=\left( \begin{matrix}{c} \alpha _{1} \\ \alpha _{2} \\ \alpha _{3} \\ \alpha _{4} \\ \alpha _{5} \\ \alpha _{6} \\ S_{11} \\ S_{12} \end{matrix} \right) ,\;\;\;\;\;\vec{b}=\left( \begin{matrix}{c} 0 \\ 0 \\ 0 \\ 0 \\ e^{i\lambda b} \\ -i\lambda e^{i\lambda b} \\ 0 \\ 0 \end{matrix} \right) . }[/math]

Note that the coefficients [math]\displaystyle{ S_{11} }[/math] and [math]\displaystyle{ S_{12} }[/math] are contained in [math]\displaystyle{ \vec{a}. }[/math]

The eigenfunctions for the wave propagating from the right [math]\displaystyle{ \Phi ^{\lt }\left( \lambda ,x\right) \lt math\gt are found similarly. Since }[/math]S_{11}</math> represents the amplitude of the reflected wave and [math]\displaystyle{ S_{12} }[/math] represents the amplitude of the transmitted wave, conservation of energy requires that [math]\displaystyle{ \left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}=1. }[/math] Similarly, since the boundary conditions are symmetric [math]\displaystyle{ S_{22}\left( \lambda \right) =S_{11}\left( \lambda \right) }[/math] and [math]\displaystyle{ S_{12}\left( \lambda \right) =S_{21}\left( \lambda \right) . }[/math]

Solution with the eigenfunctions

Equation (10) can be solved by a generalised Fourier transform based on the eigenfunctions of the operator [math]\displaystyle{ \mathcal{P} }[/math]. The eigenfunctions are orthogonal since [math]\displaystyle{ \mathcal{P} }[/math] is self-adjoint, but they must be normalised. This is accomplished by using the following identity

[math]\displaystyle{ \int_{0}^{\infty }e^{i\left( \lambda _{1}-\lambda _{2}\right) t}dt=\pi \delta \left( \lambda _{1}-\lambda _{2}\right) }[/math]

where [math]\displaystyle{ \delta }[/math] is the Dirac delta function. Therefore

[math]\displaystyle{ \begin{matrix} \left\langle \Phi ^{\gt }\left( x,\lambda _{1}\right) ,\Phi ^{\gt }\left( x,\lambda _{2}\right) \right\rangle _{\mathcal{H}} &=&\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}\left( 1+\left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}\right) \\ &&+\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}\lambda _{2}\left( 1+\left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}\right) \\ &=&4\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}, \nonumber \end{matrix} }[/math]

using the condition [math]\displaystyle{ \left| S_{11}\right| ^{2}+\left| S_{12}\right| ^{2}=1. }[/math] Similarly

[math]\displaystyle{ \left\langle \Phi ^{\lt }\left( x,\lambda _{1}\right) ,\Phi ^{\lt }\left( x,\lambda _{2}\right) \right\rangle _{\mathcal{H}}=4\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2} }[/math]

and

[math]\displaystyle{ \begin{matrix} \left\langle \Phi ^{\gt }\left( x,\lambda _{1}\right) ,\Phi ^{\lt }\left( x,\lambda _{2}\right) \right\rangle _{\mathcal{H}} &=&2\pi \delta \left( \lambda _{1}-\lambda _{2}\right) \lambda _{1}^{2}\left( S_{11}S_{21}^{\ast }+S_{12}S_{22}^{\ast }\right) \\ &=&0. \nonumber \end{matrix} }[/math]

The generalised Fourier transform which solves the evolution equation (\ref {unitary}) is

[math]\displaystyle{ \begin{matrix} U\left( x,t\right) &=&\int_{-\infty }^{\infty }\left\langle U_{0}\left( x\right) ,\frac{\Phi ^{\gt }\left( x,\lambda \right) }{4\pi \lambda ^{2}} \right\rangle _{\mathcal{H}}\Phi ^{\gt }\left( x,\lambda \right) e^{i\lambda t}d\lambda (16) \\ &&+\int_{-\infty }^{\infty }\left\langle U_{0}\left( x\right) ,\frac{\Phi ^{\lt }\left( x,\lambda \right) }{4\pi \lambda ^{2}}\right\rangle _{\mathcal{H} }\Phi ^{\lt }\left( x,\lambda \right) e^{i\lambda t}d\lambda . \nonumber \end{matrix} }[/math]

Equation (16) is the cornerstone of the approach. The integral in equation (16) can be calculated by the fast Fourier transform while the inner product can calculated by the fast Fourier transform if the initial condition [math]\displaystyle{ U_{0} }[/math] is zero underneath the plate ([math]\displaystyle{ -b\lt x\lt b). }[/math]

Numerical Calculations

The intention of this paper is to develop the solution methods rather than describe the physics of the motion and therefore only a few solutions are presented. From OhkusuISOPE we assume the plate stiffness is [math]\displaystyle{ \beta =2\times 10^{4}\lt math\gt and the plate length is }[/math]b=50</math> throughout. These values are typical for a floating runway. Figures \ref{incomingspecplot1} and \ref {incomingspecplot2} show the displacement and potential, respectively, for a pulse travelling to the right at the times [math]\displaystyle{ t=0, }[/math] 30, 60, 90, 120, 150, 180, 210, and 240. The incoming wave pulse was chosen to be a Gaussian in potential centered at [math]\displaystyle{ x=-125 }[/math] and sufficiently sharp to be negligible under the plate,

[math]\displaystyle{ U_{0}\left( x\right) =\left( \begin{matrix}{c} \phi \left( x\right) \\ i\phi ^{\prime }\left( x\right) \end{matrix} \right) }[/math]

where

[math]\displaystyle{ \phi \left( x\right) =\left\{ \begin{matrix}{c} e^{-\tfrac{(x+125)^{2}}{350}},\;\;x\lt -50, \\ 0,\;\;x\gt -50. \end{matrix} \right. }[/math]

At [math]\displaystyle{ t=0 }[/math] the plate is initially at rest and the wave is to the left of the plate propagating towards it. From [math]\displaystyle{ t=30 }[/math] the wave has reached the plate and the plate begins to undergo a complex bending motion in response to the incoming wave. The response of the plate in turn induces waves in the surrounding water which propagate away from the plate to the left and right. The final picture, [math]\displaystyle{ t=240 }[/math], shows the plate at rest with waves now propagating away from it. The majority of the wave energy has passed under the plate and continues to propagate to the right. However, the shape of the outgoing wave profile is markedly different from the incoming wave profile. Also, there is a significant reflected wave propagating away from the plate to the left.

Figures \ref{spectral1} and \ref{spectral2} show the evolution of the plate from an initial displacement in the absence of wave forcing for the times [math]\displaystyle{ t=0, }[/math] 20, 40, 60, 80, 100, 120, 140, and 160. Only the plate displacement is initially non-zero so that

[math]\displaystyle{ U_{0}\left( x\right) =\left( \begin{matrix}{c} 0 \\ i\zeta \left( x\right) \end{matrix} \right) . }[/math]

Figure \ref{spectral1} shows the motion for the symmetric initial plate displacement

[math]\displaystyle{ \zeta \left( x\right) =e^{-\tfrac{x^{2}}{350}}. }[/math]

As the plate evolves the plate vibrates, straightens, and the amplitude decays. The decay is due to the radiation of energy by the waves generated in the surrounding water. A complex wave train is produced by the plate motion and can be seen propagating away from the plate. Figure \ref {spectral2} shows the motion for the non-symmetric initial plate displacement

[math]\displaystyle{ \zeta \left( x\right) =e^{-\tfrac{\left( x-50\right) ^{2}}{350}}. }[/math]

Again as the plate evolves it straightens, vibrates, and decays and induces waves in the surrounding water.

The Lax-Phillips Scattering Solution Method.

In this section, a solution to the time-dependent motion of the plate-water system is developed using the Lax-Phillips scattering theory (\cite {laxphilips}). This solution method will only solve for an initial condition which is zero outside the plate, i.e. [math]\displaystyle{ U_{0}\left( x\right) =0 }[/math] if [math]\displaystyle{ \left| x\right| \gt b }[/math]. However, it calculates the solution by an expansion in a countable number of modes.

Lax-Phillips Scattering(17)

The Lax-Phillips scattering theory will be briefly outlined here for our specific problem. The Hilbert space [math]\displaystyle{ \mathcal{H} }[/math]\ is decomposed into three subspaces. The incoming space, denoted by [math]\displaystyle{ D_{-}, }[/math] consists of all waves travelling towards the plate, either from the left or the right, as appropriate. The outgoing subspace, denoted by [math]\displaystyle{ D_{+}, }[/math] consists of all waves travelling away from the plate, again either to the left or right, as appropriate. What remains is the scattering space, denoted by [math]\displaystyle{ K,\ }[/math] consisting of the potential and displacement under the plate.

To apply the Lax-Phillips scattering the following conditions are required: [math]\displaystyle{ D_{-}\lt math\gt and }[/math]D_{+}</math> must be orthogonal; the incoming subspace must span the entire space under temporal evolution. For our system, the first condition follows from the inner product and the second condition follows from the simple structure of the eigenfunctions of the operator [math]\displaystyle{ \mathcal{P}. }[/math] From the Lax-Phillips scattering theory, since these conditions hold, the equation of motion for the plate in the absence of incoming waves can be written

[math]\displaystyle{ \frac{1}{i}\partial _{t}U=\mathcal{B}U, (18) }[/math]

where [math]\displaystyle{ \mathcal{B} }[/math] is a non-self-adjoint operator. [math]\displaystyle{ \mathcal{B} }[/math] is related to [math]\displaystyle{ \mathcal{P} }[/math] by

[math]\displaystyle{ e^{i\mathcal{B}t}=P_{K}\left. e^{i\mathcal{P}t}\right| _{K} }[/math]

where [math]\displaystyle{ P_{K} }[/math] is the projection onto the subspace [math]\displaystyle{ K }[/math] and [math]\displaystyle{ \left. {}\right| _{K}\lt math\gt denotes a restriction to }[/math]K.[math]\displaystyle{ Therefore }[/math]e^{i\mathcal{B}t}</math> is the restricted to [math]\displaystyle{ K }[/math] of the evolution of an initial condition which is zero outside [math]\displaystyle{ K. }[/math] It should be noted that the equality in equation (\ref {nonselffirst}) is in general only true asymptotically. However the numerical results show for our case we have equality for all times.

From the Lax-Phillips scattering theory, equation (18) can be solved by finding the eigenvalues (sometimes referred to as scattering frequencies or resonances) and eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math]. The eigenvalues of [math]\displaystyle{ \mathcal{B} }[/math] occur at the singularities of the analytic extension to [math]\displaystyle{ \mathbb{C} }[/math] of the scattering matrix, [math]\displaystyle{ \mathbf{S}(\lambda ). }[/math] This is given by

[math]\displaystyle{ \mathbf{S}(\lambda )=\left( \begin{matrix}{cc} S_{11}\left( \lambda \right) & S_{12}\left( \lambda \right) \\ S_{21}\left( \lambda \right) & S_{22}\left( \lambda \right) \end{matrix} \right) (19) }[/math]

where [math]\displaystyle{ S_{11}\left( \lambda \right) }[/math], [math]\displaystyle{ S_{12}\left( \lambda \right) , }[/math] [math]\displaystyle{ S_{21}\left( \lambda \right) ,\lt math\gt and }[/math]S_{22}\left( \lambda \right) </math> are the scattered wave coefficients found from the single frequency solutions in section 12. As a consequence of the Lax-Phillips scattering theory the scattering matrix is unitary for real [math]\displaystyle{ \lambda }[/math] and the singularities must all lie in the upper complex plane ([math]\displaystyle{ {Im}\left( \lambda \right) \gt 0). }[/math] Once the singularities have been found, the eigenfunctions can be calculated. They are not orthogonal, since [math]\displaystyle{ \mathcal{B} }[/math] is a non-self-adjoint operator, but a biorthogonal system can be formed using the eigenfunctions of the adjoint operator, [math]\displaystyle{ \mathcal{B}^{\ast }. }[/math]

The eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math] are the modes of vibration for the plate-water system. These modes have a decay as well as an oscillation due to the radiation of energy into the surrounding water. The frequency of the oscillation is determined by the real part of the eigenvalue and the rate of decay is determined by the imaginary part of the eigenvalue.

While the eigenvalues of [math]\displaystyle{ \mathcal{B} }[/math] occur precisely at the singularities of the solution found by a Laplace transformation in time the Lax-Phillips scattering theory solution has three major advantages over the Laplace transform solution: the eigenvalues (singularities) can be found using the scattering matrix; the difficult equations in the Laplace space involving the initial condition do not need to be solved; the contribution of the singularity (the residue) can be found directly from the inner product of the initial condition with the corresponding eigenfunction of the adjoint operator, [math]\displaystyle{ \mathcal{B}^{\ast }. }[/math]

Finding the Singularities of the Scattering Matrix

While the analytic extension of the scattering matrix is straightforward, the linear system (15) is solved for complex [math]\displaystyle{ \lambda , }[/math] finding the singularities of the scattering matrix is non-trivial[math]\displaystyle{ . }[/math] The difficulty lies in the fact that we must search the complex plane for the singularities with no a priori knowledge about their location. We use a complex search algorithm based contour integration. The determinant of the scattering matrix is integrated around the contour of a region of the complex plane. If the value of this integral is zero, then the region is assumed to contain no singularities and the search is terminated (the possibility that the contribution of two or more singularities might cancel can be treated by considering further integrals, such as the variation of the argument around the contour). If the value of the integral is not zero, then the region must contain singularities and it is then divided into subregions and the search is repeated. Once the singularities have been located sufficiently well they are used as seeds for Newton's method and found to high accuracy.

If the eigenvalues have to be found for different parameter values then a homotopy, or continuation, method can be used, which avoids the slow complex search method. This method uses the known locations of the eigenvalues for one parameter value to determine the eigenvalues for a new parameter value by taking sufficiently small steps that Newton's method can be used with the previous solution as a seed. Unfortunately, a homotopy method requires the solution of the eigenvalues for at least one parameter value as an initial seed and this must be accomplished by a complex search algorithm.

The position of the eigenvalues for [math]\displaystyle{ \beta =2\times 10^{4} }[/math] and [math]\displaystyle{ b=50 }[/math] are shown in Figure \ref{spectrum}. They are denoted by [math]\displaystyle{ \lambda _{n}, }[/math] where [math]\displaystyle{ n\in \mathbb{Z}\lt math\gt , and ordered by increasing real part, with }[/math]n=0</math> corresponding the eigenvalue with smallest absolute real part. From the picture and on physical grounds, it seems likely that there exist asymptotics for the eigenvalues, however this theory is not developed here.

Eigenfunctions

The eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math] associated with the eigenvalue [math]\displaystyle{ \lambda _{n}\lt math\gt are denoted by }[/math]\Phi ^{+}(\lambda _{n},x),[math]\displaystyle{ and those of }[/math]\mathcal{B} ^{\ast }[math]\displaystyle{ (the adjoint of }[/math]\mathcal{B})[math]\displaystyle{ associated with the eigenvalue }[/math] \lambda _{n}^{\ast }[math]\displaystyle{ are denoted by }[/math]\hat{\Phi}^{+}\left( \lambda _{n}^{\ast },x\right) </math>. That is,

[math]\displaystyle{ \mathcal{B}\Phi ^{+}\left( \lambda _{n},x\right) =\lambda _{n}\Phi ^{+}\left( \lambda _{n},x\right) }[/math]

and

[math]\displaystyle{ \mathcal{B}^{\ast }\hat{\Phi}^{+}\left( \lambda _{n}^{\ast },x\right) =\lambda _{n}^{\ast }\hat{\Phi}^{+}\left( \lambda _{n}^{\ast },x\right) . }[/math]

The eigenfunction [math]\displaystyle{ \Phi ^{+}\left( \lambda _{n},x\right) }[/math] can be written

[math]\displaystyle{ \Phi ^{+}\left( \lambda _{n},x\right) =\left( \begin{matrix}{c} \phi ^{+}\left( \lambda _{n},x\right) \\ i\zeta ^{+}\left( \lambda _{n},x\right) \end{matrix} \right) =\left( \begin{matrix}{c} \sum\limits_{j=1}^{6}\alpha _{j}e^{\mu _{j}\left( \lambda _{n}\right) x} \\ \sum\limits_{j=1}^{6}-\frac{\alpha _{j}\mu _{j}\left( \lambda _{n}\right) ^{2}}{\lambda _{n}}e^{\mu _{j}\left( \lambda _{n}\right) x} \end{matrix} \right) }[/math]

where [math]\displaystyle{ \mu _{j}\left( \lambda \right) }[/math] are the six roots of equation (\ref {cubicinlambda}) and the coefficients [math]\displaystyle{ \alpha _{j} }[/math] are found from the boundary conditions at the end of the plate (4) and the following condition. Since the scattering matrix is singular, there are no incoming wave from either direction. We use the condition that there is no incoming wave at [math]\displaystyle{ x=-b }[/math] and that the outgoing wave is of unit amplitude, i.e.

[math]\displaystyle{ \phi ^{+}\left( \lambda _{n},-b\right) =e^{i\lambda _{n}b},\;=\ = \left. \partial _{x}\phi ^{+}\left( \lambda _{n},x\right) \right| _{x=-b}=i\lambda _{n}e^{i\lambda _{n}b}. }[/math]

We do not use the condition that there is no outgoing wave at [math]\displaystyle{ x=b }[/math] because the system will become over determined. Therefore, the coefficients [math]\displaystyle{ \alpha _{j} }[/math] satisfy the linear equation

[math]\displaystyle{ \mathbf{M}\vec{a}\mathbf{=}\vec{b}, }[/math]

where

[math]\displaystyle{ \mathbf{M}=\left( \begin{matrix}{cccccc} \mu _{1}^{4}e^{-\mu _{1}b} & \mu _{2}^{4}e^{-\mu _{2}b} & \mu _{3}^{4}e^{-\mu _{3}b} & \mu _{4}^{4}e^{-\mu _{4}b} & \mu _{5}^{4}e^{-\mu _{5}b} & \mu _{6}^{4}e^{-\mu _{6}b} \\ \mu _{1}^{5}e^{-\mu _{1}b} & \mu _{2}^{5}e^{-\mu _{2}b} & \mu _{3}^{5}e^{-\mu _{3}b} & \mu _{4}^{5}e^{-\mu _{4}b} & \mu _{5}^{5}e^{-\mu _{5}b} & \mu _{6}^{5}e^{-\mu _{6}b} \\ \mu _{1}^{4}e^{\mu _{1}b} & \mu _{2}^{4}e^{\mu _{2}b} & \mu _{3}^{4}e^{\mu _{3}b} & \mu _{4}^{4}e^{\mu _{4}b} & \mu _{5}^{4}e^{\mu _{5}b} & \mu _{6}^{4}e^{\mu _{6}b} \\ \mu _{1}^{5}e^{\mu _{1}b} & \mu _{2}^{5}e^{\mu _{2}b} & \mu _{3}^{5}e^{\mu _{3}b} & \mu _{4}^{5}e^{\mu _{4}b} & \mu _{5}^{5}e^{\mu _{5}b} & \mu _{6}^{5}e^{\mu _{6}b} \\ e^{-\mu _{1}b} & e^{-\mu _{2}b} & e^{-\mu _{3}b} & e^{-\mu _{4}b} & e^{-\mu _{5}b} & e^{-\mu _{6}b} \\ \mu _{1}e^{-\mu _{1}b} & \mu _{2}e^{-\mu _{2}b} & \mu _{3}e^{-\mu _{3}b} & \mu _{4}e^{-\mu _{4}b} & \mu _{5}e^{-\mu _{5}b} & \mu _{6}e^{-\mu _{6}b} \end{matrix} \right) }[/math]

and [math]\displaystyle{ \vec{a} }[/math] and [math]\displaystyle{ \vec{b} }[/math] are given by

[math]\displaystyle{ \vec{a}=\left( \begin{matrix}{c} \alpha _{1} \\ \alpha _{2} \\ \alpha _{3} \\ \alpha _{4} \\ \alpha _{5} \\ \alpha _{6} \end{matrix} \right) ,\;\;\;\;\;\vec{b}=\left( \begin{matrix}{c} 0 \\ 0 \\ 0 \\ 0 \\ e^{i\lambda b} \\ i\lambda e^{i\lambda b} \end{matrix} \right) . }[/math]

The eigenfunctions for the adjoint operator are found similarly.

Figure \ref{eigfunctions} shows the real and imaginary parts of the eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math] for [math]\displaystyle{ n=1, }[/math] [math]\displaystyle{ 3 }[/math], 5, and [math]\displaystyle{ 7, }[/math] again with [math]\displaystyle{ \beta =2\times 10^{4}\lt math\gt and }[/math]b=50.</math> While the eigenfunctions do not have a simple shape, increasing oscillation is apparent as [math]\displaystyle{ n }[/math] increases.

Inner products

A biorthogonal system with respect to the energy inner product (\ref {energyinnerprod}) is formed by the eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math], [math]\displaystyle{ \Phi ^{+}\left( \lambda _{n},x\right) ,\lt math\gt and the eigenfunctions of }[/math]\mathcal{B} ^{\ast },[math]\displaystyle{ }[/math]\hat{\Phi}^{+}\left( \lambda _{n},x\right) </math>. To normalise the biorthogonal system, the inner product of [math]\displaystyle{ \Phi ^{+}\left( \lambda _{n},x\right) \lt math\gt and }[/math]\hat{\Phi}^{+}\left( \lambda _{n},x\right) </math> has to be determined. From the definition of the energy inner product (\ref {energyinnerprod})

[math]\displaystyle{ \begin{matrix} \left\langle \Phi \left( \lambda _{m},x\right) ,\hat{\Phi}\left( \lambda _{n},x\right) \right\rangle _{\mathcal{H}} &=&\int_{-b}^{b}\partial _{x}\phi ^{+}\left( \lambda _{m},x\right) \left( \partial _{x}\hat{\phi}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx (20) \\ &&+\int_{-b}^{b}(1+P)i\zeta ^{+}\left( \lambda _{m}^{\ast },x\right) \left( i \hat{\zeta}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx. \nonumber \end{matrix} }[/math]

The two integrals in (20) are considered separately. The first is

[math]\displaystyle{ \begin{matrix} &&\int_{-b}^{b}\partial _{x}\phi ^{+}\left( \lambda _{m},x\right) \left( \partial _{x}\hat{\phi}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx \\ &=&\int_{-b}^{b}\left( \sum_{j=1}^{6}\mu _{j}\left( \lambda _{m}\right) \alpha _{j}e^{\mu _{j}\left( \lambda _{m}\right) x}\right) \left( \sum_{k=1}^{6}\mu _{k}\left( \lambda _{m}\right) \alpha _{k}e^{\mu _{k}\left( \lambda _{n}\right) x}\right) dx \\ &=&\sum_{j=1}^{6}\sum_{k=1}^{6}\int_{-b}^{b}-\mu _{j}\left( \lambda _{m}\right) ^{2}\alpha _{j}e^{\mu _{j}\left( \lambda _{m}\right) x}\alpha _{k}e^{\mu _{k}\left( \lambda _{n}\right) x}dx \\ &=&\sum_{j=1}^{6}\sum_{k=1}^{6}-\mu _{j}\left( \lambda _{m}\right) ^{2}\alpha _{j}\alpha _{k}\left( \frac{e^{\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b} }{\mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) } \right) \nonumber \end{matrix} }[/math]

and the second is

[math]\displaystyle{ \begin{matrix} &&\int_{-b}^{b}(1+P)i\zeta ^{+}\left( \lambda _{m}^{\ast },x\right) \left( i \hat{\zeta}^{+}\left( \lambda _{n}^{\ast },x\right) \right) ^{\ast }dx \\ &=&\int_{-b}^{b}\left( \sum_{j=1}^{6}-\frac{\alpha _{j}}{\lambda _{m}} \left( \mu _{j}\left( \lambda _{m}\right) ^{2}+\beta \mu _{j}\left( \lambda _{m}\right) ^{6}\right) e^{\mu _{j}\left( \lambda _{m}\right) x}\right) \left( \sum_{k=1}^{6}-\frac{\alpha _{k}}{\lambda _{n}}\mu _{k}\left( \lambda _{n}\right) ^{2}e^{\mu _{k}\left( \lambda _{n}\right) x}\right) dx \\ &=&\sum_{j=1}^{6}\sum_{k=1}^{6}\int_{-b}^{b}\frac{\alpha _{j}\alpha _{k}}{ \lambda _{m}\lambda _{n}}\left( \mu _{j}\left( \lambda _{m}\right) ^{2}+\beta \mu _{j}\left( \lambda _{m}\right) ^{6}\right) \mu _{k}\left( \lambda _{n}\right) ^{2}e^{\mu _{j}\left( \lambda _{m}\right) x}e^{\mu _{k}\left( \lambda _{n}\right) x}dx \\ &=&\sum_{j=1}^{6}\sum_{k=1}^{6}\frac{\alpha _{j}\alpha _{k}}{\lambda _{m}\lambda _{n}}\left( \mu _{j}\left( \lambda _{m}\right) ^{2}+\beta \mu _{j}\left( \lambda _{m}\right) ^{6}\right) \mu _{k}\left( \lambda _{n}\right) ^{2}\left( \frac{e^{\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}}{\mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) }\right) \\ &=&\sum_{j=1}^{6}\sum_{k=1}^{6}\frac{\alpha _{j}\alpha _{k}}{\lambda _{m}\lambda _{n}}\left( -\lambda _{m}^{2}\right) \mu _{k}\left( \lambda _{n}\right) ^{2}\left( \frac{e^{\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}-e^{-\left( \mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) \right) b}}{\mu _{j}\left( \lambda _{m}\right) +\mu _{k}\left( \lambda _{n}\right) }\right) . \nonumber \end{matrix} }[/math]

Therefore the calculation of the inner product in equation (\ref {innerprodall}) does not require numerical integration.

Solution

By solving (18) using the eigenfunctions of [math]\displaystyle{ \mathcal{B} }[/math] and [math]\displaystyle{ \mathcal{B}^{\ast } }[/math] we find the evolution of the plate, from an initial displacement [math]\displaystyle{ U_{0}(x), }[/math] is

[math]\displaystyle{ U\left( x,t\right) =\sum_{n=-\infty }^{\infty }e^{i\lambda _{n}t}\frac{ \left\langle U_{0}\left( x\right) ,\hat{\Phi}\left( \lambda _{n},x\right) \right\rangle _{\mathcal{H}}}{\left\langle \Phi \left( \lambda _{n},x\right) ,\hat{\Phi}\left( \lambda _{n},x\right) \right\rangle _{\mathcal{H}}}\Phi \left( \lambda _{n},x\right) . (21) }[/math]

The inner product of [math]\displaystyle{ U_{0} }[/math] with the eigenfunction [math]\displaystyle{ \hat{\Phi}\left( \lambda _{n},x\right) }[/math] is the only quantity left to compute in (\ref {evolutionB}). This inner product is written

[math]\displaystyle{ \left\langle U_{0}\left( x\right) ,\hat{\Phi}\left( \lambda _{n},x\right) \right\rangle _{\mathcal{H}}=\sum\limits_{j=1}^{6}\left( \begin{matrix}{c} \alpha _{j}\dint_{-b}^{b}-\mu _{j}\left( \lambda _{n}\right) ^{2}e^{\mu _{j}\left( \lambda _{n}\right) x}\phi _{0}\left( x\right) dx+\left. \alpha _{j}\mu _{j}\left( \lambda _{n}\right) e^{\mu _{j}\left( \lambda _{n}\right) x}\phi _{0}\left( x\right) \right| _{-b}^{b} \\ \alpha _{j}\lambda _{n}\dint_{-b}^{b}e^{\mu _{j}\left( \lambda _{n}\right) x}i\zeta _{0}\left( x\right) dx \end{matrix} \right) (22) }[/math]

and the integrals must be evaluated by numerical integration. The solutions calculated using the Lax-Phillips scattering theory are identical to those found using the self-adjoint operator method and for this reason no further figures are shown.\pagebreak

Conclusions

The spectral theory of a linear thin plate floating on shallow water has been derived. Two spectral-theory solutions have been presented which determine the time-dependent motion of the thin plate. The first method was based on self-adjoint operator theory and the second method was based on the Lax-Phillips scattering. The self-adjoint method solved both the wave forcing and the free plate problem while the Lax-Phillips method only solved for a free plate. The eigenfunctions for the self-adjoint method are orthogonal and the eigenvalues are continuous and consist of all [math]\displaystyle{ \mathbb{R} , }[/math] which makes the calculation of the eigenvalues trivial. The Lax-Phillips method has discrete eigenvalues which must be calculated numerically and the system of eigenfunctions is biorthogonal. The advantage of the Lax-Phillips method is that the modes of vibration of the plate-water system and their frequency and rate of decay are found. While the relative speeds of the two methods depends of the exact way in which they are implemented, the Lax-Phillips method should be considerably faster if the eigenvalues have been determined.

The development of a spectral theory for more complicated hydroelastic problems remains a major challenge. While this theory must be more complicated than that presented here, many features can be expected to remain. For example, ohmatsuVLFS has shown that the single frequency solutions can be used to solve certain time-dependent problems and \cite {Hazard} have shown that modes, in which the solution can be expanded, exist for other hydroelastic systems.

{\Large Acknowledgments}

\begin{acknowledgment} I would like to thank the anonymous reviewers, Dr. Kathy Ruggerio, and Prof. James Sneyd for their very helpful comments. Also, Prof. Boris Pavlov for explaining the Lax-Phillips scattering. \pagebreak \end{acknowledgment}

\bibliographystyle{jfm} \bibliography{mike,others} \pagebreak

\begin{center} {\huge Figure Captions} \end{center}

\textsc{Figure} 1. Schematic diagram of a thin plate floating on shallow water and the coordinates and dimensions of the problem.

\textsc{Figure} 2. The evolution of the displacement due to a pulse travelling to the right for the times shown. The plate occupies the region [math]\displaystyle{ -50\leq x\leq 50\lt math\gt and is shown by the bold line. }[/math]\beta =2\times 10^{4},b=50. </math>

\textsc{Figure} 3. The evolution of the potential due to a pulse travelling to the right for the times shown. The plate occupies the region [math]\displaystyle{ -50\leq x\leq 50\lt math\gt and is shown by the bold line. }[/math]\beta =2\times 10^{4},b=50.</math>

\textsc{Figure} 4. The evolution of the displacement for a plate released at [math]\displaystyle{ t=0 }[/math] for the times shown. The plate occupies the region [math]\displaystyle{ -50\leq x\leq 50 }[/math] and is shown by the bold line. [math]\displaystyle{ \beta =2\times 10^{4},b=50. }[/math]

\textsc{Figure} 5. The evolution of the displacement for a plate released at [math]\displaystyle{ t=0 }[/math] for the times shown. The plate occupies the region [math]\displaystyle{ -50\leq x\leq 50 }[/math] and is shown by the bold line. [math]\displaystyle{ \beta =2\times 10^{4},b=50. }[/math]

\textsc{Figure} 6. The location of the first 19 eigenvalues [math]\displaystyle{ \lambda _{n} }[/math] of [math]\displaystyle{ \mathcal{B} }[/math] for [math]\displaystyle{ \beta =2\times 10^{4},\;b=50. }[/math]

\textsc{Figure} 7. The real (solid) and imaginary (dashed) parts of the resonance eigenfunctions for [math]\displaystyle{ n=1, }[/math] [math]\displaystyle{ 3 }[/math], 5, and [math]\displaystyle{ 7 }[/math] as shown. [math]\displaystyle{ \beta =2\times 10^{4},\;b=50. }[/math]


\pagebreak

\FRAME{ftFU}{5.5434in}{1.6544in}{0pt}{\Qcb{}}{\Qlb{shallow}}{shallow.eps}{ \special{language "Scientific Word";type "GRAPHIC";display "ICON";valid_file "F";width 5.5434in;height 1.6544in;depth 0pt;original-width 6.9297in;original-height 1.6544in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'shallow.eps';file-properties "XNPEU";}}

Some nonsense\pagebreak

\FRAME{ftF}{4.6164in}{3.6296in}{0pt}{}{\Qlb{incomingspecplot1}}{ incomingspecplot1.eps}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.6164in;height 3.6296in;depth 0pt;original-width 6.8632in;original-height 5.3869in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'incomingspecplot1.eps';file-properties "XNPEU";}}

\FRAME{ftF}{4.5602in}{3.6288in}{0pt}{}{\Qlb{incomingspecplot2}}{ incomingspecplot2.eps}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.5602in;height 3.6288in;depth 0pt;original-width 6.7793in;original-height 5.386in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'incomingspecplot2.eps';file-properties "XNPEU";}}

\FRAME{fthFU}{4.6155in}{3.6357in}{0pt}{\Qcb{{}}}{\Qlb{spectral1}}{ spectral1.eps}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.6155in;height 3.6357in;depth 0pt;original-width 6.8493in;original-height 5.3852in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'spectral1.eps';file-properties "XNPEU";}}

\FRAME{fthFU}{4.5887in}{3.608in}{0pt}{\Qcb{{}}}{\Qlb{spectral2}}{ spectral2.eps}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.5887in;height 3.608in;depth 0pt;original-width 6.8493in;original-height 5.3852in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'spectral2.eps';file-properties "XNPEU";}}

\FRAME{ftFU}{4.5982in}{3.7853in}{0pt}{\Qcb{{}}}{\Qlb{spectrum}}{spectrum.eps }{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.5982in;height 3.7853in;depth 0pt;original-width 6.8165in;original-height 5.5486in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'spectrum.eps';file-properties "XNPEU";}}

\FRAME{ftbpFU}{4.6164in}{3.7377in}{0pt}{\Qcb{{}}}{\Qlb{eigfunctions}}{ eigfunctions.eps}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "ICON";valid_file "F";width 4.6164in;height 3.7377in;depth 0pt;original-width 6.845in;original-height 5.5365in;cropleft "0";croptop "1";cropright "1";cropbottom "0";filename 'eigfunctions.eps';file-properties "XNPEU";}}