<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
<div class="moz-cite-prefix">On 7/29/13 3:42 PM, Bishesh Khanal
wrote:<br>
</div>
<blockquote
cite="mid:CAEhex8juYKy719p90y3qkjdPejA-nwibSDXpPnBhM=Lq5PgEzQ@mail.gmail.com"
type="cite">
<meta http-equiv="Content-Type" content="text/html;
charset=ISO-8859-1">
<div dir="ltr"><br>
<div class="gmail_extra"><br>
<br>
<div class="gmail_quote">On Fri, Jul 26, 2013 at 6:05 PM, Dave
May <span dir="ltr"><<a moz-do-not-send="true"
href="mailto:dave.mayhem23@gmail.com" target="_blank">dave.mayhem23@gmail.com</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">
<div dir="ltr">
<div>Yes, the nullspace is important.<br>
</div>
<div>There is definitely a pressure null space of
constants which needs to be considered.<br>
</div>
<div>I don't believe ONLY using<br>
-fieldsplit_1_ksp_constant_null_space<br>
</div>
<div>is not sufficient. <br>
<br>
You need to inform the KSP for the outer syetem (the
coupled u-p system) that there is a null space of
constants in the pressure system. This cannot (as far
as I'm aware) be set via command line args. You need
to write a null space removal function which accepts
the complete (u,p) vector which only modifies the
pressure dofs.<br>
<br>
</div>
<div>Take care though when you write this though.
Because of the DMDA formulation you are using, you
have introduced dummy/ficticious pressure dofs on the
right/top/front faces of your mesh. Thus your null
space removal function should ignore those dofs when
your define the constant pressure constraint.<br>
<br>
</div>
</div>
</blockquote>
<div><br>
</div>
<div>I'm aware of the constant pressure null space but did
not realize that just having
-fieldsplit_1_ksp_constant_null_space would work, thanks.
Now I tried doing what you said but I get some problems!<br>
I created a null basis vector say, mNullBasis (creating a
global vector using the DMDA formulation). I set all it's
values to zero except the non-dummy pressure dofs whose
values are set to one. Now I use MatNullSpaceCreate to
create a MatNullSpace say mNullSpace using the mNullBasis
vector.<br>
</div>
<div><br>
The relevant lines in my solver function are (mKsp and mDa
are my relevant KSP and DM objects):<br>
<br>
ierr =
DMKSPSetComputeRHS(mDa,computeRHSTaras3D,this);CHKERRQ(ierr);
<br>
ierr =
DMKSPSetComputeOperators(mDa,computeMatrixTaras3D,this);CHKERRQ(ierr);<br>
ierr =
KSPSetDM(mKsp,mDa);CHKERRQ(ierr);
<br>
ierr = KSPSetNullSpace(mKsp,mNullSpace);CHKERRQ(ierr);<br>
ierr = KSPSetFromOptions(mKsp);CHKERRQ(ierr);<br>
ierr = KSPSetUp(mKsp);CHKERRQ(ierr);<br>
ierr = KSPSolve(mKsp,NULL,NULL);CHKERRQ(ierr);<br>
<br>
</div>
<div>And a corresponding addition in the function
computeMatrixTaras3D is something equivalent to:<br>
ierr =
MatNullSpaceRemove(mNullSpace,b,NULL);CHKERRQ(ierr);
//where b is the vector that got updated with the rhs
values before this line. <br>
<br>
</div>
<div>When using the -ksp_view in the option I see that a
null space now gets attached to the outer solver, while I
still need to keep -fieldsplit_1_ksp_constant_null_space
option as it is to have a null space attached to the inner
pressure block solver. However the problem is I do not get
the expected result with gamg and resulting solution blows
up.<br>
</div>
<div>Am I doing something wrong or completely missing some
point here with nullspace removal idea ? <br>
</div>
<div>(Note that I am imposing a Dirichlet boundary
condition, so in my case I think the only null space is
indeed the constant pressure and not the translation and
rigid body rotation which would be present if using for
e.g. Neumann boundary condition)<br>
<br>
</div>
</div>
</div>
</div>
</blockquote>
<meta charset="utf-8">
Bishesh,<br>
<br>
yes, pressure is the only filed that is defined up to a constant in
the incompressible case, if you properly set the velocity Dirichlet
BC.<br>
<br>
What I was saying is that you need to know the full null space of
your linear system matrix if you want to setup an efficient
algebraic multigrid preconditioner. (Near) null space is a starting
point for the tentative prolongator. I don't know how it looks like
for staggered grid, and don't know how to make PETSc use this
information (the latter should be easy I suspect). Geometric
multigrid for staggered grid is also not a straightforward option,
because you need to define and pass to PETSc a custom
interpolation-restriction operator with PCMGSetRestriction or
PCMGSetInterpolation
<meta charset="utf-8">
. Without efficient multigrid you code will be severely limited for
large-scale 3D problems (especially with variable viscosity, if this
is your goal).<br>
<br>
Which problem with variable viscosity are you solving? To my opinion
the easiest way to remove the pressure null space in your case is to
add minor compressibility. In many cases you can easily get rid of
strict incompressibility assumption, because everything is
compressible (at least slightly).<br>
<br>
If you still want to project out constant pressure, then Dave is
right, -fieldsplit_1_ksp_constant_null_space will most likely not
work, because of dummy DOFs. Another complication is that you should
define the null space only for the pressure. You probably should
retrieve KSP for the pressure Schur complement using
PCFieldSplitGetSubKSP after setup of fieldsplit preconditioner, and
try to set your custom null space directly for it. I'm not sure
however if this will work out. We have a version of staggered grid
solver without dummy pressure DOFs. For us it works both with and
without ..._ksp_constant_null_space for pressure.<br>
<br>
Anton<br>
<blockquote
cite="mid:CAEhex8juYKy719p90y3qkjdPejA-nwibSDXpPnBhM=Lq5PgEzQ@mail.gmail.com"
type="cite">
<div dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">
<div><br>
<br>
</div>
<div><br>
<br>
</div>
<div><br>
</div>
<div> <br>
<br>
</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">
<div>
<br>
</div>
<div><br>
Cheers,<br>
</div>
Dave<br>
</div>
<div class="">
<div class="h5">
<div class="gmail_extra"><br>
<br>
<div class="gmail_quote">On 26 July 2013 17:42,
Matthew Knepley <span dir="ltr"><<a
moz-do-not-send="true"
href="mailto:knepley@gmail.com"
target="_blank">knepley@gmail.com</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">
<div dir="ltr">
<div>
<div>On Fri, Jul 26, 2013 at 10:13 AM,
Bishesh Khanal <span dir="ltr"><<a
moz-do-not-send="true"
href="mailto:bisheshkh@gmail.com"
target="_blank">bisheshkh@gmail.com</a>></span>
wrote:<br>
</div>
</div>
<div class="gmail_extra">
<div class="gmail_quote">
<div>
<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"><br>
<div class="gmail_extra"><br>
<br>
<div class="gmail_quote">
On Fri, Jul 26, 2013 at 4:22
PM, Matthew Knepley <span
dir="ltr"><<a
moz-do-not-send="true"
href="mailto:knepley@gmail.com"
target="_blank">knepley@gmail.com</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">
<div dir="ltr">
<div>
<div>On Fri, Jul 26,
2013 at 9:11 AM,
Bishesh Khanal <span
dir="ltr"><<a
moz-do-not-send="true"
href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span>
wrote:<br>
</div>
</div>
<div class="gmail_extra">
<div class="gmail_quote">
<div>
<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"><br>
<div
class="gmail_extra"><br>
<br>
<div
class="gmail_quote">
On Fri, Jul
26, 2013 at
2:32 PM,
Matthew
Knepley <span
dir="ltr"><<a
moz-do-not-send="true" href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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">
<div dir="ltr">
<div>On Fri,
Jul 26, 2013
at 7:28 AM,
Bishesh Khanal
<span
dir="ltr"><<a
moz-do-not-send="true" href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span>
wrote:<br>
</div>
<div
class="gmail_extra">
<div
class="gmail_quote">
<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"><br>
<div
class="gmail_extra"><br>
<br>
<div
class="gmail_quote">
On Wed, Jul
17, 2013 at
9:48 PM, Jed
Brown <span
dir="ltr"><<a
moz-do-not-send="true" href="mailto:jedbrown@mcs.anl.gov"
target="_blank">jedbrown@mcs.anl.gov</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">
<div>Bishesh
Khanal <<a
moz-do-not-send="true" href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>>
writes:<br>
<br>
> Now, I
implemented
two different
approaches,
each for both
2D and 3D, in<br>
> MATLAB.
It works for
the smaller
sizes but I
have problems
solving it for<br>
> the
problem size I
need (250^3
grid size).<br>
> I use
staggered grid
with p on cell
centers, and
components of
v on cell<br>
> faces.
Similar split
up of K to
cell center
and faces to
account for
the<br>
> variable
viscosity
case)<br>
<br>
</div>
Okay, you're
using a
staggered-grid
finite
difference
discretization
of<br>
variable-viscosity
Stokes. This
is a common
problem and I
recommend<br>
starting with
PCFieldSplit
with Schur
complement
reduction
(make that<br>
work first,
then switch to
block
preconditioner).
</blockquote>
<div><br>
</div>
<div>Ok, I
made my 3D
problem work
with
PCFieldSplit
with Schur
complement
reduction
using the
options:<br>
-pc_fieldsplit_type
schur
-pc_fieldsplit_detect_saddle_point
-fieldsplit_1_ksp_constant_null_space<br>
</div>
<div> <br>
<br>
</div>
<blockquote
class="gmail_quote"
style="margin:0px
0px 0px
0.8ex;border-left:1px
solid
rgb(204,204,204);padding-left:1ex">You
can use PCLSC
or<br>
(probably
better for
you), assemble
a
preconditioning
matrix
containing<br>
the inverse
viscosity in
the
pressure-pressure
block. This
diagonal<br>
matrix is a
spectrally
equivalent (or
nearly so,
depending on<br>
discretization)
approximation
of the Schur
complement.
The velocity<br>
block can be
solved with
algebraic
multigrid.
Read the
PCFieldSplit<br>
docs (follow
papers as
appropriate)
and let us
know if you
get stuck.<br>
</blockquote>
</div>
<br>
</div>
<div
class="gmail_extra">Now,
I got a little
confused in
how exactly to
use command
line options
to use
multigrid for
the velocity
bock and PCLS
for the
pressure
block. After
going through
the manual I
tried the
following:<br>
</div>
</div>
</blockquote>
<div> </div>
</div>
<div>You want
Algebraic
Multigrid</div>
<div><br>
</div>
<div>-pc_type
fieldsplit
-pc_fieldsplit_detect_saddle_point</div>
<div>-pc_fieldsplit_type
schur</div>
<div>
-fieldsplit_0_pc_type
gamg</div>
<div>
<div>
-fieldsplit_1_ksp_type
fgmres</div>
<div>
-fieldsplit_1_ksp_constant_null_space</div>
<div>
-fieldsplit_1_ksp_monitor_short</div>
<div>
-fieldsplit_1_pc_type
lsc</div>
<div>-ksp_converged_reason<br>
</div>
<div><br>
</div>
</div>
</div>
</div>
</div>
</blockquote>
<div>I tried
the above set
of options but
the solution I
get seem to be
not correct.
The velocity
field I get
are quite
different than
the one I got
before without
using gamg
which were the
expected one.<br>
Note: (Also, I
had to add one
extra option
of
-fieldsplit_1_ksp_gmres_restart
100 , because
the
fieldsplit_1_ksp
residual norm
did not
converge
within default
30 iterations
before
restarting).</div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
</div>
</div>
<div>These are all
iterative solvers.
You have to make
sure everything
converges. </div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
<div>When I set restart to
100, and do -ksp_monitor, it
does converge (for the
fieldsplit_1_ksp). Are you
saying that in spite of
having -ksp_converged_reason
in the option and petsc
completing the run with the
message "Linear solve
converged due to
CONVERGED_RTOL .." not
enough to make sure that
everything is converging ?
If that is the case what
should I do for this
particular problem ?<br>
</div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
</div>
</div>
<div>If your outer iteration converges,
and you do not like the solution, there
are usually two possibilities:</div>
<div><br>
</div>
<div> 1) Your tolerance is too high,
start with it cranked down all the way
(1e-10) and slowly relax it</div>
<div><br>
</div>
<div> 2) You have a null space that you
are not accounting for</div>
<div>
<div><br>
</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">
<div class="gmail_extra">
<div class="gmail_quote">
<div>I have used the MAC scheme
with indexing as shown in:
fig 7.5, page 96 of: <a
moz-do-not-send="true"
href="http://books.google.co.uk/books?id=W83gxp165SkC&printsec=frontcover&dq=Introduction+to+Numerical+Geodynamic+Modelling&hl=en&sa=X&ei=v6TmUaP_L4PuOs3agJgE&ved=0CDIQ6AEwAA"
target="_blank">http://books.google.co.uk/books?id=W83gxp165SkC&printsec=frontcover&dq=Introduction+to+Numerical+Geodynamic+Modelling&hl=en&sa=X&ei=v6TmUaP_L4PuOs3agJgE&ved=0CDIQ6AEwAA</a><br>
<br>
Thus I have a DM with 4 dof
but there are several "ghost
values" set to 0.<br>
</div>
<div>Would this cause any
problem when using the
multigrid ? (This has worked
fine when not using the
multigrid.)</div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
</div>
<div>I don't know exactly how you have
implemented this. These should be rows
of the identity.</div>
<div>
<div> </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">
<div class="gmail_extra">
<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 dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">
<div>I do this problem
in the tutorial with</div>
<div>a constant
viscosity.</div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
<div>Which tutorial are you
referring to ? Could you
please provide me the link
please ? </div>
</div>
</div>
</div>
</blockquote>
<div><br>
</div>
</div>
<div>
There are a few on the PETSc Tutorials
page, but you can look at this</div>
<div><br>
</div>
<div> <a moz-do-not-send="true"
href="http://www.geodynamics.org/cig/community/workinggroups/short/workshops/cdm2013/presentations/SessionIV_Solvers.pdf"
target="_blank">http://www.geodynamics.org/cig/community/workinggroups/short/workshops/cdm2013/presentations/SessionIV_Solvers.pdf</a></div>
<div><br>
</div>
<div>for a step-by-step example of a
saddle-point problem at the end.</div>
<div><br>
</div>
<div> Matt</div>
<div>
<div>
<div> </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">
<div class="gmail_extra">
<div class="gmail_quote"><br>
<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">
<div class="gmail_extra">
<div class="gmail_quote">
<div><br>
</div>
<div> Matt</div>
<div>
<div>
<div> </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">
<div
class="gmail_extra">
<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 dir="ltr">
<div
class="gmail_extra">
<div
class="gmail_quote">
<div>
<div>
</div>
</div>
<div> Matt</div>
<div>
<div>
<div><br>
</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">
<div
class="gmail_extra">
<br>
</div>
<div
class="gmail_extra">but
I get the
following
errror:<br>
<br>
[0]PETSC
ERROR:
---------------------
Error Message
------------------------------------<br>
[0]PETSC
ERROR: Null
argument, when
expecting
valid pointer!<br>
[0]PETSC
ERROR: Null
Object:
Parameter # 2!<br>
[0]PETSC
ERROR:
------------------------------------------------------------------------<br>
[0]PETSC
ERROR: Petsc
Release
Version 3.4.2,
Jul, 02, 2013
<br>
[0]PETSC
ERROR: See
docs/changes/index.html
for recent
updates.<br>
[0]PETSC
ERROR: See
docs/faq.html
for hints
about trouble
shooting.<br>
[0]PETSC
ERROR: See
docs/index.html
for manual
pages.<br>
[0]PETSC
ERROR:
------------------------------------------------------------------------<br>
[0]PETSC
ERROR:
src/AdLemMain
on a
arch-linux2-cxx-debug
named edwards
by bkhanal Fri
Jul 26
14:23:40 2013<br>
[0]PETSC
ERROR:
Libraries
linked from
/home/bkhanal/Documents/softwares/petsc-3.4.2/arch-linux2-cxx-debug/lib<br>
[0]PETSC
ERROR:
Configure run
at Fri Jul 19
14:25:01 2013<br>
[0]PETSC
ERROR:
Configure
options
--with-cc=gcc
--with-fc=g77
--with-cxx=g++
--download-f-blas-lapack=1
--download-mpich=1
-with-clanguage=cxx
--download-hypre=1<br>
[0]PETSC
ERROR:
------------------------------------------------------------------------<br>
[0]PETSC
ERROR:
MatPtAP() line
8166 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/mat/interface/matrix.c<br>
[0]PETSC
ERROR:
PCSetUp_MG()
line 628 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/mg/mg.c<br>
[0]PETSC
ERROR:
PCSetUp() line
890 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/interface/precon.c<br>
[0]PETSC
ERROR:
KSPSetUp()
line 278 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c<br>
[0]PETSC
ERROR:
KSPSolve()
line 399 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c<br>
[0]PETSC
ERROR:
PCApply_FieldSplit_Schur()
line 807 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>
[0]PETSC
ERROR:
PCApply() line
442 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/interface/precon.c<br>
[0]PETSC
ERROR:
KSP_PCApply()
line 227 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/include/petsc-private/kspimpl.h<br>
[0]PETSC
ERROR:
KSPInitialResidual()
line 64 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itres.c<br>
[0]PETSC
ERROR:
KSPSolve_GMRES()
line 239 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/impls/gmres/gmres.c<br>
[0]PETSC
ERROR:
KSPSolve()
line 441 in
/home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c<br>
<br>
</div>
<div
class="gmail_extra"><br>
</div>
</div>
</blockquote>
</div>
</div>
</div>
<br>
<br
clear="all">
<span><font
color="#888888"><span><font
color="#888888">
<div>
<div><br>
</div>
-- <br>
What most
experimenters
take for
granted before
they begin
their
experiments is
infinitely
more
interesting
than any
results to
which their
experiments
lead.<br>
-- Norbert
Wiener
</div>
</font></span></font></span></div>
</div>
<span><font
color="#888888">
</font></span></blockquote>
</div>
<span><font
color="#888888"><br>
</font></span></div>
</div>
<span><font
color="#888888">
</font></span></blockquote>
</div>
</div>
</div>
<span><font
color="#888888">
<div>
<div><br>
<br clear="all">
<div><br>
</div>
-- <br>
What most
experimenters
take for granted
before they
begin their
experiments is
infinitely more
interesting than
any results to
which their
experiments
lead.<br>
-- Norbert
Wiener
</div>
</div>
</font></span></div>
</div>
</blockquote>
</div>
<br>
</div>
</div>
</blockquote>
</div>
</div>
</div>
<div>
<div><br>
<br clear="all">
<div><br>
</div>
-- <br>
What most experimenters take for granted
before they begin their experiments is
infinitely more interesting than any
results to which their experiments lead.<br>
-- Norbert Wiener
</div>
</div>
</div>
</div>
</blockquote>
</div>
<br>
</div>
</div>
</div>
</blockquote>
</div>
<br>
</div>
</div>
</blockquote>
<br>
</body>
</html>