<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On Mon, Oct 31, 2016 at 8:44 AM, Jed Brown <span dir="ltr"><<a href="mailto:jed@jedbrown.org" target="_blank">jed@jedbrown.org</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><span class="gmail-">Jeremy Theler <<a href="mailto:jeremy@seamplex.com">jeremy@seamplex.com</a>> writes:<br>
<br>
> Hi again<br>
><br>
> I have been wokring on these issues. Long story short: it is about the<br>
> ordering of the unknown fields in the vector.<br>
><br>
> Long story:<br>
> The physics is linear elastic problem, you can see it does work with LU<br>
> over a simple cube (warp the displacements to see it does represent an<br>
> elastic problem, E=200e3, nu=0.3):<br>
><br>
</span>> <a href="https://caeplex.com/demo/results.php?id=5817146bdb561" rel="noreferrer" target="_blank">https://caeplex.com/demo/<wbr>results.php?id=5817146bdb561</a><br>
<span class="gmail-">><br>
><br>
> Say my three displacements (unknowns) are u,v,w. I can define the<br>
> unknown vector as (is this called node-based ordering?)<br>
><br>
> [u1 v1 w1 u2 v2 w2 ... un vn wn]^T<br>
><br>
> Another option is (is this called unknown-based ordering?)<br>
><br>
> [u1 u2 ... un v1 v2 ... vn w1 w2 ... wn]^T<br>
><br>
><br>
> With lu/preonly the results are the same, although the stiffnes matrixes<br>
> for each case are attached as PNGs. And of course, the near-nullspace<br>
> vectors are different. So PCSetCoordinates() should work with one<br>
> ordering and not with another one, an issue I did not take into<br>
> consideration.<br>
><br>
> After understanding Matt's point about the near nullspace (and reading<br>
> some interesting comments from Jed on scicomp stackexchange) I did built<br>
> my own vectors (I had to take a look at MatNullSpaceCreateRigidBody()<br>
> because I found out by running the code the nullspace should be an<br>
> orthonormal basis, it should say so in the docs).<br>
<br>
</span>Where?<br>
<br>
"vecs   - the vectors that span the null space (excluding the constant vector); these vectors must be orthonormal."<br>
<br>
<a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatNullSpaceCreate.html" rel="noreferrer" target="_blank">https://www.mcs.anl.gov/petsc/<wbr>petsc-current/docs/<wbr>manualpages/Mat/<wbr>MatNullSpaceCreate.html</a><br>
<br>
And if you run in debug mode (default), as you always should until you<br>
are confident that your code is correct, MatNullSpaceCreate tests that<br>
your vectors are orthonormal.<br>
<span class="gmail-"><br>
> Now, there are some results I do not understand. I tried these six<br>
> combinations:<br>
><br>
> order      near-nullspace       iterations    norm<br>
> -----      --------------       ----------    ----<br>
> unknown    explicit             10            1.6e-6<br>
> unknown    PCSetCoordinates     15            1.7e-7<br>
> unknown    none                 15            2.4e-7<br>
> node       explicit             fails with error -11<br>
> node       PCSetCoordinates     fails with error -11<br>
> node       none                 13            3.8e-7<br>
<br>
</span>Did you set a block size for the "node-based" orderings?  Are you sure<br>
the above is labeled correctly?  Anyway, PCSetCoordinates uses<br>
"node-based" ordering.  Implementation performance will generally be<br>
better with node-based ordering -- it has better memory streaming and<br>
cache behavior.<br>
<br>
The AIJ matrix format will also automatically do an "inode" optimization<br>
to reduce memory bandwidth and enable block smoothing (default<br>
configuration uses SOR smoothing).  You can use -mat_no_inode to try<br>
turning that off.<br>
<span class="gmail-"><br>
> Error -11 is<br>
> PETSc's linear solver did not converge with reason<br>
> 'DIVERGED_PCSETUP_FAILED' (-11)<br>
<br>
</span>Isn't there an actual error message?<br>
<span class="gmail-"><br>
> Any explanation (for dumbs)?<br>
> Another thing to take into account: I am setting the dirichlet BCs with<br>
> MatZeroRows(), but I am not updating the columns to keep symmetry. Can<br>
> this pose a problem for GAMG?<br>
<br>
</span>Usually minor, but it is better to maintain symmetry.<br></blockquote><div><br></div><div>If the boundary values are not zero, no way to maintain symmetry unless we reduce the extra part of  the matrix. Not <span class="gmail-">updating the columns is better in this situation.<br><br></span></div><div><span class="gmail-">Fande,<br></span></div><div><br> </div></div><br></div></div>