[petsc-users] MatNest question

[Jed Brown]

Alex Fleeter <luis.saturday at gmail.com> writes:

> Hi:
>
> I have been trying MatNest. I have two questions.
>
> 1. I have to call MatAssemblyBegin/End for both sub-matrices and for the
> MatNest object. From my experiment, if I miss any of those, I will get an
> error message
> [0]PETSC ERROR: Object is in wrong state
> [0]PETSC ERROR: Not for unassembled matrix
>
> I am a bit confused because in the example
> https://www.mcs.anl.gov/petsc/petsc-current/src/snes/tutorials/ex70.c.html The
> MatAssemblyBegin/End are only called for submatrices. I assume that example
> works correctly.

The example is run in continuous integration.  You'll need to provide more details (like a minimal reproducer) if you want us to explain why assembly is needed in your usage, or if it can be relaxed.

> On the other hand, in the actual implementation,
> https://www.mcs.anl.gov/petsc/petsc-current/src/mat/impls/nest/matnest.c.html#MatAssemblyBegin_Nest,
> I can see that the MatAssemblyBegin/End are already called for each
> sub-matrix. It seems that one should only call one MatAssemblyBegin/End for
> the big nest matrix. Unfortunately, I tried that and got an error message
> shown above. I must have misunderstood something.
>
> 2. From previous threads, Jed has strongly suggested using standard Vec
> instead of VecNest, unless there are strong reasons. In my case, since
> MatNest is already used as the underlying matrix object, the local element
> assembly will generate residual vectors for separate physics. 

Please look at src/snes/tutorials/ex28.c.  Also, VecGetSubVector() is no-copy when the subvector is contiguous, so you can assemble into "separate" vectors that are actually part of the same contiguous Vec.

VecNest incurs extra overhead in all the usual (e.g., Krylov) operations.  I don't think you have a good reason to use it.  Yes, I'm that troll under the bridge trying to scare you off.

> So it is rather natural to put the right-hand side residual vector as
> a nest vector.  Of course, a consequence is that I have to put every
> vector object in the nested format, since my experiments suggest
> VecNest cannot be used with standard Vec. My question is this. How
> mature is the VecNest implementation? Is it a good design to put every
> vector in the nest format?
>
> Thanks,
>
> Alex