FrictionQPotFEM 0.23.3
|
Class that uses GMatElastoPlasticQPot to evaluate stress everywhere. More...
#include <FrictionQPotFEM/UniformSingleLayer2d.h>
Public Member Functions | |
template<class C , class E , class L , class M , class Y > | |
System (const C &coor, const E &conn, const E &dofs, const L &iip, const L &elastic_elem, const M &elastic_K, const M &elastic_G, const L &plastic_elem, const M &plastic_K, const M &plastic_G, const Y &plastic_epsy, double dt, double rho, double alpha, double eta) | |
Define the geometry, including boundary conditions and element sets. More... | |
std::string | type () const override |
Return the type of system. More... | |
double | typical_plastic_h () const |
Element height of all elements along the weak layer. More... | |
double | typical_plastic_dV () const |
Integration point volume of all elements along the weak layer. More... | |
void | initEventDrivenSimpleShear () |
Initialise event driven protocol for affine simple shear. More... | |
double | addSimpleShearToFixedStress (double target_stress, bool dry_run=false) |
Add simple shear until a target equivalent macroscopic stress has been reached. More... | |
double | addElasticSimpleShearToFixedStress (double target_stress, bool dry_run=false) |
Add simple shear until a target equivalent macroscopic stress has been reached. More... | |
double | triggerElementWithLocalSimpleShear (double deps_kick, size_t plastic_element, bool trigger_weakest=true, double amplify=1.0) |
Apply local strain to the right to a specific plastic element. More... | |
array_type::tensor< double, 2 > | plastic_ElementYieldBarrierForSimpleShear (double deps_kick=0.0, size_t iquad=0) |
Read the distance to overcome the first cusp in the element. More... | |
![]() | |
template<class C , class E , class L , class M , class Y > | |
System (const C &coor, const E &conn, const E &dofs, const L &iip, const L &elastic_elem, const M &elastic_K, const M &elastic_G, const L &plastic_elem, const M &plastic_K, const M &plastic_G, const Y &plastic_epsy, double dt, double rho, double alpha, double eta) | |
Define the geometry, including boundary conditions and element sets. More... | |
virtual size_t | N () const |
Return the linear system size (in number of blocks). More... | |
virtual std::string | type () const |
Return the type of system. More... | |
double | rho () const |
Mass density. More... | |
void | setRho (double rho) |
Overwrite the mass density to a homogeneous quantity. More... | |
template<class T > | |
void | setMassMatrix (const T &val_elem) |
Overwrite mass matrix, based on certain density that is uniform per element. More... | |
void | setEta (double eta) |
Overwrite the value of the damping at the interface. More... | |
double | eta () const |
Get the damping at the interface. More... | |
void | setAlpha (double alpha) |
Overwrite background damping density (proportional to the velocity), To use a non-homogeneous density use setDampingMatrix(). More... | |
double | alpha () const |
Background damping density. More... | |
template<class T > | |
void | setDampingMatrix (const T &val_elem) |
Overwrite damping matrix, based on certain density (taken uniform per element). More... | |
bool | isHomogeneousElastic () const |
Check if elasticity is homogeneous. More... | |
double | dt () const |
Get time step. More... | |
void | setDt (double dt) |
Overwrite time step. More... | |
template<class T > | |
void | setU (const T &u) |
Overwrite nodal displacements. More... | |
template<class T > | |
void | setV (const T &v) |
Overwrite nodal velocities. More... | |
template<class T > | |
void | setA (const T &a) |
Overwrite nodal accelerations. More... | |
virtual void | updated_u () |
Evaluate relevant forces when m_u is updated. More... | |
void | updated_v () |
Evaluate relevant forces when m_v is updated. More... | |
template<class T > | |
void | setFext (const T &fext) |
Overwrite external force. More... | |
const auto & | elastic_elem () const |
List of elastic elements. More... | |
const auto & | plastic_elem () const |
List of plastic elements. More... | |
const auto & | conn () const |
Connectivity. More... | |
const auto & | coor () const |
Nodal coordinates. More... | |
const auto & | dofs () const |
DOFs per node. More... | |
const auto & | u () const |
Nodal displacements. More... | |
const auto & | v () const |
Nodal velocities. More... | |
const auto & | a () const |
Nodal accelerations. More... | |
auto & | mass () const |
Mass matrix, see System::m_M. More... | |
auto & | damping () const |
Damping matrix, see setAlpha() and setDampingMatrix(). More... | |
const auto & | fext () |
External force. More... | |
void | applyShearStress (double sigxy) |
Modify the external force such that a shear stress is applied to the top boundary. More... | |
const auto & | fint () |
Internal force. More... | |
const auto & | fmaterial () const |
Force deriving from elasticity. More... | |
const auto & | fdamp () const |
Force deriving from damping. More... | |
double | residual () const |
Norm of the relative residual force (the external force at the top/bottom boundaries is used for normalisation). More... | |
void | quench () |
Set nodal velocities and accelerations equal to zero. More... | |
double | t () const |
Current time. More... | |
void | setT (double arg) |
Overwrite time. More... | |
size_t | inc () const |
The increment, see setInc(). More... | |
virtual void | setInc (size_t arg) |
Set increment. More... | |
const auto & | dV () const |
Integration point volume (of each integration point) More... | |
const GooseFEM::VectorPartitioned & | vector () const |
GooseFEM vector definition. More... | |
const GooseFEM::Element::Quad4::Quadrature & | quad () const |
GooseFEM quadrature definition. More... | |
GMatElastoPlasticQPot::Cartesian2d::Elastic< 2 > & | elastic () |
Elastic material mode, used for all elastic elements. More... | |
GMatElastoPlasticQPot::Cartesian2d::Cusp< 2 > & | plastic () |
Elastic material mode, used for all elastic elements. More... | |
array_type::tensor< double, 2 > | K () const |
Bulk modulus per integration point. More... | |
array_type::tensor< double, 2 > | G () const |
Shear modulus per integration point. More... | |
const array_type::tensor< double, 4 > & | Sig () |
Stress tensor of each integration point. More... | |
const array_type::tensor< double, 4 > & | Eps () |
Strain tensor of each integration point. More... | |
array_type::tensor< double, 4 > | Epsdot () const |
Strain-rate tensor (the symmetric gradient of the nodal velocities) of each integration point. More... | |
array_type::tensor< double, 4 > | Epsddot () const |
Symmetric gradient of the nodal accelerations of each integration point. More... | |
GooseFEM::MatrixPartitioned | stiffness () const |
Stiffness tensor of each integration point. More... | |
const array_type::tensor< double, 4 > & | plastic_Epsdot () |
Strain-rate tensor of integration points of plastic elements only, see System::plastic. More... | |
void | refresh () |
Recompute all forces. More... | |
template<class T > | |
double | eventDriven_setDeltaU (const T &delta_u, bool autoscale=true) |
Set purely elastic and linear response to specific boundary conditions. More... | |
const auto & | eventDriven_deltaU () const |
Get displacement perturbation used for event driven code, see eventDriven_setDeltaU(). More... | |
virtual double | eventDrivenStep (double deps, bool kick, int direction=1, bool yield_element=false, bool iterative=false) |
Add event driven step for the boundary conditions that correspond to the displacement perturbation set in eventDriven_setDeltaU(). More... | |
void | timeStep () |
Make a time-step: apply velocity-Verlet integration. More... | |
long | timeSteps (size_t n, size_t nmargin=0) |
Make a number of time steps, see timeStep(). More... | |
size_t | timeStepsUntilEvent (size_t nmargin=0, double tol=1e-5, size_t niter_tol=20, size_t max_iter=1e7) |
Perform a series of time-steps until the next plastic event, or equilibrium. More... | |
std::tuple< array_type::tensor< double, 1 >, array_type::tensor< double, 2 >, array_type::tensor< double, 2 > > | minimise_highfrequency (const array_type::tensor< size_t, 1 > &nodes, const array_type::tensor< size_t, 1 > &top, size_t t_step=1, size_t nmargin=0, double tol=1e-5, size_t niter_tol=20, size_t max_iter=1e7) |
Minimise energy while doing a high-frequency measurement of: More... | |
template<class T > | |
long | flowSteps (size_t n, const T &v, size_t nmargin=0) |
Make a number of steps with the following protocol. More... | |
long | minimise (size_t nmargin=0, double tol=1e-5, size_t niter_tol=20, size_t max_iter=1e7, bool time_activity=false, bool max_iter_is_error=true) |
Minimise energy: run System::timeStep until a mechanical equilibrium has been reached. More... | |
size_t | quasistaticActivityFirst () const |
Increment with the first plastic event. More... | |
size_t | quasistaticActivityLast () const |
Increment with the last plastic event. More... | |
template<class T > | |
size_t | minimise_truncate (const T &idx_n, size_t A_truncate=0, size_t S_truncate=0, double tol=1e-5, size_t niter_tol=20, size_t max_iter=1e7) |
Like Generic2d::minimise(), but stops when a certain number of blocks has yielded at least once. More... | |
size_t | minimise_truncate (size_t A_truncate=0, size_t S_truncate=0, double tol=1e-5, size_t niter_tol=20, size_t max_iter=1e7) |
Like Generic2d::minimise(), but stops when a certain number of blocks has yielded at least once. More... | |
array_type::tensor< double, 2 > | affineSimpleShear (double delta_gamma) const |
Get the displacement field that corresponds to an affine simple shear of a certain strain. More... | |
array_type::tensor< double, 2 > | affineSimpleShearCentered (double delta_gamma) const |
Get the displacement field that corresponds to an affine simple shear of a certain strain. More... | |
Additional Inherited Members | |
![]() | |
void | computeEpsSig () |
Compute m_Sig and m_Eps using m_u. More... | |
void | computeForceMaterial () |
Update m_fmaterial based on the current displacement field m_u. More... | |
virtual void | computeInternalExternalResidualForce () |
Compute: More... | |
![]() | |
array_type::tensor< double, 2 > | m_coor |
Nodal coordinates, see coor(). More... | |
size_t | m_N |
Number of plastic elements, alias of m_nelem_plas. More... | |
size_t | m_nelem |
Number of elements. More... | |
size_t | m_nelem_elas |
Number of elastic elements. More... | |
size_t | m_nelem_plas |
Number of plastic elements. More... | |
size_t | m_nne |
Number of nodes per element. More... | |
size_t | m_nnode |
Number of nodes. More... | |
size_t | m_ndim |
Number of spatial dimensions. More... | |
size_t | m_nip |
Number of integration points. More... | |
array_type::tensor< size_t, 1 > | m_elem_elas |
Elastic elements. More... | |
array_type::tensor< size_t, 1 > | m_elem_plas |
Plastic elements. More... | |
GooseFEM::Element::Quad4::Quadrature | m_quad |
Quadrature for all elements. More... | |
GooseFEM::Element::Quad4::Quadrature | m_quad_elas |
m_quad for elastic elements only. More... | |
GooseFEM::Element::Quad4::Quadrature | m_quad_plas |
m_quad for plastic elements only. More... | |
GooseFEM::VectorPartitioned | m_vector |
Convert vectors between 'nodevec', 'elemvec', .... More... | |
GooseFEM::VectorPartitioned | m_vector_elas |
m_vector for elastic elements only. More... | |
GooseFEM::VectorPartitioned | m_vector_plas |
m_vector for plastic elements only. More... | |
GooseFEM::MatrixDiagonalPartitioned | m_M |
Mass matrix (diagonal). More... | |
GooseFEM::MatrixDiagonal | m_D |
Damping matrix (diagonal). More... | |
GMatElastoPlasticQPot::Cartesian2d::Elastic< 2 > | m_elas |
Material for elastic el. More... | |
GMatElastoPlasticQPot::Cartesian2d::Cusp< 2 > | m_plas |
Material for plastic el. More... | |
array_type::tensor< double, 2 > | m_u |
Nodal displacements. More... | |
array_type::tensor< double, 2 > | m_v |
Nodal velocities. More... | |
array_type::tensor< double, 2 > | m_a |
Nodal accelerations. More... | |
array_type::tensor< double, 2 > | m_v_n |
Nodal velocities last time-step. More... | |
array_type::tensor< double, 2 > | m_a_n |
Nodal accelerations last time-step. More... | |
array_type::tensor< double, 3 > | m_ue |
Element vector (used for displacements). More... | |
array_type::tensor< double, 3 > | m_fe |
Element vector (used for forces). More... | |
array_type::tensor< double, 3 > | m_ue_elas |
m_ue for elastic element only More... | |
array_type::tensor< double, 3 > | m_fe_elas |
m_fe for elastic element only More... | |
array_type::tensor< double, 3 > | m_ue_plas |
m_ue for plastic element only More... | |
array_type::tensor< double, 3 > | m_fe_plas |
m_fe for plastic element only More... | |
array_type::tensor< double, 2 > | m_fmaterial |
Nodal force, deriving from elasticity. More... | |
array_type::tensor< double, 2 > | m_felas |
Nodal force from elasticity of elastic elements. More... | |
array_type::tensor< double, 2 > | m_fplas |
Nodal force from elasticity of plastic elements. More... | |
array_type::tensor< double, 2 > | m_fdamp |
Nodal force from damping. More... | |
array_type::tensor< double, 2 > | m_fvisco |
Nodal force from damping at the interface. More... | |
array_type::tensor< double, 2 > | m_ftmp |
Temporary for internal use. More... | |
array_type::tensor< double, 2 > | m_fint |
Nodal force: total internal force. More... | |
array_type::tensor< double, 2 > | m_fext |
Nodal force: total external force (reaction force) More... | |
array_type::tensor< double, 2 > | m_fres |
Nodal force: residual force. More... | |
array_type::tensor< double, 4 > | m_Eps |
Quad-point tensor: strain. More... | |
array_type::tensor< double, 4 > | m_Sig |
Quad-point tensor: stress. More... | |
array_type::tensor< double, 4 > | m_Epsdot_plas |
Quad-point tensor: strain-rate for plastic el. More... | |
GooseFEM::Matrix | m_K_elas |
Stiffness matrix for elastic elements only. More... | |
size_t | m_qs_inc_first = 0 |
First increment with plastic activity during minimisation. More... | |
size_t | m_qs_inc_last = 0 |
Last increment with plastic activity during minimisation. More... | |
size_t | m_inc |
Current increment (current time = m_dt * m_inc). More... | |
double | m_dt |
Time-step. More... | |
double | m_eta |
Damping at the interface. More... | |
double | m_rho |
Mass density (non-zero only if homogeneous). More... | |
double | m_alpha |
Background damping density (non-zero only if homogeneous). More... | |
bool | m_set_D |
Use m_D only if it is non-zero. More... | |
bool | m_set_visco |
Use m_eta only if it is non-zero. More... | |
bool | m_full_outdated |
Keep track of the need to recompute fields on full geometry. More... | |
array_type::tensor< double, 2 > | m_pert_u |
See eventDriven_setDeltaU() More... | |
array_type::tensor< double, 4 > | m_pert_Epsd_plastic |
Strain deviator for m_pert_u. More... | |
Class that uses GMatElastoPlasticQPot to evaluate stress everywhere.
Definition at line 47 of file UniformSingleLayer2d.h.
|
inlinevirtual |
Reimplemented from FrictionQPotFEM::Generic2d::System.
Definition at line 53 of file UniformSingleLayer2d.h.
|
inline |
Define the geometry, including boundary conditions and element sets.
C | Type of nodal coordinates, e.g. array_type::tensor<double, 2> |
E | Type of connectivity and DOFs, e.g. array_type::tensor<size_t, 2> |
L | Type of node/element lists, e.g. array_type::tensor<size_t, 1> |
coor | Nodal coordinates. |
conn | Connectivity. |
dofs | DOFs per node. |
iip | DOFs whose displacement is fixed. |
elastic_elem | Elastic elements. |
elastic_K | Bulk modulus per quad. point of each elastic element, see setElastic(). |
elastic_G | Shear modulus per quad. point of each elastic element, see setElastic(). |
plastic_elem | Plastic elements. |
plastic_K | Bulk modulus per quad. point of each plastic element, see Plastic(). |
plastic_G | Shear modulus per quad. point of each plastic element, see Plastic(). |
plastic_epsy | Yield strain per quad. point of each plastic element, see Plastic(). |
dt | Time step, set setDt(). |
rho | Mass density, see setMassMatrix(). |
alpha | Background damping density, see setDampingMatrix(). |
eta | Damping at the interface (homogeneous), see setEta(). |
Definition at line 79 of file UniformSingleLayer2d.h.
|
inline |
Add simple shear until a target equivalent macroscopic stress has been reached.
Depending of the target stress compared to the current equivalent macroscopic stress, the shear can be either to the left or to the right.
target_stress | Target stress (equivalent deviatoric value of the macroscopic stress tensor). |
dry_run | If true do not apply displacement, do not check. |
Throws | if yielding is triggered before the stress was reached. |
Definition at line 224 of file UniformSingleLayer2d.h.
|
inline |
Add simple shear until a target equivalent macroscopic stress has been reached.
Depending of the target stress compared to the current equivalent macroscopic stress, the shear can be either to the left or to the right.
target_stress | Target stress (equivalent deviatoric value of the macroscopic stress tensor). |
dry_run | If true do not apply displacement, do not check. |
Definition at line 173 of file UniformSingleLayer2d.h.
|
inline |
Initialise event driven protocol for affine simple shear.
Definition at line 152 of file UniformSingleLayer2d.h.
|
inline |
Read the distance to overcome the first cusp in the element.
deps_kick | Size of the local stain kick to apply. |
iquad | Which integration point to check:
|
[N, 2]
with columns [delta_eps, delta_epsxy]
. Definition at line 334 of file UniformSingleLayer2d.h.
|
inline |
Apply local strain to the right to a specific plastic element.
This 'triggers' one element while keeping the boundary conditions unchanged. Note that by applying shear to the element, yielding can also be triggered in the surrounding elements.
deps_kick | Size of the local stain kick to apply. |
plastic_element | Which plastic element to trigger: System::plastic_elem()(plastic_element). |
trigger_weakest | If true , trigger the weakest integration point. If false , trigger the strongest. |
amplify | Amplify the strain kick with a certain factor. |
Definition at line 255 of file UniformSingleLayer2d.h.
|
inlineoverridevirtual |
Return the type of system.
Reimplemented from FrictionQPotFEM::Generic2d::System.
Reimplemented in FrictionQPotFEM::UniformSingleLayerThermal2d::System.
Definition at line 115 of file UniformSingleLayer2d.h.
|
inline |
Integration point volume of all elements along the weak layer.
Definition at line 140 of file UniformSingleLayer2d.h.
|
inline |
Element height of all elements along the weak layer.
Definition at line 125 of file UniformSingleLayer2d.h.