<html>
<head>
<meta http-equiv="content-type" content="text/html; charset=UTF-8">
</head>
<body text="#000000" bgcolor="#FFFFFF">
Dear pdelab users,
<div class="moz-forward-container">
<p>I am currently experiencing a rather strange problem during
parallel solution of my finite volume code. I have written a
short outline of my code below for reference. <br>
</p>
<p>At some point during computation, if dune throws an error, the
code catches this error, resets the solution vector to the old
value, halves the time step size, and tries to redo the
calculation (osm.apply()). <br>
</p>
<p>However, if I get the error "FMatrixError: matrix is singular",
the solver seems to freeze. Even the initial defect is not
shown! (See the terminal output below.) I am not sure why this
is so, and I have not experienced this issue before.</p>
<p>I will be very thankful if someone can help me figure out a way
around this problem.</p>
<p>Thanks, and warm wishes, Shubhangi</p>
<p><br>
</p>
<font size="-1"><b>// code layout</b><br>
</font><br>
<font size="-1"> ...UG grid, generated using gmsh, GV, ...</font><br>
<br>
<font size="-1"> typedef
Dune::PDELab::QkDGLocalFiniteElementMap<GV::Grid::ctype,
double, 0, dim, Dune::PDELab::QkDGBasisPolynomial::lagrange>
FEMP0;</font><br>
<font size="-1"> FEMP0 femp0;</font><br>
<font size="-1"> typedef
Dune::PDELab::GridFunctionSpace<GV,FEMP0,Dune::PDELab::P0ParallelConstraints,Dune::PDELab::ISTL::VectorBackend<>>
GFS0;</font><br>
<font size="-1"> GFS0 gfs0(gv,femp0);</font><br>
<font size="-1"> typedef
Dune::PDELab::PowerGridFunctionSpace< GFS0,num_of_vars,
Dune::PDELab::ISTL::VectorBackend<Dune::PDELab::ISTL::Blocking::fixed>,
Dune::PDELab::EntityBlockedOrderingTag> GFS_TCH;</font><font
size="-1"><br>
</font><br>
<font size="-1"> ... LocalOperator LOP lop, TimeLocalOperator
TOP top, GridOperator GO go, InstationaryGridOperator IGO igo,
... </font><font size="-1"><br>
</font><br>
<font size="-1"> typedef
Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR<IGO> LS;</font><br>
<font size="-1"> LS ls(gfs,50,1,false,true);</font><br>
<font size="-1"> typedef Dune::PDELab::Newton< IGO, LS, U
> PDESOLVER;</font><br>
<font size="-1"> PDESOLVER pdesolver( igo, ls );</font><br>
<font size="-1">
Dune::PDELab::ImplicitEulerParameter<double> method;</font><br>
<br>
<font size="-1"> Dune::PDELab::OneStepMethod< double, IGO,
PDESOLVER, U, U > osm( method, igo, pdesolver );</font><font
size="-1"><br>
</font><br>
<font size="-1"> //TIME-LOOP</font><br>
<font size="-1"> while( time < t_END - 1e-8){</font><br>
<font size="-1"> try{</font><br>
<font size="-1"> //PDE-SOLVE</font><br>
<font size="-1"> osm.apply( time, dt, uold, unew );</font><br>
<font size="-1"> exceptionCaught = false;</font><br>
<font size="-1"> }catch ( Dune::Exception &e ) {</font><br>
<font size="-1"> //RESET</font><br>
<font size="-1"> exceptionCaught = true;</font><br>
<font size="-1"> std::cout << "Catched Error,
Dune reported error: " << e << std::endl;</font><br>
<font size="-1"> unew = uold;</font><br>
<font size="-1"> dt *= 0.5;</font><br>
<font size="-1">
osm.getPDESolver().discardMatrix();</font><br>
<font size="-1"> continue;</font><br>
<font size="-1"> }</font><br>
<font size="-1"> uold = unew;</font><br>
<font size="-1"> time += dt;</font><br>
<font size="-1"> }</font><br>
<p><font size="-1"><br>
</font></p>
<p><font size="-1"><b>// terminal output showing FMatrixError...</b><br>
</font></p>
<p><font size="-1"> <br>
<font color="#3333ff"> time = 162.632 , time+dt = 164.603 ,
opTime = 180 , dt : 1.97044<br>
<br>
READY FOR NEXT ITERATION. <br>
_____________________________________________________<br>
current opcount = 2<br>
****************************<br>
TCH HYDRATE:<br>
****************************<br>
TIME STEP [implicit Euler] 89 time (from): 1.6263e+02
dt: 1.9704e+00 time (to): 1.6460e+02<br>
STAGE 1 time (to): 1.6460e+02.<br>
Initial defect: 2.1649e-01<br>
Using a direct coarse solver (SuperLU)<br>
Building hierarchy of 2 levels (inclusive coarse solver)
took 0.2195 seconds.<br>
=== BiCGSTABSolver<br>
12.5 6.599e-11<br>
=== rate=0.1733, T=1.152, TIT=0.09217, IT=12.5<br>
Newton iteration 1. New defect: 3.4239e-02. Reduction
(this): 1.5816e-01. Reduction (total): 1.5816e-01<br>
Using a direct coarse solver (SuperLU)<br>
Building hierarchy of 2 levels (inclusive coarse solver)
took 0.195 seconds.<br>
=== BiCGSTABSolver<br>
17 2.402e-11<br>
=== rate=0.2894, T=1.655, TIT=0.09738, IT=17<br>
Newton iteration 2. New defect: 3.9906e+00. Reduction
(this): 1.1655e+02. Reduction (total): 1.8434e+01<br>
Using a direct coarse solver (SuperLU)<br>
Building hierarchy of 2 levels (inclusive coarse solver)
took 0.8697 seconds.<br>
=== BiCGSTABSolver<br>
Catched Error, Dune reported error: FMatrixError
[luDecomposition:/home/sgupta/dune_2_6/source/dune/dune-common/dune/common/densematrix.hh:909]:
matrix is singular<br>
_____________________________________________________<br>
current opcount = 2<br>
****************************<br>
TCH HYDRATE:<br>
****************************<br>
TIME STEP [implicit Euler] 89 time (from): 1.6263e+02
dt: 9.8522e-01 time (to): 1.6362e+02<br>
STAGE 1 time (to): 1.6362e+02.<br>
<br>
</font></font></p>
<p><font color="#3333ff"><b><font size="-1" color="#cc0000">...
nothing happens here... the terminal appears to freeze...</font></b><br>
</font></p>
<p><br>
</p>
<p><br>
</p>
<pre class="moz-signature" cols="72">--
Dr. Shubhangi Gupta
Marine Geosystems
GEOMAR Helmholtz Center for Ocean Research
Wischhofstraße 1-3,
D-24148 Kiel
Room: 12-206
Phone: +49 431 600-1402
Email: <a class="moz-txt-link-abbreviated" href="mailto:sgupta@geomar.de" moz-do-not-send="true">sgupta@geomar.de</a></pre>
</div>
</body>
</html>