Difference between revisions of "Infinite Array Green Function"

From WikiWaves
Jump to navigationJump to search
m (Reverted edits by DronvArroc (DronvArroc); changed back to last version by Meylan)
 
(31 intermediate revisions by 4 users not shown)
Line 62: Line 62:
 
  \phi(\mathbf{x}) = \phi^{\rm in}(\mathbf{x})
 
  \phi(\mathbf{x}) = \phi^{\rm in}(\mathbf{x})
 
  +\int_{\Delta_{0}}  
 
  +\int_{\Delta_{0}}  
\sum_{m=-\infty}^{\infty} \left(G^{\mathbf{P}} {n_\xi}(\mathbf{x},\xi+(0,ml,0))e^{im\sigma l} \phi(\xi)  
+
\sum_{m=-\infty}^{\infty} \left(G^{\mathbf{P}}_{n_\xi}(\mathbf{x},\xi+(0,ml,0))e^{im\sigma l} \phi(\xi)  
 
- G^{\mathbf{P}} (\mathbf{x},\xi)e^{im\sigma l} \phi_{n_\xi}(\xi) \right) d\xi
 
- G^{\mathbf{P}} (\mathbf{x},\xi)e^{im\sigma l} \phi_{n_\xi}(\xi) \right) d\xi
 
</math></center>
 
</math></center>
Line 75: Line 75:
  
 
The spatial representation of the periodic Green function given by
 
The spatial representation of the periodic Green function given by
equation ((G_periodic)) is slowly convergent
+
equation is slowly convergent
 
and in the far field the terms decay in  
 
and in the far field the terms decay in  
 
magnitude like <math>O(n^{-1/2})</math>. In this section we
 
magnitude like <math>O(n^{-1/2})</math>. In this section we
 
show how to accelerate the convergence. We begin with the asymptotic
 
show how to accelerate the convergence. We begin with the asymptotic
approximation of the three-dimensional Green function ((green))
+
approximation of the Three-dimensional [[Free-Surface Green Function]]
 
far from the source point,
 
far from the source point,
 
<center><math>
 
<center><math>
  G(\mathbf{x},0;\mbox{\boldmath<math>\xi</math>}) \sim -\frac{ik}{2}
+
  G(\mathbf{x},\xi) \sim -\frac{ik}{2}
  \,H_{0}( k |\mathbf{x}-\mbox{\boldmath<math>\xi</math>}|),  
+
  \,H_{0}( k |\mathbf{x}-\xi|),  
  \qquad \mbox{as <math>|\mathbf{x}-\mbox{\boldmath</math>\xi<math>}| \to \infty</math>},
+
  |\mathbf{x}-\xi| \to \infty
  (asyG)
 
 
</math></center>
 
</math></center>
[[Weh_Lait]] where <math>H_0 \equiv H_{0}^{(1)}</math> is the Hankel function  
+
[[Wehausen and Laitone 1960]] where <math>H_0 \equiv H_{0}^{(1)}</math> is the Hankel function  
of the first kind of order zero [[abr_ste]]. In Linton [[Linton98]]  
+
of the first kind of order zero [[Abramowitz and Stegun 1964]]. In Linton [[Linton 1998]]  
 
various methods were described in which the convergence of the periodic  
 
various methods were described in which the convergence of the periodic  
 
Green functions was improved. One such method, which suits the particular  
 
Green functions was improved. One such method, which suits the particular  
 
problem being considered here, involves writing the periodic
 
problem being considered here, involves writing the periodic
 
Green function as
 
Green function as
<center><math>\begin{matrix}
+
<center><math>
  G_{\mathbf{P}} (\mathbf{x};\mbox{\boldmath<math>\xi</math>})  
+
  G_{\mathbf{P}} (\mathbf{x};\xi)  
  &\! = &\!\!\!\! \displaystyle{\sum_{m=-\infty}^{\infty}} \!
+
  = \sum_{m=-\infty}^{\infty}
\left[ G\left(\mathbf{x};\mbox{\boldmath<math>\xi</math>}+(0,ml)\right)
+
\left[
+ \frac{ik}{2} H_{0} \Big(k\sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big)
+
G\left(\mathbf{x};\xi)+(0,ml)\right)  
\right] \exp^{im\sigma l}
+
+ \frac{ik}{2} H_{0}  
\\
+
\Big(k\sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big) e^{im\sigma l}
  &\! &\! -\displaystyle{\sum_{m=-\infty}^{\infty}}
+
\right]
 +
  -\sum_{m=-\infty}^{\infty}
 
  \frac{ik}{2}H_{0} \Big(k\sqrt{ (X+cl)^2 + Y_{m}^2 }\Big)
 
  \frac{ik}{2}H_{0} \Big(k\sqrt{ (X+cl)^2 + Y_{m}^2 }\Big)
  \exp^{im\sigma l}
+
  e^{im\sigma l}
  (near_accelerated)
+
</math></center>
\end{matrix}</math></center>
 
 
where <math>c</math> is a numerical smoothing parameter, introduced to avoid the  
 
where <math>c</math> is a numerical smoothing parameter, introduced to avoid the  
singularity at <math>\mathbf{x} = \mbox{\boldmath</math>\xi<math>}</math> in the Hankel  
+
singularity at <math>\mathbf{x} = \xi</math> in the Hankel  
 
function and
 
function and
 
<center><math>
 
<center><math>
X = x-\xi,\quad \mbox{and} \quad
+
X = x-\xi,\quad \mathrm{and} \quad
 
Y_{m} = (y-\eta)-ml.
 
Y_{m} = (y-\eta)-ml.
 
</math></center>
 
</math></center>
Furthermore we use the fact that second slowly convergent sum in ((near_accelerated)) can be transformed to
+
Furthermore we use the fact that second slowly convergent sum can be transformed to
<center><math> (help_accelerated)
+
<center><math>  
  -\frac{i}{l} \sum_{m=-\infty}^{\infty}
+
-\sum_{m=-\infty}^{\infty}
  \frac{\exp^{ik \mu_{m} |X+c| }\,\exp^{i \sigma_{m} Y_0}}{\mu_{m}}
+
\frac{ik}{2}H_{0} \Big(k\sqrt{ (X+cl)^2 + Y_{m}^2 }\Big)
 +
  e^{im\sigma l}
 +
-\frac{i}{l} \sum_{m=-\infty}^{\infty}
 +
  \frac{e^{ik \mu_{m} |X+c| }\,e^{i \sigma_{m} Y_0}}{\mu_{m}}
 
