# Finite-difference software for solving custom equations

Are there any good, easy to use, software for simulating the evolution of systems of generic differential equations? I know there are custom programs for various specific circumstances (such as electrodynamics). But say I have some equations and I quickly want to simulate the evolution of some initial state; without reinventing the wheel each time, I would love some plug-and-play program where I can just state the equations, initial conditions, boundary conditions, discretisation criteria etc. and simulate the system. Maybe even an option for easily switching between different solver methods.

Is there anything like this? Thanks :)

• This looks interesting: fenicsproject.org
– Tom
Aug 25 '20 at 10:19
• PDE's or ODE's? Aug 25 '20 at 14:55
• PDEs. Specifically things like Gross-Pitaevskii equations for simulating the evolution of the macroscopic wavefunction of exciton-polariton condensates.
– Tom
Aug 26 '20 at 17:56
• That does not look too bad; but the equation is for a complex-valued function. I am pretty sure all those standard PDE packages assume real-valued equations, so you'd have to split your equation into real and imaginary parts. Aug 26 '20 at 21:20
• The julia library should be able to deal with complex numbers just fine (including things like arbitrary precision complex numbers and matrices with arbitrary element types). Aug 26 '20 at 23:45

I'm going to assume since you mention electrodynamics that you're interested in PDEs. You've already mentioned FEniCS in your comment. FEniCS offers a domain-specific language (DSL) called UFL. This language makes it simple to express the weak form of a PDE for discretization via the finite element method. There's another package called Firedrake*, which also uses UFL to express PDEs. Firedrake has built-in support for extruded meshes, where you take a 2D domain and lift it into 3D. I develop applications in the geophysical sciences, and this feature is absolutely indispensable for me. On the other hand, the FEniCS developers are much more heavily invested in biological problems and consequently have worked a lot on solving mixed-dimensional problems (see this paper).

You mentioned wanting to use the finite difference method, so you should also be aware of Devito. Like FEniCS and Firedrake, Devito also uses an embedded DSL to describe PDEs. For an example, you can see this sample code for a seismic wave propagation problem.

With each of these libraries, you should be freed from the burden of writing the low-level kernel code to form the linear or nonlinear system corresponding to the discretization of your problem. This is possible because they either call out to skeletonized low-level kernels and inject code for your particular problem, or they completely generate the low-level kernels on the spot. In other words, the library automates the translation from your PDE to a discretized linear or nonlinear system.

The library does not automate away the process of deciding how to solve that linear or nonlinear system once it has been formed. Firedrake, for example, is very tightly integrated with PETSc. You can make choices about methods by passing options to PETSc's linear solvers, for example to use LU, GMRES, or CG to solve linear systems.

Opinions of whether these libraries are A Good Thing among computational scientists vary. On the one hand, they provide you with a convenient way to set up and solve physics problems without having to write a mess of low-level code, so they are great for prototyping. Some have argued that they are ultimately a hindrance if you want to do real work because (1) there are insufficient escape hatches for when you need to do weird or unusual things and (2) you give up all control over how your problem is formed. Personally, I am happy to give up this control and I avoid doing weird or unusual things, but this is a choice you have to make for yourself.

*I'm an occasional contributor to Firedrake and I have grant funding to develop software built on it.

• It looks like Firedrake is finite element, not finite difference. Devito is definitely finite-difference based, and I've not looked at FEnicCS (but the name suggests finite element too...). If the OP is really interested in only finite-difference code then devito might be the only real option amongst these. (I don't want to enter a FD vs FE war here... it just seems the apps support different things) Aug 26 '20 at 0:08
• Yes both Firedrake and FEniCS are more in the finite element mind but I included them because the question also has the "finite element" tag. I was assuming Tom was most interested in tools that would help with rapid prototyping and application development with a preference for FDM but that FEM was ok too. Also not interested in arguing about which approach is best, just pointing out which tools are out there. Aug 26 '20 at 0:30
• Thanks a lot for the suggestions. I will check everything out. I currently don't understand all the subtle differences between the various methods. I just want something that solves PDEs as a function of time and space in very simple meshes. I have developed my own C++ code in the past based on the Runga-Kutta method, but I'd much rather rely on well-developed open source code. When working with non-linear equations it's very easy for bugs in the code to hide and pretend they are giving real physics!
– Tom
Aug 27 '20 at 10:10
• @Tom that's for sure! Many years ago I was writing some code to simulate a non-Newtonian fluid, where viscosity $\sim$ strain rate${}^{-2/3}$. Like a doofus I put 2 / 3 into my code, forgetting that integer division gives integer results, so this value is equal to 0. I spent months simulating a Newtonian viscous fluid before I thought to test against an analytical solution and found out that this was totally wrong, but in the meantime the results were more or less reasonable -- it still walked and talked like a fluid! Aug 27 '20 at 16:34

You might want to check out DifferentialEquations.jl. It supports ODEs, PDEs, stochastic equations, delay equations, and basically everything else. It also has really good automatic sparsity detection and ability to work on GPU for big systems. DiffEqOperators.jl is the submodule which has the automated finite difference operators (with lazy stencil application) for easy PDE discretization.

The BOUT++ project http://boutproject.github.io offers a set of tools for finite-difference solution of systems of PDEs, primarily targeting fluid dynamics and plasma physics but not limited to those application areas. Equations in human-readable form are automatically discretized in space and integrated in time, using efficient parallel solvers.