Difference between revisions of "Variable Depth Shallow Water Wave Equation"

From WikiWaves
Jump to navigationJump to search
Line 83: Line 83:
 
Taking a separable solution <math>\ w (x,t) = \Tau (t) \hat{w} (x)</math>  gives the eigenvalue problem
 
Taking a separable solution <math>\ w (x,t) = \Tau (t) \hat{w} (x)</math>  gives the eigenvalue problem
 
<center>
 
<center>
<math>\partial_x \left( c(x)^2 \partial_x\hat{w} \right) = \lambda\hat{w} \quad (2)</math>
+
<math>\partial_x \left( h(x) \partial_x\hat{w} \right) = -\kappa^{2}\hat{w} \quad (1)</math>
 
</center>
 
</center>
  
Given boundary conditions <math>\hat{w} (0) = a</math> and <math>\hat{w} (1) = b</math> we can take
+
Given boundary conditions <math>\hat{w} (0) = a</math> and <math>\hat{w} (1) = b</math> we can take <math>\hat{w} = (b-a)x + a + u </math> With <math> u  </math>  satisfying <math> u (0) = u (1) = 0</math>
 +
 
 +
Substituting this form into (1) gives
 +
<center>
 +
<math> (b-a)\partial_xh(x)+\partial_x(h(x)\partial_xu) = -\kappa^{2}\left((b-a)x+a+u\right)</math>
 +
</center>
 +
Or, on rearranging
 +
<center>
 +
<math> \partial_x(h(x)\partial_xu)+\kappa^{2}u = -(b-a)\partial_xh(x)-\kappa^{2}\left((b-a)x+a\right) \quad (2)</math>
 +
</center>
 +
Now consider the homogenous Sturm-Liouville problem for u
 +
<center>
 +
<math> \partial_x(h(x)\partial_xu)+\lambda u = 0\quad u(0)=u(1)=0 \quad (3)</math>
 +
</center>
 +
By Sturm-Liouville theory this has an infinite set of eigenvalues <math> \lambda_k </math> with corresponding eigenfunctions <math> u_k </math>.  Also since <math> u_k(0)=u_k(1)=0\quad \forall k </math> Each <math> u_k</math> can be expanded as a fourier series in terms of sine functions.
 +
<center>
 +
<math> u_k = \sum_{n=1}^{\infty} a_{n,k}\sin(n\pi x)</math>
 +
</center>
 +
Transforming (3) into the equivalent variational problem gives
 
<center>
 
<center>
<math>\hat{w} = (b-a)x + a + u \quad (3)</math>
+
<math> J[u] = \int_{0}^{1}\,hu'^{2}-\lambda u^{2} \, dx \quad (4) </math>
 
</center>
 
</center>
 
+
Substituting the fourier expansion into (4) implies J must be stationary at <math> \frac{\partial J}{\partial a_{n}}=0 \quad \forall n </math>
With <math> u = \sum_{n=1}^{N} a_n \sin (n \pi x)</math> a series solution satisfying <math> u (0) = u (1) = 0</math>
 
 
 
We wish to solve for <math> a_n </math> given <math> a </math> and <math> b </math>. Equation (2) can be transformed into the Sturm-Liouville problem of minimizing the functional
 
 
<center>
 
<center>
<math> J [u] = \int_{0}^{1} c^2 \left(\frac{d \hat{w}}{d x}\right)^2 + \lambda \hat{w}^2 \,\mathrm{d}x</math>
+
<math> \frac{\partial J}{\partial a_{n}}=\int_{0}^{1}\,hn\pi \cos(n\pi x)\sum_{m=1}^{\infty} a_{m}\cos(m\pi x)\,dx-\frac{\lambda} {2}a_{n}=0</math>
 
</center>
 
</center>
which is minimal exactly when <math>\partial_{a_n} J = 0 \quad \forall a_n </math>
+
By defining a vector <math> \textbf{a} = \left(a_{n}\right)</math> and a matrix <math>M_{(n,m)} = 2\int_{0}^{1}\,hnm\pi^{2} \cos(n\pi x)\cos(m\pi x)\,dx </math> we have the linear system <math> M\textbf{a} = \lambda\textbf{a} </math> which returns the eigenvalues and eigenvectors of equation (3), with eigenvectors <math> \textbf{a} </math> representing coefficient vectors of the fourier expansions of eigenfunctions.
 
 
Substituting in <math>\hat{w} = (b-a)x + a + \sum_{n=1}^{N} a_n \sin (n \pi x)</math> and <math>\frac {\mathrm{d} \hat{w}}{\mathrm{d} x} = (b-a) + \sum_{n=1}^{N} a_n n \pi \cos (n \pi x)</math>
 
 
 
gives
 
  
 +
If we now construct <math> u = \sum_{1}^{\infty} b_k u_k </math> and substitute this into equation (2) we get
 