</math></center>
 
</math></center>
[[Linton98,Jorgenson90,Singh90]] where
+
[[Linton 1998]] where
 
<math>\sigma_{m} = \sigma + 2 m \pi/l</math> and
 
<math>\sigma_{m} = \sigma + 2 m \pi/l</math> and
 
<center><math>
 
<center><math>
Line 125: Line 127:
 
where the positive real or positive imaginary part of  
 
where the positive real or positive imaginary part of  
 
the square root is taken.
 
the square root is taken.
Combining equations ((near_accelerated)) and ((help_accelerated))
+
Combining these equations  
 
we obtain the accelerated version of the periodic Green function  
 
we obtain the accelerated version of the periodic Green function  
<center><math>\begin{matrix}
+
<center><math>
  G_{\mathbf{P}} (\mathbf{x};\mbox{\boldmath<math>\xi</math>})  
+
  G_{\mathbf{P}} (\mathbf{x};\xi)  
  &\! = &\!\!\!\! \displaystyle{\sum_{m=-\infty}^{\infty}} \!
+
  = \sum_{m=-\infty}^{\infty}
  \left[ G\left(\mathbf{x};\mbox{\boldmath<math>\xi</math>}+(0,ml)\right)   
+
  \left[ G\left(\mathbf{x};\xi+(0,ml)\right)   
 
  + \frac{ik}{2} H_{0} \Big(k \sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big)
 
  + \frac{ik}{2} H_{0} \Big(k \sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big)
  \right] \exp^{im\sigma l}
+
  \right] e^{im\sigma l}
\\
+
  -\frac{i}{l}\sum_{m=-\infty}^{\infty}
&\! &\! -\frac{i}{l}\displaystyle{\sum_{m=-\infty}^{\infty}
+
  \frac{e^{ik \mu_{m}|X+cl| }e^{i\sigma_{m}Y_0}}{\mu_{m}}.
  \frac{\exp^{ik \mu_{m}|X+cl| }\exp^{i\sigma_{m}Y_0}}{\mu_{m}}}.
+
</math></center>
  (Gp_fast)
+
The convergence of the two sums depends on the value
\end{matrix}</math></center>
 
The convergence of the two sums in ((Gp_fast)) depends on the value
 
 
of <math>c</math>. For small <math>c</math> the first sum converges rapidly while the second converges
 
of <math>c</math>. For small <math>c</math> the first sum converges rapidly while the second converges
 
slowly. For large <math>c</math> the second sum converges rapidly while the first converges
 
slowly. For large <math>c</math> the second sum converges rapidly while the first converges
Line 144: Line 144:
 
The smoothing parameter <math>c</math> must be carefully chosen to balance these two
 
The smoothing parameter <math>c</math> must be carefully chosen to balance these two
 
effects. Of course, the convergence also depends strongly on how close together the  
 
effects. Of course, the convergence also depends strongly on how close together the  
points <math>\mathbf{x}</math> and <math>\mbox{\boldmath</math>\xi<math>}</math> are.
+
points <math>\mathbf{x}</math> and <math>\xi</math> are.
  
 
Note that some special combinations of wavelength <math>\lambda</math> and angle
 
Note that some special combinations of wavelength <math>\lambda</math> and angle
 
of incidence <math>\theta</math> cause the periodic Green function to diverge
 
of incidence <math>\theta</math> cause the periodic Green function to diverge
( [[Scott98]],[[Jorgenson90]]).  This singularity is closely
+
( [[Scott 1998]]).  This singularity is closely
 
related to the diffracted waves and will be explained shortly.
 
related to the diffracted waves and will be explained shortly.
  
Line 154: Line 154:
  
 
We begin with the accelerated periodic Green function, equation
 
We begin with the accelerated periodic Green function, equation
((Gp_fast)) setting <math>c=0</math> and considering the case when <math>X</math> is large
+
setting <math>c=0</math> and considering the case when <math>X</math> is large
 
(positive or negative). We also note that for <math>m</math> sufficiently small
 
(positive or negative). We also note that for <math>m</math> sufficiently small
 
or large <math>i\mu_m</math> will be negative and the corresponding terms will
 
or large <math>i\mu_m</math> will be negative and the corresponding terms will
 
decay. Therefore
 
decay. Therefore
<center><math> (G_p_accelerate)
+
<center><math>
  G_{\mathbf{P}} (\mathbf{x};\mbox{\boldmath<math>\xi</math>})
+
  G_{\mathbf{P}} (\mathbf{x};\xi)
  \sim - \frac{i}{l} \sum_{m=-M}^{N} \frac{\exp^{ik\mu_{m}|X|}\, \exp^{i\sigma_{m}Y_0}}
+
  \sim - \frac{i}{l} \sum_{m=-M}^{N} \frac{e^{ik\mu_{m}|X|}\, e^{i\sigma_{m}Y_0}}
  {\mu_{m}}, \qquad \mbox{as <math>X \to \pm \infty</math>}
+
  {\mu_{m}}, X \to \pm \infty
  (new_G_p_trunc)
 
 
</math></center>
 
</math></center>
 
where the integers <math>M</math> and <math>N</math>  satisfy the following
 
where the integers <math>M</math> and <math>N</math>  satisfy the following
Line 169: Line 168:
 
\left.
 
\left.
 
\begin{matrix}
 
\begin{matrix}
[c]{c}
 
 
\sigma_{-M-1}<-k<\sigma_{-M},\\
 
\sigma_{-M-1}<-k<\sigma_{-M},\\
 
\sigma_{N}<k<\sigma_{N+1}.
 
\sigma_{N}<k<\sigma_{N+1}.
 
\end{matrix}
 
\end{matrix}
\right\} (cut_off)
+
\right\}  
 
</math></center>
 
</math></center>
Equations ((M_N)) can be written as  
+
These equations can be written as  
 
<center><math>
 
<center><math>
 
\frac{l}{2\pi}\left(\sigma+k-2\pi \right) < M < \frac{l}{2\pi}\left(
 
\frac{l}{2\pi}\left(\sigma+k-2\pi \right) < M < \frac{l}{2\pi}\left(
 
\sigma+k \right),
 
\sigma+k \right),
(diffracted_M)
 
 
</math></center>
 
</math></center>
 
and
 
