Hi!

I have developed a Mathematica code to solve the Boltzmann equation governing the evolution of a physical system. It is based on the so-called Direct Simulation Monte Carlo approach: one starts with a huge number of particles and simulates their collisions. Specifically, I build a table with each row representing a particle (with, say 10^7 rows), split it into n >> 1 sub-tables ("cells"), and launch the interaction kernels for sub-tables in parallel to simulate the interactions within each cell, isolated. Each interaction may modify the row (elastic scatterings), add new rows (backward annihilation), or remove rows (annihilation).

I heavily used built-in compilation features, with CompilationTarget->"C". However, I have reached the limits of my expertise in tuning its performance and overall design. Some of the most annoying points are that at the beginning of each timestep, I have to shuffle the whole table (mimic particles' migration from cell to cell), i.e., use RandomSample[table], which takes 20% of all the computational time; also, the initialization of the parallel routine governing dynamics within sub-tables (I do this using ParallelMap[compiledcode[#],subtables] - it is much faster than ParallelTable[compiledcode[tab],{tab,subtables}]) takes 40% of all the time just because of redistributing the large table over kernels.

For some reason, I prefer the main environment to simulate the dynamics to be Mathematica. While it may be, in principle, possible to make a C++ library and load it using LibraryLink (potentially preserving the C++ performance), I am wondering if it may be possible to improve the code such that it is maximally efficient within Mathematica only, to maintain flexibility. The C++ version (not implemented in Mathematica, though) already exists, having a much faster shuffling of the table and distributing the sub-table over kernels.

I seek advice regarding these points, and would be grateful for any help. However, given the large code, it may be unsuitable to copy-paste pieces of code here. Is there an option for how to proceed?