v0.8.23
Classes | Public Member Functions | Public Attributes | List of all members
SurfaceSlidingConstrains Struct Reference

Shape preserving constrains, i.e. nodes sliding on body surface. More...

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

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

Classes

struct  DriverElementOrientation
 Class implemented by user to detect face orientation. More...
 
struct  MyTriangleFE
 
struct  OpJacobian
 

Public Member Functions

MyTriangleFEgetLoopFeRhs ()
 
MyTriangleFEgetLoopFeLhs ()
 
MoFEMErrorCode getOptions ()
 
 SurfaceSlidingConstrains (MoFEM::Interface &m_field, DriverElementOrientation &orientation)
 
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, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
 

Public Attributes

MoFEM::InterfacemField
 
boost::shared_ptr< MyTriangleFEfeRhsPtr
 
boost::shared_ptr< MyTriangleFEfeLhsPtr
 
MyTriangleFEfeRhs
 
MyTriangleFEfeLhs
 
DriverElementOrientationcrackFrontOrientation
 
double aLpha
 

Detailed Description

Shape preserving constrains, i.e. nodes sliding on body surface.

Derivation and implementation of constrains preserving body surface, i.e. body shape and volume.

The idea starts form observation that body shape can be globally characterized by constant calculated as volume over its area

\[ \frac{V}{A} = C \]

Above equation expressed in integral form is

\[ \int_\Omega \textrm{d}V = C \int_\Gamma \textrm{d}S \]

where notting that,

\[ \frac{1}{3} \int_\Omega \textrm{div}[\mathbf{X}] \textrm{d}\Omega = C \int_\Gamma \textrm{d}S \]

and applying Gauss theorem we get

\[ \int_\Gamma \mathbf{X}\cdot \frac{\mathbf{N}}{\|\mathbf{N}\|} \textrm{d}\Gamma = 3C \int_\Gamma \textrm{d}S. \]

Drooping integrals on both sides, and linearizing equation, we get

\[ \frac{\mathbf{N}}{\|\mathbf{N}\|} \cdot \delta \mathbf{X} = 3C - \frac{\mathbf{N}}{\|\mathbf{N}\|}\cdot \mathbf{X} \]

where \(\delta \mathbf{X}\) is displacement sub-inctrement. Above equation is a constrain if satisfied in body shape and volume is conserved. Final form of constrain equation is

\[ \mathcal{r} = \frac{\mathbf{N}}{\|\mathbf{N}\|}\cdot \mathbf{X} - \frac{\mathbf{N_0}}{\|\mathbf{N_0}\|}\cdot \mathbf{X}_0 = \frac{\mathbf{N}}{\|\mathbf{N}\|}\cdot (\mathbf{X}-\mathbf{X}_0) \]

In the spirit of finite element method the constrain equation is multiplied by shape functions and enforce using Lagrange multiplier method

\[ \int_\Gamma \mathbf{N}^\mathsf{T}_\lambda \left( \frac{\mathbf{N}}{\|\mathbf{N}\|}\mathbf{N}_\mathbf{X}\cdot (\overline{\mathbf{X}}-\overline{\mathbf{X}}_0) \right) \|\mathbf{N}\| \textrm{d}\Gamma = \mathbf{0}. \]

Above equation is nonlinear, applying to it Taylor expansion, we can get form which can be used with Newton interactive method

\[ \begin{split} &\int_\Gamma \mathbf{N}^\mathsf{T}_\lambda \left\{ \mathbf{N}\mathbf{N}_\mathbf{X} + \left(\mathbf{X}-\mathbf{X}_0\right) \cdot \left( \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\xi}\right]\cdot\mathbf{B}_\eta - \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\eta}\right]\cdot\mathbf{B}_\xi \right) \right\} \textrm{d}\Gamma \cdot \delta \overline{\mathbf{X}}\\ = &\int_\Gamma \mathbf{N}^\mathsf{T}_\lambda \mathbf{N}\cdot(\mathbf{X}-\mathbf{X}_0) \textrm{d}\Gamma \end{split}. \]