and
Line 185: Line 182:
 
  \frac{l}{2\pi}\left( k - \sigma \right) > N > \frac{l}{2\pi}
 
  \frac{l}{2\pi}\left( k - \sigma \right) > N > \frac{l}{2\pi}
 
  \left( k-\sigma - 2\pi \right)
 
  \left( k-\sigma - 2\pi \right)
  (diffracted_N)
 
 
</math></center>
 
</math></center>
[[Linton98]].
+
[[Linton 1998]].
 
It is obvious that <math>G_{\mathbf{P}}</math> will diverge if <math>\sigma_m = \pm k</math>;
 
It is obvious that <math>G_{\mathbf{P}}</math> will diverge if <math>\sigma_m = \pm k</math>;
 
these values correspond to cut-off frequencies which are an expected
 
these values correspond to cut-off frequencies which are an expected
Line 196: Line 192:
 
The diffracted waves are the plane waves which are observed as <math>x \to \pm
 
The diffracted waves are the plane waves which are observed as <math>x \to \pm
 
\infty</math>. Their amplitude and form are obtained by substituting the limit
 
\infty</math>. Their amplitude and form are obtained by substituting the limit
of the periodic Green function ((new_G_p_trunc)) as  <math>x\to\pm\infty</math>
+
of the periodic Green function as  <math>x\to\pm\infty</math>
into the boundary integral equation for the potential ((bem_eq_2)).
+
into the boundary integral equation for the potential.
 
This gives us
 
This gives us
 
<center><math>
 
<center><math>
 
  \lim_{x\to\pm\infty}
 
  \lim_{x\to\pm\infty}
 
  \phi^{s} ( \mathbf{x},0 )  = - \frac{i}{l} \sum_{m=-M}^{N} \int_{\Delta_0}
 
  \phi^{s} ( \mathbf{x},0 )  = - \frac{i}{l} \sum_{m=-M}^{N} \int_{\Delta_0}
  \frac{\exp^{ik\mu_{m} |X| } \exp^{i\sigma_{m}Y_0}}{\mu_{m}}
+
  \frac{e^{ik\mu_{m} |X| } e^{i\sigma_{m}Y_0}}{\mu_{m}}
  \left[ k\phi(\mbox{\boldmath<math>\xi</math>},0)
+
  \left[ k\phi(\xi,0)
  - w(\mbox{\boldmath<math>\xi</math>}) \right]  
+
  - w(\xi) \right]  
  d\mbox{\boldmath<math>\xi</math>},
+
  d\xi,
  (phi_s)
 
 
</math></center>
 
</math></center>
 
where <math>\phi^{s} = \phi-\phi^{\rm in}</math> is the scattered wave which
 
where <math>\phi^{s} = \phi-\phi^{\rm in}</math> is the scattered wave which
Line 212: Line 207:
 
<center><math>
 
<center><math>
 
  \lim_{x\to-\infty}\phi^{s}  
 
  \lim_{x\to-\infty}\phi^{s}  
  (\mathbf{x},0) = A_{m}^{-}\,\exp^{ik\mu_{m}x}\exp^{i\sigma_{m}y},
+
  (\mathbf{x},0) = A_{m}^{-}\,e^{ik\mu_{m}x}e^{i\sigma_{m}y},
  (phi_m_min)
 
 
</math></center>
 
</math></center>
where the amplitudes <math>A_{m}^{-}</math> are identified from ((phi_s))
+
where the amplitudes <math>A_{m}^{-}</math> are  
as
 
 
<center><math>
 
<center><math>
 
  A_{m}^{-} = -\frac{i}{\mu_{m}l} \int_{\Delta_0}
 
  A_{m}^{-} = -\frac{i}{\mu_{m}l} \int_{\Delta_0}
  \exp^{ik\mu_{m}\xi } \exp^{-i\sigma_{m}\eta}
+
  e^{ik\mu_{m}\xi } e^{-i\sigma_{m}\eta}
  \left[  k\phi\left( \mbox{\boldmath<math>\xi</math>}\right)
+
  \left[  k\phi\left( \xi\right)
  - w (\mbox{\boldmath<math>\xi</math>}) \right]  
+
  - w (\xi) \right]  
  d\mbox{\boldmath<math>\xi</math>}.
+
  d\xi.
  (Ad_m)
 
 
</math></center>
 
</math></center>
 
Likewise as <math>x \to \infty</math> the scattered wave is given by
 
Likewise as <math>x \to \infty</math> the scattered wave is given by
 
<center><math>
 
<center><math>
 
\lim_{x\to\infty}\phi^{s} (\mathbf{x},0) =
 
\lim_{x\to\infty}\phi^{s} (\mathbf{x},0) =
  A_{m}^{+} \exp^{-ik\mu_{m}x} \exp^{i\sigma_{m}y},
+
  A_{m}^{+} e^{-ik\mu_{m}x} e^{i\sigma_{m}y},
  (phi_m_plus)
 
 
</math></center>
 
</math></center>
 
where <math>A_{m}^{+}</math> are
 
where <math>A_{m}^{+}</math> are
 
<center><math>
 
<center><math>
 
  A_{m}^{+} = -\frac{i}{\mu_{m}l}\int_{\Delta_0}
 
  A_{m}^{+} = -\frac{i}{\mu_{m}l}\int_{\Delta_0}
  \exp^{-ik\mu_{m}\xi }\exp^{-i\sigma_{m}\eta}
+
  e^{-ik\mu_{m}\xi }e^{-i\sigma_{m}\eta}
  \left[  k\phi (\mbox{\boldmath<math>\xi</math>},0)
+
  \left[  k\phi (\xi,0)
  - w(\mbox{\boldmath<math>\xi</math>}) \right]
+
  - w(\xi) \right]
  d\mbox{\boldmath<math>\xi</math>}.
+
  d\xi.
  (Ad_p)
 
 
</math></center>
 
</math></center>
 
The diffracted waves propagate at various angles with respect to the  
 
The diffracted waves propagate at various angles with respect to the  
Line 268: Line 258:
 
<center><math>
 
<center><math>
 
R = A_{0}^{-}
 
R = A_{0}^{-}
  = -\frac{i}{\mu_{0}l}\int_{\Delta_0}\exp^{ik (\xi\cos\theta
+
  = -\frac{i}{\mu_{0}l}\int_{\Delta_0}e^{ik (\xi\cos\theta
  -\eta\sin\theta)}\left[ k\phi(\mbox{\boldmath<math>\xi</math>},0)
+
  -\eta\sin\theta)}\left[ k\phi(\xi,0)
  - w(\mbox{\boldmath<math>\xi</math>})\right]
+
  - w(\xi)\right]
  d\mbox{\boldmath<math>\xi</math>}.
