How to work around lack of dynamic profiles in porting from Epetra to Tpetra?

[jjramsey]

As I mentioned in an earlier question, I've started an attempt to port an old code from Epetra to Tpetra, since Epetra is being deprecated.

There are a couple uses of Epetra in the old code. In one, I'm doing the usual routing of assembling a global matrix for a linear finite element calculation: looping over elements "owned" by a processing element and summing their stiffness matrices into the global matrix. I do this in two steps, first assembling an Epetra_FECrsGraph by looping over the elements almost as if I were assembling the matrix proper (except not calculating the stiffness matrices), and second by constructing an Epetra_FECrsMatrix from the graph and then looping over the elements again and summing their stiffness matrices into this matrix. (I think this was done in an example code from Didasko.)

The other use was to generate an adjacency list from a list of elements and the nodes that belonged to each of them. Again, I was looping over elements "owned" by a processing element, and the process went much like the aforementioned assembly of the graph.

In both of these uses, I relied on Epetra_FECrsGraph having a dynamic profile. It was unclear to me how many nonzeros per row I needed, so I let the process of finite element assembly do it for me. Now, with the constructors of Tpetra's CrsGraph and FECrsGraph all requiring upfront knowledge of the number of nonzeros (or something that entails that knowledge), the lack of a dynamic profile in Tpetra presents a problem.

Are there known workarounds for this issue, or some sane straightforward way of finding the number of nonzeros that I hadn't thought of?

Tagging @trilinos/tpetra

[cgcgcg]

@trilinos/tpetra

[cgcgcg]

I think the most reasonable option is to loop over the elements and construct an upper bound of the entries per row. Cause you know whether a dof is shared between elements or interior. Tpetra does not know that. I believe the decision for not supporting the dynamic profile anymore was prompted by performance considerations.

[jjramsey]

@cgcgcg:

I think the most reasonable option is to loop over the elements and construct an upper bound of the entries per row. Cause you know whether a dof is shared between elements or interior.

When you were talking about looping over elements, did you have in mind looping over the elements one owns or owned+ghost elements?

[cgcgcg]

I have never used it, but here is my intuition for how this works. You want numEntPerRow in the construction of Tpetra::FECrsGraph to account for all your local inserts into the overlapping graph. So you need to loop over owned+ghosted and either determine the rank-local max per row or a row-wise upper bound. Once you call globalAssemble on the FECrsGraph Tpetra will move the matrix to the non-overlapping row map.

[csiefer2]

IiRC, you can also loop over only your owned elements and the graph will expand itself during endFill()

[jjramsey]

IiRC, you can also loop over only your owned elements and the graph will expand itself during endFill()

@csiefer2, I thought that the lack of a dynamic profile meant that the graph couldn't expand itself. Doesn't maxNumEntriesPerRow mean the max number of entries per row when taking account contributions to the graph from all processes, rather than just the contributions from the owned process itself?

[csiefer2]

The graph can't expand itself in the insert*Entries() phase. That rule doesn't quite apply to fillComplete / endFill. I suspect we don't actually document that fact.

[jjramsey]

The graph can't expand itself in the insert*Entries() phase. That rule doesn't quite apply to fillComplete / endFill.

So if I understand correctly, maxNumEntriesPerRow means the max number of entries per row that are made locally (and probably is used to size some internal arrays in the FECrsGraph object that cannot be resized during the the insert*Entries() phase). It should not require communication across processes to estimate.

I am a bit concerned about the qualification "quite" in "That rule doesn't quite apply to fillComplete / endFill," which makes me wonder if I am right about it not requiring communication across processes to estimate maxNumEntriesPerRow.

I suspect we don't actually document that fact.

Yes, that is not documented, or at least isn't documented clearly. Indeed, documentation for the newer parts of Trilinos (and even many of the older parts) is fairly sparse, and some of it is out of date. For example, there's a Tpetra tutorial describes the second integer argument to the CrsMatrix constructor as a size hint, while the documentation for that constructor calls it a "strict upper bound". But that's a bit off topic for this discussion.

[csiefer2]

I believe the qualification was unnecessary on my part.