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

Constructor & Destructor Documentation

◆ SurfaceSlidingConstrains()

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

Definition at line 408 of file SurfaceSlidingConstrains.hpp.

410  : mField(m_field), feRhsPtr(new MyTriangleFE(m_field)),
411  feLhsPtr(new MyTriangleFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
412  crackFrontOrientation(orientation), aLpha(1) {
413 
414  ierr = getOptions();
415  CHKERRABORT(PETSC_COMM_WORLD, ierr);
416  }
DriverElementOrientation & crackFrontOrientation
boost::shared_ptr< MyTriangleFE > feLhsPtr
boost::shared_ptr< MyTriangleFE > feRhsPtr

Member Function Documentation

◆ getLoopFeLhs()

MyTriangleFE& SurfaceSlidingConstrains::getLoopFeLhs ( )

Definition at line 390 of file SurfaceSlidingConstrains.hpp.

390 { return feLhs; }

◆ getLoopFeRhs()

MyTriangleFE& SurfaceSlidingConstrains::getLoopFeRhs ( )

Definition at line 388 of file SurfaceSlidingConstrains.hpp.

388 { return feRhs; }

◆ getOptions()

MoFEMErrorCode SurfaceSlidingConstrains::getOptions ( )

Definition at line 395 of file SurfaceSlidingConstrains.hpp.

395  {
397  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
398  "Get surface sliding constrains element scaling",
399  "none");
400  CHKERRQ(ierr);
401  CHKERR PetscOptionsScalar("-surface_sliding_alpha", "scaling parameter", "",
402  aLpha, &aLpha, PETSC_NULL);
403  ierr = PetscOptionsEnd();
404  CHKERRQ(ierr);
406  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
CHKERRQ(ierr)
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

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

569  {
571 
572  if (alpha != nullptr) {
573  aLpha = *alpha;
574  }
575 
576  boost::shared_ptr<VectorDouble> active_variables_ptr(
577  new VectorDouble(3 + 9));
578  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
579  boost::shared_ptr<MatrixDouble> jacobian_ptr(
580  new MatrixDouble(3 + 9, 3 + 9));
581 
582  feRhs.getOpPtrVector().clear();
583  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
584  lagrange_multipliers_field_name, active_variables_ptr));
585  feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
586  material_field_name, active_variables_ptr));
587  feRhs.getOpPtrVector().push_back(new OpJacobian(
588  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
589  jacobian_ptr, crackFrontOrientation, false, aLpha));
590  feRhs.getOpPtrVector().push_back(
591  new OpAssembleRhs<3, 9>(lagrange_multipliers_field_name, results_ptr));
592  feRhs.getOpPtrVector().push_back(
593  new OpAssembleRhs<3, 9>(material_field_name, results_ptr));
594 
595  // Adding operators to calculate the left hand side
596  feLhs.getOpPtrVector().clear();
597  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
598  lagrange_multipliers_field_name, active_variables_ptr));
599  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
600  material_field_name, active_variables_ptr));
601  feLhs.getOpPtrVector().push_back(new OpJacobian(
602  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
603  jacobian_ptr, crackFrontOrientation, true, aLpha));
604  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
605  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
606  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
607  material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
608  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
609  material_field_name, material_field_name, jacobian_ptr));
610 
612  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:75
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
DriverElementOrientation & crackFrontOrientation
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:413
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:73

◆ setOperatorsConstrainOnly()

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

Definition at line 615 of file SurfaceSlidingConstrains.hpp.

617  {
619 
620  boost::shared_ptr<VectorDouble> active_variables_ptr(
621  new VectorDouble(3 + 9));
622  boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
623  boost::shared_ptr<MatrixDouble> jacobian_ptr(
624  new MatrixDouble(3 + 9, 3 + 9));
625 
626  // Adding operators to calculate the left hand side
627  feLhs.getOpPtrVector().clear();
628  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
629  lagrange_multipliers_field_name, active_variables_ptr));
630  feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
631  material_field_name, active_variables_ptr));
632  feLhs.getOpPtrVector().push_back(new OpJacobian(
633  tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
634  jacobian_ptr, crackFrontOrientation, true, aLpha));
635  feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
636  lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
637 
639  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:75
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
DriverElementOrientation & crackFrontOrientation
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:413
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:73

Member Data Documentation

◆ aLpha

double SurfaceSlidingConstrains::aLpha

Definition at line 394 of file SurfaceSlidingConstrains.hpp.

◆ crackFrontOrientation

DriverElementOrientation& SurfaceSlidingConstrains::crackFrontOrientation

Definition at line 392 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyTriangleFE& SurfaceSlidingConstrains::feLhs

Definition at line 389 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr<MyTriangleFE> SurfaceSlidingConstrains::feLhsPtr

Definition at line 385 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyTriangleFE& SurfaceSlidingConstrains::feRhs

Definition at line 387 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr<MyTriangleFE> SurfaceSlidingConstrains::feRhsPtr

Definition at line 385 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface& SurfaceSlidingConstrains::mField

Definition at line 356 of file SurfaceSlidingConstrains.hpp.


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