[petsc-dev] Jed!!!
Eric Chamberland
Eric.Chamberland at giref.ulaval.ca
Tue Apr 29 12:32:31 CDT 2014
On 04/27/2014 12:32 PM, Jed Brown wrote:
> Barry Smith <bsmith at mcs.anl.gov> writes:
>
>> I am totally confused and likely users are to. Please point me to an
>> example that actually uses Matnest the way it is “supposed to be
>> used”.
>
> src/snes/examples/tutorials/ex28.c
>
> Note that this does not contain the word "nest". That is as intended.
>
>> I cannot even find a MatSetValues_Nest() so I do not see how one can
>> even interchange code between AIJ and Nest.
>
> You're supposed to use MatGetLocalSubMatrix(). Translating from global
> indices to local indices is a disaster that we want to avoid. So we go
> the other way. Speak the language of "split local" spaces during
> assembly and the data structure itself can live wherever is most
> efficient for the solver.
We have a question about this:
Where is the "disaster"? In our mind, we think the MatSetValues_Nest
should:
1) do the assembly for local lines in the proper submatrices
2) communicate the non-local lines at MatAssemblyEnd() so that each
process receives some lines which will be assembled by method 1)
So to do 1) (the assembly of local lines), you "just have to" convert
column indices to:
a) retrieve which sub-matrices it belongs to
b) for each sub-matrix, do the assembly in that sub-matrix which
involves converting lines/columns indices to -1 for elementary values
not in that matrix.
We coded something (in sequential) in our code to do the work, and
observed surprisingly that it is faster to do the assembly in a nested
matrix:
3.97s vs 4.24s for a MatNest vs 357 369^2 CSR Matrix
(we suppose it could be faster because of the dichotomous search done in
MatSetValues_SeqAIJ)
However, to do the correct work in parallel we think it should be
implemented directly into something like MatSetValues_Nest... Have we
missed something?
The use case scenario is the following: we have a P2-hierarchical (for
example, a velocity field) which is made of a linear part and a
quadratic "correction" part. For a specific preconditioner
(http://onlinelibrary.wiley.com/doi/10.1002/nla.757/abstract)
we want to do the "normal" assembly then work on the submatrices which
are split by the linear vs quadratic parts of the velocity field. The
assembly is done "normally" in our code, meaning we don't want to
duplicate the computations done by the formulation (physical
properties) but instead we want to have the elementary matrix to be
split by the assembly functionality.
Thanks for your insights!
Eric
>
>> The users manual has the word nest several times but no indication
>> of how to use the damn thing.
>
> From the manual:
>
> The key to format-independent assembly is the function
>
> MatGetLocalSubMatrix(Mat A,IS isrow,IS iscol,Mat *submat);
>
More information about the petsc-dev
mailing list