A quick question, as I read the comment from @ivanpribec
He mentioned many old f77 fortran ode solver code are not thread-safe.
Uhm, I guess that is probably an issue for openMP, since those old code uses some common block perhaps?
But I wonder, if I use MPI, does the ‘thread-safe’ issue still matter?
It seems in MPI, each core or thread actually run the program independently, and there will be no communication between the cores or threads if I do not let them communicate.
So, in MPI, is there still thread-safe issue?
Answering your specific question, as you say, each MPI “rank” is an independent process (not a thread) so there’s no issue of thread safety. But adding MPI to an application takes significant work, especially if the code uses global variables for communication. (Coarrays might be a better bet, but still requires effort.) You would be better served to look for solvers that were designed to support parallelism.
I may be misunderstanding things, but MPI is a library and if you are building your own library on top of MPI, you will have no control over the execution environment. People could be calling your library routine from multiple threads and if any internal data structures are shared, results can be unexpected.
It’s not only common blocks, it’s also procedure variables with save attribute.
If you need to solve multiple systems of parametrized ODE’s (e.g. for different combinations of parameters or many different initial conditions) then OpenMPI or co-arrays could be an option. Here’s a table (original source) which highlights some differences between SPMD and co-arrays:
If you need to solve a single large system of ODEs in parallel I’d recommend taking a look at Sundials. Here’s the list of NVector types available in Sundials for different accelerator/parallel programming paradigms:
I see. I usually wrote MPI code and run on supercomputer and never have the thread-safe issue because each rank run independently in MPI.
Roughly speaking, what I did is like solving ODEs from 1000 different initial conditions then average the results. In this case instead of solving the ODEs on a single core for 1000 loops, I run MPI on 1000 cores each with a different initial condition and then collect the results on rank 0 and average the results. In this case, MPI works just fine and it does not matter if the ODE solver is thread-safe or not, because each rank actually run independently. I am not seeking to solve ODEs parallelly.
But I guess I got what you mean. If I were to really solve the ODEs parallelly, then it is better to use those parallel ODE solver. Otherwise modifying those non-parallel ODE solver to make them parallel is not an easy task.