[petsc-users] Computing residual norm in KSPFGMRESCycle()

[Jed Brown]

Dave, have you considered using GCR instead of FGMRES?  It's a flexible
method that is equivalent in many circumstances, but provides the
residual and solution at each iteration without needing to "build" it.

Dave Lee via petsc-users <petsc-users at mcs.anl.gov> writes:

> Hi Matt and Barry,
>
> thanks for the good ideas.
>
> I wasn't aware of the function KSPSetConvergenceTest(), thanks for alerting
> me to this Matt. Unfortunately I really do need to exclude a couple of
> specific degrees of freedom from the norm, and as such I want to make sure
> I can replicate the norm as computed within UpdateHessenberg() before I
> replace this.
>
> Apologies Barry, there was a stupid bug in that code, it should have read:
>     KSPFGMRESUpdateHessenberg(ksp,loc_it,hapend,&res_norm);
>     KSPFGMRESBuildSoln(fgmres->nrs,ksp->vec_sol,tmp_sol,ksp,loc_it+1);
>     VecAXPY(tmp_sol,-1.0,ksp->vec_rhs);
>     VecNorm(tmp_sol,NORM_2,&res_tmp);
>
> Which does not give the same result as the norm for UpdateHessenberg().
> Unfortunately I want to avoid using
> KSPBuildResidual(ksp,NULL,NULL,&resid);
> as this calls KSP_MatMult(), which will in turn call my snes residual
> assembly routine. This routine in turn calls a Navier Stokes solver to run
> for a fixed (long) time period, and so is expensive.
>
> Replacing VecNorm() with my own version is indeed what I'm doing elsewhere,
> I'm just not sure of how to replace the UpdateHessenberg() norm right now.
>
> Cheers, Dave.
>
> On Fri, Jul 5, 2019 at 12:59 PM Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
>
>>
>>
>> > On Jul 4, 2019, at 7:29 PM, Dave Lee via petsc-users <
>> petsc-users at mcs.anl.gov> wrote:
>> >
>> > Hi PETSc,
>> >
>> > I have a problem for which I need to exclude certain degrees of freedom
>> from the evaluation of the norm as used to monitor the residual in FGMRES
>> (don't ask!). Therefore I would like to replace the existing norm
>> evaluation with my own version.
>> >
>> > Currently this is done within:
>> > KSPFGMRESUpdateHessenberg(ksp,loc_it,hapend,&res_norm)
>> > as:
>> > *res = PetscAbsScalar(*RS(it+1));
>> >
>> > Firstly, I would like to make sure I can replicate the existing residual
>> norm. I tried to do this as:
>> > KSPFGMRESUpdateHessenberg(ksp,loc_it,hapend,&res_norm);
>> > KSPFGMRESBuildSoln(fgmres->nrs,ksp >vec_sol,tmp_sol,ksp,loc_it);
>>
>>    I don't understand the lines of code below. I don't see how they could
>> compute the norm of the residual.
>>
>> From KSPMonitorTrueResidualNorm() you can compute the unpreconditioned
>> residual with
>>
>>   ierr = KSPBuildResidual(ksp,NULL,NULL,&resid);CHKERRQ(ierr);
>>   ierr = VecNorm(resid,NORM_2,&truenorm);CHKERRQ(ierr);
>>   ierr = VecDestroy(&resid);CHKERRQ(ierr);
>>
>> since FGMRES always uses right preconditioning (and hence the
>> unpreconditioned residual) this should match the value computed internally
>> by FGMRES (subject to roundoff differences)
>>
>> To exclude certain values in the residual from the norm calculation you
>> could write a custom "VecNorm" norm that skipped those entries.
>>
>>   Barry
>>
>>
>>
>>
>>
>> > VecNorm(tmp_sol,NORM_2,&res_tmp);
>> > VecNorm(ksp->vec_rhs,NORM_2,&res_rhs);
>> > res_norm = fabs(res_tmp-res_rhs);
>> >
>> > But this doesn't match the norm that comes out of UpdateHessenberg.
>> >
>> > Any ideas on how I can re-create the norm derived from UpdateHessenberg
>> as an expansion over the Kyrlov vectors in FGMRES?
>> >
>> > Cheers, Dave.
>>
>>