v0.15.0
Loading...
Searching...
No Matches
SurfaceSlidingConstrains Struct Reference

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

#include "tools/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)
 
virtual ~SurfaceSlidingConstrains ()=default
 
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)
 
MyTriangleFEgetLoopFeRhs ()
 
MyTriangleFEgetLoopFeLhs ()
 
MoFEMErrorCode getOptions ()
 
 SurfaceSlidingConstrains (MoFEM::Interface &m_field, DriverElementOrientation &orientation)
 
virtual ~SurfaceSlidingConstrains ()=default
 
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 Member Functions inherited from GenericSliding
 GenericSliding ()=default
 
virtual ~GenericSliding ()=default
 
 GenericSliding ()=default
 
virtual ~GenericSliding ()=default
 

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

Constructor & Destructor Documentation

◆ SurfaceSlidingConstrains() [1/2]

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

Definition at line 403 of file SurfaceSlidingConstrains.hpp.

405 : mField(m_field), feRhsPtr(new MyTriangleFE(m_field)),
406 feLhsPtr(new MyTriangleFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
407 crackFrontOrientation(orientation), aLpha(1) {
408
409 ierr = getOptions();
410 CHKERRABORT(PETSC_COMM_WORLD, ierr);
411 }
static PetscErrorCode ierr
boost::shared_ptr< MyTriangleFE > feRhsPtr
DriverElementOrientation & crackFrontOrientation
boost::shared_ptr< MyTriangleFE > feLhsPtr

◆ ~SurfaceSlidingConstrains() [1/2]

virtual SurfaceSlidingConstrains::~SurfaceSlidingConstrains ( )
virtualdefault

◆ SurfaceSlidingConstrains() [2/2]

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

Definition at line 402 of file SurfaceSlidingConstrains.hpp.

404 : mField(m_field), feRhsPtr(new MyTriangleFE(m_field)),
405 feLhsPtr(new MyTriangleFE(m_field)), feRhs(*feRhsPtr), feLhs(*feLhsPtr),
406 crackFrontOrientation(orientation), aLpha(1) {
407
408 ierr = getOptions();
409 CHKERRABORT(PETSC_COMM_WORLD, ierr);
410 }

◆ ~SurfaceSlidingConstrains() [2/2]

virtual SurfaceSlidingConstrains::~SurfaceSlidingConstrains ( )
virtualdefault

Member Function Documentation

◆ getLoopFeLhs() [1/2]

MyTriangleFE & SurfaceSlidingConstrains::getLoopFeLhs ( )
inline

Definition at line 387 of file SurfaceSlidingConstrains.hpp.

387{ return feLhs; }

◆ getLoopFeLhs() [2/2]

MyTriangleFE & SurfaceSlidingConstrains::getLoopFeLhs ( )
inline

Definition at line 387 of file SurfaceSlidingConstrains.hpp.

387{ return feLhs; }

◆ getLoopFeRhs() [1/2]

MyTriangleFE & SurfaceSlidingConstrains::getLoopFeRhs ( )
inline

Definition at line 385 of file SurfaceSlidingConstrains.hpp.

385{ return feRhs; }

◆ getLoopFeRhs() [2/2]

MyTriangleFE & SurfaceSlidingConstrains::getLoopFeRhs ( )
inline

Definition at line 385 of file SurfaceSlidingConstrains.hpp.

385{ return feRhs; }

◆ getOptions() [1/2]

MoFEMErrorCode SurfaceSlidingConstrains::getOptions ( )
inline

Definition at line 392 of file SurfaceSlidingConstrains.hpp.

392 {
394 PetscOptionsBegin(PETSC_COMM_WORLD, "",
395 "Get surface sliding constrains element scaling",
396 "none");
397 CHKERR PetscOptionsScalar("-surface_sliding_alpha", "scaling parameter", "",
398 aLpha, &aLpha, PETSC_NULLPTR);
399 PetscOptionsEnd();
401 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.

◆ getOptions() [2/2]

MoFEMErrorCode SurfaceSlidingConstrains::getOptions ( )
inline

Definition at line 392 of file SurfaceSlidingConstrains.hpp.

392 {
394 PetscOptionsBegin(PETSC_COMM_WORLD, "",
395 "Get surface sliding constrains element scaling", "none");
396 CHKERR PetscOptionsScalar("-surface_sliding_alpha", "scaling parameter", "",
397 aLpha, &aLpha, PETSC_NULLPTR);
398 PetscOptionsEnd();
400 }

◆ setOperators() [1/2]

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

Definition at line 561 of file SurfaceSlidingConstrains.hpp.

564 {
566
567 if (alpha != nullptr) {
568 aLpha = *alpha;
569 }
570
571 boost::shared_ptr<VectorDouble> active_variables_ptr(
572 new VectorDouble(3 + 9));
573 boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
574 boost::shared_ptr<MatrixDouble> jacobian_ptr(
575 new MatrixDouble(3 + 9, 3 + 9));
576
577 feRhs.getOpPtrVector().clear();
578 feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
579 lagrange_multipliers_field_name, active_variables_ptr));
580 feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
581 material_field_name, active_variables_ptr));
582 feRhs.getOpPtrVector().push_back(new OpJacobian(
583 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
584 jacobian_ptr, crackFrontOrientation, false, aLpha));
585 feRhs.getOpPtrVector().push_back(
586 new OpAssembleRhs<3, 9>(lagrange_multipliers_field_name, results_ptr));
587 feRhs.getOpPtrVector().push_back(
588 new OpAssembleRhs<3, 9>(material_field_name, results_ptr));
589
590 // Adding operators to calculate the left hand side
591 feLhs.getOpPtrVector().clear();
592 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
593 lagrange_multipliers_field_name, active_variables_ptr));
594 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
595 material_field_name, active_variables_ptr));
596 feLhs.getOpPtrVector().push_back(new OpJacobian(
597 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
598 jacobian_ptr, crackFrontOrientation, true, aLpha));
599 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
600 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
601 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
602 material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
603 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
604 material_field_name, material_field_name, jacobian_ptr));
605
607 }
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.

