9

can you give me a link for a good and simple explanation of the Shortley-Weller finite-difference scheme? I tried to google it but all I get is (inaccessible) academic publications. I also tried reading the dedicated chapter (4.8) in the Wolfgang Hackbusch's book "Elliptic Differential Equations" but I found it rather difficult. Thank you

To Christian Clason: thank you for you answer , but there is a thing I still don't understand: what do I have to do to apply this method be to arbitrary boundaries, for example an asymmetric airfoil in a flow?

Community
  • 1
Michael
  • 91
  • 2

1 Answers1

9

As far as I can tell, this scheme just consists in replacing the uniform finite difference stencil near the boundary with a non-uniform stencil (with at least one point shifted to lie on the boundary). Basically, you take your arbitrarily shaped domain, put it in a box, discretize the box with a uniform grid, throw away all grid points that do not have at least one neighbor inside the domain, and shift the remaining grid points outside the domain either horizontally or vertically (whichever is shortest) so that they lie on the boundary. (The actual implementation is much more tedious, of course.)

To obtain the non-uniform stencil at one of the nodes next to a boundary node, one proceeds similarly to (one of) the derivations of the uniform stencil: Interpolate the (unknown) function by a quadratic polynomial in the nodes and take the second derivative. It suffices to consider the one-dimensional case with the nodes $x_1=x-h_1,x_2=x,x_3=x-h_2$. Then

$$D^2_h u(x) \approx u(x-h_1)\ell''_1(x) + u(x) \ell_2''(x) + u(x+h_2)\ell''_3(x),$$

where $\ell_j=\Pi_{i\neq j}(x-x_i)/(x_j-x_i)$ are the Lagrange polynomials corresponding to the nodes. Computing the derivatives yields

$$D^2_h u(x) = \frac{2}{h_1(h_1+h_2)} u(x-h_1) - \frac{2}{h_1h_2} u(x) + \frac{2}{h_2(h_1+h_2)} u(x+h_2)$$

as claimed. (You can also use the Newton form of the interpolating polynomial, which simplifies computing the derivatives, especially for higher orders.) Doing the same in $y$ and summing the stencils gives equation (4.8.7).

You can find more detailed examples in Randy LeVeque's Finite Difference Methods for Ordinary and Partial Differential Equations (e.g., page 9), or on this blog post (which also contains NumPy code for computing the coefficients given arbitrary $h_1$ and $h_2$). This is also treated in detail in Morton and Mayers, Numerical solution of Partial Differential Equations, section 3.4.

How you treat the boundary nodes depends on your boundary conditions. For Dirichlet conditions, you proceed as you would for a uniform mesh. For Neumann conditions, you use the above approach (non-uniform interpolation -- now simultaneously in $x$ and $y$ -- and differentiation) to approximate the normal derivative at the boundary node to get a local stencil; see Morton and Mayers, page 75ff.

Christian Clason
  • 12,301
  • 3
  • 48
  • 68
  • If we replace $h_1, h_2$ with just $h$ then $D_h = \Delta_h$. But $\Delta_h$ is derived using finite differences, i.e. Taylor series. Can you please clarify how the Lagrange interpolation can coincide with Taylor series / finite differences? – sequence Feb 22 '19 at 13:31
  • 2
    @sequence, one can derive finite difference approximations as derivatives of the interpolation polynomial constructed using the nodes of the stencil. So, you take the points of the stencil which you want to use, construct the Lagrange polynomial and differentiate it to derive the formula for approxation of the derivative you want. – VorKir Feb 23 '19 at 02:37
  • @VorKir It is remarkable, nevertheless, that these approximations will coincide with Taylor series approximations. – sequence Feb 23 '19 at 14:59
  • 1
    @sequence Not really, since a Taylor approximation is a linear (quadratic etc.) polynomial approximation to a function at a point, and so is the interpolating polynomial. By the fundamental theorem of algebra, these have to be the same. (If you are still sceptical, remember that the points for the stencil are not chosen randomly.) – Christian Clason Feb 24 '19 at 09:32
  • @ChristianClason Thanks for clarifying. Indeed, I had not paid significant attention to the fact that both approximations are polynomials, which must be unique, and hence about the applicability of the fundamental theorem of algebra. – sequence Feb 24 '19 at 18:38