v0.9.0
Classes | Public Member Functions | Public Attributes | List of all members
EdgeSlidingConstrains Struct Reference

#include <users_modules/basic_finite_elements/src/SurfaceSlidingConstrains.hpp>

Inheritance diagram for EdgeSlidingConstrains:
[legend]
Collaboration diagram for EdgeSlidingConstrains:
[legend]

Classes

struct  CalculateEdgeBase
 
struct  MyEdgeFE
 
struct  OpJacobian
 

Public Member Functions

MyEdgeFEgetLoopFeRhs ()
 
MyEdgeFEgetLoopFeLhs ()
 
MoFEMErrorCode getOptions ()
 
 EdgeSlidingConstrains (MoFEM::Interface &m_field)
 
MoFEMErrorCode setOperators (int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
 
MoFEMErrorCode setOperators (int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name, const double *alpha=nullptr)
 
MoFEMErrorCode setOperatorsConstrainOnly (int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
 
MoFEMErrorCode setOperatorsConstrainOnly (int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
 

Public Attributes

MoFEM::InterfacemField
 
boost::shared_ptr< MyEdgeFEfeRhsPtr
 
boost::shared_ptr< MyEdgeFEfeLhsPtr
 
MyEdgeFEfeRhs
 
MyEdgeFEfeLhs
 
double aLpha
 

Detailed Description

Examples
mesh_smoothing.cpp.

Definition at line 660 of file SurfaceSlidingConstrains.hpp.

Constructor & Destructor Documentation

◆ EdgeSlidingConstrains()

EdgeSlidingConstrains::EdgeSlidingConstrains ( MoFEM::Interface m_field)

Definition at line 952 of file SurfaceSlidingConstrains.hpp.

953  : mField(m_field), feRhsPtr(new MyEdgeFE(m_field)),
954  feLhsPtr(new MyEdgeFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
955  aLpha(1) {
956  ierr = getOptions();
957  CHKERRABORT(PETSC_COMM_WORLD, ierr);
958  }
boost::shared_ptr< MyEdgeFE > feRhsPtr
boost::shared_ptr< MyEdgeFE > feLhsPtr

Member Function Documentation

◆ getLoopFeLhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeLhs ( )

Definition at line 935 of file SurfaceSlidingConstrains.hpp.

935 { return feLhs; }

◆ getLoopFeRhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeRhs ( )

Definition at line 933 of file SurfaceSlidingConstrains.hpp.

933 { return feRhs; }

◆ getOptions()

MoFEMErrorCode EdgeSlidingConstrains::getOptions ( )

Definition at line 939 of file SurfaceSlidingConstrains.hpp.

939  {
941  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
942  "Get edge sliding constrains element scaling",
943  "none");
944  CHKERRQ(ierr);
945  CHKERR PetscOptionsScalar("-edge_sliding_alpha", "scaling parameter", "",
946  aLpha, &aLpha, PETSC_NULL);
947  ierr = PetscOptionsEnd();
948  CHKERRQ(ierr);
950  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
CHKERRQ(ierr)
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ setOperators() [1/2]

MoFEMErrorCode EdgeSlidingConstrains::setOperators ( int  tag,
Range  edges,
Range  faces,
const std::string  lagrange_multipliers_field_name,
const std::string  material_field_name 
)

Definition at line 1119 of file SurfaceSlidingConstrains.hpp.

1121  {
1124  edges, faces);
1125  CHKERR setOperators(tag, lagrange_multipliers_field_name,
1126  material_field_name);
1128  }
MoFEMErrorCode setOperators(int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
static MoFEMErrorCode setTags(moab::Interface &moab, Range edges, Range tris, bool number_pathes=true)

◆ setOperators() [2/2]

MoFEMErrorCode EdgeSlidingConstrains::setOperators ( int  tag,
const std::string  lagrange_multipliers_field_name,
const std::string  material_field_name,
const double alpha = nullptr 
)

Definition at line 1130 of file SurfaceSlidingConstrains.hpp.

1133  {
1135 
1136  if (alpha) {
1137  aLpha = *alpha;
1138  }
1139 
1140  boost::shared_ptr<VectorDouble> active_variables_ptr(
1141  new VectorDouble(4 + 6));
1142  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1143  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1144  new MatrixDouble(4 + 6, 4 + 6));
1145 
1146  feRhs.getOpPtrVector().clear();
1147  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1148  lagrange_multipliers_field_name, active_variables_ptr));
1149  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1150  material_field_name, active_variables_ptr));
1151  feRhs.getOpPtrVector().push_back(new OpJacobian(
1152  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1153  jacobian_ptr, false, aLpha));
1154  feRhs.getOpPtrVector().push_back(
1155  new OpAssembleRhs<4, 6>(lagrange_multipliers_field_name, results_ptr));
1156  feRhs.getOpPtrVector().push_back(
1157  new OpAssembleRhs<4, 6>(material_field_name, results_ptr));
1158 
1159  // Adding operators to calculate the left hand side
1160  feLhs.getOpPtrVector().clear();
1161  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1162  lagrange_multipliers_field_name, active_variables_ptr));
1163  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1164  material_field_name, active_variables_ptr));
1165  feLhs.getOpPtrVector().push_back(new OpJacobian(
1166  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1167  jacobian_ptr, true, aLpha));
1168  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1169  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1170  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1171  material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
1172  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1173  material_field_name, material_field_name, jacobian_ptr));
1174 
1176  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72

◆ setOperatorsConstrainOnly() [1/2]

MoFEMErrorCode EdgeSlidingConstrains::setOperatorsConstrainOnly ( int  tag,
Range  edges,
Range  faces,
const std::string  lagrange_multipliers_field_name,
const std::string  material_field_name 
)

Definition at line 1179 of file SurfaceSlidingConstrains.hpp.

1181  {
1183 
1185  edges, faces);
1186  CHKERR setOperatorsConstrainOnly(tag, lagrange_multipliers_field_name,
1187  material_field_name);
1189  }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
MoFEMErrorCode setOperatorsConstrainOnly(int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
static MoFEMErrorCode setTags(moab::Interface &moab, Range edges, Range tris, bool number_pathes=true)

◆ setOperatorsConstrainOnly() [2/2]

MoFEMErrorCode EdgeSlidingConstrains::setOperatorsConstrainOnly ( int  tag,
const std::string  lagrange_multipliers_field_name,
const std::string  material_field_name 
)

Definition at line 1192 of file SurfaceSlidingConstrains.hpp.

1194  {
1196 
1197  boost::shared_ptr<VectorDouble> active_variables_ptr(
1198  new VectorDouble(4 + 6));
1199  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1200  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1201  new MatrixDouble(4 + 6, 4 + 6));
1202 
1203  // Adding operators to calculate the left hand side
1204  feLhs.getOpPtrVector().clear();
1205  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1206  lagrange_multipliers_field_name, active_variables_ptr));
1207  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1208  material_field_name, active_variables_ptr));
1209  feLhs.getOpPtrVector().push_back(new OpJacobian(
1210  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1211  jacobian_ptr, true, aLpha));
1212  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1213  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1214 
1216  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72

Member Data Documentation

◆ aLpha

double EdgeSlidingConstrains::aLpha

Definition at line 937 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyEdgeFE& EdgeSlidingConstrains::feLhs

Definition at line 934 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feLhsPtr

Definition at line 930 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyEdgeFE& EdgeSlidingConstrains::feRhs

Definition at line 932 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feRhsPtr

Definition at line 930 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface& EdgeSlidingConstrains::mField

Definition at line 882 of file SurfaceSlidingConstrains.hpp.


The documentation for this struct was generated from the following file: