PDE : Can not solve Helmholtz equation

[sompongt]

PDE : Can not solve Helmholtz equation

(This is not a homework. I doing my research on numerical boundary integral. I need the analytical solution to compare the results with my computer program. I try to solve this equation, but it not success. I need urgent help.)

I working on anti-plane elasticity problem. The physical problem is described on a rectangular domain,
0 \leq x \leq a and 0 \leq y \leq b.

Let u_z = u(x,y;t) is the displacement function in z direction.

u_z = u(x,y;t) is the displacement function in z direction.
The boundary condition are specified by:

(1) u = 0 on the lower edge (y=0)

(2) \text{Traction} = 0 on the left edge (x=0)
(3) \text{Traction} = 0 on the right edge (x=a)
(4) \text{Traction} = P e^{i \omega t} on the top edge (y=b) in the z direction

This problem can be formulate in the frequency domain to obtain the Helmholtz equation, as below.

The problem seem not difficult. But I can not solve it. Could anybody help me?
Any suggestion are welcome.

Mathematical derivation start here. The problem statement describe as,

\begin{equation}
\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial x^2} + k^2 u = 0,
0 \leq x \leq a, 0 \leq y \leq b \text{...(1)}
\end{equation}

where
\begin{equation}
k^2 = \frac{\omega^2}{c^2}
\end{equation}

\omega = angular frequency of the excited force
c = speed of wave = function of material properties

subjected by the following boundary conditions,

\begin{align}
\frac{\partial u}{\partial x}(x=0, y) &= 0 \text{...(2)}\label{bc1}\\
\frac{\partial u}{\partial x}(x=a, y) &= 0 \text{...(3)}\label{bc2} \\
u(x, y = 0) &= 0 \text{...(4)}\label{bc3}\\
\frac{\partial u}{\partial y}(x, y=b) &= \text{constant}\label{bc4} = P \text{...(5)}
\end{align}

Solution:

By separation of variables, let
\begin{equation}
u(x,y) = X(x)Y(y)
\end{equation}

Then, working on the standard process,

\begin{align}
X'' Y + X Y'' + k^2 X Y = 0 \\
\frac{X''}{X} + \frac{Y''}{Y} + k^2 = 0
\end{align}

let

\begin{align}
\frac{X''}{X} &= -\alpha^2 &\rightarrow X''+ \alpha^2 X &= 0 \\
\frac{Y''}{Y} &= -\beta^2 &\rightarrow Y''+ \beta^2 Y &= 0
\end{align}

Relation between \alpha and \beta,
\begin{equation}
\alpha ^2 + \beta^2 = k^2 \label{abrelation} \text{...(6)}
\end{equation}

Solve for X(x) and Y(y)

\begin{align}
X(x) &= C_1 \cos \alpha x + C_2 \sin \alpha x \\
Y(y) &= C_3 \cos \beta y + C_4 \sin \beta y
\end{align}

So that,

\begin{align}
u(x,y) &= X(x)Y(y) \\
&= (C_1 \cos \alpha x + C_2 \sin \alpha x)(C_3 \cos \beta y + C_4 \sin \beta y) \label{reduce1}
\end{align}

Determine coefficients C_1, C_2, C_3, C_4 by using boundary conditions (BC) eq(2-5),

First, use boundary condition in eq(4):
\begin{align}
u(x,y=0) &= 0 \\
(C_1 \cos \alpha x + C_2 \sin \alpha x)(C_3 \cdot 1 + C_4 \cdot 0) &= 0 \\
C_3(C_1 \cos \alpha x + C_2 \sin \alpha x) &= 0
\end{align}

We obtain,
\begin{equation}
C_3 = 0
\end{equation}

So that, u(x,y) reduce to
\begin{align}
u(x,y) &= C_4 \sin \beta y (C_1 \cos \alpha x + C_2 \sin \alpha x) \\
&= \sin \beta y (C_4C_1 \cos \alpha x + C_4C_2 \sin \alpha x) \\
&= \sin \beta y (C_5 \cos \alpha x + C_6 \sin \alpha x) \label{reduce2} \\
\frac{\partial u}{\partial x} (x,y) &= \alpha \sin \beta y (-C_5 \sin \alpha x + C_6 \cos \alpha x)
\end{align}

Second, use boundary condition in eq(2):
\begin{align}
\frac{\partial u}{\partial x}(x=0,y) &= 0 \\
\alpha \sin \beta y (-C_5 \cdot 0 + C_6 \cdot 1) &= 0 \\
C_6 \cdot \alpha \sin \beta y &= 0
\end{align}

We obtain,

\begin{equation}
C_6 = 0
\end{equation}

So that, u(x,y) reduce to
\begin{align}
u(x,y) &= C_5 \cos \alpha x \sin \beta y \\
\frac{\partial u}{\partial x} (x,y) &= - \alpha C_5 \sin \alpha x \sin \beta y
\end{align}

Third step, use boundary condition in eq(3):