+
  d\xi.
  (R)
 
 
</math></center>
 
</math></center>
 
The coefficient, <math>T</math>, for the fundamental transmitted wave for  
 
The coefficient, <math>T</math>, for the fundamental transmitted wave for  
Line 278: Line 267:
 
<center><math>
 
<center><math>
 
  T = 1 + A_{0}^{+}
 
  T = 1 + A_{0}^{+}
  = 1 - \frac{i}{\mu_{0}l} \int_{\Delta_0} \exp^{-ik(\xi\cos\theta
+
  = 1 - \frac{i}{\mu_{0}l} \int_{\Delta_0} e^{-ik(\xi\cos\theta
  +\eta\sin\theta)}\left[ k\phi(\mbox{\boldmath<math>\xi</math>},0)
+
  +\eta\sin\theta)}\left[ k\phi(\xi,0)
  - w(\mbox{\boldmath<math>\xi</math>}) \right]
+
  - w(\xi) \right]
  d\mbox{\boldmath<math>\xi</math>}.
+
  d\xi.
  (T)
 
 
</math></center>
 
</math></center>
  
Line 293: Line 281:
 
<center><math>
 
<center><math>
 
\cos\theta = \left( |R|^2+|T|^2 \right) \cos\theta
 
\cos\theta = \left( |R|^2+|T|^2 \right) \cos\theta
  + \sum_{\stackrel{\scriptstyle{m=-M} }{m ~ \neq 0}}^{N}
+
  + \sum_{m=-M,\,m \neq 0}^{N}
 
  \left( |A_{m}^{-}|^2 \cos\psi_{m}^{-} +|A_{m}^{+}|^2  
 
  \left( |A_{m}^{-}|^2 \cos\psi_{m}^{-} +|A_{m}^{+}|^2  
 
  \cos\psi_{m}^{+} \right).
 
  \cos\psi_{m}^{+} \right).
  (NRGbal)
 
 
</math></center>
 
</math></center>
The energy balance equation ((NRGbal)) can be used as an accuracy
+
The energy balance equation can be used as an accuracy
 
check on the numerical results.
 
check on the numerical results.
 
=Results=
 
 
We tested the convergence of our accelerated version of the Green function
 
and we use <math>c = 0.05</math> and <math>44</math> terms in the first sum (the spatial term) and <math>46</math> terms in the second sum (the spectral term) of ((Gp_fast))).
 
These values were determined by a convergence study, the details of which
 
can be found in Wang [[wangphd04]].
 
These values will be used in all our subsequent
 
calculations. We consider four geometries for the plates which are shown
 
in Figure (floes4jfs).
 
 
 
==Scattering from a Dock==
 
 
Aside from the energy balance equation or wide spacing, it is difficult
 
to compare our results to establish their validity. However,
 
there is one case in which we can make comparisons. If we
 
consider the case when we have the dock boundary condition
 
under the plate (so that <math>w=0</math>) and the plates are square and joined
 
then the problem reduces to a two dimensional dock problem
 
which is discussed extensively in \cite[Chapter 2]{Linton_mciver01}.
 
To impose the condition of
 
a dock we simply solve equation ((phi_plate)) setting
 
<math>\vec{w}</math> to zero (we do not require equation ((comp_eqn))),
 
choosing a square plate (geometry 1) and setting the plate separation to <math>l=4</math>. 
 
 
Figure (fig_RnT4stiff) shows the reflection and
 
transmission coefficients for a square plate (geometry 1) with
 
the plate separation <math>l=4</math> and the dock boundary condition
 
(crosses) and the solution to the two-dimensional dock problem using
 
the method of [[Linton_mciver01]] (solid and dashed lines).
 
As expected the results agree.
 
 
Table (table_stiff) shows the values of the coefficient
 
<math>A_{m}^{\pm}</math> for a dock of geometry 1 with <math>{\lambda=4}</math>,
 
<math>{l=6}</math> and <math>{\theta=\pi/6}</math>. These results are given
 
to assist in numerical comparisons. 
 
\begin{table}
 
\begin{tabular}{ccc}
 
\hline
 
<math>m</math> & <math>A_{m}^{-}</math> & <math>A_{m}^{+}</math> \\ \hline
 
<math>-2</math> & <math>-0.214-0.042i</math> & <math>0.232+0.023i</math> \\
 
<math>-1</math> & <math>0.266-0.268i</math> & <math>-0.185+0.349i</math> \\
 
<math>0</math> & <math>0.631-0.210i</math> & <math>-0.702-0.141i</math> \\
 
\end{tabular}
 
\caption[]{The coefficients <math>A_{m}^{\pm}</math> for the
 
case of a dock of geometry 1 with <math>{\lambda=4}</math>,
 
<math>{l=6}</math> and <math>{\theta=\pi/6}</math>.} (table_stiff)
 
\end{table}
 
 
==Scattering from Elastic Plates==
 
 
We begin with a short table of numerical results. Table (flexible_table)
 
is equivalent to
 
Table (table_stiff) except that the plate is now elastic
 
with <math>\beta = 0.1</math> and <math>\gamma = 0</math>. As expected the reflected
 
energy is less because the waves can propagate under the elastic
 
plates.
 
\begin{table}
 
\begin{tabular}{ccc}
 
\hline
 
<math>m</math> & <math>A_{m}^{-}</math> & <math>A_{m}^{+} </math> \\ \hline
 
<math>-2</math> & <math>0.001+0.014i</math> & <math>-0.040-0.016i</math> \\
 
<math>-1</math> & <math>-0.016-0.008i</math> & <math>-0.070-0.099i</math> \\
 
<math>0</math> & <math>-0.058-0.072i</math> & <math>-0.209-0.582i</math> \\
 
\end{tabular}
 
\caption[]{The coefficients <math>A_{m}^{\pm}</math> for the
 
case of an elastic plate of geometry 1 with  <math>{\beta=0.1}</math>,
 
<math>{\gamma=0}</math>, <math>{\lambda=4}</math>,
 
<math>{l=6}</math> and <math>{\theta=\pi/6}</math>.} (flexible_table)
 
\end{table}
 
Figures (fig_square_vartheta_beta0p1) and
 
(fig_triangle_vartheta_beta0p1) show the amplitudes of
 
