On Thu, Jun 14, 2012 at 10:07 PM, Thomas Witkowski <span dir="ltr"><<a href="mailto:thomas.witkowski@tu-dresden.de" target="_blank">thomas.witkowski@tu-dresden.de</a>></span> wrote:<br><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">
Am 14.06.2012 15:06, schrieb Matthew Knepley:
<blockquote type="cite">On Thu, Jun 14, 2012 at 9:01 PM, Thomas Witkowski <span dir="ltr"><<a href="mailto:thomas.witkowski@tu-dresden.de" target="_blank">thomas.witkowski@tu-dresden.de</a>></span>
wrote:<br>
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
I try to implement an efficient solver for my FEM based
unsteady Navier Stokes code. The scenario I consider is realy
simple: 2D FEM with Taylor-Hood element, equi-spaced grid,
simple channel flow with prescribed inflow boundaries for the
velocity, no boundary conditions for the pressure. Because I
want to</blockquote>
<div><br>
</div>
<div>With no pressure BC, you need to specify the null space.
Look at SNES ex62.</div>
</div>
</blockquote>
As so often, I forgotten about this :) As I have free outflow
boundary for the velocity (so no Dirichlet boundary), I think, the
nullspace must be somehow different than just setting a constant on
all pressure nodes. So I want back to the even more simple case, the
driven cavity scenario. Here I have dirichlet BC on all velocity
components and boundaries and no pressure BC. There is no influence
on the solver, whether I set the null space or not (I checked for it
with -ksp_view that there is a null space set). Only when going from
fgmres to gmres for the outer solver, setting the null space has
some (positive) change in the solver convergence. </div></blockquote><div><br></div><div>It has a null space, so whatever you observe in a single run is an outlier. Try removing the null space in ex62. It will fail to converge.</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"> <br><blockquote type="cite"><div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">solve
unsteady flows, at each timepoint an Ossen problem is solved
(implicit Euler time discretization, "trivial" linearization).
For using PETSc, I created a fieldsplit preconditioner, that
is configured in the following way:<br>
<br>
mpirun -np 2 ./navier_stokes init/channel_parallel.dat.2d \<br>
-ksp_type fgmres \<br>
-pc_type fieldsplit \<br>
-pc_fieldsplit_type SCHUR \<br>
-pc_fieldsplit_schur_fact_type FULL \<br>
-ksp_converged_reason \<br>
-ksp_monitor_true_residual \<br>
-fieldsplit_velocity_ksp_type preonly \<br>
-fieldsplit_velocity_pc_type lu \<br>
-fieldsplit_velocity_pc_factor_mat_solver_package mumps \<br>
-fieldsplit_pressure_ksp_type gmres \<br>
-fieldsplit_pressure_ksp_max_it 10 \<br>
-fieldsplit_pressure_ksp_rtol 1.0e-2 \<br>
-fieldsplit_pressure_pc_type lsc \<br>
</blockquote>
<div><br>
</div>
<div>This makes no sense unless you specify an auxiliary
operator. Just leave it at jacobi. When you use LU</div>
<div>for velocity, it will converge in 1 iterate. Since it
doesn't, it means you have a problem with the null space.</div>
</div>
</blockquote>
So even for the driven cavity, the null space is set, and again I
use the same solver, it needs something like 10 or 15 iterations. I
don't understand your argument, that using LU for the velocity part,
this solver should converge in 1 iterations. It still has an inexact
solver for the Schur complement and the solver for the composite
matrix in LSC is also inexact.<br></div></blockquote><div><br></div><div>1) Start with -fieldsplit_pressure_ksp_rtol 1e-10. It will converge in 1 iterate. Back off the Schur tolerance</div><div> to balance that cost with the cost of extra iterates. This is the game.</div>
<div><br></div><div>2) Again, LSC is not doing anything for you right now, unless you are indeed giving an auxiliary operator. Use jacobi.</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">
Thomas<br>
<br>
<blockquote type="cite">
<div class="gmail_quote">
<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">
-fieldsplit_pressure_lsc_ksp_type bcgs \<br>
-fieldsplit_pressure_lsc_ksp_max_it 10 \<br>
-fieldsplit_pressure_lsc_ksp_rtol 1.0e-2 \<br>
-fieldsplit_pressure_lsc_pc_type hypre<br>
<br>
The direct solver for the velocity part is just for debugging
and should be replaced when everything else works fine. The
point is, that I found this constellation to be not very
efficient. It takes around 20 to 30 iterations, which takes
around 30 seconds on a very small problem size (around 20000
global unknows for each velocity component and 5000 global
unknowns for the pressure) on a very fast desktop CPU (some
new Intel Xeno with 4 core). Any hints for improvements?<span><font color="#888888"><br>
<br>
Thomas<br>
</font></span></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<br>
</font></span></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<br>