v0.10.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 Member Functions inherited from GenericSliding
 GenericSliding ()=default
 
virtual ~GenericSliding ()=default
 

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 647 of file SurfaceSlidingConstrains.hpp.

Constructor & Destructor Documentation

◆ EdgeSlidingConstrains()

EdgeSlidingConstrains::EdgeSlidingConstrains ( MoFEM::Interface m_field)

Definition at line 920 of file SurfaceSlidingConstrains.hpp.

921  : mField(m_field), feRhsPtr(new MyEdgeFE(m_field)),
922  feLhsPtr(new MyEdgeFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
923  aLpha(1) {
924  ierr = getOptions();
925  CHKERRABORT(PETSC_COMM_WORLD, ierr);
926  }
boost::shared_ptr< MyEdgeFE > feRhsPtr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
boost::shared_ptr< MyEdgeFE > feLhsPtr

Member Function Documentation

◆ getLoopFeLhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeLhs ( )

Definition at line 903 of file SurfaceSlidingConstrains.hpp.

903 { return feLhs; }

◆ getLoopFeRhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeRhs ( )

Definition at line 901 of file SurfaceSlidingConstrains.hpp.

901 { return feRhs; }

◆ getOptions()

MoFEMErrorCode EdgeSlidingConstrains::getOptions ( )

Definition at line 907 of file SurfaceSlidingConstrains.hpp.

907  {
909  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
910  "Get edge sliding constrains element scaling",
911  "none");
912  CHKERRQ(ierr);
913  CHKERR PetscOptionsScalar("-edge_sliding_alpha", "scaling parameter", "",
914  aLpha, &aLpha, PETSC_NULL);
915  ierr = PetscOptionsEnd();
916  CHKERRQ(ierr);
918  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ 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 1087 of file SurfaceSlidingConstrains.hpp.

1089  {
1092  edges, faces);
1093  CHKERR setOperators(tag, lagrange_multipliers_field_name,
1094  material_field_name);
1096  }
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:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
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 1098 of file SurfaceSlidingConstrains.hpp.

1101  {
1103 
1104  if (alpha) {
1105  aLpha = *alpha;
1106  }
1107 
1108  boost::shared_ptr<VectorDouble> active_variables_ptr(
1109  new VectorDouble(4 + 6));
1110  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1111  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1112  new MatrixDouble(4 + 6, 4 + 6));
1113 
1114  feRhs.getOpPtrVector().clear();
1115  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1116  lagrange_multipliers_field_name, active_variables_ptr));
1117  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1118  material_field_name, active_variables_ptr));
1119  feRhs.getOpPtrVector().push_back(new OpJacobian(
1120  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1121  jacobian_ptr, false, aLpha));
1122  feRhs.getOpPtrVector().push_back(
1123  new OpAssembleRhs<4, 6>(lagrange_multipliers_field_name, results_ptr));
1124  feRhs.getOpPtrVector().push_back(
1125  new OpAssembleRhs<4, 6>(material_field_name, results_ptr));
1126 
1127  // Adding operators to calculate the left hand side
1128  feLhs.getOpPtrVector().clear();
1129  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1130  lagrange_multipliers_field_name, active_variables_ptr));
1131  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1132  material_field_name, active_variables_ptr));
1133  feLhs.getOpPtrVector().push_back(new OpJacobian(
1134  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1135  jacobian_ptr, true, aLpha));
1136  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1137  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1138  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1139  material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
1140  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1141  material_field_name, material_field_name, jacobian_ptr));
1142 
1144  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:415
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74

◆ 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 1147 of file SurfaceSlidingConstrains.hpp.

1149  {
1151 
1153  edges, faces);
1154  CHKERR setOperatorsConstrainOnly(tag, lagrange_multipliers_field_name,
1155  material_field_name);
1157  }
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:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
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:415
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 1160 of file SurfaceSlidingConstrains.hpp.

1162  {
1164 
1165  boost::shared_ptr<VectorDouble> active_variables_ptr(
1166  new VectorDouble(4 + 6));
1167  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1168  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1169  new MatrixDouble(4 + 6, 4 + 6));
1170 
1171  // Adding operators to calculate the left hand side
1172  feLhs.getOpPtrVector().clear();
1173  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1174  lagrange_multipliers_field_name, active_variables_ptr));
1175  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1176  material_field_name, active_variables_ptr));
1177  feLhs.getOpPtrVector().push_back(new OpJacobian(
1178  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1179  jacobian_ptr, true, aLpha));
1180  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1181  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1182 
1184  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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:415
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74

Member Data Documentation

◆ aLpha

double EdgeSlidingConstrains::aLpha

Definition at line 905 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyEdgeFE& EdgeSlidingConstrains::feLhs

Definition at line 902 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feLhsPtr

Definition at line 898 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyEdgeFE& EdgeSlidingConstrains::feRhs

Definition at line 900 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feRhsPtr

Definition at line 898 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface& EdgeSlidingConstrains::mField

Definition at line 869 of file SurfaceSlidingConstrains.hpp.


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