<br><br><div class="gmail_quote">On Sat, May 26, 2012 at 9:07 PM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

<div class="gmail_quote"><div class="im">On Sat, May 26, 2012 at 8:04 PM, Dmitry Karpeev <span dir="ltr"><<a href="mailto:karpeev@mcs.anl.gov" target="_blank">karpeev@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">


<div><br><br>On Saturday, May 26, 2012, Jed Brown <<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>> wrote:<br>> Dmitry, I have problems with this commit.<br>> <a href="http://petsc.cs.iit.edu/petsc/petsc-dev/rev/440ec9350d51" target="_blank">http://petsc.cs.iit.edu/petsc/petsc-dev/rev/440ec9350d51</a><br>



><br>> If I try running -pc_type asm with a vector problem defined on a DM, PCASM calls DMCreateDecomposition which splits fields on each process. This is of very limited value since usually every field is coupled to every other, so after increasing the overlap by one, you end up solving the entire subdomain twice. However, at this point, you are not using overlap so instead, we get subdomain and field blocked Jacobi. (Equivalent to BJacobi with additive fieldsplit inside.) This is certainly unexpected behavior.<br>



> It also crashes with Seq(S)BAIJ because MatGetSubMatrices_SeqBAIJ can only extract when the block size matches.<br>> And it does the wrong thing when used with MPI(S)BAIJ because MatGetSubMatrices_MPIBAIJ calls ISCompressIndicesGeneral which silently "fills in" the blocks, leading to a Mat and Vec size mismatch. Hong, how should we enforce that the original index set actually addresses blocks? Why is it not okay to sort? (Don't we rely on those indices being sorted somewhere else anyway?)<br>



<br></div>The problem is that this API mixes field splits and (spatial) domain decompositions.<br>Initially I thought the two kinds could be served up via the same API because they use the same syntax, but the semantics are different. In particular,  domain decompositions don't change the block size, while field splits typically do. Then it's usually downwards -- when splitting up a block of fields pointwise, but sometimes upwards -- in a primal-dual split of a KKT system the primal subsystem might have a higher natural block size than the full system.<br>



<br>I think the solution is to split up the API: DMCreateFieldSplit() and DMCreateDomainDecomposition(). I'll look into it later tonight, when I get to my computer.</blockquote><div><br></div></div><div>If the arguments are definitely going to be the same, there could be an enum to designate which form is being requested.</div>

</div></blockquote><div>I think that might make it more cumbersome for impls to provide only DD or only FS decompositions. </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

<div class="gmail_quote">
<div><br></div><div><br></div><div>What are we going to do now? Leave the code there deactivated like it is now for the release and discuss later?</div></div></blockquote><div>I want to fix and reactivate  the code for the release so I can use it from libMesh.</div>

<div>I'll push roughly this fix by tomorrow: DMCreateDecomposition(DM)()--> DMCreateFieldSplit(DM)()/DMCreateDomainDecomposition(DM)() used, respectively, by PCFIELDSPLIT/PC(G)ASM.</div><div>That should address the problems you were encountering -- unintended forwarding of field splits to PCASM.</div>

<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="gmail_quote"><div class="im"><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">


<div><br><br>> Our test coverage must be pretty sucky if we went two months with ASM this broken for systems.<br>
><br>> As a lower priority, -pc_asm_decomposition should provide a list of valid decompositions rather than be a magic string. Also, this form cannot be called outside of a PetscOptionsBegin/End pair, otherwise it doesn't print -help output. It should be moved to PCSetUp_ASM.<br>



> ierr = PetscOptionsString("-pc_asm_decomposition", "Name of the DM defining the composition", "PCSetDM", ddm_name, ddm_name,1024,&flg); CHKERRQ(ierr);<br><br></div>The problem here is that the PC object doesn't really know what decompositions are available -- the DM does. We could have DMGetSubdomainDecompositionNames() advertise what'd available. Another problem is that what decompositions make sense is DM-impl-dependent.  It is more stark for fieldsplits: I have a fork of libMesh that can generate decompositions based on the system's variables, and essentially any groupings of variable names are valid splits ( e.g., v1,v2;v3,v4. I think it makes sense to maintain that flexibility).</blockquote>


<div><br></div></div><div>Cool, what about deeper nestings like [v1,[[v2,v3],v4]].</div></div></blockquote><div>Nested splits are to be handled recursively by the individual splits' DMs.  And each split should also be domain-decomposable, so that geometrically-defined subdomains (e.g., Exodus II blocks) can be used with PC(G)ASM.</div>

<div>I am testing that functionality in libMesh/Moose right now. </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="gmail_quote"><div class="im"><div> </div>

<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"> Because of that I think the best help line would only exhibit the syntax, rather than enumerate all possible splits.</blockquote>


<div><br></div></div><div>How about having a DM-specific function for that designation?</div></div></blockquote><div>DMXXXGetFieldSplits()/DMXXXGetDomainDecompositions(), for example?</div><div>I think that would work, if the implementation defined only a small number of decompositions.</div>

<div>If instead there is a short grammar defining well-formed strings naming such decompositions, I'm not sure how to </div><div>communicate that to the user, other than in documentation.  I have DMLibMeshGetVars() to return all vars.  FieldSplits</div>

<div>can be defined as arbitrary groupings of those.  Each split will then advertise only the vars in the groupings that defines the split.</div><div><br></div><div>Dmitry.</div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

<div class="gmail_quote"><div class="im"><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"> Then DMSetUp_Xxx should handle that.<span><font color="#888888"><br>



<br>Dmitry.</font></span><div><div><br>><br>><br>> I have disabled that code for now, but we have to decide whether it's conceptually right (considering the effect on overlap) and if so, how to fix it.<br>
> <a href="http://petsc.cs.iit.edu/petsc/petsc-dev/rev/d9a0e3f9dacd" target="_blank">http://petsc.cs.iit.edu/petsc/petsc-dev/rev/d9a0e3f9dacd</a>
</div></div></blockquote></div></div><br>
</blockquote></div><br>