<html><head><meta http-equiv="Content-Type" content="text/html; charset=us-ascii"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Oct 7, 2020, at 8:11 AM, Mark Adams <<a href="mailto:mfadams@lbl.gov" class="">mfadams@lbl.gov</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class=""><div dir="ltr" class=""><div class=""><br class=""></div></div><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, Oct 7, 2020 at 8:27 AM Victoria Hamtiaux <<a href="mailto:victoria.hamtiaux@uclouvain.be" target="_blank" class="">victoria.hamtiaux@uclouvain.be</a>> wrote:<br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div class=""><p class="">Thanks for all the answers,</p><p class=""><br class="">
</p><p class="">How can I do the "semi-coarsening"? I don't really know how those
preconditionners work so I don't how how to change them or so..</p><p class=""><br class=""></p></div></blockquote>You have to write custom code to do semi-coarsening. PETSc does not provide that and you would not want to do it yourself, most likely.</div></div></div></blockquote><div><br class=""></div> We do not provide it directly but if you are using PCMG and DMDA it is relatively straight-forward. You create a coarse DM and then refine it but each refinement you only do in the directions you want set each time with DMDASetRefinementFactor(). Once you have the collections of refined DM's you provide them to PCMG.</div><div><br class=""></div><div> Barry</div><div><br class=""><blockquote type="cite" class=""><div class=""><div dir="ltr" class=""><div class="gmail_quote"><div class=""></div><div class=""> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><div class="">
<br class="webkit-block-placeholder"></div><p class="">I have a question because you both seem to say that my matrix is
supposed to be symmetric which is not the case. \</p></div></blockquote><div class="">You said "my matrix is symmetric." </div><div class=""><br class=""></div><div class="">Then you said " I suspect that by stretching the grid, my matrix is not symmetric anymore and that it might cause a problem."</div><div class=""><br class=""></div><div class="">We are saying that by stretchin the grid the matrix is still symmetric even if the grid has lost a symmetry. I don't know of a mechanism for stretching the grid to make the matrix asymmetric. So we are suggesting that you verify your suspicion that the matrix is symmetric.</div><div class=""><br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><p class="">And in fact, I
don't get how it can be symmetric. Because you will have something
close to symmetric. For example when you are at the center of your
domain it will be symmetric, but when your at a point at the
boundaries I don't get how you can be symmetric, you won't have
something at the left and the right of your main diagonal... (I
don't know if my explanations are understandable)</p></div></blockquote><div class="">You can make a discretization that is not symmetric because of boundary conditions but I assume that is not the case because you said your matrix is symmetric.</div><div class=""> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><p class="">Best regards, <br class="">
</p><p class=""><br class="">
</p><p class="">Victoria<br class="">
</p><p class=""><br class="">
</p><p class=""><br class="">
</p>
<div class="">On 7/10/20 14:20, Mark Adams wrote:<br class="">
</div>
<blockquote type="cite" class="">
<div dir="ltr" class="">GMG (geometric MG) is stronger as Matt said, but it
is affected by stretched grids in a similar way. A way to fix
this in GMG is semi-coarsening, which AMG _can_ do
automatically.<br class="">
</div>
<br class="">
<div class="gmail_quote">
<div dir="ltr" class="gmail_attr">On Wed, Oct 7, 2020 at 8:02 AM
Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank" class="">knepley@gmail.com</a>> wrote:<br class="">
</div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr" class="">
<div dir="ltr" class="">On Wed, Oct 7, 2020 at 7:07 AM Victoria
Hamtiaux <<a href="mailto:victoria.hamtiaux@uclouvain.be" target="_blank" class="">victoria.hamtiaux@uclouvain.be</a>>
wrote:<br class="">
</div>
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div class=""><p class="">Hello Matt, <br class="">
</p><p class=""><br class="">
</p><p class="">I just checked the symmetry of my matrix and it is
not symmetric. But it is not symmetric either when I
use a uniform grid.</p><p class="">The domain is 3D and I'm using finite differences,
so I guess it is normal that at multiple places
(when I deal with points near the boundaries), the
matrix is not symmetric.<br class="">
</p><p class="">So I was wrong, the problem doesn't come from the
fact that the matrix is not symmetric. I don't know
where it comes from, but when I switch from uniform
to stretched grid, the solver stops working
properly. Could it be from the preconditionner of
the solver that I use?<br class="">
</p><p class="">Do you have any other idea ? </p>
</div>
</blockquote>
<div class="">I would consider using GMG. As Mark says, AMG is very
fiddly with stretched grids. For Poisson, GMG works
great and you seem to have regular grids.</div>
<div class=""><br class="">
</div>
<div class=""> Thanks,</div>
<div class=""><br class="">
</div>
<div class=""> Matt </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div class=""><p class="">Thanks for your help, <br class="">
</p><p class=""><br class="">
</p><p class="">Victoria<br class="">
</p><p class=""><br class="">
</p>
<div class="">On 7/10/20 12:48, Matthew Knepley wrote:<br class="">
</div>
<blockquote type="cite" class="">
<div dir="ltr" class="">
<div dir="ltr" class="">On Wed, Oct 7, 2020 at 6:40 AM
Victoria Hamtiaux <<a href="mailto:victoria.hamtiaux@uclouvain.be" target="_blank" class="">victoria.hamtiaux@uclouvain.be</a>>
wrote:<br class="">
</div>
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Dear all,<br class="">
<br class="">
<br class="">
After the discretization of a poisson equation
with purely Neumann (or <br class="">
periodic) boundary conditions, I get a matrix
which is singular.<br class="">
<br class="">
<br class="">
The way I am handling this is by using a
NullSpace with the following <br class="">
code :<br class="">
<br class="">
MatNullSpace nullspace;<br class="">
MatNullSpaceCreate(PETSC_COMM_WORLD,
PETSC_TRUE, 0, 0, &nullspace);<br class="">
MatSetNullSpace(p_solverp->A, nullspace);<br class="">
MatSetTransposeNullSpace(p_solverp->A,
nullspace);<br class="">
MatNullSpaceDestroy(&nullspace);<br class="">
<br class="">
<br class="">
Note that I am using the hypre preconditionner
BOOMERANG and the default <br class="">
solver GMRES.<br class="">
<br class="">
<br class="">
KSPCreate(PETSC_COMM_WORLD,&p_solverp->ksp);<br class="">
KSPSetOperators(p_solverp->ksp,
p_solverp->A, p_solverp->A);<br class="">
PC prec;<br class="">
KSPGetPC(p_solverp->ksp, &prec);<br class="">
PCSetType(prec,PCHYPRE);//PCHYPRE seems
the best<br class="">
PCHYPRESetType(prec,"boomeramg");
//boomeramg is the best<br class="">
KSPSetInitialGuessNonzero(p_solverp->ksp,PETSC_TRUE);<br class="">
KSPSetFromOptions(p_solverp->ksp);<br class="">
KSPSetTolerances(p_solverp->ksp,
1.e-10, 1e-10, PETSC_DEFAULT, <br class="">
PETSC_DEFAULT);<br class="">
KSPSetReusePreconditioner(p_solverp->ksp,PETSC_TRUE);<br class="">
KSPSetUseFischerGuess(p_solverp->ksp,1,5);<br class="">
KSPGMRESSetPreAllocateVectors(p_solverp->ksp);<br class="">
KSPSetUp(p_solverp->ksp);<br class="">
<br class="">
<br class="">
<br class="">
And this works fine when my grid is uniform,
so that my matrix is <br class="">
symmetric.<br class="">
<br class="">
<br class="">
But when I stretch the grid near the boundary
(my grid is then <br class="">
non-uniform), it doesn't work properly
anymore. I suspect that by <br class="">
stretching the grid, my matrix is not
symmetric anymore and that it <br class="">
might cause a problem.<br class="">
</blockquote>
<div class=""><br class="">
</div>
<div class="">Symmetry is a property of the operator, so
you should be symmetric on your</div>
<div class="">stretched grid. If not, I think you have
the discretization wrong. You can check</div>
<div class="">symmetry using</div>
<div class=""><br class="">
</div>
<div class=""> <a href="https://eur03.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.mcs.anl.gov%2Fpetsc%2Fpetsc-current%2Fdocs%2Fmanualpages%2FMat%2FMatIsSymmetric.html&data=02%7C01%7Cvictoria.hamtiaux%40uclouvain.be%7C494a0f05bd214b5974f008d86abb6e02%7C7ab090d4fa2e4ecfbc7c4127b4d582ec%7C0%7C0%7C637376700525925749&sdata=nWt0sejio7o1PzMoc7tPu7JOvcNqofRuMQ91ynW54r4%3D&reserved=0" target="_blank" class="">https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatIsSymmetric.html</a></div>
<div class=""><br class="">
</div>
<div class="">Also, if you suspect your discretization,
you should probably do an MMS test to</div>
<div class="">verify that you discretization converges at
the correct rate.</div>
<div class=""><br class="">
</div>
<div class=""> Thanks,</div>
<div class=""><br class="">
</div>
<div class=""> Matt</div>
<div class=""> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"> I tried
fixing the solution at an arbitrary point, but
sometimes doing <br class="">
this, I get errors near that fixed point. I
've seen on the petsc-users <br class="">
forum that you usually don't recommend to fix
a point, but I don't <br class="">
really know how to proceed differently.<br class="">
<br class="">
<br class="">
What would you recommend to solve this
problem?<br class="">
<br class="">
<br class="">
Thanks for your help,<br class="">
<br class="">
<br class="">
Best regards,<br class="">
<br class="">
<br class="">
Victoria<br class="">
<br class="">
<br class="">
</blockquote>
</div>
<br clear="all" class="">
<div class=""><br class="">
</div>
-- <br class="">
<div dir="ltr" class="">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">What most experimenters take for
granted before they begin their
experiments is infinitely more
interesting than any results to
which their experiments lead.<br class="">
-- Norbert Wiener</div>
<div class=""><br class="">
</div>
<div class=""><a href="https://eur03.safelinks.protection.outlook.com/?url=http:%2F%2Fwww.cse.buffalo.edu%2F~knepley%2F&data=02%7C01%7Cvictoria.hamtiaux%40uclouvain.be%7C494a0f05bd214b5974f008d86abb6e02%7C7ab090d4fa2e4ecfbc7c4127b4d582ec%7C0%7C0%7C637376700525925749&sdata=LHARUv3BxSwWnxN2LUJnX3vr2ZJ9f50EMQzw44Hy%2FqY%3D&reserved=0" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class="">
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</blockquote>
</div>
</blockquote>
</div>
<br clear="all" class="">
<div class=""><br class="">
</div>
-- <br class="">
<div dir="ltr" class="">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">What most experimenters take for granted
before they begin their experiments is
infinitely more interesting than any results
to which their experiments lead.<br class="">
-- Norbert Wiener</div>
<div class=""><br class="">
</div>
<div class=""><a href="https://eur03.safelinks.protection.outlook.com/?url=http:%2F%2Fwww.cse.buffalo.edu%2F~knepley%2F&data=02%7C01%7Cvictoria.hamtiaux%40uclouvain.be%7C494a0f05bd214b5974f008d86abb6e02%7C7ab090d4fa2e4ecfbc7c4127b4d582ec%7C0%7C0%7C637376700525935751&sdata=kOKe%2FLj7pvAdzldpTNlRfC7BS6Vv4S5mU6Cb8pPpmrE%3D&reserved=0" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class="">
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</blockquote>
</div>
</blockquote>
</div>
</blockquote></div>
</div>
</div></blockquote></div><br class=""></body></html>