[Trilinos-Users] [EXTERNAL] Re: Trouble Ensuring x (as Guess) is used in Belos Solver

Corey A. Henderson cahenderson at wisc.edu
Sat Mar 16 13:28:15 EDT 2019


Heidi, thank you for your comprehensive explanation. I added the parameters
your suggested to my setup and immediately saw the results I expected.

Thank you very much!

Corey

On Wed, Feb 13, 2019 at 2:30 PM Thornquist, Heidi K <hkthorn at sandia.gov>
wrote:

> Hi Corey,
>
>
>
> The tiny hiccup in the pseudocode is that it doesn’t seem to account for
> the fact that the GMRES solvers in Belos use relative residuals to
> determine convergence.  By default, the relative scaling is the norm of the
> initial residual vector.  If you recursively solve the same linear system,
> without resetting the solution vector to zero, the linear solver will
> continue to try and drop the relative residual by the tolerance specified.
> That means that after the second solve, if you leave b unchanged, the
> relative residual with respect to the norm of b will be twice around the
> accuracy specified by the original tolerance.
>
>
>
> As an example, take test_resolve_gmres_hb.cpp in the directory
> packages/belos/epetra/test/BlockGmres.  That test illustrates multiple ways
> to reset a solver object when a user is looking to reuse a solver object.
> If I refrain from resetting the solution vector after the first solve (line
> 198), the second solve will resolve the linear problem using the solution
> vector from the first solve.  That looks like this:
>
>
>
> First solve:
>
>
>
> Belos::StatusTestGeneralOutput: Passed
>
>   (Num calls,Mod test,State test): (211, 1, Passed)
>
>    Passed.......OR Combination ->
>
>      OK...........Number of Iterations = 210 < 1029
>
>      Converged....(2-Norm Res Vec) / (2-Norm Prec Res0)
>
>                   residual [ 0 ] = 9.48466e-06 < 1e-05
>
>
>
> Second solve:
>
>
>
> Belos::StatusTestGeneralOutput: Passed
>
>   (Num calls,Mod test,State test): (401, 1, Passed)
>
>    Passed.......OR Combination ->
>
>      OK...........Number of Iterations = 400 < 1029
>
>      Converged....(2-Norm Res Vec) / (2-Norm Prec Res0)
>
>                   residual [ 0 ] = 9.63506e-06 < 1e-05
>
>
>
> However, after the first solve my residual with respect to the norm of the
> initial residual vector is 9.48466e-06 and after the second solve it is
> 9.13853e-11.  The best way to handle this is have the solver use a scaling
> different from the norm of the initial residual vector.  The norm of the b
> vector is a good choice here, so add these parameters to the Belos
> parameter list:
>
>
>
> belosList.set( "Implicit Residual Scaling", "Norm of RHS" ); // Implicit
> residual scaling for convergence
>
> belosList.set( "Explicit Residual Scaling", "Norm of RHS" ); // Explicit
> residual scaling for convergence
>
> Changing this option in the test_resolve_gmres_hb.cpp code, results in the
> second solve not performing any iterations because convergence was already
> achieved.  This is the desired result:
>
>
>
> Belos::StatusTestGeneralOutput: Passed
>
>   (Num calls,Mod test,State test): (1, 1, Passed)
>
>    Passed.......OR Combination ->
>
>      OK...........Number of Iterations = 0 < 1029
>
>      Converged....(2-Norm Res Vec) / (2-Norm RHS )
>
>                   residual [ 0 ] = 9.48466e-06 < 1e-05
>
>
>
> Hope this helps.
>
>
>
> Thanks,
>
> Heidi
>
>
>
> --
>
>
>
> Heidi K. Thornquist
>
> Electrical Models & Simulation
>
> Sandia National Laboratories
>
>
>
>
>
>
>
> *From: *Trilinos-Users <trilinos-users-bounces at trilinos.org> on behalf of
> "Corey A. Henderson" <cahenderson at wisc.edu>
> *Date: *Wednesday, February 13, 2019 at 10:36 AM
> *To: *"Heroux, Mike" <MHeroux at csbsju.edu>
> *Cc: *"trilinos-users at trilinos.org" <trilinos-users at trilinos.org>
> *Subject: *[EXTERNAL] Re: [Trilinos-Users] Trouble Ensuring x (as Guess)
> is used in Belos Solver
>
>
>
> I did not. I was going to ping the thread to be sure it had been seen.
>
>
>
> I am still not able to convince myself that x is being used as a guess.
>
>
>
> On Sat, Feb 9, 2019, 8:19 PM Heroux, Mike <MHeroux at csbsju.edu wrote:
>
> Corey,
>
>
>
> Did you get past this issue?  Messages being sent to
> trilinos-users at trilinos.org are not being reliably received by everyone
> on the list.
>
>
>
> Thanks.
>
>
>
> Mike
>
>
>
>
>
> *From: *Trilinos-Users <trilinos-users-bounces at trilinos.org> on behalf of
> "Corey A. Henderson" <cahenderson at wisc.edu>
> *Date: *Wednesday, December 26, 2018 at 1:30 PM
> *To: *"trilinos-users at trilinos.org" <trilinos-users at trilinos.org>
> *Subject: *[Trilinos-Users] Trouble Ensuring x (as Guess) is used in
> Belos Solver
>
>
>
> Hello,
>
>
>
> I am trying to set up an iterative call to the Pseudoblock GMRES solver in
> Belos where I expect x to vary only slightly between calls (because b
> varies only slightly). When I inspect the number of iterations after
> subsequent calls, even if b is not changed, I get no reduction in
> iterations as I would expect from using the previous x as a guess. I
> believe I am missing something in how Belos is reset().
>
>
>
> I need some help in learning how to ensure the previously-solved x is
> reused as a guess for subsequent iterations.
>
>
>
> Pseudocode:
>
>
>
> Set up A, x=0, preconditioner M, initial b
>
> Create problem from M, A, x, b
>
> Create solver ("GMRES")
>
> Solver.setProblem()
>
> loop
>
>     solve()
>
>     print numIters
>
>     copy off of x
>
>     Solver.reset(Problem)  <-- Where I think my mistake is
>
>     update b
>
> end loop
>
>
>
> Even if b is unchanged, I get the same number of iterations (and duration)
> on subsequent solves, where I would expect them to drop significantly since
> x is near (or just is) the solution. If I omit the solver.reset() call I
> get an error on the next solve(). I believe x is not being used as a guess
> the way I'm setting up my code, but I don't understand why.
>
>
>
> Can someone explain how to ensure Belos uses x as a guess? I have been
> over the documentation and can't figure it what I'm doing wrong.
>
>
>
> Thanks!
>
>
>
> Corey
>
> --
>
> Corey A. Henderson
>
> PhD Candidate and NSF Graduate Fellow
>
> Dept. of Engineering Physics
>
> Univ. of Wisconsin - Madison
>
>

-- 
Corey A. Henderson
PhD Candidate and NSF Graduate Fellow
Dept. of Engineering Physics
Univ. of Wisconsin - Madison
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://trilinos.org/pipermail/trilinos-users/attachments/20190316/17c2ecda/attachment.html>


More information about the Trilinos-Users mailing list