the diffracted waves
 
due to the array as a function of the incident angle
 
for plates of geometry one and two respectively with <math>\beta=0.1</math>,
 
<math>\gamma=0</math>, <math>k=\pi/2</math>, and <math>l=6</math>. There are 3 pairs of diffracted
 
waves (including the reflected-transmitted pair) for any angle.
 
<math>G_{\bf P}</math> diverges if <math>\sigma_n = \pm k</math> which for our values
 
of <math>l</math> and <math>k</math> means that <math>\theta = \pm 0.3398</math>. As <math>\theta</math> moves
 
across these points one of the diffracted waves disappears
 
(at <math>\pm\pi/2</math>) and an other appears (at <math>\mp\pi/2</math>). In the plots
 
we have plotted <math>A^{\pm}_{-2}</math> and <math>A^{\pm}_{1}</math> with the same
 
line style and also <math>A^{\pm}_{-1}</math> and <math>A^{\pm}_{2}</math> since they represent
 
diffracted waves which appear and disappear together. Interestingly
 
the result of doing this is to produce smooth curves for
 
<math>-\pi/2<\theta<\pi/2</math>.  Figures (fig_square_vartheta_beta0p1) and
 
(fig_triangle_vartheta_beta0p1) show that there is a very
 
strong dependence on the amount of reflected energy as a function
 
of incident angle with small angles giving the smallest reflection
 
and large angles giving the greatest reflection. 
 
 
Figures (fig_disp_p_square) to (fig_disp_p_trapezoid)
 
show the real part of the displacement for
 
five plates (<math>\Delta_{j}</math>, <math>j=-2,-1,0,1,2 </math>)
 
of the array for plates of geometry one to four
 
respectively, with <math>\beta=0.1</math>, <math>\gamma=0</math>, and <math>l=6</math>. The angle
 
of incidence is <math>\theta=\pi/6</math>.
 
We consider two values of the wavenumber, <math>k=\pi/2</math> (a) and
 
<math>k=\pi/4</math> (b). The complex response of the elastic plates is
 
apparent in these figures as is the coupling between
 
the water and the plate. It is also clear that there is a great
 
deal of difference in the individual behaviour of plates of
 
different geometries. This has practical implications for
 
experiments which might be performed on individual ice floes.
 
For example, from these figures it appears that it will be very difficult
 
to make
 
measurements of  wave spectra using an accelerometer
 
deployed on an ice floe if the floe size is comparable to the
 
wavelength.
 
 
In Figure (energy_vs_k) we consider the total reflected energy
 
<math>E_R</math> given by
 
<center><math>
 
E_R =  |R|^2 \cos\theta
 
+ \sum_{\stackrel{\scriptstyle{m=-M} }{m ~ \neq 0}}^{N}
 
  |A_{m}^{-}|^2 \cos\psi_{m}^{-}
 
</math></center>
 
as a function of <math>k</math> for <math>l= 6</math>, <math>\beta=0.1</math>, and <math>\gamma=0</math>
 
for floes of geometry 1 and 2 and for angles of <math>\theta = 0</math>, <math>\pi/6</math>,
 
and <math>\pi/3</math>. In Figure (energy_vs_k) we have
 
divided <math>E_R</math> by <math>\cos\theta</math> to normalise the curves,
 
since if all the energy is reflected
 
<math>E_R/\cos\theta=1</math>.
 
This figure shows the kind of important results which we
 
can produce from our model even with some simple calculations.
 
The results show that for <math>k<1</math> there is no scattering of energy
 
at all. While it is to be expected that long wavelength waves
 
will not be strongly scattered this figure gives us
 
a quantitative value for the <math>k</math> below which there
 
is negligible scattering. As <math>k</math> increases the reflection
 
increases and this increase is much more marked for waves incident
 
at a large angle. This is to be expected since waves which
 
are normally incident can be expected to pass through the plates
 
even for short wavelengths. However, the effect of angle appears to
 
be much stronger than might be expected on a simple geometric
 
argument (which is also what was found in Figures
 
(fig_square_vartheta_beta0p1) and (fig_triangle_vartheta_beta0p1)).
 
Interestingly, the effect of geometry is much less
 
significant than either <math>\theta</math> or <math>k</math>. This
 
implies that for practical purposes it
 
may be sufficient to take one or two representative geometries.
 
This seems surprising given the differences in the responses of
 
plates of different geometries shown in Figures  (fig_disp_p_square)
 
to (fig_disp_p_trapezoid). 
 
Figure (energy_vs_k) can be regarded as a preliminary figure,
 
showing the kinds of calculations which are required to develop
 
a model for wave scattering in the marginal ice zone from the present
 
results.
 
 
=Concluding Remarks=
 
 
Motivated by the problem of modelling wave propagation in the
 
marginal ice zone and by the general problem of scattering by
 
arrays of bodies we have presented a solution to the problem
 
of wave scattering by an infinite array of floating elastic
 
plates. The model is based on the linear theory and assumes that
 
the floe submergence is negligible. For this reason it is applicable
 
to low to moderate wave heights and to floes whose size is large compared
 
to the floe thickness.
 
The solution method is similar to that used to solve
 
for a single plate except that the periodic Green function must
 
be used. The periodic Green function is obtained by summing the free-surface
 
Green function for all <math>y</math>. However the periodic Green function is slowly
 
convergent and therefore a method, which is commonly used in Optics, is
 
devised to accelerate the convergence. The acceleration method involves
 
expressing the infinite sum of the Green function in its far-field form.
 
From the far-field representation of the periodic Green Function, we
 
calculate the diffracted wave
 
far from the array and the cut-off frequencies. We have checked our
 
numerical calculations for energy balance and against the
 
limiting case when the plates are rigid and joined where the
 
solution reduces to that of a rigid dock. We have also presented
 
solutions for a range of elastic plate geometries.
 
 
The solution method could be extended to water of finite depth using
 
a similar periodic Green function, but in this case based on the free
 
surface Green function for finite depth water. The same problems of slow
 
convergence would arise and a similar method for convergence of the
 
series would be required. We believe that a method similar to the one
 
presented here could be used to accelerate the convergence of the periodic
 
Green function for water of finite depth. Another extension would be to consider
 
a doubly periodic lattice in which a doubly periodic Green function
 
would be required. This should also be able to be computed quickly by
 
methods similar to the ones presented here.
 
 
The most important extension of this work concerns using it to construct
 
