Currently, the calculations are performed using a single array for each computational size. For multi-core calculations, multithreading is an option. Several frameworks exist that, with varying emphases, allow for the parallelization of NumPy computations. Numba is a strong contender, but it doesn’t currently support masked arrays. Dask is another fascinating project, though its focus is on parallelizing the processing of massive datasets across multiple machines and handling ‘out-of-memory’ computations. As a result, testing with Dask for this specific scenario didn’t yield any significant performance gains.
I chose to handle parallelization by implementing my own subdivision. The model area is automatically split into multiple tiles based on its size. Then, these tiles are processed in parallel, depending on the number of available CPU cores, using Python’s ThreadPoolExecutor. In this approach, tiles for each region are recreated from the main array each time, and after computation, the results are written back to the main tile.