v0.9.1
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 642 of file SurfaceSlidingConstrains.hpp.

Constructor & Destructor Documentation

◆ EdgeSlidingConstrains()

EdgeSlidingConstrains::EdgeSlidingConstrains ( MoFEM::Interface m_field)

Definition at line 915 of file SurfaceSlidingConstrains.hpp.

916  : mField(m_field), feRhsPtr(new MyEdgeFE(m_field)),
917  feLhsPtr(new MyEdgeFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
918  aLpha(1) {
919  ierr = getOptions();
920  CHKERRABORT(PETSC_COMM_WORLD, ierr);
921  }
boost::shared_ptr< MyEdgeFE > feRhsPtr
boost::shared_ptr< MyEdgeFE > feLhsPtr

Member Function Documentation

◆ getLoopFeLhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeLhs ( )

Definition at line 898 of file SurfaceSlidingConstrains.hpp.

898 { return feLhs; }

◆ getLoopFeRhs()

MyEdgeFE& EdgeSlidingConstrains::getLoopFeRhs ( )

Definition at line 896 of file SurfaceSlidingConstrains.hpp.

896 { return feRhs; }

◆ getOptions()

MoFEMErrorCode EdgeSlidingConstrains::getOptions ( )

Definition at line 902 of file SurfaceSlidingConstrains.hpp.

902  {
904  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
905  "Get edge sliding constrains element scaling",
906  "none");
907  CHKERRQ(ierr);
908  CHKERR PetscOptionsScalar("-edge_sliding_alpha", "scaling parameter", "",
909  aLpha, &aLpha, PETSC_NULL);
910  ierr = PetscOptionsEnd();
911  CHKERRQ(ierr);
913  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
CHKERRQ(ierr)
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

1084  {
1087  edges, faces);
1088  CHKERR setOperators(tag, lagrange_multipliers_field_name,
1089  material_field_name);
1091  }
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:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
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 1093 of file SurfaceSlidingConstrains.hpp.

1096  {
1098 
1099  if (alpha) {
1100  aLpha = *alpha;
1101  }
1102 
1103  boost::shared_ptr<VectorDouble> active_variables_ptr(
1104  new VectorDouble(4 + 6));
1105  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1106  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1107  new MatrixDouble(4 + 6, 4 + 6));
1108 
1109  feRhs.getOpPtrVector().clear();
1110  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1111  lagrange_multipliers_field_name, active_variables_ptr));
1112  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1113  material_field_name, active_variables_ptr));
1114  feRhs.getOpPtrVector().push_back(new OpJacobian(
1115  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1116  jacobian_ptr, false, aLpha));
1117  feRhs.getOpPtrVector().push_back(
1118  new OpAssembleRhs<4, 6>(lagrange_multipliers_field_name, results_ptr));
1119  feRhs.getOpPtrVector().push_back(
1120  new OpAssembleRhs<4, 6>(material_field_name, results_ptr));
1121 
1122  // Adding operators to calculate the left hand side
1123  feLhs.getOpPtrVector().clear();
1124  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1125  lagrange_multipliers_field_name, active_variables_ptr));
1126  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1127  material_field_name, active_variables_ptr));
1128  feLhs.getOpPtrVector().push_back(new OpJacobian(
1129  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1130  jacobian_ptr, true, aLpha));
1131  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1132  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1133  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1134  material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
1135  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1136  material_field_name, material_field_name, jacobian_ptr));
1137 
1139  }
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:482
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:412
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 1142 of file SurfaceSlidingConstrains.hpp.

1144  {
1146 
1148  edges, faces);
1149  CHKERR setOperatorsConstrainOnly(tag, lagrange_multipliers_field_name,
1150  material_field_name);
1152  }
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:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
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:412
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 1155 of file SurfaceSlidingConstrains.hpp.

1157  {
1159 
1160  boost::shared_ptr<VectorDouble> active_variables_ptr(
1161  new VectorDouble(4 + 6));
1162  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(4 + 6));
1163  boost::shared_ptr<MatrixDouble> jacobian_ptr(
1164  new MatrixDouble(4 + 6, 4 + 6));
1165 
1166  // Adding operators to calculate the left hand side
1167  feLhs.getOpPtrVector().clear();
1168  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
1169  lagrange_multipliers_field_name, active_variables_ptr));
1170  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<4>(
1171  material_field_name, active_variables_ptr));
1172  feLhs.getOpPtrVector().push_back(new OpJacobian(
1173  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1174  jacobian_ptr, true, aLpha));
1175  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<4, 6>(
1176  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1177 
1179  }
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:482
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:412
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72

Member Data Documentation

◆ aLpha

double EdgeSlidingConstrains::aLpha

Definition at line 900 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyEdgeFE& EdgeSlidingConstrains::feLhs

Definition at line 897 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feLhsPtr

Definition at line 893 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyEdgeFE& EdgeSlidingConstrains::feRhs

Definition at line 895 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr<MyEdgeFE> EdgeSlidingConstrains::feRhsPtr

Definition at line 893 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface& EdgeSlidingConstrains::mField

Definition at line 864 of file SurfaceSlidingConstrains.hpp.


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