Best Methodologies for Managing a Mesh in Parallel Finite Element Computation?

Question

I am currently developing a domain decomposition method for the solution of the scattering problem. Basically I am solving a system of Helmholtz BVPs iteratively. I discretize the equations using finite element method over triangular or tetrahedral meshes. I'm developing the code towards my Phd thesis. I am aware of some of the existing finite element libraries out there such as deal.ii or DUNE and though I think they are great, with inspirational design and API, for learning purposes I wanted to develop my own little application from scratch.

I am at a point where I have my serial versions running and now I want to parallelize them. After all, it is one of the strengths of the domain decomposition framework to formulate algorithms that are easy to parallelize, at least in principle. In practice however, there are many details one must consider. Mesh management is one of them. If the applications is to achieve high resolution while scaling well to many CPUs the replication of an entire mesh on every CPU is inefficient.

I wanted to ask those developers who work on similar applications in high performance computing environments how they deal with this issue.

There is p4est library for distributed mesh management. I do not need AMR so it might be an overkill since I'm only interested in using uniform meshes and I'm not sure if it can refine triangular meshes. I could also simply create a uniform mesh then feed it into one of the mesh partitioners and do some post processing of the output.

The simplest approach seems to create a separate file for each partition containing mesh information relevant only to that particular partition. This file would be read by a single CPU which would be responsible for assembly of the discrete system on that portion of the mesh. Of course, some global partition connectivity/neighborhood information would also need to be stored in a file read by all CPUs for inter process communication.

What other approaches are out there? If some of you could share, what are some of the commonly used methodologies in the industry, or government research institutions related to handling this issue? I am quite new to programming a parallel finite element solver and I wanted to get a feel for whether or not I'm thinking about this problem correctly and how others are approaching it. Any advice or pointers to relevant research articles would be greatly appreciated!

Thanks in advance!

Answer 1

If you are not using AMR and do not want to scale beyond 1K-4K cores then simply do this.

1. Rank 0 reads the entire mesh and partitions it using METIS/Scotch etc. (Note: This is a serial operation).

2. Rank 0 broadcasts the element/node partitioning info to all other ranks and frees the memory (used to store the mesh)

3. All ranks read the nodes/elements they own (including ghost nodes) from the same input file (Note: 2000 ranks accessing the same input file might sound slow but is not in practice, though it may be bad for the file system but then we are doing it only once).

4. All ranks need to create the local to global node/element/dof mappings for application of BCs and assembling of matrices and renumber the nodes.

After everything is said and done all data on a rank will be local so you should be able to scale well (memory wise). I do all this in about 100 lines (see lines 35-132 here) in a small code of mine.

Now if your mesh is too large (e.g., >100-250 million elements) that you cannot partition it using METIS on a single node and need ParMETIS/PT-Scotch then you have the additional work of partitioning it in parallel before all cores/ranks can read it. In such a scenario it might be easier to keep the partitioning phase separate from the main code for logistical reasons.

Btw AMR libs usually dont do tets. Also PETSc is good choice for parallelization of your code.

Edit: Also see here and here.

Answer 2

This may not come as a surprise to you given that I develop deal.II, but here's my perspective: When I talk to students, I typically tell them to develop their own prototype in the beginning so they can see how it's done. But then, once they've got something small running, I make them use a library that allows them to go so much further because they don't have to reinvent the wheel with basically every step they take.

In your case, you've already seen how to implement a simple Helmholtz solver. But you'll spend the next 6 months writing the code necessary to do it in parallel, you'll spend another 3 months if you want to use more complicated geometries. You'll then spend 6 more months if you want an efficient solver. And all of this time you're writing code that's already been written by someone else and that, in a sense, doesn't get you any closer to what you actually need to do for your PhD: develop something new that hasn't been done before. If you go down this road, you'll spend 2-3 years of your PhD time re-doing what others have done, and maybe 1 year doing something new.

The alternative is that you now spend 6 months learning one of the existing libraries, but after that you'll have 2-3 years where you really do new stuff, things where every other week you can walk into your adviser's office and show him/her something that's truly new, that runs on massively large scales, or is just very cool in other regards. I think you probably see where I'm going with this by now.

Answer 3

This is not a complete answer.

For implementation of parallel domain decomposition methods, I encountered some complications. First, one can use many processors for one subdomain or feed one processor with many subdomains and one might want to implement both paradigms. Second, substructured form of domain decomposition methods requires separating out subdomains' (not element's) faces, edges, vertices. I do not think these complications are readily included in the parallel mesh management. The situation becomes simpler if you consider one processor for one subdomain and using the overlapping RAS/RASHO method. Even in this case, you would better manage your parallel layout yourself, because when your subdomains have overlap the usual parallel mesh management may require communications for the entire overlapping regions while the necessary data for RAS/RASHO are only on the boundaries.