Solving vector field PDE with DiffEqOperators.jl

Hello all,

I am wondering if there is a workaround to solving a vector field for lets say Maxwell’s equations such as

\frac{\partial \mathbf{E}}{\partial t} = \frac{1}{\varepsilon_0} \times \mathbf{H}

\frac{\partial \mathbf{H}}{\partial t} = \frac{1}{\mu_0} \times \mathbf{E}

where \mathbf{E} and \mathbf{H} are vectors, with DiffEqOperators.

From my last post, I figured out how to calculate scalar fields, so essentially now I would just need to calculate three scalar fields at once per vector field (if in 3D) since in practice, one would just solve for each component of the vector field at a time and iterate through them.

This is similar to the example from the https://diffeq.sciml.ai/stable/tutorials/ode_example/ where a system of ODEs is solved using


function parameterized_lorenz!(du,u,p,t)

  x,y,z = u

  σ,ρ,β = p

  du[1] = dx = σ*(y-x)

  du[2] = dy = x*(ρ-z) - y

  du[3] = dz = x*y - β*z

end

The main problem with extending it with lets say


function maxwell!(du,u,p,t)

  e_x, e_y, e_z, h_x, h_y, h_z = u

  eps, mu = p

  du[1] = de_x = (1/eps) .* ((D_y*Q*h_z) - (D_z*Q*h_y))

  du[2] = de_y = (1/eps) .* ((D_z*Q*h_x) - (D_x*Q*h_z))

  du[3] = de_z = (1/eps) .* ((D_x*Q*h_y) - (D_y*Q*h_x))

  du[4] = dh_x = (1/mu) .* ((D_z*Q*e_y) - (D_y*Q*e_z))

  du[5] = dh_y = (1/mu) .* ((D_x*Q*e_z) - (D_z*Q*e_x))

  du[6] = dh_z = (1/mu) .* ((D_y*Q*e_x) - (D_x*Q*e_y))

end

(mind that this is a very rough example) is that D*Q*u will already give a vector (or a matrix I’m not exactly sure), so we cannot index du in this case in the same way that is possible for the ODE system of equations. I’m wondering, how would someone set up this system using the ODE function form given that several fields must be solved concurrently and share values between each other. Would this approach of defining everything in one function be best, or instead would it be better to define each component as a separate function and pass the values between them. Likewise, is it possible to solve this using the ODE solve framework present in DifferentialEquations.jl or would I have to create my own custom time stepper to iterate each piece. Thanks.

Take a look at RecursiveArrayTools.jl, which exports a VectorOfArray type that should work well here.

1 Like