Limits of the wave equation with piecewise constant propagation speed

approximation-theoryordinary differential equationsreference-request

Consider a wave equation
$$\frac{\partial^2 u}{\partial t^2} = c(x)^2 \frac{\partial^2 u}{\partial x^2} \tag{1}$$
In frequency domain this becomes an ODE:
$$-\omega^2 u = c(x)^2 \frac{\partial^2 u}{\partial x^2} \tag{2}$$
We can solve (2) analytically when $c(x)$ is constant, then the solution is $\exp\left(\pm\frac{i \omega}{c} x\right)$. We can also solve it analytically if $c(x)$ is piecewise constant: write in every region-of-constant-$c$ a linear combination of leftward and rightward propagating waves, impose the appropriate continuity conditions on the interfaces between regions with different values of $c$, solve for the coefficients of the linear combinations.

Suppose now that $c(x)$ is some arbitrary function. We can still write it as a limit of piecewise constant functions:
$$c(x)=\lim_{\Delta\rightarrow 0^+} \sum_{i=-\infty}^{\infty} \left\{\begin{matrix}c(i\Delta) & x\in [i\Delta,(i+1)\Delta[ \\ 0 & \text{otherwise}\end{matrix}\right.$$

Let us call these piecewise constant approximations to $c$, $c_{\Delta}$
$$c_{\Delta}(x)=\sum_{i=-\infty}^{\infty} \left\{\begin{matrix}c(i\Delta) & x\in [i\Delta,(i+1)\Delta[ \\ 0 & \text{otherwise}\end{matrix}\right.$$
We can analytically solve, for any strictly positive $\Delta$,
$$-\omega^2 u_{\Delta}= c_{\Delta}(x)^2 \frac{\partial^2 u_{\Delta}}{\partial x^2}$$

Question

Is $\lim_{\Delta\rightarrow 0^+} u_{\Delta} = u$?

In other words, can we approximate general solutions of (2) by substituting in a piecewise constant approximation of $c(x)$, and then solving analytically?

I expect the answer to be "no".

  • because $\lim_{\Delta\rightarrow 0^+} c_{\Delta} = c$, but $\lim_{\Delta\rightarrow 0^+} \frac{\partial}{\partial x}c_{\Delta} \neq \frac{\partial}{\partial x} c$, in fact the limit doesn't even exist
  • In literature on numerical solutions of PDEs, I have never seen this suggested as a viable way of approximating solutions of wave equations.

Best Answer

New:

For the initial conditions of $u_\Delta$ choose $u_\Delta(0) = u(0)$ and $u'_\Delta(0) = u'(0)$. Define $v = u - u_\Delta$ and $\delta c = c - c_\Delta$. So $v(0) = v'(0) = 0$. Suppose $c$ is sufficiently well behaved so $u$ is also. $u_\Delta$ should be well behaved excluding a discontinuous second derivative and if ever $c_\Delta = 0$ then $u_\Delta$ and $u''_\Delta$ jumps to zero.

Define $v = u - u_\Delta$ and $\delta c = c - c_\Delta$. Subtracting $$ -\omega^2 v = c_\Delta^2 v'' - (2c\delta c - \delta c^2) u'' $$ Suppose for some $x_0 \geq 0$ we want to know $v(x_0)$. Note if $\Delta \rightarrow 0$ then the sup norm $||\delta c||_{\infty, [0,x_0]} \rightarrow 0$. Since $u''$ and $c$ are well behaved we're arguing $||(2c\delta c - \delta c^2) u''||_{\infty,[0,x_0]} \rightarrow 0$. If this term on the RHS were zero, by the uniqueness theorem $v(x)=0$ on $[0,x_0]$. The fact it approaches zero simply means it perturbs $v(x)$ from $0$ an arbitrarily small amount.

So $u_\Delta \rightarrow u$ as $\Delta \rightarrow 0$ if $c$ is sufficiently nice. I haven't tried to address what a necessary and sufficient condition may be.


Original:

It's hard for me to picture the limit does not hold for a large class of well behaved functions. On the other hand surely we can contrive functions which will break the limit. By work such as Lebesgue's I feel finding a necessary and sufficient condition on functions for the limit to hold may not be easy.

Thoughts for proof: If we can prove the limit holds when $c(x)$ is a polynomial perhaps we can take a limit on the degree to argue the limit holds for many Taylor series.

For convenience scale $c(x)$ by $\omega$. Suppose $c(x)$ is a polynomial and denote $c(x)^2 = c_0 + c_1 x + \cdots + c_n x^{2n}$. If we were able to solve for $u$ and $u_\Delta$ we'd be able to take the limit and see if they match. Applying the Laplace transform which is generally a simplification

\begin{align} \mathcal{L}\{-u\} = -U(s) =& \sum_{0\leq k \leq 2n} \mathcal{L}\{c_k x^k u''\} \\ =& \sum_{0\leq k \leq 2n} (-1)^k c_k \frac{d^k}{ds^k}\bigg[s^2U-su(0)-u'(0)\bigg] \\ =& -c_0 (su(0)+ u'(0)) + c_1 u(0) + \sum_{0\leq k\leq 2n} (-1)^k (c_k s^2-2c_{k+1}s+2c_{k+2}) U^{(k)} \end{align} where $c_{2n+1},c_{2n+2} = 0$.

\begin{align} \mathcal{L}\{-u_\Delta\}=-U_\Delta(s)= & \sum_{0\leq k} \mathcal{L}\Big\{\big(H(x-k\Delta)-H(x-(k+1)\Delta)\big)c(k\Delta)^2 u_\Delta''(x)\Big\} \\ =& \sum_{0\leq k} c(k\Delta)^2 e^{-sk\Delta}\Big[\mathcal{L}\{u''_\Delta(x+k\Delta)\} - e^{-s\Delta}\mathcal{L}\{u''_\Delta(x+(k+1)\Delta)\}\Big] \end{align}

using $\mathcal{L}\{H(x-k\Delta)u''_\Delta(x)\}=\int_0^\infty e^{-sx}H(x-k\Delta)u''_\Delta(x) dx = e^{-sk\Delta}\mathcal{L}\{u''_\Delta(x+k\Delta)\}$ where $H$ is the heaviside function. Unfortunately $u''_\Delta(x)$ is not of the right form. Yet going back and directly solving the original equations for any $\Delta$ seems recursive without a generating formula. Perhaps there's less explicit theorems which can be used.

Related Question