\begin{align}
\frac{\partial u}{\partial x}(x=a,y) &= 0 \\
-\alpha C_5 \sin \alpha a \sin \beta y = 0
\end{align}

Which can be conclude that,

\begin{align}
\sin \alpha a &= 0 \\
\alpha a &= n \pi \\
\alpha_n &= \frac{n \pi}{a}
\end{align}

After we find \alpha, we can determine \beta from eq(6)

\begin{align}
\beta_n^2 = k^2 - \alpha_n^2
\end{align}

At this point, we can represent u(x,y) as infinite series by using principle of superposition

\begin{align}
u(x,y) &= \sum_{n = 0}^{\infty} C_n \cos \alpha_n x \sin \beta_n y \\
\frac{\partial u}{\partial y} (x,y) &= \sum_{n = 0}^{\infty} \beta_n C_n \cos \alpha_n x \cos \beta_n y
\end{align}

We can determine the unknowns C_n by using the last boundary condition in eq(5)

\begin{align}
P &= \frac{\partial u}{\partial y} (x,y = b) \\
&= \sum_{n = 0}^{\infty} \beta_n C_n \cos \alpha_n x \cos \beta_n b \\
&= \sum_{n = 0}^{\infty} \underbrace{\left(C_n\beta_n \cos \beta_n b \right)}_{\text{constant} = \bar{C_n}} \cos \alpha_n x \\
&= \sum_{n = 0}^{\infty} \bar{C_n}\cos \alpha_n x \text{...(7)}
\end{align}

Consider \bar{C_n} as coefficient of Fourier cosine series,

\begin{align}
P &= \frac{c_0}{2} + \sum_{ n = 1 }^{\infty} c_n \cos \frac{n\pi x}{a} \\
c_n &= \frac{2}{a} \int_{0}^{a} P \cdot \cos \frac{n\pi x}{a} dx
\end{align}

It is not so hard to find that,
\begin{align}
c_0 &= 2P \\
c_n &= 0
\end{align}

For me, the problem arise here. How can I find \bar{C_n} in eq(7)? Could anybody help me? What point that I am wrong? I feel very headache. Please, please, please help me.

 

[hunt_mat]

Take a step back, you separated variables correctly but you must say that:
<br /> \frac{X&#039;&#039;}{X}+\frac{Y&#039;&#039;}{Y}+k^{2}=0\rightarrow \frac{X&#039;&#039;}{X}=\mu\rightarrow \frac{Y&#039;&#039;}{Y}+k^{2}=-\mu<br />
Start from this point.

 

[rajb245]

You basically got it. Your equation (7) is setting a constant P equal to a Fourier cosine series. The lowest order term in such a series is a constant, and you even correctly identified the value c_0=2P; if you take half that, this IS your \bar C_0 value. Working that back to your original C_0, you get C_0 = \frac{\bar C_0}{\beta_0 \cos(\beta_0 b)}=\frac{P}{k\cos(k b)}, and all the other terms in the series are zero. This gives your final solution, u(x,y)=\frac{P}{k\cos(k b)}\sin(k y).

Sanity check this to make sure it solves your problem:
1. Solves the PDE? Yes.
2. Satisfies the \frac{\partial}{\partial x} boundary conditions? Yes, it doesn't depend on x, so the x partials are zero everywhere.
3. Zero on the bottom boundary? Yes, the sine makes sure that happens.
4. Has the right derivative value at the top boundary? Yes, the constants are such that the function outward normal derivative is equal to P at the top.

Its unsurprising that the solution is independent of x. There is no forcing in that direction, since the top boundary is moving up and down in unison at the frequency \omega. This launches plane waves down the sheet. Any variation in that top boundary condition will launch more interesting waves that move side to side and interact with the side boundaries.

 

[sompongt]

You basically got it. Your equation (7) is setting a constant P equal to a Fourier cosine series. The lowest order term in such a series is a constant, and you even correctly identified the value c_0=2P; if you take half that, this IS your \bar C_0 value. Working that back to your original C_0, you get C_0 = \frac{\bar C_0}{\beta_0 \cos(\beta_0 b)}=\frac{P}{k\cos(k b)}, and all the other terms in the series are zero. This gives your final solution, u(x,y)=\frac{P}{k\cos(k b)}\sin(k y).

Sanity check this to make sure it solves your problem:
1. Solves the PDE? Yes.
2. Satisfies the \frac{\partial}{\partial x} boundary conditions? Yes, it doesn't depend on x, so the x partials are zero everywhere.
3. Zero on the bottom boundary? Yes, the sine makes sure that happens.
4. Has the right derivative value at the top boundary? Yes, the constants are such that the function outward normal derivative is equal to P at the top.

Its unsurprising that the solution is independent of x. There is no forcing in that direction, since the top boundary is moving up and down in unison at the frequency \omega. This launches plane waves down the sheet. Any variation in that top boundary condition will launch more interesting waves that move side to side and interact with the side boundaries.

rajb245, thank you very much.