<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
<font face="Ubuntu">Matt,<br>
<br>
thank you.<br>
The matrix arises from discretization of the Poisson equation in
incompressible flow calculations.<br>
<br>
Michele<br>
<br>
</font>
<div class="moz-cite-prefix">On 08/13/2013 06:34 PM, Matthew Knepley
wrote:<br>
</div>
<blockquote
cite="mid:CAMYG4GmNG1kQffmtVj2PajMfJkGfS8_x115K7e_kq1-vdTv44Q@mail.gmail.com"
type="cite">
<div dir="ltr">On Tue, Aug 13, 2013 at 8:25 PM, Michele Rosso <span
dir="ltr"><<a moz-do-not-send="true"
href="mailto:mrosso@uci.edu" target="_blank">mrosso@uci.edu</a>></span>
wrote:<br>
<div class="gmail_extra">
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0
.8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000"> <font
face="Ubuntu">Matt,<br>
<br>
thank you! I will try to reduce the number of levels
and see how it goes. <br>
I asked about the speed since CG + Block Jacobi with
ICC in each block runs faster then CG + MG, so I
thought I was missing something.<br>
</font></div>
</blockquote>
<div><br>
</div>
<div>What is the operator? If its the mass matrix, we
expect this.</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0
.8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000"><font face="Ubuntu">
Could you please tell me how to get rid of Chebichev?<br>
</font></div>
</blockquote>
<div><br>
</div>
<div>I believe its</div>
<div><br>
</div>
<div> -mg_levels_ksp_type richardson</div>
<div><br>
</div>
<div>but you can check.</div>
<div><br>
</div>
<div> matt</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0
.8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000"><font face="Ubuntu">
Best,<br>
Michele<br>
<br>
</font>
<div>On 08/13/2013 05:51 M, Matthew Knepley wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr">On Tue, Aug 13, 2013 at 7:05 PM,
Michele Rosso <span dir="ltr"><<a
moz-do-not-send="true"
href="mailto:mrosso@uci.edu" target="_blank">mrosso@uci.edu</a>></span>
wrote:<br>
<div class="gmail_extra">
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0
0 0 .8ex;border-left:1px #ccc
solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000"> <font
face="Ubuntu">Hi Matt,<br>
<br>
I attached the output of the commands you
suggested.<br>
The options I used are:<br>
<br>
-log_summary -ksp_monitor -ksp_view
-ksp_converged_reason -pc_type mg
-pc_mg_galerkin -pc_mg_levels 5
-options_left<br>
</font></div>
</blockquote>
<div><br>
</div>
<div>The convergence is great. I notice that
your coarse solve takes no time. You could
probably use fewer levels for</div>
<div>this problem. For this problem there is no
easy things left I think. We are currently
debating how you can squeeze</div>
<div>something extra out of the smoother. Here
you could probably get rid of Chebychev and
use only SOR.</div>
<div><br>
</div>
<div> Matt</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0
0 0 .8ex;border-left:1px #ccc
solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000"><font
face="Ubuntu"> and here are the lines of
codes where I setup the solution process:<br>
<br>
call DMDACreate3d( PETSC_COMM_WORLD
, &<br>
&
DMDA_BOUNDARY_PERIODIC ,
DMDA_BOUNDARY_PERIODIC, &<br>
&
DMDA_BOUNDARY_PERIODIC ,
DMDA_STENCIL_STAR, &<br>
& N_Z , N_Y , N_X
, N_B3 , N_B2 , 1_ip, 1_ip , 1_ip , &<br>
& NNZ ,NNY , NNX,
da , ierr)<br>
<br>
<br>
! Create Global Vectors <br>
call DMCreateGlobalVector(da,b,ierr)<br>
call VecDuplicate(b,x,ierr)<br>
<br>
! Set initial guess for first use of
the module to 0<br>
call VecSet(x,0.0_rp,ierr) <br>
<br>
! Create matrix <br>
call DMCreateMatrix(da,MATAIJ,A,ierr)<br>
<br>
! Create solver<br>
call
KSPCreate(PETSC_COMM_WORLD,ksp,ierr) <br>
call KSPSetDM(ksp,da,ierr)<br>
call
KSPSetDMActive(ksp,PETSC_FALSE,ierr)<br>
call
KSPSetOperators(ksp,A,A,SAME_NONZERO_PATTERN,ierr)<br>
call KSPSetType(ksp,KSPCG,ierr)<br>
call
KSPSetNormType(ksp,KSP_NORM_UNPRECONDITIONED,ierr)<br>
call
KSPSetInitialGuessNonzero(ksp,PETSC_TRUE,ierr)<br>
call KSPSetTolerances(ksp, tol
,PETSC_DEFAULT_DOUBLE_PRECISION,&<br>
&
PETSC_DEFAULT_DOUBLE_PRECISION,PETSC_DEFAULT_INTEGER,ierr)
<br>
<br>
! Nullspace removal <br>
call MatNullSpaceCreate(
PETSC_COMM_WORLD,PETSC_TRUE,PETSC_NULL_INTEGER,&<br>
&
PETSC_NULL_INTEGER,nullspace,ierr)<br>
call
KSPSetNullspace(ksp,nullspace,ierr)<br>
call
MatNullSpaceDestroy(nullspace,ierr) <br>
<br>
! To allow using option from command
line<br>
call KSPSetFromOptions(ksp,ierr)
<br>
<br>
<br>
Hope I did not omit anything useful.<br>
Thank you for your time.<br>
<br>
Best,<br>
Michele<br>
<br>
<br>
<br>
<br>
</font>
<div>On 08/13/2013 04:26 PM, Matthew Knepley
wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr">On Tue, Aug 13, 2013 at
6:09 PM, Michele Rosso <span dir="ltr"><<a
moz-do-not-send="true"
href="mailto:mrosso@uci.edu"
target="_blank">mrosso@uci.edu</a>></span>
wrote:<br>
<div class="gmail_extra">
<div class="gmail_quote">
<blockquote class="gmail_quote"
style="margin:0 0 0
.8ex;border-left:1px #ccc
solid;padding-left:1ex">
<div bgcolor="#FFFFFF"
text="#000000"> <font
face="Ubuntu">Hi Karli,<br>
<br>
thank you for your hint: now
it works. <br>
Now I would like to speed up
the solution: I was counting
on increasing the number of
levels/the number of
processors used, but now I see
I cannot do that.<br>
Do you have any hint to
achieve better speed?<br>
Thanks!<br>
</font></div>
</blockquote>
<div><br>
</div>
<div>"Better speed" is not very
helpful for us, and thus we cannot
offer much help. You could</div>
<div><br>
</div>
<div> 1) Send the output of
-log_summary -ksp_monitor
-ksp_view</div>
<div><br>
</div>
<div> 2) Describe the operator
succintly</div>
<div><br>
</div>
<div> Matt</div>
<div> </div>
<blockquote class="gmail_quote"
style="margin:0 0 0
.8ex;border-left:1px #ccc
solid;padding-left:1ex">
<div bgcolor="#FFFFFF"
text="#000000"> <font
face="Ubuntu"> Best,<br>
Michele<br>
</font>
<blockquote type="cite">
<blockquote type="cite">
<blockquote type="cite">
<blockquote type="cite">
<blockquote type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite">
<blockquote
type="cite"><br>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</blockquote>
</div>
</blockquote>
</div>
</div>
</div>
</blockquote>
<br>
</div>
</blockquote>
</div>
<br>
<br clear="all">
<span class="HOEnZb"><font color="#888888">
<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 </font></span></div>
</div>
</blockquote>
<br>
</div>
</blockquote>
</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>
</blockquote>
<br>
</body>
</html>