Equation expressing forces on shape as result of constrains, as result Lagrange multiplier method have following form

\[ \begin{split} &\int_\Gamma \mathbf{N}^\mathsf{T}_\mathbf{X} \cdot \mathbf{N} \mathbf{N}_\lambda \textrm{d}\Gamma \cdot \delta\overline{\lambda}\\ + &\int_\Gamma \lambda \mathbf{N}^\mathsf{T}_\mathbf{X} \left( \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\xi}\right]\cdot\mathbf{B}_\eta - \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\eta}\right]\cdot\mathbf{B}_\xi \right) \textrm{d}\Gamma \delta\overline{\mathbf{X}}\\ = &\int_\Gamma \lambda \mathbf{N}^\mathsf{T}_\mathbf{X} \cdot \mathbf{N} \textrm{d}\Gamma \end{split} \]

Above equations are assembled into global system of equations as following

\[ \left[ \begin{array}{cc} \mathbf{K} + \mathbf{B} & \mathbf{C}^\mathsf{T} \\ \mathbf{C} + \mathbf{A} & 0 \end{array} \right] \left\{ \begin{array}{c} \delta \overline{\mathbf{X}} \\ \delta \overline{\lambda} \end{array} \right\}= \left[ \begin{array}{c} \mathbf{f} - \mathbf{C}^\mathsf{T}\overline{\lambda} \\ \overline{\mathbf{r}} \end{array} \right] \]

where

\[ \mathbf{C}= \int_\Gamma \mathbf{N}_\lambda^\mathsf{T} \mathbf{N} \cdot \mathbf{N}_\mathbf{X} \textrm{d}\Gamma, \]

\[ \mathbf{B}= \int_\Gamma \lambda (\mathbf{X}-\mathbf{X}_0)\cdot \left( \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\xi}\right]\cdot\mathbf{B}_\eta - \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\eta}\right]\cdot\mathbf{B}_\xi \right) \textrm{d}\Gamma \]

and

\[ \mathbf{A}= \int_\Gamma \mathbf{N}^\mathsf{T}_\lambda \left(\mathbf{X}-\mathbf{X}_0\right) \cdot \left( \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\xi}\right]\cdot\mathbf{B}_\eta - \textrm{Spin}\left[\frac{\partial\mathbf{X}}{\partial\eta}\right]\cdot\mathbf{B}_\xi \right). \]

Examples
mesh_smoothing.cpp.

Definition at line 328 of file SurfaceSlidingConstrains.hpp.

Constructor & Destructor Documentation

◆ SurfaceSlidingConstrains()

SurfaceSlidingConstrains::SurfaceSlidingConstrains ( MoFEM::Interface m_field,
DriverElementOrientation orientation 
)

Definition at line 423 of file SurfaceSlidingConstrains.hpp.

425  : mField(m_field), feRhsPtr(new MyTriangleFE(m_field)),
426  feLhsPtr(new MyTriangleFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
427  crackFrontOrientation(orientation), aLpha(1) {
428 
429  ierr = getOptions();
430  CHKERRABORT(PETSC_COMM_WORLD, ierr);
431  }
DriverElementOrientation & crackFrontOrientation
boost::shared_ptr< MyTriangleFE > feLhsPtr
boost::shared_ptr< MyTriangleFE > feRhsPtr

Member Function Documentation

◆ getLoopFeLhs()

MyTriangleFE& SurfaceSlidingConstrains::getLoopFeLhs ( )

Definition at line 404 of file SurfaceSlidingConstrains.hpp.

404 { return feLhs; }

◆ getLoopFeRhs()

MyTriangleFE& SurfaceSlidingConstrains::getLoopFeRhs ( )

Definition at line 402 of file SurfaceSlidingConstrains.hpp.

402 { return feRhs; }

◆ getOptions()

MoFEMErrorCode SurfaceSlidingConstrains::getOptions ( )

Definition at line 410 of file SurfaceSlidingConstrains.hpp.

410  {
412  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
413  "Get surface sliding constrains element scaling",
414  "none");
415  CHKERRQ(ierr);
416  CHKERR PetscOptionsScalar("-surface_sliding_alpha", "scaling parameter", "",
417  aLpha, &aLpha, PETSC_NULL);
418  ierr = PetscOptionsEnd();
419  CHKERRQ(ierr);
421  }
#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()

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