a scattering model for the marginal ice zone. Such a model would be based
 
on the solutions presented here but would be a far from trivial extension
 
of it. For a model to be realistic the effects which have been induced by
 
the assumption of periodicity would have to be removed. One method to do
 
this would be to consider random spacings of the ice floes and to average the
 
result over these. This should lead to a scattered wave field over
 
all directions rather than discrete scattering angles as well as a reflected
 
and transmitted wave. These solutions would then have to be combined, again with
 
some kind of averaging and wide spacing approximation, to compute the
 
scattering by multiple rows of ice floes.
 
 
 
 
\bibliography{/home/groups/seaice/bibdata/mike,/home/groups/seaice/bibdata/others,/home/meylan/Papers/Periodic_Cynthia/ALMOSTALL,/home/meylan/Papers/Periodic_Cynthia/Cynthia}
 
 
 
 
\pagebreak
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7]
 
{array}
 
\caption{A schematic diagram of the periodic array of floating elastic
 
plates.}
 
(fig_array)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}[ptb]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{floes4jfs}
 
\caption{Diagram of the four plate geometries for which we will
 
calculate solutions. }
 
(floes4jfs)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7]{stiff_compare}
 
\caption{The reflection coefficient <math>R</math> (solid line)
 
and the transmission coefficient <math>T</math> (dashed line) as a function
 
of <math>k</math> for a two-dimensional dock of length 4 for the incident
 
angles shown. The crosses are the same problem solved using
 
the three-dimensional array code with the dock boundary
 
condition and using plates of geometry 1 with <math>l=4</math>.}
 
(fig_RnT4stiff)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{diffracted_square}
 
\caption{The diffracted waves <math>A_m^\pm</math> for a periodic array
 
of geometry one plates with <math>k=\pi/2</math>, <math>l=6</math>, <math>\beta=0.1,</math> <math>\gamma=0</math>,
 
and <math>l=6</math>. The solid line
 
is <math>A_0^\pm</math>, the chained line is <math>A_{-2}^\pm</math> and <math>A_1^\pm</math> and the
 
dashed line is <math>A_{-1}^\pm</math> and <math>A_2^\pm</math>.}
 
(fig_square_vartheta_beta0p1)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{diffracted_triangle}
 
\caption{As for figure (fig_square_vartheta_beta0p1) except that
 
the pate has geometry 2. }
 
(fig_triangle_vartheta_beta0p1)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{disp_P_square}
 
\caption{The real part of the displacement <math>w</math> for five plates of geometry one
 
which are part of a periodic
 
array, <math>l= 6</math>, <math>\theta=\pi/6</math>, <math>\beta=0.1</math>, <math>\gamma=0</math> and
 
(a) <math>k=\pi/2</math> and (b) <math>k=\pi/4=8</math>.}
 
(fig_disp_p_square)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{disp_P_triangle}
 
\caption{As for figure (fig_disp_p_square) except that the
 
plate is of geometry two.}
 
(fig_disp_p_triangle)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{disp_P_parallelogram}
 
\caption{As for figure (fig_disp_p_square) except that the
 
plate is of geometry three.}
 
(fig_disp_p_parallelogram)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{disp_P_trapezoid}
 
\caption{As for figure (fig_disp_p_square) except that the
 
plate is of geometry four.}
 
(fig_disp_p_trapezoid)
 
\end{center}
 
\end{figure}
 
 
\begin{figure}
 
[p]
 
\begin{center}
 
\includegraphics[scale = 0.7
 
 
 
]
 
{energy_vs_k}
 
\caption{The total reflected energy <math>E_R</math> divided by
 
<math>\cos\theta</math> as a function of <math>k</math> for floes of geometry 1 (a)
 
and 2 (b). <math>l= 6</math>, <math>\beta=0.1</math>, <math>\gamma=0</math>, and
 
<math>\theta = 0</math> (solid line), <math>\pi/6</math> (dashed line), and <math>\pi/3</math>
 
(chained line).}
 
(energy_vs_k)
 
\end{center}
 
\end{figure}
 
 
 
  
  
 
[[Category:Infinite Array]]
 
[[Category:Infinite Array]]

Latest revision as of 09:11, 9 January 2009

Introduction

We present here the solution to the Infinite Array based on an infinite image system of Free-Surface Green Functions

Problem Formulation

We begin by formulating the problem. Cartesian coordinates [math]\displaystyle{ (x,y,z) }[/math] are chosen with [math]\displaystyle{ z }[/math] vertically upwards such that [math]\displaystyle{ z=0 }[/math] coincides with the mean free surface of the water. An infinite array of identical bodies are periodically spaced along the [math]\displaystyle{ y }[/math]-axis with uniform separation [math]\displaystyle{ l }[/math]. The problem is to determine the motion of the water and the bodies when plane waves are obliquely-incident from [math]\displaystyle{ x=-\infty }[/math] upon the periodic array of bodies.

The bodies occupy [math]\displaystyle{ \Delta_m }[/math], [math]\displaystyle{ -\infty \lt m \lt \infty }[/math]. Periodicity implies that if [math]\displaystyle{ (x,y) \in \Delta_0 }[/math], then [math]\displaystyle{ (x,y+ml) \in \Delta_m }[/math], [math]\displaystyle{ -\infty \lt m \lt \infty }[/math].

We assume that we have the Standard Linear Wave Scattering Problem. The incident wave potential given by

[math]\displaystyle{ \phi^{{\rm in}} = \frac{A}{k} e^{ik (x\cos\theta+y\sin\theta)}\,e^{kz}, }[/math]

