# Solving a complex ODE with large number of variables (>1e6 variables) - best practise?

I have to solve a non-linear ODE of the shape $$\partial_zA=f(A)$$ with $$f$$ a non-linear function and $$A$$ a matrix/vector with >1e6 variables (i.e. $$A$$ is a matrix with >1000x1000 entries). For each calculation of the right hand side of the system (i.e. $$f(A)$$) I have to calculate roughly five matrix-matrix multiplications $$B\cdot A$$ with $$B$$ a constant matrix and roughly five FFTs over each row of $$A$$. $$A$$ is using complex numbers.

To solve this system I intended to use a solver such as odeint from boost. I also intended to run the calculation on a GPU, to use the increased speed of matrix-matrix-multiplications and batch FFTs. I had to find out, though, that odeint does not support thrust::complex as value type. Therefore, my current approach is to copy $$A$$ to the GPU in each step, and return $$f(A)$$ afterwards to the solver. This is rather costly, and therefore I would like to avoid it.

Thus, are there other strategies (apart from writing my own solver) which could help solving this PDE in reasonable time? Should I still try to go for the GPU, or rather try things like MPI (or a hybrid approach)? Moreover, I noticed that for several sets of parameters the equation behaves rather stiff. Are there options for solving the system in that case?

• Good answers are likely to depend on properties of $f$. What is the original equation you're trying to solve? – David Ketcheson May 6 '20 at 9:57
• @DavidKetcheson: It's a rather long equation, which can be condensed into that function $f(A)$. Moreover, depending on the parameters for the starting conditions I can either add or remove terms from $f(A)$, therefore I intended to keep the question more general. Is an answer still possible without that information? – arc_lupus May 8 '20 at 10:54
• In my opinion, no, it's not possible to give a very useful answer with the information you have provided. – David Ketcheson May 9 '20 at 9:33