I have two fields, A (complex scalar field) and B (real vector field), which are sampled on staggered grids. This means that (in, say, three dimensions) A is an complex array with N^3 elements, while B is a real array with 3*(N-1)^3 elements. I also have corresponding equations of motion with dA = fA(A,B) and dB = fB(A,B). I’m rolling my own integrator, but I would prefer to switch to DifferentialEquations.jl.
I can see a trivial, even if inconvenient, solution to the difference of array shapes by flattening and stacking both into a vector and then view-ing reshape-ing them into correct form inside the problem definition.
My bigger concern is the complex/real split between the arrays. This type difference has physical relevance and is enforced through the code base, so I would prefer to not cast the real field into a complex-valued field just for sake of this numerical step. There is also space, where in practical situations, N will be large (easily N=10^4, could increase magnitude or two) and I’m not too keen on almost doubling the memory footprint just to store billions of zeros.
What is the most graceful way of solving a system of equations like this?