<div class="gmail_quote">On Thu, Jan 19, 2012 at 02:42, Klaij, Christiaan <span dir="ltr"><<a href="mailto:C.Klaij@marin.nl">C.Klaij@marin.nl</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="im">> > I have two DMs which I add to a DMComposite and then use MATNEST when<br>
> > getting the corresponding matrix. This gives me block (0,0) and block<br>
> > (1,1). How do I set/get blocks (0,1) and (1,0)? Looking at ex28 I tried<br>
> > MatGetLocalSubMatrix but it gives a null arg...<br>
> ><br>
><br>
> So the problem is that we have no way of knowing what preallocation<br>
> (nonzero pattern) _should_ go in the off-diagonal part. Unfortunately, the<br>
> current preallocation mechanism (DMCompositeSetCoupling()) is a difficult<br>
> thing to implement and the mechanism does not directly apply to MatNest. If<br>
> you have ideas for a good preallocation API, I would like to hear it. I<br>
> need to get back to the preallocation issue because it's an obvious wart in<br>
> the multiphysics support (as long as we don't have fast dynamic<br>
> preallocation, which is a somewhat viable alternative). What I would like<br>
> is for the user to call MatGetLocalSubMatrix() for any blocks that they<br>
> want allocated and set preallocation in terms of the local ordering.<br>
><br>
> The current (unfortunate) solution for MatNest with off-diagonal parts is<br>
> to create the submatrices after DMGetMatrix(), preallocate as you like, and<br>
> copy the ISLocalToGlobalMappings over.<br>
<br>
</div>I see the problem, no ideas for a good general preallocation<br>
mechanism, sorry. Preallocation would depend on how the user<br>
discretizes the cross terms, wouldn't it? So how can we expect<br>
PETSc to "deduce" it from the diagonal blocks? As a user I would<br>
be happy to preallocate myself.<br>
<br>
Out of curiosity, how does it happen in ex28.c: I only see two<br>
matrices being defined (line 347 and 355), I don't see<br>
DMCompositeSetCoupling() at all, yet on line 271 block (0,1) is<br>
available...<br></blockquote><div><br></div><div>It only assembles the block diagonal when you use MatNest.</div><div><br></div><div> if (!Buk) PetscFunctionReturn(0); /* Not assembling this block */</div><div><br></div>
<div>It builds the whole matrix when you use AIJ, but preallocation isn't correct.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
What about a less general (but important) case: saddle point<br>
problems arising from incompressible Stokes, Oseen and<br>
Navier-Stokes eqs. with Schur type preconditioning. In 2D with N<br>
cells and co-located variables arranged<br>
as (u1,...,uN,v1,...,vN,p1,...pN) the matrix would have the form<br>
[Q G, D 0] with Q a 2N-by-2N matrix, G a 2N-by-N matrix and D a<br>
N-by-2N matrix. Since the variables are co-located, they share<br>
the same partitioning but could have different stencils. How to use<br>
the "split local space", DMComposite and MATNEST in this case?<br></blockquote><div><br></div><div>If you order this way, then you don't need DMComposite or MatNest (although you can still make a MatNest that operates in this ordering, we just don't have a way to make it automatically).</div>
</div>