I'm trying to solve a nonlinear partial differential equation \begin{equation} L(u_{xxtt},u_{xx}u_{tt},u_{xt}^2,u_{xt},u_{tt})=0 \end{equation} using finite difference methods. In order to remove the second derivatives with respect to time, I substitute \begin{equation} v := u_t \end{equation} such that it becomes \begin{equation} L(v_{xxt},u_{xx}v_{t},v_{x}^2,v_{x},v_{t}) = 0 \end{equation} So in my algorithm, I first implicitly calculate $v$ at the next time step using central difference approximations (in order to ensure that local truncation error is second order both in time and space), i.e. \begin{align} (v_i{}^m)_{xx} & = (v^m_{i+1} - 2 v_i{}^m + v^m_{i-1})/\Delta x^2 \\ (v_i{}^m)_x & = (v^m_{i+1}-v^m_{i-1}) / 2\Delta x \\ (v_i{}^m)_t & = (v_i^{m+1} - v_i{}^{m-1})/2 \Delta t \end{align} Once I've calculated $v_i^{m+1}$, I calculate $u_i^{m+1}$ as follows \begin{equation} v_i{}^m = (u_i{}^m)_t = \frac{u^{m+1}_i - u_i^{m-1}}{2 \Delta t} \implies u^{m+1} = 2 \Delta t v_i{}^m + u_i^{m-1} \tag{1} \end{equation} I would like to determine the stability of this program, and any tips on how to do that are greatly appreciated.

I first wanted to determine the stability using von Neumann's method, but I don't think this is possible (but please correct me if I'm wrong). Let me explain:

Based on the answer given to my previous question, I've decided to "freeze" $u_{xx}$ and $v_x$ in order to linearize the PDE \begin{equation} L(v_{xxt},v_{x},v_{t}) = 0 \end{equation} However, if I substitute equation $(1)$ into my finite difference operator, then I will get terms involving $u^{m+2}, u^{m+1}, u^m, u^{m-1}, u^{m-2}$. Thus, if I now substitute in the following (as per usual for the von Neumann stability analysis) \begin{equation} u_i{}^m = A^m \exp(\hat{i} k i \Delta x) \; , \; \hat{i}:= \sqrt{-1} \end{equation} then I will end up with a fifth order polynomial in $A$ \begin{equation} aA^4 + bA^3 + cA^2 + dA + e = 0 \end{equation} This means it will not be possible to solve for $A$ and thus I cannot determine the amplification factor. Hence, unless I've missed something, it's not possible to use the von Neumann analysis to determine stability of the scheme.

  • 1
    $\begingroup$ Of course you can solve for $A$. You just can't write down an explicit formula for it in terms of radicals. $\endgroup$ – David Ketcheson Mar 8 '16 at 8:00

You can always try to find the roots numerically.

But in terms of alternate ways of analysis. If you are still referring to the same paper, then it shows that there is an analytic solution to the equation. What you can try to do is substitute this solution in your numerical scheme and see whether you can do something with it.

  • $\begingroup$ No, I'm currently working on my own project solving a nonlinear PDE, but I use that paper as inspiration. $\endgroup$ – Hunter Mar 7 '16 at 17:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.