v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | Private Attributes | List of all members
SphericalArcLengthControl Struct Reference

Implementation of spherical arc-length method. More...

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

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

Public Member Functions

DEPRECATED SphericalArcLengthControl (ArcLengthCtx *arc_ptr_raw)
 
 SphericalArcLengthControl (boost::shared_ptr< ArcLengthCtx > &arc_ptr)
 
virtual ~SphericalArcLengthControl ()
 
MoFEMErrorCode preProcess ()
 
MoFEMErrorCode operator() ()
 
MoFEMErrorCode postProcess ()
 
virtual double calculateLambdaInt ()
 Calculate f_lambda(dx,lambda) More...
 
virtual MoFEMErrorCode calculateDb ()
 Calculate db. More...
 
virtual MoFEMErrorCode calculateDxAndDlambda (Vec x)
 
virtual MoFEMErrorCode calculateInitDlambda (double *dlambda)
 
virtual MoFEMErrorCode setDlambdaToX (Vec x, double dlambda)
 

Public Attributes

ArcLengthCtxarcPtrRaw
 

Private Attributes

boost::shared_ptr< ArcLengthCtxarcPtr
 

Detailed Description

Implementation of spherical arc-length method.

\[ \alpha \| \Delta\mathbf{x} \|^2 + \Delta\lambda^2 \beta^2 \| \mathbf{F}_\lambda \|^2 = s^2 \]

This is particular implementation of ArcLength control, i.e. spherical arc length control. If beta is set to 0 and alpha is non-zero it is cylindrical arc-length control. Works well with general problem with non-linear geometry. It not guarantee dissipative loading path in case of physical nonlinearities.

Definition at line 374 of file ArcLengthTools.hpp.

Constructor & Destructor Documentation

◆ SphericalArcLengthControl() [1/2]

SphericalArcLengthControl::SphericalArcLengthControl ( ArcLengthCtx arc_ptr_raw)
Deprecated:
use constructor with shared_ptr

Definition at line 661 of file ArcLengthTools.cpp.

662 : FEMethod(), arcPtrRaw(arc_ptr_raw) {}
structure for User Loop Methods on finite elements

◆ SphericalArcLengthControl() [2/2]

SphericalArcLengthControl::SphericalArcLengthControl ( boost::shared_ptr< ArcLengthCtx > &  arc_ptr)

Definition at line 664 of file ArcLengthTools.cpp.

666 : FEMethod(), arcPtrRaw(arc_ptr.get()), arcPtr(arc_ptr) {}
boost::shared_ptr< ArcLengthCtx > arcPtr

◆ ~SphericalArcLengthControl()

SphericalArcLengthControl::~SphericalArcLengthControl ( )
virtual

Definition at line 668 of file ArcLengthTools.cpp.

668{}

Member Function Documentation

◆ calculateDb()

MoFEMErrorCode SphericalArcLengthControl::calculateDb ( )
virtual

Calculate db.

\[ \textrm{d}\mathbf{b} = 2 \alpha \Delta\mathbf{x} \]

Reimplemented in ShellArcLength.

Definition at line 691 of file ArcLengthTools.cpp.

691 {
693 CHKERR VecCopy(arcPtrRaw->dx, arcPtrRaw->db);
694 CHKERR VecScale(arcPtrRaw->db, 2 * arcPtrRaw->alpha);
696}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
SmartPetscObj< Vec > db
db derivative of f(dx*dx), i.e. db = d[ f(dx*dx) ]/dx
SmartPetscObj< Vec > dx
dx = x-x0
double alpha
displacement scaling factor

◆ calculateDxAndDlambda()

MoFEMErrorCode SphericalArcLengthControl::calculateDxAndDlambda ( Vec  x)
virtual

Definition at line 742 of file ArcLengthTools.cpp.

742 {
744 // dx
745 CHKERR VecCopy(x, arcPtrRaw->dx);
746 CHKERR VecAXPY(arcPtrRaw->dx, -1, arcPtrRaw->x0);
747 CHKERR VecGhostUpdateBegin(arcPtrRaw->dx, INSERT_VALUES, SCATTER_FORWARD);
748 CHKERR VecGhostUpdateEnd(arcPtrRaw->dx, INSERT_VALUES, SCATTER_FORWARD);
749 // dlambda
750 if (arcPtrRaw->getPetscLocalDofIdx() != -1) {
751 double *array;
752 CHKERR VecGetArray(arcPtrRaw->dx, &array);
754 array[arcPtrRaw->getPetscLocalDofIdx()] = 0;
755 CHKERR VecRestoreArray(arcPtrRaw->dx, &array);
756 }
757 CHKERR VecGhostUpdateBegin(arcPtrRaw->ghosTdLambda, INSERT_VALUES,
758 SCATTER_FORWARD);
759 CHKERR VecGhostUpdateEnd(arcPtrRaw->ghosTdLambda, INSERT_VALUES,
760 SCATTER_FORWARD);
761 // dx2
763 MOFEM_LOG_C("WORLD", Sev::verbose, "\tdlambda = %6.4e dx2 = %6.4e\n",
766}
#define ArcFunctionBegin
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
DofIdx getPetscLocalDofIdx()
Get local index of load factor.
double dx2
inner_prod(dX,dX)
SmartPetscObj< Vec > x0
displacement vector at beginning of step
double dLambda
increment of load factor
SmartPetscObj< Vec > ghosTdLambda

◆ calculateInitDlambda()

MoFEMErrorCode SphericalArcLengthControl::calculateInitDlambda ( double dlambda)
virtual

Definition at line 769 of file ArcLengthTools.cpp.

769 {
771 *dlambda = std::sqrt(pow(arcPtrRaw->s, 2) /
772 (pow(arcPtrRaw->beta, 2) * arcPtrRaw->F_lambda2));
773 if (!(*dlambda == *dlambda)) {
774 MOFEM_LOG("WORLD", Sev::error)
775 << "s " << arcPtrRaw->s << " " << arcPtrRaw->beta << " "
777 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE,
778 "Increment of lambda is not a number");
779 }
781}
@ MOFEM_IMPOSSIBLE_CASE
Definition: definitions.h:35
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
double s
arc length radius
double F_lambda2
inner_prod(F_lambda,F_lambda);
double beta
force scaling factor

