[petsc-users] Context behind SNESComputeFunction call

Kong, Fande fande.kong at inl.gov
Fri Jan 26 17:34:40 CST 2018


Hi Barry,

I made minor changes on src/snes/examples/tutorials/ex2.c to demonstrate
this issue.  Please see the attachment.

./ex2 -snes_monitor -ksp_monitor -snes_mf_operator 1




























*atol=1e-50, rtol=1e-08, stol=1e-08, maxit=50, maxf=10000 FormFunction is
called   0 SNES Function norm 5.414682427127e+00     0 KSP Residual norm
9.559082033938e-01  FormFunction is called  FormFunction is called     1
KSP Residual norm 1.703870633386e-09  FormFunction is called  FormFunction
is called   1 SNES Function norm 2.952582481151e-01     0 KSP Residual norm
2.672054855433e-02  FormFunction is called  FormFunction is called     1
KSP Residual norm 1.519298012177e-10  FormFunction is called  FormFunction
is called   2 SNES Function norm 4.502289047587e-04     0 KSP Residual norm
4.722075651268e-05  FormFunction is called  FormFunction is called     1
KSP Residual norm 3.834927363659e-14  FormFunction is called  FormFunction
is called   3 SNES Function norm 1.390073376131e-09 number of SNES
iterations = 3Norm of error 1.49795e-10, Iterations 3*

"FormFunction" is called TWICE at "0 KSP".

If we comment out MatMFFDSetFunction:

*/* ierr =
MatMFFDSetFunction(Jacobian,FormFunction_MFFD,(void*)snes);CHKERRQ(ierr);
*/*


./ex2 -snes_monitor -ksp_monitor -snes_mf_operator 1
























*atol=1e-50, rtol=1e-08, stol=1e-08, maxit=50, maxf=10000 FormFunction is
called   0 SNES Function norm 5.414682427127e+00     0 KSP Residual norm
9.559082033938e-01  FormFunction is called     1 KSP Residual norm
1.703870633386e-09  FormFunction is called  FormFunction is called   1 SNES
Function norm 2.952582481151e-01     0 KSP Residual norm 2.672054855433e-02
 FormFunction is called     1 KSP Residual norm 1.519298012177e-10
 FormFunction is called  FormFunction is called   2 SNES Function norm
4.502289047587e-04     0 KSP Residual norm 4.722075651268e-05  FormFunction
is called     1 KSP Residual norm 3.834927363659e-14  FormFunction is
called  FormFunction is called   3 SNES Function norm 1.390073376131e-09
number of SNES iterations = 3Norm of error 1.49795e-10, Iterations 3*

"FormFunction" is called ONCE at "0 KSP".

Hopefully, this example makes the point clear.


Fande,

On Fri, Jan 26, 2018 at 3:50 PM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:

>
>
>   So you are doing something non-standard? Are you not just using -snes_mf
> or -snes_mf_operator? Can you send me a sample code that has the extra
> function evaluations? Because if you run through regular usage with the
> debugger you will see there is no extra evaluation.
>
>    Barry
>
>
> > On Jan 26, 2018, at 4:32 PM, Kong, Fande <fande.kong at inl.gov> wrote:
> >
> >
> >
> > On Fri, Jan 26, 2018 at 3:10 PM, Smith, Barry F. <bsmith at mcs.anl.gov>
> wrote:
> >
> >
> > > On Jan 26, 2018, at 2:15 PM, Kong, Fande <fande.kong at inl.gov> wrote:
> > >
> > >
> > >
> > > On Mon, Jan 8, 2018 at 2:15 PM, Smith, Barry F. <bsmith at mcs.anl.gov>
> wrote:
> > >
> > >
> > > > On Jan 8, 2018, at 2:59 PM, Alexander Lindsay <
> alexlindsay239 at gmail.com> wrote:
> > > >
> > > > Is there any elegant way to tell whether SNESComputeFunction is
> being called under different conceptual contexts?
> > > >
> > > > E.g. non-linear residual evaluation vs. Jacobian formation from
> finite differencing vs. Jacobian-vector products from finite differencing?
> > >
> > >   Under normal usage with the options database no.
> > >
> > > Hi Barry,
> > >
> > > How difficult to provide an API? Is it possible?
> > >
> > >
> > >
> > >   If you have some reason to know you could write three functions and
> provide them to SNESSetFunction(), MatMFFDSetFunction(), and
> MatFDColoringSetFunction(). Note that these functions have slightly
> different calling sequences but you can have all of them call the same
> underlying common function if you are only interested in, for example, how
> many times the function is used for each purpose.
> > >
> > > If we use this way for the Jacobian-free Newton, the function
> evaluation will be called twice at the first linear iteration because the
> computed residual vector at the nonlinear step  is not reused. Any way to
> reuse the function residual of the Newton step instead of recomputing a new
> residual at the first linear iteration?
> >
> >    It does reuse the function evaluation. Why do you think it does not?
> If you look at MatMult_MFFD() you will see the lines of code
> >
> >   /* compute func(U) as base for differencing; only needed first time in
> and not when provided by user */
> >   if (ctx->ncurrenth == 1 && ctx->current_f_allocated) {
> >     ierr = (*ctx->func)(ctx->funcctx,U,F);CHKERRQ(ierr);
> >   }
> >
> > since the if is satisfied it does not compute the function at the base
> location.  To double check I ran src/snes/examples/tutorials/ex19 with
> -snes_mf in the debugger and verified that the "extra" function evaluations
> are not done.
> >
> > In MatAssemblyEnd_SNESMF,
> >
> >   if (j->func == (PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction)
> {
> >     ierr = SNESGetFunction(snes,&f,NULL,NULL);CHKERRQ(ierr);
> >     ierr = MatMFFDSetBase_MFFD(J,u,f);CHKERRQ(ierr);
> >   } else {
> >     /* f value known by SNES is not correct for other differencing
> function */
> >     ierr = MatMFFDSetBase_MFFD(J,u,NULL);CHKERRQ(ierr);
> >   }
> >
> >
> > Will hit ierr = MatMFFDSetBase_MFFD(J,u,NULL);CHKERRQ(ierr), because
> SNES and MAT have different function pointers.
> >
> > In MatMFFDSetBase_MFFD(Mat J,Vec U,Vec F),
> >
> >   if (F) {
> >     if (ctx->current_f_allocated) {ierr = VecDestroy(&ctx->current_f);
> CHKERRQ(ierr);}
> >     ctx->current_f           = F;
> >     ctx->current_f_allocated = PETSC_FALSE;
> >   } else if (!ctx->current_f_allocated) {
> >     ierr = MatCreateVecs(J,NULL,&ctx->current_f);CHKERRQ(ierr);
> >
> >     ctx->current_f_allocated = PETSC_TRUE;
> >   }
> >
> > Because F=NULL, we then have ctx->current_f_allocated = PETSC_TRUE.
> >
> > Then, the following if statement is true:
> >
> >   if (ctx->ncurrenth == 1 && ctx->current_f_allocated) {
> >     ierr = (*ctx->func)(ctx->funcctx,U,F);CHKERRQ(ierr);
> >   }
> >
> >
> > Fande,
> >
> >
> >
> >   Barry
> >
> >
> > >
> > > Fande,
> > >
> > >
> > >
> > >    Barry
> > >
> > >
> > >
> > > >
> > > > Alex
> > >
> > >
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20180126/6e5e4cd5/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: ex2.c
Type: text/x-csrc
Size: 12999 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20180126/6e5e4cd5/attachment.bin>


More information about the petsc-users mailing list