Second order 2D differentiation DiffEqOperators.jl

Specific Domains Numerics
2

Hi,

What you want to do is:

dx = 0.01
axis = collect(-1:dx:1);
umesh = [(i, j) for i in axis, j in axis];

sines((x,y)) = sin(pi*x)*sin(pi*y)
u = sines.(umesh)

Qx, Qy = Dirichlet0BC(Float64, size(u))

Dxx = CenteredDifference{1}(2,2,dx,N)
Dyy = CenteredDifference{2}(2,2,dx,N)

A = (Dxx + Dyy)*compose(Qx, Qy)

laplacian_of_u = A*u

Doing things in this way improves performance, as DiffEqOperators fuses the stencils of the operators, fusing loops and leveraging GPU convolutions from NNlib.

For the extension to 3D you do the obvious thing, note that the BC constructors return a variable number of Q operators depending on ndims(u), and compose takes any number of arguments.

I’m working on DiffEqOperators and am happy to support if you have further questions, let me know how you get on!

1 Like