I've been working on a finite element code on unstructured methods, which I've parallelized using the Schur complement method. Here's a summary of how I did it:
- Assign each triangle of the mesh to a domain
- For each node, determine which domain it is in or whether it is on the artificial boundary
- Send each processor its part of the mesh -- the triangles it owns as well as any nodes it owns, even if it shares those nodes with some other domain
- Build the matrices, in parallel, on each processor
- Solve some stuff, using some way of exchanging node data between processors
I'm solving a quasilinear system of elliptic PDE, so each linear solve is just an iteration of Picard's method. In particular, I need to change the entries of the stiffness matrix at every step, since it depends on the gradient of the solution. (The non-zero structure of the matrix doesn't change.)
This all works just fine. However, it disagrees with the usual approach I see taken.
My code divides up the triangles between processes, and some of them have to share nodes. This means a extra storage, but the entries of the stiffness matrix can be filled completely in parallel. On the other hand, libraries like PETSc and METIS divide up all the nodes into disjoint sets and send them, along with their corresponding rows in the stiffness matrix, to each process. But, you have to communicate to fill the stiffness matrix.
So: why do the big scientific libraries have a marked preference for dividing up the nodes instead of the elements? Am I missing something?