◆ setOperators() [2/2]

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

Definition at line 560 of file SurfaceSlidingConstrains.hpp.

563 {
565
566 if (alpha != nullptr) {
567 aLpha = *alpha;
568 }
569
570 boost::shared_ptr<VectorDouble> active_variables_ptr(
571 new VectorDouble(3 + 9));
572 boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
573 boost::shared_ptr<MatrixDouble> jacobian_ptr(
574 new MatrixDouble(3 + 9, 3 + 9));
575
576 feRhs.getOpPtrVector().clear();
577 feRhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
578 lagrange_multipliers_field_name, active_variables_ptr));
579 feRhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
580 material_field_name, active_variables_ptr));
581 feRhs.getOpPtrVector().push_back(new OpJacobian(
582 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
583 jacobian_ptr, crackFrontOrientation, false, aLpha));
584 feRhs.getOpPtrVector().push_back(
585 new OpAssembleRhs<3, 9>(lagrange_multipliers_field_name, results_ptr));
586 feRhs.getOpPtrVector().push_back(
587 new OpAssembleRhs<3, 9>(material_field_name, results_ptr));
588
589 // Adding operators to calculate the left hand side
590 feLhs.getOpPtrVector().clear();
591 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
592 lagrange_multipliers_field_name, active_variables_ptr));
593 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
594 material_field_name, active_variables_ptr));
595 feLhs.getOpPtrVector().push_back(new OpJacobian(
596 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
597 jacobian_ptr, crackFrontOrientation, true, aLpha));
598 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
599 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
600 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
601 material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
602 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
603 material_field_name, material_field_name, jacobian_ptr));
604
606 }

◆ setOperatorsConstrainOnly() [1/2]

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

Definition at line 610 of file SurfaceSlidingConstrains.hpp.

612 {
614
615 boost::shared_ptr<VectorDouble> active_variables_ptr(
616 new VectorDouble(3 + 9));
617 boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
618 boost::shared_ptr<MatrixDouble> jacobian_ptr(
619 new MatrixDouble(3 + 9, 3 + 9));
620
621 // Adding operators to calculate the left hand side
622 feLhs.getOpPtrVector().clear();
623 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
624 lagrange_multipliers_field_name, active_variables_ptr));
625 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
626 material_field_name, active_variables_ptr));
627 feLhs.getOpPtrVector().push_back(new OpJacobian(
628 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
629 jacobian_ptr, crackFrontOrientation, true, aLpha));
630 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
631 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
632
634 }

◆ setOperatorsConstrainOnly() [2/2]

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

Definition at line 609 of file SurfaceSlidingConstrains.hpp.

611 {
613
614 boost::shared_ptr<VectorDouble> active_variables_ptr(
615 new VectorDouble(3 + 9));
616 boost::shared_ptr<VectorDouble> results_ptr(new VectorDouble(3 + 9));
617 boost::shared_ptr<MatrixDouble> jacobian_ptr(
618 new MatrixDouble(3 + 9, 3 + 9));
619
620 // Adding operators to calculate the left hand side
621 feLhs.getOpPtrVector().clear();
622 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsLambda(
623 lagrange_multipliers_field_name, active_variables_ptr));
624 feLhs.getOpPtrVector().push_back(new OpGetActiveDofsPositions<3>(
625 material_field_name, active_variables_ptr));
626 feLhs.getOpPtrVector().push_back(new OpJacobian(
627 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
628 jacobian_ptr, crackFrontOrientation, true, aLpha));
629 feLhs.getOpPtrVector().push_back(new OpAssembleLhs<3, 9>(
630 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
631
633 }

Member Data Documentation

◆ aLpha

double SurfaceSlidingConstrains::aLpha

Definition at line 391 of file SurfaceSlidingConstrains.hpp.

◆ crackFrontOrientation

DriverElementOrientation & SurfaceSlidingConstrains::crackFrontOrientation

Definition at line 389 of file SurfaceSlidingConstrains.hpp.

◆ feLhs

MyTriangleFE & SurfaceSlidingConstrains::feLhs

Definition at line 386 of file SurfaceSlidingConstrains.hpp.

◆ feLhsPtr

boost::shared_ptr< MyTriangleFE > SurfaceSlidingConstrains::feLhsPtr

Definition at line 382 of file SurfaceSlidingConstrains.hpp.

◆ feRhs

MyTriangleFE & SurfaceSlidingConstrains::feRhs

Definition at line 384 of file SurfaceSlidingConstrains.hpp.

◆ feRhsPtr

boost::shared_ptr< MyTriangleFE > SurfaceSlidingConstrains::feRhsPtr

Definition at line 382 of file SurfaceSlidingConstrains.hpp.

◆ mField

MoFEM::Interface & SurfaceSlidingConstrains::mField

Definition at line 353 of file SurfaceSlidingConstrains.hpp.


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