<center>
 
<center>
<math> \partial_{a_n} J = \int_{0}^{1} c^2 \left( n \pi \cos (n \pi x) \left((b-a) + \sum_{m=1}^{N} a_m m \pi \cos (m \pi x) \right) \right) \, \mathrm{d}x  + \int_{0}^{1} \lambda \left( \sin (n \pi x) \left((b-a)x + a + \sum_{m=1}^{N} a_m \sin (m \pi x) \right) \right) \, \mathrm{d}x = 0 \quad (4) </math>
+
<math> \sum_{1}^{\infty} (\kappa^{2}-\lambda_k) b_k u_k = -(b-a)\partial_xh(x)-\kappa^{2}\left((b-a)x+a\right) \quad (5) </math>
 
</center>
 
</center>
 
+
And defining the RHS of equation (5) as <math> f(x) </math>, a known function, we can retrieve the coefficients <math> b_k </math> by integrating against <math> u_k </math>
Since <math> \sin(n \pi x) </math> and <math> \sin(m \pi x) </math> are orthogonal  the second integral in (4) can be calculated.
 
 
 
<math> \int_{0}^{1} \lambda \left( \sin (n \pi x) \left((b-a)x + a + \sum_{m=1}^{N} a_m \sin (m \pi x) \right) \right) \, \mathrm{d}x =\lambda \int_{0}^{1} \left( (b-a)x+a \right) \sin (n \pi x) \, \mathrm{d}x + \frac{\lambda a_n}{2}</math>
 
 
 
and
 
 
 
<math> \int_{0}^{1} \lambda \left( \sin (n \pi x) \left((b-a)x + a + \sum_{m=1}^{N} a_m \sin (m \pi x) \right) \right) \, \mathrm{d}x </math>
 
<math>= \lambda \left((b-a) \int_{0}^{1} x \sin(n \pi x) \, dx + a \int_{0}^{1} \sin(n \pi x) \, \mathrm{d}x \right) </math>
 
<math>= \frac{\lambda (b-a)}{(n \pi)^2} \Big[ \sin(n \pi x) -n \pi x \cos( n \pi x) \Big]_{0}^{1} - \frac{a \lambda}{n \pi} \Big[\cos(n \pi x) \Big]_{0}^{1} </math>
 
<math> = \frac{\lambda (b-a)}{(n \pi)^2} \left(n \pi (-1)^{n+1} \right) + \frac{a \lambda}{n \pi} \left(1 - (-1)^n \right) </math>
 
<math> = \frac{\lambda}{n \pi} \left( b (-1)^{n+1}+a \right) </math>
 
 
 
So equation (4) can be written as
 
<math>\int_{0}^{1} c^2 \left( n \pi \cos (n \pi x) \left((b-a) + \sum_{m=1}^{N} a_m m \pi \cos (m \pi x) \right) \right) \, \mathrm{d}x + \frac{\lambda}{n \pi} \left( b (-1)^{n+1}+a \right) + \frac{\lambda a_n}{2}= 0 </math>
 
 
 
The remaining integral can be split into two parts and with the sum taken outside the integral we obtain
 
<math>\int_{0}^{1} c^2 (b-a) n \pi \cos (n \pi x)  + \sum_{m=1}^{N} \left( \int_{0}^{1} c^2 n m \pi^2 \cos(n \pi x) \cos (m \pi x) \, \mathrm{d}x \right) a_m + \frac{\lambda}{n \pi} \left( b (-1)^{n+1}+a \right) + \frac{\lambda a_n}{2}= 0 </math>
 
 
 
or, on rearranging
 
<math> \sum_{m=1}^{N} \left( \int_{0}^{1} c^2 n m \pi^2 \cos(n \pi x) \cos (m \pi x) \, \mathrm{d}x \right) a_m+ \frac{\lambda a_n}{2}= -\int_{0}^{1} c^2 (b-a) n \pi \cos (n \pi x) - \frac{\lambda}{n \pi} \left( b (-1)^{n+1}+a \right) </math>
 
 
 
Now if we form the vector <math> \mathbf{a} = \begin{pmatrix} a_1 \\ \vdots \\ a_N \end{pmatrix} </math> and recalling that we have the above expression for all n, we can write the above as a matrix multiplication of <math> \mathbf{a} </math>
 
 
<center>
 
<center>
<math> \mathrm{M} \mathbf{a} = \mathbf{f} </math>
+
<math> b_k = \frac{\int_{0}^{1}\,f u_k\,dx}{(\kappa^{2}-\lambda_k) \int_{0}^{1}\, u_{k}^{2}\,dx} </math>
 
</center>
 
</center>
  
