First Order Schemes

[jettrich]

Hi, in case I want to use findiff for educational purpose, showing the differences in accuracy for first, second and higher order schemes, how can I enforce the usage of pure forward/backward differencing? Furthermore, e.g. thinking of advection, is it possible to realize upwind or TVD schemes with findiff? Any help is greatly appreciated, kind regards!

[maroba]

Hi,

if you want to enforce a specific difference scheme, you can use the Stencil class in findiff.stencils. Assuming you are using the currently latest version, 0.9.2, an example in 1D would read:

import numpy as np
from findiff.stencils import Stencil

# some example data:
x = np.linspace(0, 1, 101)
dx = x[1] - x[0]
f = np.sin(x)

# define an explicit forward stencil for the first derivative along axis 0:
stencil = Stencil(offsets=[0, 1, 2], partials:{(1,): 1}, spacings=dx)

Note that the partials keyword is kind of awkward. It allows you to define a linear combination of partial derivative with constant coefficients with one single dict. The keys are tuples and each item of the tuple specifies the degree of the derivative along the axis. The values are the constant coefficients. That is, if you have a 2D grid, a partial derivative derivative along the axis 1 would read partials={(0, 1): 1} and a mixed one (1, 1): 1. Anyways, continuing with the 1D example:

stencil.accuracy
Out: 2

and

stencil.values
Out: {(0,): -150.0, (1,): 200.0, (2,): -50.0}

Apply the stencil at a given point, e.g. at x[5]:

stencil(f, at=5)
Out: 0.9987835384105308

Doing this point-wise is inefficient, though. But you can also apply the stencil to slices and masks instead of single points by using the on argument instead of at. This is much faster, of course.

I guess one could easily realize an upwind or TVD schemes using the Stencil class.

By the way, I'm about to release a new major version (1.0.0) in the next two weeks or so. It will contain symbolic capabilities which allows one to derive discretization schemes and analyze their stability behavior analytically. Maybe that's interesting for you.