# Stability Analysis

The partial differential equation, \begin{align} \dfrac{\partial f}{\partial x} + a(x)\dfrac{\partial f}{\partial y} = 0 \qquad & f(0,y) = f(L_1,y) = c_0e^{-y} \\ & f(x,0) = c_0 \;,\; f(x,L_2) = c_0e^{-L_2} \end{align} for $$f=f(x,y)$$ has the nontrivial solution, $$$$f(x,y) = Ae^{(-y-h(x))} \;,\quad a(x)=-\dfrac{dh(x)}{dx} \;.$$$$ If I discretize the equation by box integration over $$[y_{i},y_{i+1}]$$, I get: $$$$\dfrac{\Delta y}{2} \left[\dfrac{\partial f(x,y_{i})}{\partial x} + \dfrac{\partial f(x,y_{i+1})}{\partial x}\right] + a(x)\big[f(x,y_{i+1}) - f(x,y_i)\big] = 0$$$$ which converges to the nontrivial solution as $$\Delta y\rightarrow 0$$. However, I have noticed such discretizations are numerically unstable for $$a(x)\neq 0$$ and result in oscillations when solved on a computer. How can I investigate the above equation in terms of numerical stability? Does this equation belong to a specific class of PDEs, and are there techniques to improve its stability?

• Are you making sure that (i) you are doing backward integration if $a>0$, and (ii) that you are satisfying the CFL (stability) condition? – Wolfgang Bangerth Oct 5 '18 at 22:34

This is a first-order PDE but you are trying to impose Dirichlet BC on the whole boundary. This is not well-posed in general but maybe your data is such that you expect it to have a unique solution. You can write this PDE as $$v \cdot \nabla f = \nabla\cdot (v f) = 0$$ where $$v(x,y) = (1, a(x))$$ This is a hyperbolic PDE. If $$n$$ is outward normal to your domain boundary, then you can specify Dirichlet BC only on those portions of the boundary where $$v \cdot n < 0$$. So you can specify BC on the left boundary $$(x=0)$$ but not the right boundary. On the bottom and top boundaries, you have to check the sign of $$a(x)$$ to decide this.
Draw the characteristic curves $$dy/dx = a(x)$$ to see how BC will flow into your domain. This can identify a time-like direction along which you can march to find the solution.
For example, if $$a(x)=x$$ and your domain is $$[0,L_1]\times[0,L_2]$$ so that $$a(x) \ge 0$$ in your domain, then you can specify BC on the left and bottom sides. You can march along $$x$$. As Wolfgang wrote, for $$\partial f/\partial y$$ you can use backward-difference type scheme. An FD scheme may look like this $$\frac{f_{i+1,j}-f_{i,j}}{\Delta x} + a_i \frac{f_{i,j} - f_{i,j-1}}{\Delta y} = 0, \qquad a_i \ge 0$$ Note that $$f_{0,j}$$ and $$f_{i,0}$$ are the given boundary values and you determine the remaining values using the scheme above.