with
+
Also to find the coefficients <math> c_k </math> of the fourier expansion of u are just <math> \sum_{k=1}^{\infty}a_{n,k}b_{k} </math> with <math> a_{n.k} </math> being the <math>n</math>th coefficient of the <math>k</math>th eigenfunction of the Sturm-Liouville problem.
<math> \mathrm{M}_{(n,m)} = \int_{0}^{1} c^2  n m \pi^2 \cos(n \pi x) \cos (m \pi x) \, \mathrm{d}x + \delta_{nm} \frac{\lambda}{2} \quad \delta_{nm} = \left\{ \begin{matrix} 1 &  \mathrm{if} \quad n = m  \\ 0 & \mathrm{if} \quad n \neq m  \end{matrix} \right. </math>
 
 
 
and
 
<math> \mathbf{f} = \begin{pmatrix} f_1 \\ \vdots \\ f_N \end{pmatrix} </math> with <math> f_n =  -\int_{0}^{1} c^2 (b-a) n \pi \cos (n \pi x) + \frac{\lambda}{n \pi} \left( b (-1)^n-a \right) </math>
 
 
 
So, given a suitable <math> c(x) </math> and boundary conditions <math>\hat{w} (0) = a</math> and <math>\hat{w} (1) = b</math> we have a system of linear equations that can be solved to give the coefficients <math> a_n </math> which in turn define the function <math>\hat{w} </math>.
 
  
 
[[Category:Shallow Depth]]
 
[[Category:Shallow Depth]]

Revision as of 22:49, 7 January 2009

Introduction

We consider here the problem of waves reflected by a region of variable depth in an otherwise uniform depth region assuming the equations of Shallow Depth.

Equations

We begin with the shallow depth equation

[math]\displaystyle{ \rho(x)\partial_t^2 \zeta = \partial_x \left(h(x) \partial_x \zeta \right). }[/math]

subject to the initial conditions

[math]\displaystyle{ \zeta_{t=0} = \zeta_0(x)\,\,\,{\rm and}\,\,\, \partial_t\zeta_{t=0} = \partial_t\zeta_0(x) }[/math]

where [math]\displaystyle{ \zeta }[/math] is the displacement, [math]\displaystyle{ \rho }[/math] is the string density and [math]\displaystyle{ h(x) }[/math] is the variable depth (note that we are unifying the variable density string and the wave equation in variable depth because the mathematical treatment is identical).

Waves in a finite basin

We consider the problem of waves in a finite basin [math]\displaystyle{ 0\lt x\lt 1 }[/math]. At the edge of the basin the boundary conditions are

[math]\displaystyle{ \left.\partial_x \zeta\right|_{x=0} = \left.\partial_x \zeta\right|_{x=1} =0 }[/math]

.

We solve the equations by expanding in the modes for the basin which satisfy

[math]\displaystyle{ \partial_x \left(h(x) \partial_x \zeta_n \right) = -\lambda_n, }[/math]

normalised so that

[math]\displaystyle{ \int_0^1 \zeta_n \zeta_m = \delta_{mn}. }[/math]

The solution is then given by

[math]\displaystyle{ \zeta(x,t) = \sum_{n=0}^{\infty} \left(\int_0^1 \zeta_n(x^\prime) \zeta_0 (x^\prime) dx^\prime \right) \zeta_n(x) \cos(\sqrt{\lambda_n} t ) }[/math]
[math]\displaystyle{ + \sum_{n=1} ^{\infty} \left(\int_0^1 \zeta_n(x^\prime) \partial_t\zeta_0 (x^\prime) dx^\prime \right) \zeta_n(x) \frac{\sin(\sqrt{\lambda_n} t )}{\sqrt{\lambda_n}} }[/math]

where we have assumed that [math]\displaystyle{ \lambda_0 = 0 }[/math].

Calculation of [math]\displaystyle{ \zeta_n }[/math]

We can calculate the eigenfunctions [math]\displaystyle{ \zeta_n }[/math] by an expansion in the modes for the case of uniform depth. We use the Rayleigh-Ritz method. The eigenfunctions are local minimums of

[math]\displaystyle{ J[\zeta] = \int_0^1 \frac{1}{2}\left\{ \left(h(x) \partial_x \zeta\right)^2 - \lambda \zeta^2 \right\} }[/math]

subject to the boundary conditions that the normal derivative vanishes (where [math]\displaystyle{ \lambda }[/math] is the eigenvalue).

We expand the displacement in the eigenfunctions for constant depth [math]\displaystyle{ h=1 }[/math]

[math]\displaystyle{ \zeta = \sum_{n=1}^{N} a_n \psi_n(x) }[/math]

where

[math]\displaystyle{ \psi_n = \sqrt{2} \cos( (n-1) \pi x),\,\,n\ne 1 }[/math]

[math]\displaystyle{ \psi_0 = 1,\, }[/math]

and substitute this expansion into the variational equation we obtain

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

where the elements of the matrix M are

[math]\displaystyle{ m_{mn} = \int_0^1 \left\{ \left(\partial_x \psi_m h(x) \partial_x \psi_n\right) \right\} }[/math]

Matlab code

Waves in an infinite basin

We assume that the depth is constant and equal to one outside the region [math]\displaystyle{ 0\lt x\lt 1 }[/math]. We can therefore write the wave as

Solution using Separation of Variables

Taking a separable solution [math]\displaystyle{ \ w (x,t) = \Tau (t) \hat{w} (x) }[/math] gives the eigenvalue problem

[math]\displaystyle{ \partial_x \left( h(x) \partial_x\hat{w} \right) = -\kappa^{2}\hat{w} \quad (1) }[/math]

Given boundary conditions [math]\displaystyle{ \hat{w} (0) = a }[/math] and [math]\displaystyle{ \hat{w} (1) = b }[/math] we can take [math]\displaystyle{ \hat{w} = (b-a)x + a + u }[/math] With [math]\displaystyle{ u }[/math] satisfying [math]\displaystyle{ u (0) = u (1) = 0 }[/math]

Substituting this form into (1) gives

[math]\displaystyle{ (b-a)\partial_xh(x)+\partial_x(h(x)\partial_xu) = -\kappa^{2}\left((b-a)x+a+u\right) }[/math]

Or, on rearranging

[math]\displaystyle{ \partial_x(h(x)\partial_xu)+\kappa^{2}u = -(b-a)\partial_xh(x)-\kappa^{2}\left((b-a)x+a\right) \quad (2) }[/math]

Now consider the homogenous Sturm-Liouville problem for u

[math]\displaystyle{ \partial_x(h(x)\partial_xu)+\lambda u = 0\quad u(0)=u(1)=0 \quad (3) }[/math]

By Sturm-Liouville theory this has an infinite set of eigenvalues [math]\displaystyle{ \lambda_k }[/math] with corresponding eigenfunctions [math]\displaystyle{ u_k }[/math]. Also since [math]\displaystyle{ u_k(0)=u_k(1)=0\quad \forall k }[/math] Each [math]\displaystyle{ u_k }[/math] can be expanded as a fourier series in terms of sine functions.

[math]\displaystyle{ u_k = \sum_{n=1}^{\infty} a_{n,k}\sin(n\pi x) }[/math]

Transforming (3) into the equivalent variational problem gives

[math]\displaystyle{ J[u] = \int_{0}^{1}\,hu'^{2}-\lambda u^{2} \, dx \quad (4) }[/math]

Substituting the fourier expansion into (4) implies J must be stationary at [math]\displaystyle{ \frac{\partial J}{\partial a_{n}}=0 \quad \forall n }[/math]

[math]\displaystyle{ \frac{\partial J}{\partial a_{n}}=\int_{0}^{1}\,hn\pi \cos(n\pi x)\sum_{m=1}^{\infty} a_{m}\cos(m\pi x)\,dx-\frac{\lambda} {2}a_{n}=0 }[/math]

By defining a vector [math]\displaystyle{ \textbf{a} = \left(a_{n}\right) }[/math] and a matrix [math]\displaystyle{ M_{(n,m)} = 2\int_{0}^{1}\,hnm\pi^{2} \cos(n\pi x)\cos(m\pi x)\,dx }[/math] we have the linear system [math]\displaystyle{ M\textbf{a} = \lambda\textbf{a} }[/math] which returns the eigenvalues and eigenvectors of equation (3), with eigenvectors [math]\displaystyle{ \textbf{a} }[/math] representing coefficient vectors of the fourier expansions of eigenfunctions.

If we now construct [math]\displaystyle{ u = \sum_{1}^{\infty} b_k u_k }[/math] and substitute this into equation (2) we get

[math]\displaystyle{ \sum_{1}^{\infty} (\kappa^{2}-\lambda_k) b_k u_k = -(b-a)\partial_xh(x)-\kappa^{2}\left((b-a)x+a\right) \quad (5) }[/math]

And defining the RHS of equation (5) as [math]\displaystyle{ f(x) }[/math], a known function, we can retrieve the coefficients [math]\displaystyle{ b_k }[/math] by integrating against [math]\displaystyle{ u_k }[/math]

[math]\displaystyle{ b_k = \frac{\int_{0}^{1}\,f u_k\,dx}{(\kappa^{2}-\lambda_k) \int_{0}^{1}\, u_{k}^{2}\,dx} }[/math]

Also to find the coefficients [math]\displaystyle{ c_k }[/math] of the fourier expansion of u are just [math]\displaystyle{ \sum_{k=1}^{\infty}a_{n,k}b_{k} }[/math] with [math]\displaystyle{ a_{n.k} }[/math] being the [math]\displaystyle{ n }[/math]th coefficient of the [math]\displaystyle{ k }[/math]th eigenfunction of the Sturm-Liouville problem.