where [math]\displaystyle{ A }[/math] is the dimensionless amplitude and [math]\displaystyle{ \theta }[/math] is the direction of propagation of the wave (with [math]\displaystyle{ \theta = 0 }[/math] corresponding to normal incidence.

Transformation to an Integral Equation

We now Floquet's theorem (Scott 1998) (also called the assumption of periodicity in the water wave context) which states the displacement from adjacent plates differ only by a phase factor. If the potential under the central plate [math]\displaystyle{ \Delta_{0} }[/math] is given by [math]\displaystyle{ \phi( \mathbf{x}_{0},0) }[/math], [math]\displaystyle{ \mathbf{x}_{0}\in\Delta_{0} }[/math], then by Floquet's theorem the potential satisfies

[math]\displaystyle{ \phi(\mathbf{x}_{m},0) = \phi(\mathbf{x}_{0},0) e^{im\sigma l}, }[/math]

and the displacement of the plate [math]\displaystyle{ \Delta_{m} }[/math] satisfies

[math]\displaystyle{ w(\mathbf{x}_{m}) = w(\mathbf{x}_{0}) e^{im\sigma l}, }[/math]

where [math]\displaystyle{ \mathbf{x}_{m} \in \Delta_{m} }[/math], [math]\displaystyle{ -\infty \lt m \lt \infty }[/math] and the phase difference is [math]\displaystyle{ \sigma = k\sin\theta }[/math] (see, for example, Linton 1998).

A standard approach to the solution of the equations of motion for the water is the Green Function Solution Method in which we transform the equations into a boundary integral equation using the Free-Surface Green Function. In doing so we obtain

[math]\displaystyle{ \phi(\mathbf{x}) = \phi^{\rm in} (\mathbf{x},0) +\sum_{m=-\infty}^{\infty} \int_{\Delta_{m}} \left(G_{n_\xi}(\mathbf{x},\xi) \phi(\xi) - G(\mathbf{x},\xi) \phi_{n_\xi}(\xi) \right) d\xi }[/math]

[math]\displaystyle{ G(\mathbf{x},\xi) }[/math] is the Free-Surface Green Function This can be written alternatively as

[math]\displaystyle{ \phi(\mathbf{x}) = \phi^{\rm in}(\mathbf{x}) +\int_{\Delta_{0}} \sum_{m=-\infty}^{\infty} \left(G^{\mathbf{P}}_{n_\xi}(\mathbf{x},\xi+(0,ml,0))e^{im\sigma l} \phi(\xi) - G^{\mathbf{P}} (\mathbf{x},\xi)e^{im\sigma l} \phi_{n_\xi}(\xi) \right) d\xi }[/math]

where the kernel [math]\displaystyle{ G_{\mathbf{P}} }[/math] (referred to as the periodic Green function) is given by

[math]\displaystyle{ G^{\mathbf{P}} (\mathbf{x};\xi) = \sum_{m=-\infty}^{\infty} G(\mathbf{x},\xi+(0,ml,0))e^{im\sigma l}. }[/math]

Accelerating the Convergence of the Periodic Green Function

The spatial representation of the periodic Green function given by equation is slowly convergent and in the far field the terms decay in magnitude like [math]\displaystyle{ O(n^{-1/2}) }[/math]. In this section we show how to accelerate the convergence. We begin with the asymptotic approximation of the Three-dimensional Free-Surface Green Function far from the source point,

[math]\displaystyle{ G(\mathbf{x},\xi) \sim -\frac{ik}{2} \,H_{0}( k |\mathbf{x}-\xi|), |\mathbf{x}-\xi| \to \infty }[/math]

Wehausen and Laitone 1960 where [math]\displaystyle{ H_0 \equiv H_{0}^{(1)} }[/math] is the Hankel function of the first kind of order zero Abramowitz and Stegun 1964. In Linton Linton 1998 various methods were described in which the convergence of the periodic Green functions was improved. One such method, which suits the particular problem being considered here, involves writing the periodic Green function as

[math]\displaystyle{ G_{\mathbf{P}} (\mathbf{x};\xi) = \sum_{m=-\infty}^{\infty} \left[ G\left(\mathbf{x};\xi)+(0,ml)\right) + \frac{ik}{2} H_{0} \Big(k\sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big) e^{im\sigma l} \right] -\sum_{m=-\infty}^{\infty} \frac{ik}{2}H_{0} \Big(k\sqrt{ (X+cl)^2 + Y_{m}^2 }\Big) e^{im\sigma l} }[/math]

where [math]\displaystyle{ c }[/math] is a numerical smoothing parameter, introduced to avoid the singularity at [math]\displaystyle{ \mathbf{x} = \xi }[/math] in the Hankel function and

[math]\displaystyle{ X = x-\xi,\quad \mathrm{and} \quad Y_{m} = (y-\eta)-ml. }[/math]

Furthermore we use the fact that second slowly convergent sum can be transformed to

[math]\displaystyle{ -\sum_{m=-\infty}^{\infty} \frac{ik}{2}H_{0} \Big(k\sqrt{ (X+cl)^2 + Y_{m}^2 }\Big) e^{im\sigma l} -\frac{i}{l} \sum_{m=-\infty}^{\infty} \frac{e^{ik \mu_{m} |X+c| }\,e^{i \sigma_{m} Y_0}}{\mu_{m}} }[/math]

Linton 1998 where [math]\displaystyle{ \sigma_{m} = \sigma + 2 m \pi/l }[/math] and

[math]\displaystyle{ \mu_m = \left[ 1-\left(\frac{\sigma_{m}}{k} \right)^{2} \right]^{\frac{1}{2}}, }[/math]

where the positive real or positive imaginary part of the square root is taken. Combining these equations we obtain the accelerated version of the periodic Green function

[math]\displaystyle{ G_{\mathbf{P}} (\mathbf{x};\xi) = \sum_{m=-\infty}^{\infty} \left[ G\left(\mathbf{x};\xi+(0,ml)\right) + \frac{ik}{2} H_{0} \Big(k \sqrt{\left( X+cl\right)^2 + Y_{m}^2}\Big) \right] e^{im\sigma l} -\frac{i}{l}\sum_{m=-\infty}^{\infty} \frac{e^{ik \mu_{m}|X+cl| }e^{i\sigma_{m}Y_0}}{\mu_{m}}. }[/math]

The convergence of the two sums depends on the value of [math]\displaystyle{ c }[/math]. For small [math]\displaystyle{ c }[/math] the first sum converges rapidly while the second converges slowly. For large [math]\displaystyle{ c }[/math] the second sum converges rapidly while the first converges slowly. The smoothing parameter [math]\displaystyle{ c }[/math] must be carefully chosen to balance these two effects. Of course, the convergence also depends strongly on how close together the points [math]\displaystyle{ \mathbf{x} }[/math] and [math]\displaystyle{ \xi }[/math] are.

Note that some special combinations of wavelength [math]\displaystyle{ \lambda }[/math] and angle of incidence [math]\displaystyle{ \theta }[/math] cause the periodic Green function to diverge ( Scott 1998). This singularity is closely related to the diffracted waves and will be explained shortly.

The scattered waves (modes)

