I am trying to solve a reaction-diffusion in Julia.
I’ve worked with ModelingToolkit and the SciML ecosystem before, but it seems that with PDEs there is not a single best way to solve them.
The simplest way I found to solve the PDE is using MethodOfLines.jl, but the discretization is quite slow for slightly large 2D problems.
I’ve also seen ParallelStencil.jl which looks promising and has been suggested (Issue solving 2D wave equation with ModelingToolkit.jl - #10 by ChrisRackauckas), but I was not able to figure it out how to use with ModelingToolkit.jl.
I am in no need of complex geometries or anything, just want to work with simple rectangular toroidal domains or use Neumann boundary conditions at the edges (du/dt = 0).
Is there a simple way to solve a system like that? Any reccomendations are appreciated.
There is not a single best way to solve such systems right now. The best way performance is still to write it out by hand, and there are some tutorials that show how to do this the best way:
ParallelStencil does not make sense for these types of equations as the overhead of the reaction term makes it so that way the stenciling loop shouldn’t be optimized and instead its more optimal to fold it in with the reaction loop, meaning that you do a single naive pass. This is only not true if the reactions are linear, and thus you have a fully linear system, and then you would still merge the two components into a loop but now you’d have a fully linear stencil and use ParallelStencil. But having a linear reaction essentially never happens so that’s just an oddity to mention in passing.
MethodOfLines.jl does this codegen to the fused form but is currently not loop preserving. We’re working on getting that into JuliaSimCompiler and via array equation support, but currently it won’t do as the by-hand solution here.
So the best solution is to keep using DifferentialEquations.jl and writing the equations as functions by hand instead of creating ModelingToolkit.jl models? Or is there a way to do the same with ModelingToolkit.jl?
Is there a way to use ModelingToolkit.jl models/equations in Gridap.jl? Or interoperability with SciML packages? I know about the package and would be awesome if that’s the case.
We recently set up an example for connecting VoronoiFVM.jl with Catalyst.jl, which is bases on ModelingToolkit.jl: Coupling with Catalyst.jl · VoronoiFVM.jl This is an example for hererogeneous reactions. May be I will make another example just for reaction-diffusion.