Definition at line 581 of file SurfaceSlidingConstrains.hpp.

584  {
586 
587  if (alpha != nullptr) {
588  aLpha = *alpha;
589  }
590 
591  boost::shared_ptr<VectorDouble> active_variables_ptr(new VectorDouble(3+9));
592  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3+9));
593  boost::shared_ptr<MatrixDouble> jacobian_ptr(new MatrixDouble(3+9,3+9));
594 
595  feRhs.getOpPtrVector().clear();
596  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
597  lagrange_multipliers_field_name, active_variables_ptr));
598  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
599  material_field_name, active_variables_ptr));
600  feRhs.getOpPtrVector().push_back(new OpJacobian(
601  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
602  jacobian_ptr, crackFrontOrientation, false, aLpha));
603  feRhs.getOpPtrVector().push_back(
604  new OpAssembleRhs<3, 9>(lagrange_multipliers_field_name, results_ptr));
605  feRhs.getOpPtrVector().push_back(
606  new OpAssembleRhs<3, 9>(material_field_name, results_ptr));
607 
608  // Adding operators to calculate the left hand side
609  feLhs.getOpPtrVector().clear();
610  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
611  lagrange_multipliers_field_name, active_variables_ptr));
612  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
613  material_field_name, active_variables_ptr));
614  feLhs.getOpPtrVector().push_back(new OpJacobian(
615  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
616  jacobian_ptr, crackFrontOrientation, true, aLpha));
617  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
618  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
619  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
620  material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
621  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
622  material_field_name, material_field_name, jacobian_ptr));
623 
625  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
DriverElementOrientation & crackFrontOrientation
#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:76

◆ setOperatorsConstrainOnly()

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

Definition at line 628 of file SurfaceSlidingConstrains.hpp.

630  {
632 
633  boost::shared_ptr<VectorDouble> active_variables_ptr(
634  new VectorDouble(3 + 9));
635  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
636  boost::shared_ptr<MatrixDouble> jacobian_ptr(
637  new MatrixDouble(3 + 9, 3 + 9));
638 
639  // Adding operators to calculate the left hand side
640  feLhs.getOpPtrVector().clear();
641  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
642  lagrange_multipliers_field_name, active_variables_ptr));
643  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
644  material_field_name, active_variables_ptr));
645  feLhs.getOpPtrVector().push_back(new OpJacobian(
646  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
647  jacobian_ptr, crackFrontOrientation, true, aLpha));
648  feLhs.getOpPtrVector().push_back(
649  new OpAssembleLhs<3,9>(lagrange_multipliers_field_name, material_field_name,
650  jacobian_ptr));
651 
653  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
DriverElementOrientation & crackFrontOrientation
#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:76

Member Data Documentation

◆ aLpha

double SurfaceSlidingConstrains::aLpha

Definition at line 409 of file SurfaceSlidingConstrains.hpp.

◆ crackFrontOrientation

DriverElementOrientation& SurfaceSlidingConstrains::crackFrontOrientation

Definition at line 407 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyTriangleFE& SurfaceSlidingConstrains::feLhs

Definition at line 403 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr<MyTriangleFE> SurfaceSlidingConstrains::feLhsPtr

Definition at line 399 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyTriangleFE& SurfaceSlidingConstrains::feRhs

Definition at line 401 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr<MyTriangleFE> SurfaceSlidingConstrains::feRhsPtr

Definition at line 399 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface& SurfaceSlidingConstrains::mField

Definition at line 352 of file SurfaceSlidingConstrains.hpp.


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