|
| v0.14.0
|
Go to the documentation of this file.
28 : data_ctx(
PetscData::CtxSetNone),
f(PETSC_NULL),
A(PETSC_NULL),
29 B(PETSC_NULL), x(PETSC_NULL), x_t(PETSC_NULL), x_tt(PETSC_NULL) {}
33 this->
f = petsc_data.
f;
34 this->
A = petsc_data.
A;
35 this->
B = petsc_data.
B;
36 this->
x = petsc_data.
x;
37 this->
x = petsc_data.
x;
38 this->
x_t = petsc_data.
x_t;
52 : ksp_ctx(CTX_KSPNONE), ksp(PETSC_NULL), ksp_f(
PetscData::
f),
79 this->snes =
snes.snes;
87 *iface =
const_cast<TSMethod *
>(
this);
92 :
ts_ctx(CTX_TSNONE), ts_step(-1), ts_a(0), ts_t(0), ts_u(
PetscData::x),
111 loopSize(0), rAnk(-1), sIze(-1), refinedEntitiesPtr(nullptr),
112 refinedFiniteElementsPtr(nullptr), problemPtr(nullptr),
113 fieldsPtr(nullptr), entitiesPtr(nullptr), dofsPtr(nullptr),
114 finiteElementsPtr(nullptr), finiteElementsEntitiesPtr(nullptr),
115 adjacenciesPtr(nullptr) {}
145 "should be implemented by user in derived class (preProcess)");
156 "should be implemented by user in derived class (postProcess)");
167 "should be implemented by user in derived class (operator)");
174 const bool reset_dofs) {
179 auto get_nodes_field_data =
180 [&](boost::shared_ptr<FEDofEntity_multiIndex> &&dofs,
188 bit_number, get_id_for_min_type<MBVERTEX>()));
191 bit_number, get_id_for_max_type<MBVERTEX>()));
194 auto &first_dof = **dit;
196 const int nb_dof_idx = first_dof.getNbOfCoeffs();
197 const int max_nb_dofs = nb_dof_idx * num_nodes;
198 nodes_data.resize(max_nb_dofs,
false);
201 for (; dit != hi_dit;) {
202 const auto &dof_ptr = *dit;
203 const auto &dof = *dof_ptr;
204 const auto &sn = *dof.getSideNumberPtr();
205 const int side_number = sn.side_number;
207 int pos = side_number * nb_dof_idx;
208 auto ent_filed_data_vec = dof.getEntFieldData();
209 for (
int ii = 0; ii != nb_dof_idx; ++ii) {
210 nodes_data[pos] = ent_filed_data_vec[ii];
216 }
else if (reset_dofs) {
217 nodes_data.resize(0,
false);
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
std::bitset< 8 > Switches
const RefEntity_multiIndex * refinedEntitiesPtr
container of mofem dof entities
PetscReal ts_a
shift for U_t (see PETSc Time Solver)
int nInTheLoop
number currently of processed method
const EntFiniteElement_multiIndex * finiteElementsEntitiesPtr
static constexpr Switches CtxSetB
auto getDataDofsPtr() const
static constexpr Switches CtxSetA
MoFEMErrorCode copyPetscData(const PetscData &petsc_data)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
data structure for snes (nonlinear solver) context
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
const Field_multiIndex * fieldsPtr
raw pointer to fields container
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
const FieldEntityEntFiniteElementAdjacencyMap_multiIndex * adjacenciesPtr
int sIze
number of processors in communicator
#define CHKERR
Inline error check.
const FieldEntity_multiIndex * entitiesPtr
raw pointer to container of field entities
implementation of Data Operators for Forces and Sources
data structure for ksp (linear solver) context
auto getNumberOfNodes() const
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
static constexpr Switches CtxSetX_T
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
const RefElement_multiIndex * refinedFiniteElementsPtr
container of mofem finite element entities
static constexpr Switches CtxSetTime
@ MOFEM_OPERATION_UNSUCCESSFUL
static constexpr Switches CtxSetNone
Data structure to exchange data between mofem and User Loop Methods.
int loopSize
local number oe methods to process
constexpr auto field_name
static constexpr Switches CtxSetX_TT
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
static constexpr Switches CtxSetX
PetscInt ts_step
time step number
base class for all interface classes
const FiniteElement_multiIndex * finiteElementsPtr
raw pointer to container finite elements
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
unsigned int getFieldBitNumber(std::string field_name) const
const DofEntity_multiIndex * dofsPtr
raw pointer container of dofs
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
data structure for TS (time stepping) context
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
boost::weak_ptr< CacheTuple > cacheWeakPtr
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
const Problem * problemPtr
raw pointer to problem
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
static constexpr Switches CtxSetF
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscReal ts_dt
time step size
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
virtual MoFEMErrorCode operator()()
function is run for every finite element
KSPContext ksp_ctx
Context.