: set sqrt rel error in function (default %g)\n",p,mfctx->error_rel);CHKERRQ(ierr);
ierr = (*PetscHelpPrintf)(mfctx->snes->comm," %ssnes_mf_period : how often h is recomputed (default 1, everytime)\n",p);CHKERRQ(ierr);
if (mfctx->ops->printhelp) {
ierr = (*mfctx->ops->printhelp)(mfctx);CHKERRQ(ierr);
}
}
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFGetH"
/*@
MatSNESMFGetH - Gets the last value that was used as the differencing
parameter.
Not Collective
Input Parameters:
. mat - the matrix obtained with MatCreateSNESMF()
Output Paramter:
. h - the differencing step size
Level: advanced
.keywords: SNES, matrix-free, parameters
.seealso: MatCreateSNESMF(),MatSNESMFSetHHistory(),
MatSNESMFResetHHistory(),MatSNESMFKSPMonitor()
@*/
int MatSNESMFGetH(Mat mat,Scalar *h)
{
MatSNESMFCtx ctx;
int ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
if (ctx) {
*h = ctx->currenth;
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFKSPMonitor"
/*
MatSNESMFKSPMonitor - A KSP monitor for use with the default PETSc
SNES matrix free routines. Prints the differencing parameter used at
each step.
*/
int MatSNESMFKSPMonitor(KSP ksp,int n,PetscReal rnorm,void *dummy)
{
PC pc;
MatSNESMFCtx ctx;
int ierr;
Mat mat;
MPI_Comm comm;
PetscTruth nonzeroinitialguess;
PetscFunctionBegin;
ierr = PetscObjectGetComm((PetscObject)ksp,&comm);CHKERRQ(ierr);
ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
ierr = KSPGetInitialGuessNonzero(ksp,&nonzeroinitialguess);CHKERRQ(ierr);
ierr = PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
if (!ctx) {
SETERRQ(1,1,"Matrix is not a matrix free shell matrix");
}
if (n > 0 || nonzeroinitialguess) {
#if defined(PETSC_USE_COMPLEX)
ierr = PetscPrintf(comm,"%d KSP Residual norm %14.12e h %g + %g i\n",n,rnorm,
PetscRealPart(ctx->currenth),PetscImaginaryPart(ctx->currenth));CHKERRQ(ierr);
#else
ierr = PetscPrintf(comm,"%d KSP Residual norm %14.12e h %g \n",n,rnorm,ctx->currenth);CHKERRQ(ierr);
#endif
} else {
ierr = PetscPrintf(comm,"%d KSP Residual norm %14.12e\n",n,rnorm);CHKERRQ(ierr);
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFSetFunction"
/*@C
MatSNESMFSetFunction - Sets the function used in applying the matrix free.
Collective on Mat
Input Parameters:
+ mat - the matrix free matrix created via MatCreateSNESMF()
. v - workspace vector
. func - the function to use
- funcctx - optional function context passed to function
Level: advanced
Notes:
If you use this you MUST call MatAssemblyBegin()/MatAssemblyEnd() on the matrix free
matrix inside your compute Jacobian routine
If this is not set then it will use the function set with SNESSetFunction()
.keywords: SNES, matrix-free, function
.seealso: MatCreateSNESMF(),MatSNESMFGetH(),
MatSNESMFSetHHistory(), MatSNESMFResetHHistory(),
MatSNESMFKSPMonitor(), SNESetFunction()
@*/
int MatSNESMFSetFunction(Mat mat,Vec v,int (*func)(SNES,Vec,Vec,void *),void *funcctx)
{
MatSNESMFCtx ctx;
int ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
if (ctx) {
ctx->func = func;
ctx->funcctx = funcctx;
ctx->funcvec = v;
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFSetPeriod"
/*@
MatSNESMFSetPeriod - Sets how often h is recomputed, by default it is everytime
Collective on Mat
Input Parameters:
+ mat - the matrix free matrix created via MatCreateSNESMF()
- period - 1 for everytime, 2 for every second etc
Options Database Keys:
+ -snes_mf_period
Level: advanced
.keywords: SNES, matrix-free, parameters
.seealso: MatCreateSNESMF(),MatSNESMFGetH(),
MatSNESMFSetHHistory(), MatSNESMFResetHHistory(),
MatSNESMFKSPMonitor()
@*/
int MatSNESMFSetPeriod(Mat mat,int period)
{
MatSNESMFCtx ctx;
int ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
if (ctx) {
ctx->recomputeperiod = period;
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFSetFunctionError"
/*@
MatSNESMFSetFunctionError - Sets the error_rel for the approximation of
matrix-vector products using finite differences.
Collective on Mat
Input Parameters:
+ mat - the matrix free matrix created via MatCreateSNESMF()
- error_rel - relative error (should be set to the square root of
the relative error in the function evaluations)
Options Database Keys:
+ -snes_mf_err - Sets error_rel
Level: advanced
Notes:
The default matrix-free matrix-vector product routine computes
.vb
F'(u)*a = [F(u+h*a) - F(u)]/h where
h = error_rel*u'a/||a||^2 if |u'a| > umin*||a||_{1}
= error_rel*umin*sign(u'a)*||a||_{1}/||a||^2 else
.ve
.keywords: SNES, matrix-free, parameters
.seealso: MatCreateSNESMF(),MatSNESMFGetH(),
MatSNESMFSetHHistory(), MatSNESMFResetHHistory(),
MatSNESMFKSPMonitor()
@*/
int MatSNESMFSetFunctionError(Mat mat,PetscReal error)
{
MatSNESMFCtx ctx;
int ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
if (ctx) {
if (error != PETSC_DEFAULT) ctx->error_rel = error;
}
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFAddNullSpace"
/*@
MatSNESMFAddNullSpace - Provides a null space that an operator is
supposed to have. Since roundoff will create a small component in
the null space, if you know the null space you may have it
automatically removed.
Collective on Mat
Input Parameters:
+ J - the matrix-free matrix context
- nullsp - object created with MatNullSpaceCreate()
Level: advanced
.keywords: SNES, matrix-free, null space
.seealso: MatNullSpaceCreate(), MatSNESMFGetH(), MatCreateSNESMF(),
MatSNESMFSetHHistory(), MatSNESMFResetHHistory(),
MatSNESMFKSPMonitor(), MatSNESMFErrorRel()
@*/
int MatSNESMFAddNullSpace(Mat J,MatNullSpace nullsp)
{
int ierr;
MatSNESMFCtx ctx;
MPI_Comm comm;
PetscFunctionBegin;
ierr = PetscObjectGetComm((PetscObject)J,&comm);CHKERRQ(ierr);
ierr = MatShellGetContext(J,(void **)&ctx);CHKERRQ(ierr);
/* no context indicates that it is not the "matrix free" matrix type */
if (!ctx) PetscFunctionReturn(0);
ctx->sp = nullsp;
ierr = PetscObjectReference((PetscObject)nullsp);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFSetHHistory"
/*@
MatSNESMFSetHHistory - Sets an array to collect a history of the
differencing values (h) computed for the matrix-free product.
Collective on Mat
Input Parameters:
+ J - the matrix-free matrix context
. histroy - space to hold the history
- nhistory - number of entries in history, if more entries are generated than
nhistory, then the later ones are discarded
Level: advanced
Notes:
Use MatSNESMFResetHHistory() to reset the history counter and collect
a new batch of differencing parameters, h.
.keywords: SNES, matrix-free, h history, differencing history
.seealso: MatSNESMFGetH(), MatCreateSNESMF(),
MatSNESMFResetHHistory(),
MatSNESMFKSPMonitor(), MatSNESMFSetFunctionError()
@*/
int MatSNESMFSetHHistory(Mat J,Scalar *history,int nhistory)
{
int ierr;
MatSNESMFCtx ctx;
PetscFunctionBegin;
ierr = MatShellGetContext(J,(void **)&ctx);CHKERRQ(ierr);
/* no context indicates that it is not the "matrix free" matrix type */
if (!ctx) PetscFunctionReturn(0);
ctx->historyh = history;
ctx->maxcurrenth = nhistory;
ctx->currenth = 0;
PetscFunctionReturn(0);
}
#undef __FUNC__
#define __FUNC__ /**/"MatSNESMFResetHHistory"
/*@
MatSNESMFResetHHistory - Resets the counter to zero to begin
collecting a new set of differencing histories.
Collective on Mat
Input Parameters:
. J - the matrix-free matrix context
Level: advanced
Notes:
Use MatSNESMFSetHHistory() to create the original history counter.
.keywords: SNES, matrix-free, h history, differencing history
.seealso: MatSNESMFGetH(), MatCreateSNESMF(),
MatSNESMFSetHHistory(),
MatSNESMFKSPMonitor(), MatSNESMFSetFunctionError()
@*/
int MatSNESMFResetHHistory(Mat J)
{
int ierr;
MatSNESMFCtx ctx;
PetscFunctionBegin;
ierr = MatShellGetContext(J,(void **)&ctx);CHKERRQ(ierr);
/* no context indicates that it is not the "matrix free" matrix type */
if (!ctx) PetscFunctionReturn(0);
ctx->ncurrenth = 0;
PetscFunctionReturn(0);
}