[Dune] [#928] breakdown of BiCGSTAB
Bernd Flemisch
bernd at iws.uni-stuttgart.de
Thu Jun 9 19:37:54 CEST 2011
Sorry to be so noisy today. The only difference in the
algorithm between ISTL and Octave seems to be that ISTL
does not waste any space to include an intermediate vector
"s." I did not pursue that further.
There's a crucial difference in checking for convergence,
though: Octave only checks whether two succeeding residual
norms are the same, while ISTL checks whether some numbers
are too small. If I simply comment these three tests in
lines 760ff and 790f, everything runs smoothly and the
solver terminates successfully.
I am not sure whether this can be considered a solution...
Kind regards
Bernd
On Thu, 09 Jun 2011 19:06:30 +0200
Bernd Flemisch <bernd at iws.uni-stuttgart.de> wrote:
> Dear all,
>
> sorry, Markus, I misread your first answer and saw "b"
>instead of "x."
>
> Yes, it is fair to assume that one is clever enough to
>insert known values in the initial guess. I am just
>surprised that the BiCGSTAB algorithm is not able to
>handle that. So I guess the lesson learned is "put your
>Dirichlet BCs in the initial guess."
>
> Matlab claims to start with an all zero initial guess by
>default. It can return the residual norms. They look the
>same as ISTL's in the beginning and then slowly start to
>deviate. I have no idea what the reason for this can be.
>I attach the norms below.
>
> So I looked at Octave. Octave also solves it. And the
>bicgstab routine is textbook style and GPL'd, so I attach
>it. I will try to detect differences.
>
> Thank you again for your replies and efforts.
>
> Kind regards
> Bernd
>
> Matlab:
> 547.722557505166e+003
> 70.6686272173148e+003
> 16.0919594944600e+003
> 13.6905293585520e+000
> 7.97989495493456e+000
> 7.96042781743378e+000
> 6.27748702554648e+000
> 6.13752800810622e+000
> 4.38501049646480e+000
> 4.34518363816910e+000
> 3.86625793790769e+000
> 3.25040879545531e+000
> 2.87194723734970e+000
> 2.99085129702484e+000
> 2.31632828961268e+000
> 1.97561111997182e+000
> 1.74999000800660e+000
> 1.67753388249665e+000
> 339.915287594688e-003
> 105.666640663521e-003
> 87.0123511764460e-003
> 377.111960192396e-003
> 113.631126810502e-003
> 87.1530655241023e-003
> 20.8750454417466e-003
> 12.9096724926136e-003
> 5.72985162518586e-003
> 4.39704303378432e-003
> 1.65452338886661e-003
> 1.32933748983879e-003
> 1.05134740698569e-003
> 983.690065855372e-006
> 730.930761387681e-006
> 559.097067469436e-006
> 481.176318592935e-006
> 2.11860499559708e-003
> 292.479733029574e-006
> 334.437457236489e-006
> 197.114212603315e-006
> 177.406543554174e-006
> 147.862682181941e-006
>
>
> ISTL:
> Iter Defect Rate
> 0 547723
> 0.5 70668.6 0.129023
> 1 16092 0.22771
> 1.5 13.6905 0.000850768
> 2 7.97989 0.582877
> 2.5 7.9812 1.00016
> 3 6.28596 0.787596
> 3.5 6.10793 0.971679
> 4 4.38601 0.718085
> 4.5 4.08159 0.930592
> 5 3.29403 0.807045
> 5.5 2.75974 0.837801
> 6 2.55524 0.925898
> 6.5 3.29023 1.28764
> 7 1.67157 0.508041
> 7.5 1.67186 1.00017
> 8 0.436379 0.261014
> 8.5 0.176398 0.404231
> 9 0.112041 0.635158
> 9.5 0.227358 2.02924
> 10 0.0535411 0.235493
> 10.5 0.0378322 0.706601
> 11 0.0126924 0.335492
> 11.5 0.00882755 0.695499
> 12 0.00671734 0.760951
> 12.5 0.00856232 1.27466
> 13 0.00698652 0.815961
> 13.5 0.00436435 0.624682
> 14 0.00386916 0.886537
> 14.5 0.00348031 0.8995
> 15 0.0021925 0.629972
> 15.5 0.0222249 10.1368
> 16 0.0013511 0.0607922
> 16.5 0.00116566 0.862745
> 17 0.000935422 0.802485
> 17.5 0.000591171 0.631983
> 18 0.000375963 0.635964
> 18.5 0.000333767 0.887764
> 19 0.000286355 0.857951
> 19.5 0.000169482 0.591859
> 20 0.000139933 0.825651
> 20.5 0.000125247 0.895053
> 21 0.000104286 0.832641
> 21.5 9.08415e-05 0.87108
> 22 8.19629e-05 0.902263
> Dune reported error: ISTLError
>[apply:../../../dune/istl/solvers.hh:791]: h=0 in
>BiCGSTAB
>
>
>
>
>
>
> On Thu, 9 Jun 2011 14:58:10 +0200
> Andreas Lauser <andreas.lauser at iws.uni-stuttgart.de>
>wrote:
>> Hi All,
>>
>> I might be wrong (In this case I appologize for the
>>noise), but I always thought that BiCG(Stab) reduces to
>>CG in the symmetric positive definite case. If this is
>>true, then the breakdown of BiCGStab but not CG indicates
>>a problem with the ISTL implementation of BiCGStab,
>>right?
>>
>> cheers
>> Andreas
>>
>> Am Donnerstag, 9. Juni 2011, 14:44:52 schrieb Markus
>>Blatt:
>>> Hi,
>>>
>>> maybe my initial answer was not clear enough.
>>>
>>> On Thu, Jun 09, 2011 at 02:14:05PM +0200, Bernd Flemisch
>>>wrote:
>>> > thank you for your answer. Although you apparently
>>>think that my
>>> > answer is not interesting since you already closed
>>>this task, I
>>> > provide it here anyway.
>>>
>>> I closed it because with my suggested initial guess it
>>>works!
>>>
>>> > This matrix and right hand side are coming from a to
>>>my knowledge
>>> > correct CVFE discretization of a simple diffusion
>>>equation on the
>>> > unit square discretized with 8x4 rectangular elements.
>>>Dirichlet BCs
>>> > on top and bottom, Neumann zero left and right. The
>>>solution of the
>>> > equation system is the physically correct one, meaning
>>>constant in x
>>> > and linear in y direction.
>>> >
>>> > You are somehow right about suggesting a different
>>>right hand side.
>>> > The unusual form is coming from the fact that our
>>>Dirichlet BCs are
>>> > implemented a bit differently leading to 1's and 2's
>>>on the
>>> > diagonals in the original example. I changed that in
>>>the newly
>>> > attached example with only 1's on the Dirichlet
>>>diagonals and your
>>> > suggested rhs. BiCGSTAB still fails.
>>>
>>> Actually I suggested to change the initial guess by
>>>incorporating the
>>> already known Dirichlet values into it. This seems fair.
>>>(I changed x
>>> and not b and used your original matrix!) This approach
>>>worked for me
>>> and somehow makes sense, doesn't it?
>>>
>>> > Matlab's bicgstab works.
>>>
>>> Using the same linear system, right hand side, initial
>>>guess and
>>> preconditioner? (Or are the Dirichlet Boundaries not
>>>represented as
>>> dofs?) Then they maybe use same look ahead techniques to
>>> prevent the breakdown. As I do not have Matlab I can
>>>only guess this.
>>>
>>> Of course incorporating such techniques would be nice
>>>for ISTL, too?
>>> Maybe filing another feature request is a good idea.
>>>
>>>
>>> Cheers,
>>>
>>> Markus
>>
>>
>>
>> --
>> Andreas Lauser
>> Lehrstuhl für Hydromechanik und Hydrosystemmodellierung
>> Universität Stuttgart - Institut für Wasserbau
>> Pfaffenwaldring 61
>> D-70569 Stuttgart
>> Tel: (+49) 0711/ 685-64719
>>Fax: (+49) 0711/ 685-60430
>> e-mail: andreas.lauser at iws.uni-stuttgart.de
>>
>> _______________________________________________
>> Dune mailing list
>> Dune at dune-project.org
>> http://lists.dune-project.org/mailman/listinfo/dune
>
> ___________________________________________________________
>
> Bernd Flemisch phone: +49 711 685
>69162
> IWS, Universitaet Stuttgart fax: +49 711 685
>67020
> Pfaffenwaldring 61 email:
>bernd at iws.uni-stuttgart.de
> D-70569 Stuttgart url:
>www.hydrosys.uni-stuttgart.de
> ___________________________________________________________
___________________________________________________________
Bernd Flemisch phone: +49 711 685
69162
IWS, Universitaet Stuttgart fax: +49 711 685
67020
Pfaffenwaldring 61 email:
bernd at iws.uni-stuttgart.de
D-70569 Stuttgart url:
www.hydrosys.uni-stuttgart.de
___________________________________________________________
More information about the Dune
mailing list