Large-Amplitude Elastic Free-Surface Waves: Geometric Nonlinearity and Peakons

Abstract

An instantaneous sub-surface disturbance in a two-dimensional elastic half-space is considered. The disturbance propagates through the elastic material until it reaches the free surface, after which it propagates out along the surface. In conventional theory, the free-surface conditions on the unknown surface are projected onto the flat plane $y = 0$, so that a linear model may be used. Here, however, we present a formulation that takes explicit account of the fact that the location of the free surface is unknown a priori, and we show how to solve this more difficult problem numerically. This reveals that, while conventional linearized theory gives an accurate account of the decaying waves that travel outwards along the surface, it can under-estimate the strength of the elastic rebound above the location of the disturbance. In some circumstances, the non-linear solution fails in finite time, due to the formation of a “peakon” at the free surface. We suggest that brittle failure of the elastic material might in practice be initiated at those times and locations.

Appendix

In this Appendix, we derive the stability condition for the explicit finite-difference scheme (17) presented in Sect. 3.1. To do this, we employ von Neumann’s Fourier stability method [31, Sect. 2-2].

Since the difference scheme (17) is linear, any error introduced into the scheme satisfies the same equations, and the purpose of a stability analysis is to determine whether that error grows or decays as the index $k$ (corresponding to time) increases. We suppose that the errors have a Fourier form

$$\begin{aligned} \textrm{Error in } X_{i,j}^{k} = & A_{k} \exp \left ( \sqrt{-1} \left [ i p \Delta x + j q \Delta y \right ] \right ) \\ \textrm{Error in } Y_{i,j}^{k} = & B_{k} \exp \left ( \sqrt{-1} \left [ i p \Delta x + j q \Delta y \right ] \right ) , \end{aligned}$$

(58)

in which the two constants $p$ and $q$ depend on the details of the numerical grids in the variables $x$ and $y$. These forms (58) are substituted into the numerical scheme, and after some analysis it is found that the two growth factors $A_{k}$ and $B_{k}$ satisfy the coupled linear difference scheme

$$\begin{aligned} A_{k+1} = & 2 \lambda _{A} A_{k} - \gamma _{AB} B_{k} - A_{k-1} , \\ B_{k+1} = & 2 \lambda _{B} B_{k} - \gamma _{AB} A_{k} - B_{k-1} . \end{aligned}$$

(59)

In this scheme, we have defined three additional constants for convenience. They are

$$\begin{aligned} \lambda _{A} = & 1 - 2 \alpha ^{2} r^{2} P^{2} - 2 s^{2} Q^{2} , \\ \lambda _{B} = & 1 - 2 \alpha ^{2} s^{2} Q^{2} - 2 r^{2} P^{2} , \\ \gamma _{AB} = & 4 \left ( \alpha ^{2} - 1 \right ) r s P Q \sqrt{1 - P^{2}} \sqrt{1 - Q^{2}} . \end{aligned}$$

(60)

Here, the two grid-based constants $r$ and $s$ are

$$ r = \frac{\Delta t}{\Delta x} \quad \textrm{;} \quad s = \frac{\Delta t}{\Delta y} $$

(61)

and the two remaining constants $P$ and $Q$ involve the Fourier modes $p$ and $q$ in the assumed form (58) for the errors in the scheme; they are

$$ P = \sin \left ( \frac{1}{2} p \Delta x \right ) \quad \textrm{;} \quad Q = \sin \left ( \frac{1}{2} q \Delta y \right ) . $$

(62)

The coupled difference scheme (59) can be written as the single fourth-order linear difference equation

$$\begin{aligned} & A_{k+2} - 2 \left ( \lambda _{A} + \lambda _{B} \right ) A_{k+1} + \left ( 2 + 4 \lambda _{A} \lambda _{B} -\gamma _{AB}^{2} \right ) A_{k} \\ & \quad - 2 \left ( \lambda _{A} + \lambda _{B} \right ) A_{k-1} + A_{k-2} = 0 , \end{aligned}$$

and since this has constant coefficients, its solution may be sought by assuming $A_{k} \sim \xi ^{k}$. It follows at once that $\xi$ must solve the quartic polynomial equation

