20 #ifndef PISM_SNESPROBLEM_H
21 #define PISM_SNESPROBLEM_H
23 #include "pism/util/Grid.hh"
24 #include "pism/util/Vector2d.hh"
25 #include "pism/util/petscwrappers/SNES.hh"
26 #include "pism/util/Logger.hh"
38 virtual const std::string&
name();
74 template<
int DOF,
class U>
79 cb->
solver->compute_local_function(info,x,f);
81 MPI_Comm com = MPI_COMM_SELF;
82 PetscErrorCode ierr = PetscObjectGetComm((PetscObject)cb->
da, &com); CHKERRQ(ierr);
84 SETERRQ(com, 1,
"A PISM callback failed");
89 template<
int DOF,
class U>
94 cb->
solver->compute_local_jacobian(info, x,
J);
96 MPI_Comm com = MPI_COMM_SELF;
97 PetscErrorCode ierr = PetscObjectGetComm((PetscObject)cb->
da, &com); CHKERRQ(ierr);
99 SETERRQ(com, 1,
"A PISM callback failed");
104 template<
int DOF,
class U>
114 PISM_CHK(ierr,
"DMCreateGlobalVector");
123 ierr = DMDASNESSetFunctionLocal(*
m_DA, INSERT_VALUES,
126 PISM_CHK(ierr,
"DMDASNESSetFunctionLocal");
130 PISM_CHK(ierr,
"DMDASNESSetJacobianLocal");
132 ierr = DMSetMatType(*
m_DA,
"baij");
136 PISM_CHK(ierr,
"DMSetApplicationContext");
141 ierr = SNESSetFromOptions(
m_snes);
142 PISM_CHK(ierr,
"SNESSetFromOptions");
145 template<
int DOF,
class U>
150 template<
int DOF,
class U>
152 return "UnnamedProblem";
155 template<
int DOF,
class U>
160 ierr = SNESSolve(m_snes,NULL,m_X);
PISM_CHK(ierr,
"SNESSolve");
163 SNESConvergedReason reason;
164 ierr = SNESGetConvergedReason(m_snes, &reason);
PISM_CHK(ierr,
"SNESGetConvergedReason");
167 name().c_str(), SNESConvergedReasons[reason]);
170 m_grid->ctx()->log()->message(1,
"SNESProblem %s converged (SNES reason %s)\n",
171 name().c_str(), SNESConvergedReasons[reason]);
static RuntimeError formatted(const ErrorLocation &location, const char format[],...) __attribute__((format(printf
build a RuntimeError with a formatted message
SNESProblem(std::shared_ptr< const Grid > g)
CallbackData m_callbackData
virtual const std::string & name()
virtual void compute_local_function(DMDALocalInfo *info, const U **xg, U **yg)=0
static PetscErrorCode jacobian_callback(DMDALocalInfo *info, const U **x, Mat B, CallbackData *)
virtual void compute_local_jacobian(DMDALocalInfo *info, const U **x, Mat B)=0
std::shared_ptr< const Grid > m_grid
static PetscErrorCode function_callback(DMDALocalInfo *info, const U **x, U **f, CallbackData *)
std::shared_ptr< Wrapper > Ptr
#define PISM_CHK(errcode, name)
#define PISM_ERROR_LOCATION
SNESProblem< 1, double > SNESScalarProblem
void handle_fatal_errors(MPI_Comm com)
SNESProblem< 2, Vector2d > SNESVectorProblem
SNESProblem< DOF, U > * solver