◆ calculateLambdaInt()

double SphericalArcLengthControl::calculateLambdaInt ( )
virtual

Calculate f_lambda(dx,lambda)

\[ f_\lambda(\Delta\mathbf{x},\lambda) = \alpha \| \Delta\mathbf{x} \|^2 + \Delta\lambda^2 \beta^2 \| \mathbf{F}_\lambda \|^2 \]

Reimplemented in ShellArcLength.

Definition at line 685 of file ArcLengthTools.cpp.

685 {
686 return arcPtrRaw->alpha * arcPtrRaw->dx2 + pow(arcPtrRaw->dLambda, 2) *
687 pow(arcPtrRaw->beta, 2) *
689}

◆ operator()()

MoFEMErrorCode SphericalArcLengthControl::operator() ( )

Definition at line 698 of file ArcLengthTools.cpp.

698 {
700 switch (snes_ctx) {
701 case CTX_SNESSETFUNCTION: {
703 CHKERR VecSetValue(snes_f, arcPtrRaw->getPetscGlobalDofIdx(),
704 arcPtrRaw->res_lambda, ADD_VALUES);
705 MOFEM_LOG_C("WORLD", Sev::verbose, "\tres_lambda = %6.4e\n",
707 } break;
708 case CTX_SNESSETJACOBIAN: {
709 arcPtrRaw->dIag =
711 CHKERR MatSetValue(snes_B, arcPtrRaw->getPetscGlobalDofIdx(),
712 arcPtrRaw->getPetscGlobalDofIdx(), 1, ADD_VALUES);
713 } break;
714 default:
715 break;
716 }
718}
double res_lambda
f_lambda - s
double dIag
diagonal value
DofIdx getPetscGlobalDofIdx()
Get global index of load factor.
virtual double calculateLambdaInt()
Calculate f_lambda(dx,lambda)

◆ postProcess()

MoFEMErrorCode SphericalArcLengthControl::postProcess ( )

Definition at line 720 of file ArcLengthTools.cpp.

720 {
722 switch (snes_ctx) {
723 case CTX_SNESSETFUNCTION: {
724 MOFEM_LOG_C("WORLD", Sev::verbose, "\tlambda = %6.4e\n",
726 } break;
727 case CTX_SNESSETJACOBIAN: {
728 CHKERR VecGhostUpdateBegin(arcPtrRaw->ghostDiag, INSERT_VALUES,
729 SCATTER_FORWARD);
730 CHKERR VecGhostUpdateEnd(arcPtrRaw->ghostDiag, INSERT_VALUES,
731 SCATTER_FORWARD);
732 CHKERR MatAssemblyBegin(snes_B, MAT_FLUSH_ASSEMBLY);
733 CHKERR MatAssemblyEnd(snes_B, MAT_FLUSH_ASSEMBLY);
734 MOFEM_LOG_C("WORLD", Sev::verbose, "\tdiag = %6.4e", arcPtrRaw->dIag);
735 } break;
736 default:
737 break;
738 }
740}
FieldData & getFieldData()
Get value of load factor.
SmartPetscObj< Vec > ghostDiag

◆ preProcess()

MoFEMErrorCode SphericalArcLengthControl::preProcess ( )

Definition at line 670 of file ArcLengthTools.cpp.

670 {
672 switch (snes_ctx) {
673 case CTX_SNESSETFUNCTION: {
676 } break;
677 case CTX_SNESSETJACOBIAN: {
678 } break;
679 default:
680 break;
681 }
683}
virtual MoFEMErrorCode calculateDxAndDlambda(Vec x)
virtual MoFEMErrorCode calculateDb()
Calculate db.

◆ setDlambdaToX()

MoFEMErrorCode SphericalArcLengthControl::setDlambdaToX ( Vec  x,
double  dlambda 
)
virtual

Definition at line 783 of file ArcLengthTools.cpp.

783 {
785 // check if local dof idx is non zero, i.e. that lambda is accessible from
786 // this processor
787 if (arcPtrRaw->getPetscLocalDofIdx() != -1) {
788 double *array;
789 CHKERR VecGetArray(x, &array);
790 double lambda_old = array[arcPtrRaw->getPetscLocalDofIdx()];
791 if (!(dlambda == dlambda)) {
792 MOFEM_LOG("WORLD", Sev::error)
793 << "s " << arcPtrRaw->s << " " << arcPtrRaw->beta << " "
795 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE,
796 "Increment of lambda is not a number");
797 }
798 array[arcPtrRaw->getPetscLocalDofIdx()] = lambda_old + dlambda;
799 MOFEM_LOG_C("WORLD", Sev::verbose, "\tlambda = %6.4e, %6.4e (%6.4e)\n",
800 lambda_old, array[arcPtrRaw->getPetscLocalDofIdx()], dlambda);
801 CHKERR VecRestoreArray(x, &array);
802 }
804}

Member Data Documentation

◆ arcPtr

boost::shared_ptr<ArcLengthCtx> SphericalArcLengthControl::arcPtr
private

Definition at line 412 of file ArcLengthTools.hpp.

◆ arcPtrRaw

ArcLengthCtx* SphericalArcLengthControl::arcPtrRaw

Definition at line 376 of file ArcLengthTools.hpp.


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