3#ifndef __SMOOTHER_HPP__
4#define __SMOOTHER_HPP__
7#error "MoFEM need to be compiled with ADOL-C"
23 CHKERRABORT(PETSC_COMM_SELF,
ierr);
28 PetscOptionsBegin(PETSC_COMM_WORLD,
"",
29 "Get stabilisation element options",
"none");
31 PetscBool smoothing_on =
sTabilised ? PETSC_TRUE : PETSC_FALSE;
32 CHKERR PetscOptionsBool(
"-smoothing_stabilise",
33 "all nodes controlled by smoothing element",
"",
34 smoothing_on, &smoothing_on, PETSC_NULLPTR);
35 sTabilised = (smoothing_on == PETSC_TRUE) ?
true :
false;
43 CHKERRABORT(PETSC_COMM_WORLD,
ierr);
45 CHKERRABORT(PETSC_COMM_WORLD,
ierr);
51 std::map<int, NonlinearElasticElement::BlockData>
setOfBlocks;
65 CHKERR VolumeElementForcesAndSourcesCore::preProcess();
67 if (A != PETSC_NULLPTR)
70 if (
F != PETSC_NULLPTR)
112 CHKERR VolumeElementForcesAndSourcesCore::postProcess();
130 struct OpJacobianSmoother
136 int tag,
bool jacobian)
138 field_name, data, common_data, tag, jacobian, false, false) {}
147 for (
int dd1 = 0; dd1 < 3; dd1++) {
148 for (
int dd2 = 0; dd2 < 3; dd2++) {
172 MoFEMErrorCode
aSemble(
int row_side, EntityType row_type,
173 EntitiesFieldData::EntData &row_data) {
176 int nb_dofs = row_data.getIndices().size();
177 int *indices_ptr = &row_data.getIndices()[0];
180 iNdices.resize(nb_dofs,
false);
181 noalias(
iNdices) = row_data.getIndices();
187 VectorDofs &dofs = row_data.getFieldDofs();
188 VectorDofs::iterator dit = dofs.begin();
189 for (
int ii = 0; dit != dofs.end(); dit++, ii++) {
219 const std::string vel_field,
const std::string
field_name,
224 crack_area_tangent_constrains
233 MoFEMErrorCode
aSemble(
int row_side,
int col_side, EntityType row_type,
235 EntitiesFieldData::EntData &row_data,
236 EntitiesFieldData::EntData &col_data) {
239 int nb_row = row_data.getIndices().size();
240 int nb_col = col_data.getIndices().size();
241 int *row_indices_ptr = &row_data.getIndices()[0];
242 int *col_indices_ptr = &col_data.getIndices()[0];
252 VectorDofs &dofs = row_data.getFieldDofs();
253 VectorDofs::iterator dit = dofs.begin();
254 for (
int ii = 0; dit != dofs.end(); dit++, ii++) {
265 nb_col, col_indices_ptr, &
k(0, 0), ADD_VALUES);
269 const auto bit_number_for_crack_area_tangent_constrain =
271 const auto bit_number_for_mesh_position =
275 double *f_tangent_front_mesh_array;
277 &f_tangent_front_mesh_array);
282 for (
int nn = 0; nn < 4; nn++) {
285 auto dit = row_dofs->get<Unique_mi_tag>().lower_bound(
286 DofEntity::getLoFieldEntityUId(
287 bit_number_for_crack_area_tangent_constrain,
getConn()[nn]));
288 auto hi_dit = row_dofs->get<Unique_mi_tag>().upper_bound(
289 DofEntity::getHiFieldEntityUId(
290 bit_number_for_crack_area_tangent_constrain,
getConn()[nn]));
293 if (std::distance(dit, hi_dit) > 0) {
296 auto diit = row_dofs->get<Unique_mi_tag>().lower_bound(
297 DofEntity::getLoFieldEntityUId(bit_number_for_mesh_position,
300 auto hi_diit = row_dofs->get<Unique_mi_tag>().upper_bound(
301 DofEntity::getHiFieldEntityUId(bit_number_for_mesh_position,
305 for (; diit != hi_diit; diit++) {
307 for (
int ddd = 0; ddd < nb_col; ddd++) {
310 diit->get()->getPetscGlobalDofIdx()) {
313 "data inconsistency %d != %d",
315 diit->get()->getPetscGlobalDofIdx());
318 if (diit->get()->getPetscLocalDofIdx() == -1) {
320 "data inconsistency");
323 f_tangent_front_mesh_array[diit->get()
324 ->getPetscLocalDofIdx()] *
325 k(3 * nn + diit->get()->getDofCoeffIdx(), ddd);
326 int lambda_idx = dit->get()->getPetscGlobalDofIdx();
328 &col_indices_ptr[ddd], &
g, ADD_VALUES);
334 &f_tangent_front_mesh_array);
static PetscErrorCode ierr
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
constexpr auto field_name
unsigned int getFieldBitNumber(std::string field_name) const
Deprecated interface functions.
auto getRowDofsPtr() const
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Mat & snes_B
preconditioner of jacobian matrix
const EntityHandle * getConn()
get element connectivity
data for calculation heat conductivity and heat capacity elements
Range forcesOnlyOnEntitiesRow
boost::shared_ptr< FunctionsToCalculatePiolaKirchhoffI< adouble > > materialAdoublePtr
common data used by volume elements
std::vector< MatrixDouble3by3 > sTress
definition of volume element
Operator performs automatic differentiation.
OpJacobianPiolaKirchhoffStress(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool jacobian, bool ale, bool field_disp)
Construct operator to calculate Piola-Kirchhoff stress or its derivatives over gradient deformation.
ublas::vector< int > rowIndices
OpLhsPiolaKirchhoff_dx(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
ublas::vector< int > iNdices
OpRhsPiolaKirchhoff(const std::string field_name, BlockData &data, CommonData &common_data)
structure grouping operators and data used for calculation of nonlinear elastic element
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MyVolumeFE(MoFEM::Interface &m_field, SmootherBlockData &smoother_data)
MoFEMErrorCode postProcess()
function is run at the end of loop
SmootherBlockData & smootherData
OpJacobianSmoother(const std::string field_name, NonlinearElasticElement::BlockData &data, NonlinearElasticElement::CommonData &common_data, int tag, bool jacobian)
MoFEMErrorCode calculateStress(const int gg)
Calculate Paola-Kirchhoff I stress.
ublas::vector< int > rowFrontIndices
OpLhsSmoother(const std::string vel_field, const std::string field_name, NonlinearElasticElement::BlockData &data, NonlinearElasticElement::CommonData &common_data, SmootherBlockData &smoother_data, const std::string crack_area_tangent_constrains)
const std::string fieldCrackAreaTangentConstrains
SmootherBlockData & smootherData
MoFEMErrorCode aSemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
ublas::vector< int > frontIndices
SmootherBlockData & smootherData
OpRhsSmoother(const std::string field_name, NonlinearElasticElement::BlockData &data, NonlinearElasticElement::CommonData &common_data, SmootherBlockData &smoother_data)
MoFEMErrorCode aSemble(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
virtual ~SmootherBlockData()
MoFEMErrorCode getOptions()
std::map< int, NonlinearElasticElement::BlockData > setOfBlocks
SmootherBlockData smootherData
MyVolumeFE & getLoopFeLhs()
get lhs volume element
Smoother(MoFEM::Interface &m_field)
NonlinearElasticElement::CommonData commonData
boost::shared_ptr< MyVolumeFE > feLhsPtr
MyVolumeFE & getLoopFeRhs()
get rhs volume element
boost::shared_ptr< MyVolumeFE > feRhsPtr
MyVolumeFE & feRhs
calculate right hand side for tetrahedral elements