$$\begin{aligned} & \xi ^{4} - 2 \left ( \lambda _{A} + \lambda _{B} \right ) \xi ^{3} + \left ( 2 + 4 \lambda _{A} \lambda _{B} -\gamma _{AB}^{2} \right ) \xi ^{2} \\ & \quad - 2 \left ( \lambda _{A} + \lambda _{B} \right ) \xi + 1 = 0 . \end{aligned}$$

(63)

Stability of the original system (17) is now determined by the four roots $\xi$ of this equation, since if $| \xi | < 1$ then $A_{k}$ and $B_{k}$ will remain bounded as $k \rightarrow \infty$ and in this case, errors will not grow and the system will therefore be stable.

Remarkably, the four roots $\xi$ of (63) can be written in fairly simple closed form.¹ They are

$$ \xi = \rho \pm \sqrt{ \rho ^{2} - 1} $$

(64)

in which

$$ \rho = \frac{1}{2} \left [ \lambda _{A} + \lambda _{B} \pm \sqrt{ \gamma _{AB}^{2} + \left ( \lambda _{A} - \lambda _{B} \right )^{2} } \right ] . $$

(65)

This leads at once to the following result:

Lemma 1

The system is stable if $| \rho | < 1$ and it is unstable if $| \rho | > 1$.

Proof

If $| \rho | < 1$, then from (64) all roots are complex and

$$ \xi = \rho \pm \sqrt{-1} \sqrt{ 1 -\rho ^{2} } . $$

(66)

In this case, the absolute value is $| \xi | = 1$. Consequently the error neither grows nor decays, so that the system is non-dissipative and stable. Conversely, if $| \rho | > 1$, then the two options for $\xi$ in (64) are both real and $| \xi | > 1$ for at least one of them. The error amplitude $A_{k} \sim \xi ^{k}$ thus grows without bound and the system is therefore unstable. □

We now state our main result.

Theorem 1

If

$$ r^{2} + s^{2} < \frac{1}{1 + \alpha ^{2}} $$

(67)

the explicit finite-difference scheme (17) is stable.

Proof

If condition (67) holds, then after some algebra, (60) demonstrates that

$$\begin{aligned} \lambda _{A} + \lambda _{B} > & 0 , \\ \gamma _{AB}^{2} < & 4 \left [ 1 - \left ( \lambda _{A} + \lambda _{B} \right ) + \lambda _{A} \lambda _{B} \right ] \end{aligned}$$

both apply. These two inequalities then show that

$$\begin{aligned} \gamma _{AB}^{2} < & 4 \left [ 1 - \lambda _{A} - \lambda _{B} + \lambda _{A} \lambda _{B} \right ] \\ < & 4 \left [ 1 + \lambda _{A} + \lambda _{B} + \lambda _{A} \lambda _{B} \right ] . \end{aligned}$$

These can be manipulated further to give

$$\begin{aligned} - 2 - \lambda _{A} - \lambda _{B} < & \pm \sqrt{ \gamma _{AB}^{2} + \left ( \lambda _{A} - \lambda _{B} \right )^{2} } \\ < & 2 - \lambda _{A} - \lambda _{B} . \end{aligned}$$

From this result, it follows that

$$ -2 < 2\rho < 2 , $$

(68)

using the definition (65). Thus $| \rho | < 1$ and by Lemma 1 the system is stable. □

In addition to finding the condition (67) that is necessary for stability, it is also possible to determine when the scheme is guaranteed to be unstable. This is summarized in the following result:

Theorem 2

If

$$ r^{2} + s^{2} > 1 $$

(69)

(and $\alpha > 1$), the explicit finite-difference scheme (17) is unstable.

Proof

In the application of the numerical scheme (17), all the Fourier modes (58) will be present, including some for which $P = 1$ and $Q = 1$ in (62). For those modes, (65) shows that

$$ \rho = \left \{ \textstyle\begin{array}{c} 1 - 2 \alpha ^{2} r^{2} - 2 s^{2} \\ 1 - 2 r^{2} - 2 \alpha ^{2} s^{2} \end{array}\displaystyle \right . . $$

(70)

Since $\alpha ^{2} > 1$, it follows that, for these modes,

$$ \rho < 1 - 2 r^{2} - 2 s^{2} . $$

(71)

Now if the condition (69) holds, then these modes will have the property that $\rho < - 1$. In that case, Lemma 1 shows that the difference scheme (17) is unstable. □

Theorems 1 and 2 can be regarded as a generalization of the Courant-Friedrichs-Lewy condition in the numerical solution of hyperbolic partial differential equations [31].