I’m trying to solve the Laplace equation in 3d on orthogonal grids that are partially trimmed, i.e. some of the internal grid points are removed from the domain. These grids come from micro-ct tomography images and correspond to the void space inside a 3d porous object.
I only need to find the effective diffusivity of the tomography images and for that, I simply need to solve the Laplace equation with two Dirichlet boundary conditions imposed on two opposed faces of the simulation box.
The grid size is about 500-cubed, but roughly half the internal grid points correspond to the solid phase and thus are removed from the original domain. So, I’m left with 50 to 70 million unknowns to solve. Needless to say, the resulting sparse matrix is banded (7 diagonals on average), symmetric and positive-definite.
I have access to a workstation with up to 500 GB of memory and up to 100 CPU cores. I’m looking for the quickest way possible in Julia (that harnesses the 100 CPU cores I have access to) to solve this system since I need to do this for several tomography images.