#352 allow higher-order central difference; add MOL_2D test with Dirichlet

[valentinsulzer]

Uses CartesianIndices to apply higher-order central difference method.
Needs special handling for cases close to edges (for example index(2) - 2 gives Index(0) which does exist, so instead we need to use indices 1-5 to get 4th order). Fornberg does this for you but I'm not 100% sure it comes out as being the right order near the boundaries, since forward differences need more points to get the same accuracy https://en.wikipedia.org/wiki/Finite_difference_coefficient

[ChrisRackauckas]

There's some theory on this IIRC. With ghost point boundaries at order 1 you can prove you still get a consistent order 2 discretization because of rate at which the boundaries fall off because of their lower dimensionality. Because of how the proof works, I'd be pretty sure it extends to higher order.

[ChrisRackauckas]

You could take non-centered stencils near the boundary, which is what the derivative operators do, but that seems like too much work with little gain.

[valentinsulzer]

This actually does the non-centered stencils near the boundary (turns out this is what the x0 in fornberg is for, too) - I was just wondering whether we'd need an extra point at the boundary to get the same accuracy but agree the gain is minimal

DiffEqOperators.jl/src/MOLFiniteDifference/MOL_discretization.jl

Lines 118 to 122 in 45f79f2

	I1 = oneunit(first(indices))
	Imin = first(indices) + I1 * (order÷2)
	Imax = last(indices) - I1 * (order÷2)
	# Use max and min to apply buffers
	central_neighbor_idxs(II,j) = stencil(j) .+ max(Imin,min(II,Imax))