FrictionQPotFEM 0.23.3
|
Compared to UniformSingleLayer2d::System() this class adds thermal fluctuations. More...
#include <FrictionQPotFEM/UniformSingleLayerThermal2d.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, double temperature_dinc, size_t temperature_seed, double temperature) | |
Define the geometry, including boundary conditions and element sets. More... | |
std::string | type () const override |
Return the type of system. More... | |
const array_type::tensor< double, 2 > & | fthermal () const |
The force vector that represents the effect of temperature (on the weak layer only). More... | |
double | temperature () const |
Return the target temperature. More... | |
void | setInc (size_t arg) override |
Set increment. 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... | |
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... | |
Protected Member Functions | |
void | updateCache (int64_t index) |
Update cache of thermal forces on the weak layer. More... | |
void | computeThermalForce (bool force=false) |
Update the thermal force vector. More... | |
void | computeInternalExternalResidualForce () override |
Compute: More... | |
![]() | |
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... | |
Protected Attributes | |
size_t | m_computed |
Increment at which the thermal stress tensor was generated. More... | |
size_t | m_dinc |
Duration to keep the same random thermal stress tensor. More... | |
double | m_T |
Definition of temperature (units of equivalent stress). More... | |
double | m_std |
Standard deviation for signed equivalent thermal stress. More... | |
array_type::tensor< double, 2 > | m_fthermal |
Nodal force from temperature. More... | |
array_type::tensor< double, 4 > | m_cache |
Cache [ncache, m_elem_plas.size(), 3, 2]. More... | |
int64_t | m_cache_start |
Start index of the cache. More... | |
int64_t | m_ncache |
Number of cached items. More... | |
prrng::pcg32 | m_gen |
Random generator for thermal forces. 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... | |
Compared to UniformSingleLayer2d::System() this class adds thermal fluctuations.
The fluctuations are implemented as a dipolar force on each element edge. Both its components are drawn from a normal distribution with zero mean and 'temperature' standard deviation. The temperature specification is thereby in units of the equivalent stress, see GMatElastoPlasticQPot::Cartesian2d::Sigd().
Definition at line 52 of file UniformSingleLayerThermal2d.h.
|
inlinevirtual |
Reimplemented from FrictionQPotFEM::UniformSingleLayer2d::System.
Definition at line 58 of file UniformSingleLayerThermal2d.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(). |
temperature_dinc | Duration to keep the same random thermal stress tensor. |
temperature_seed | Seed random generator thermal stress. |
temperature | 'Temperature': in units of the the equivalent stress. |
Definition at line 87 of file UniformSingleLayerThermal2d.h.
|
inlineoverrideprotectedvirtual |
Compute:
Internal rule: all relevant forces are expected to be updated before this function is called.
Reimplemented from FrictionQPotFEM::Generic2d::System.
Definition at line 266 of file UniformSingleLayerThermal2d.h.
|
inlineprotected |
Update the thermal force vector.
force | If false updating is skipped if already computed for the current increment. |
Definition at line 215 of file UniformSingleLayerThermal2d.h.
|
inline |
The force vector that represents the effect of temperature (on the weak layer only).
Definition at line 170 of file UniformSingleLayerThermal2d.h.
|
inlineoverridevirtual |
Set increment.
arg | size_t. |
Reimplemented from FrictionQPotFEM::Generic2d::System.
Definition at line 204 of file UniformSingleLayerThermal2d.h.
|
inline |
Return the target temperature.
Definition at line 179 of file UniformSingleLayerThermal2d.h.
|
inlineoverridevirtual |
Return the type of system.
Reimplemented from FrictionQPotFEM::UniformSingleLayer2d::System.
Definition at line 161 of file UniformSingleLayerThermal2d.h.
|
inlineprotected |
Update cache of thermal forces on the weak layer.
Definition at line 188 of file UniformSingleLayerThermal2d.h.
|
protected |
Cache [ncache, m_elem_plas.size(), 3, 2].
Definition at line 155 of file UniformSingleLayerThermal2d.h.
|
protected |
Start index of the cache.
Definition at line 156 of file UniformSingleLayerThermal2d.h.
|
protected |
Increment at which the thermal stress tensor was generated.
Definition at line 150 of file UniformSingleLayerThermal2d.h.
|
protected |
Duration to keep the same random thermal stress tensor.
Definition at line 151 of file UniformSingleLayerThermal2d.h.
|
protected |
Nodal force from temperature.
Definition at line 154 of file UniformSingleLayerThermal2d.h.
|
protected |
Random generator for thermal forces.
Definition at line 158 of file UniformSingleLayerThermal2d.h.
|
protected |
Number of cached items.
Definition at line 157 of file UniformSingleLayerThermal2d.h.
|
protected |
Standard deviation for signed equivalent thermal stress.
Definition at line 153 of file UniformSingleLayerThermal2d.h.
|
protected |
Definition of temperature (units of equivalent stress).
Definition at line 152 of file UniformSingleLayerThermal2d.h.