Mesh partitioning for Multidomain problem

[lmolin3]

Hi all,
I'm working on a simple multidomain problem using domain decomposition.
I managed to extend the multidomain miniapp and solve the coupling between the two domains using Neumann-Dirichlet and a relaxation step.
I wanted to test on a geometry with 3 domains and developed the mesh for this in Gmsh.
The code seems to work in serial, but the problem I'm having in parallel seems related to the partitioning of the mesh.
Each domain has its own solver (advection-reaction-diffusion) and bchandler, which rely on the domain's ParSubMesh.
The code fails due to one or more ranks not having any elements of the specific domain.
Any suggestions in these regards?

Best,
Leo

[pazner]

In principle, MFEM should work even when there are ranks with no elements, but there are some parts of the code that don't currently support this (see e.g. PR #4260 for some improvements to this).

If you are running into problems with empty ranks, it should probably be considered a bug in MFEM and should be fixed.

[jandrej]

Agree. If you can provide a minimal failing example and file an issue, that would be great.

[lmolin3]

Thank you so much for the prompt reply!

I tracked down the specific segfault when running on 4 cores and it was actually related to

bool tensor = UsesTensorBasis(H1FESpace);

which seems to be fixed by PR #4260, since for empty partition GetFE(0) was failing.

Running on 2 cores there seems to be no issue with empty partitioning, but the CG solver still fails.

I'll come up with a minimal failing example, as at the current stage the code relies on the Solver, Operator, BCHandler classes I wrote (https://github.com/lmolin3/Cardiac-Ablation-ECM2.git, multidomain-ecm2 branch).

Once it's done I'll open an issue for this.

Best,
Leo

[jandrej]

CG should work for empty ranks

[lmolin3]

My bad, the failing test was with setting order>1 of the H1 fec or with serial/parallel refinements of the above mesh, but doesn't have to do with the empty ranks issue.

Besides that segfault fixed by the PR, the issue seems to be on the first assert in FindPointsGSLIB::DistributePointInfoToOwningMPIRanks due to the rank not having pts on the specific interface. For the Neumann-Dirichlet coupling I'm using the SubMesh::Transfer for Dirichlet, and GSLIB custom interpolation (as suggested in #4258 ) for transferring the heat-flux across interfaces.

[lmolin3]

I think #4461 and #4260 solved it for me. I was having some issues with the initialization of the meshgen variable for ParSubMesh on empty rank but i opened a separate issue for it.
Thank you again!