We begin with the accelerated periodic Green function, equation setting [math]\displaystyle{ c=0 }[/math] and considering the case when [math]\displaystyle{ X }[/math] is large (positive or negative). We also note that for [math]\displaystyle{ m }[/math] sufficiently small or large [math]\displaystyle{ i\mu_m }[/math] will be negative and the corresponding terms will decay. Therefore

[math]\displaystyle{ G_{\mathbf{P}} (\mathbf{x};\xi) \sim - \frac{i}{l} \sum_{m=-M}^{N} \frac{e^{ik\mu_{m}|X|}\, e^{i\sigma_{m}Y_0}} {\mu_{m}}, X \to \pm \infty }[/math]

where the integers [math]\displaystyle{ M }[/math] and [math]\displaystyle{ N }[/math] satisfy the following inequalities

[math]\displaystyle{ (M_N) \left. \begin{matrix} \sigma_{-M-1}\lt -k\lt \sigma_{-M},\\ \sigma_{N}\lt k\lt \sigma_{N+1}. \end{matrix} \right\} }[/math]

These equations can be written as

[math]\displaystyle{ \frac{l}{2\pi}\left(\sigma+k-2\pi \right) \lt M \lt \frac{l}{2\pi}\left( \sigma+k \right), }[/math]

and

[math]\displaystyle{ \frac{l}{2\pi}\left( k - \sigma \right) \gt N \gt \frac{l}{2\pi} \left( k-\sigma - 2\pi \right) }[/math]

Linton 1998. It is obvious that [math]\displaystyle{ G_{\mathbf{P}} }[/math] will diverge if [math]\displaystyle{ \sigma_m = \pm k }[/math]; these values correspond to cut-off frequencies which are an expected feature of periodic structures.

The diffracted waves

The diffracted waves are the plane waves which are observed as [math]\displaystyle{ x \to \pm \infty }[/math]. Their amplitude and form are obtained by substituting the limit of the periodic Green function as [math]\displaystyle{ x\to\pm\infty }[/math] into the boundary integral equation for the potential. This gives us

[math]\displaystyle{ \lim_{x\to\pm\infty} \phi^{s} ( \mathbf{x},0 ) = - \frac{i}{l} \sum_{m=-M}^{N} \int_{\Delta_0} \frac{e^{ik\mu_{m} |X| } e^{i\sigma_{m}Y_0}}{\mu_{m}} \left[ k\phi(\xi,0) - w(\xi) \right] d\xi, }[/math]

where [math]\displaystyle{ \phi^{s} = \phi-\phi^{\rm in} }[/math] is the scattered wave which is composed of a finite number of plane waves. For [math]\displaystyle{ x \to -\infty }[/math] the scattered wave is given by

[math]\displaystyle{ \lim_{x\to-\infty}\phi^{s} (\mathbf{x},0) = A_{m}^{-}\,e^{ik\mu_{m}x}e^{i\sigma_{m}y}, }[/math]

where the amplitudes [math]\displaystyle{ A_{m}^{-} }[/math] are

[math]\displaystyle{ A_{m}^{-} = -\frac{i}{\mu_{m}l} \int_{\Delta_0} e^{ik\mu_{m}\xi } e^{-i\sigma_{m}\eta} \left[ k\phi\left( \xi\right) - w (\xi) \right] d\xi. }[/math]

Likewise as [math]\displaystyle{ x \to \infty }[/math] the scattered wave is given by

[math]\displaystyle{ \lim_{x\to\infty}\phi^{s} (\mathbf{x},0) = A_{m}^{+} e^{-ik\mu_{m}x} e^{i\sigma_{m}y}, }[/math]

where [math]\displaystyle{ A_{m}^{+} }[/math] are

[math]\displaystyle{ A_{m}^{+} = -\frac{i}{\mu_{m}l}\int_{\Delta_0} e^{-ik\mu_{m}\xi }e^{-i\sigma_{m}\eta} \left[ k\phi (\xi,0) - w(\xi) \right] d\xi. }[/math]

The diffracted waves propagate at various angles with respect to the normal direction of the array. The angles of diffraction, [math]\displaystyle{ \psi_{m}^{\pm} }[/math], are given by

[math]\displaystyle{ \psi_{m}^{\pm} = \tan^{-1}\left( \frac{\sigma_{m}}{\pm k\mu_{m}}\right). (psi_m) }[/math]

Notice that for [math]\displaystyle{ m=0 }[/math] we have

[math]\displaystyle{ \psi_{0}^{\pm}=\pm\theta, (psi_0) }[/math]

where [math]\displaystyle{ \theta }[/math] is the incident angle. This is exactly as expected since we should always have a transmitted wave which travels in the same direction as the incident wave and a reflected wave which travels in the negative incident angle direction.

The fundamental reflected and transmitted waves

We need to be precise when we determine the wave of order zero at [math]\displaystyle{ x\to\infty }[/math] because we have to include the incident wave. There is always at least one set of propagating waves corresponding to [math]\displaystyle{ m=0 }[/math] which correspond to simple reflection and transmission. The coefficient, [math]\displaystyle{ R }[/math], for the fundamental reflected wave for the [math]\displaystyle{ m=0 }[/math] mode is given by

[math]\displaystyle{ R = A_{0}^{-} = -\frac{i}{\mu_{0}l}\int_{\Delta_0}e^{ik (\xi\cos\theta -\eta\sin\theta)}\left[ k\phi(\xi,0) - w(\xi)\right] d\xi. }[/math]

The coefficient, [math]\displaystyle{ T }[/math], for the fundamental transmitted wave for the [math]\displaystyle{ m=0 }[/math] mode is given by

[math]\displaystyle{ T = 1 + A_{0}^{+} = 1 - \frac{i}{\mu_{0}l} \int_{\Delta_0} e^{-ik(\xi\cos\theta +\eta\sin\theta)}\left[ k\phi(\xi,0) - w(\xi) \right] d\xi. }[/math]

Conservation of energy

The diffracted wave, taking into account the correction for [math]\displaystyle{ T }[/math], must satisfy the energy flux equation. This simply says that the energy of the incoming wave must be equal to the energy of the outgoing waves. This gives us

[math]\displaystyle{ \cos\theta = \left( |R|^2+|T|^2 \right) \cos\theta + \sum_{m=-M,\,m \neq 0}^{N} \left( |A_{m}^{-}|^2 \cos\psi_{m}^{-} +|A_{m}^{+}|^2 \cos\psi_{m}^{+} \right). }[/math]

The energy balance equation can be used as an accuracy check on the numerical results.