Seismo-VLAB  1.3
An Open-Source Finite Element Software for Meso-Scale Simulations
lin2DTria6.hpp
Go to the documentation of this file.
1 //==============================================================================
2 //
3 // Seismo Virtual Laboratory
4 // Module for Serial and Parallel Analysis of seismic
5 // wave propagation and soil-structure interaction simulation
6 // Copyright (C) 2018-2021, The California Institute of Technology
7 // All Rights Reserved.
8 //
9 // Commercial use of this program without express permission of the California
10 // Institute of Technology, is strictly prohibited. See file "COPYRIGHT" in
11 // main directory for information on usage and redistribution, and for a
12 // DISCLAIMER OF ALL WARRANTIES.
13 //
14 //==============================================================================
15 //
16 // Written by:
17 // Danilo S. Kusanovic (dkusanov@caltech.edu)
18 // Elnaz E. Seylabi (elnaze@unr.edu)
19 //
20 // Supervised by:
21 // Domniki M. Asimaki (domniki@caltech.edu)
22 //
23 // References :
24 // [1] Finite Element Procedures, Bathe, K.J., Chapter 5: pages 371-376.
25 // Prentice-Hall, 1996.
26 //
27 // Description:
30 //------------------------------------------------------------------------------
31 
32 #ifndef _LIN2DTRIA6_HPP_
33 #define _LIN2DTRIA6_HPP_
34 
35 #include <map>
36 #include <memory>
37 #include <string>
38 #include <memory>
39 #include <Eigen/Dense>
40 
41 #include "Node.hpp"
42 #include "Load.hpp"
43 #include "Material.hpp"
44 #include "Element.hpp"
45 #include "Damping.hpp"
46 #include "QuadratureRule.hpp"
47 
55 class lin2DTria6 : public Element{
56 
57  public:
66  lin2DTria6(const std::vector<unsigned int> nodes, std::unique_ptr<Material> &material, const double th, const std::string quadrature="GAUSS", const unsigned int nGauss=4);
67 
69  ~lin2DTria6();
70 
73  void CommitState();
74 
77  void ReverseState();
78 
81  void InitialState();
82 
85  void UpdateState();
86 
91  void SetDomain(std::map<unsigned int, std::shared_ptr<Node> > &nodes);
92 
96  void SetDamping(const std::shared_ptr<Damping> &damping);
97 
100  std::vector<unsigned int> GetTotalDegreeOfFreedom() const;
101 
105  Eigen::MatrixXd GetStrain() const;
106 
110  Eigen::MatrixXd GetStress() const;
111 
115  Eigen::MatrixXd GetStrainRate() const;
116 
122  Eigen::MatrixXd GetStrainAt(double x3, double x2) const;
123 
129  Eigen::MatrixXd GetStressAt(double x3, double x2) const;
130 
135  Eigen::VectorXd GetVTKResponse(std::string response) const;
136 
139  double ComputeEnergy();
140 
145  Eigen::MatrixXd ComputeMassMatrix();
146 
151  Eigen::MatrixXd ComputeStiffnessMatrix();
152 
157  Eigen::MatrixXd ComputeDampingMatrix();
158 
163  Eigen::MatrixXd ComputePMLMatrix();
164 
169  Eigen::VectorXd ComputeInternalForces();
170 
175  Eigen::VectorXd ComputeInternalDynamicForces();
176 
183  Eigen::VectorXd ComputeSurfaceForces(const std::shared_ptr<Load> &surface, unsigned int face);
184 
191  Eigen::VectorXd ComputeBodyForces(const std::shared_ptr<Load> &body, unsigned int k=0);
192 
199  Eigen::VectorXd ComputeDomainReductionForces(const std::shared_ptr<Load> &drm, unsigned int k);
200 
201  private:
203  double t;
204 
206  std::shared_ptr<Damping> theDamping;
207 
209  std::vector<std::shared_ptr<Node> > theNodes;
210 
212  std::vector<std::unique_ptr<Material> > theMaterial;
213 
215  std::unique_ptr<QuadratureRule> QuadraturePoints;
216 
220  Eigen::VectorXd ComputeStrain(const Eigen::MatrixXd &Bij) const;
221 
225  Eigen::VectorXd ComputeStrainRate(const Eigen::MatrixXd &Bij) const;
226 
231  Eigen::MatrixXd ComputeJacobianMatrix(const double ri, const double si) const;
232 
237  Eigen::MatrixXd ComputeShapeFunctionMatrix(const double ri, const double si) const;
238 
244  Eigen::MatrixXd ComputeStrainDisplacementMatrix(const double ri, const double si, const Eigen::MatrixXd &Jij) const;
245 
249  Eigen::MatrixXd ComputeInitialStiffnessMatrix() const;
250 };
251 
252 #endif
void InitialState()
Brings the material/section state to its initial state in this element.
Eigen::MatrixXd ComputeShapeFunctionMatrix(const double ri, const double si) const
Evaluates the shape function matrix at a given Gauss point.
Eigen::MatrixXd GetStressAt(double x3, double x2) const
Gets the material stress in section at coordinate (x3,x2).
This file contains the "Load" class declarations, which defines a load that can act in a node or an e...
double t
Element thickness.
Definition: lin2DTria6.hpp:203
Eigen::MatrixXd ComputeStrainDisplacementMatrix(const double ri, const double si, const Eigen::MatrixXd &Jij) const
Evaluates the strain-displacement matrix at a given Gauss point.
Eigen::MatrixXd ComputeDampingMatrix()
Compute the damping matrix of the element using gauss-integration.
Eigen::MatrixXd GetStrainRate() const
Gets the material/section (generalised) strain-rate.
This file contains the abstract "Material object" declarations, which computes the strain...
Eigen::VectorXd ComputeInternalForces()
Compute the internal (elastic) forces acting on the element.
std::vector< std::shared_ptr< Node > > theNodes
The Element Nodes.
Definition: lin2DTria6.hpp:209
void CommitState()
Save the material states in the element.
Eigen::VectorXd ComputeStrainRate(const Eigen::MatrixXd &Bij) const
Update strain rate in the element.
double ComputeEnergy()
Computes the element energy for a given deformation.
void SetDamping(const std::shared_ptr< Damping > &damping)
Sets the damping model.
std::vector< std::unique_ptr< Material > > theMaterial
Element material.
Definition: lin2DTria6.hpp:212
Eigen::VectorXd ComputeSurfaceForces(const std::shared_ptr< Load > &surface, unsigned int face)
Compute the surface forces acting on the element.
void ReverseState()
Reverse the material/section states to previous converged state in this element.
std::shared_ptr< Damping > theDamping
The Damping model.
Definition: lin2DTria6.hpp:206
std::unique_ptr< QuadratureRule > QuadraturePoints
Coordinate of Gauss points.
Definition: lin2DTria6.hpp:215
Eigen::MatrixXd ComputeJacobianMatrix(const double ri, const double si) const
Computes the jacobian of the transformation.
~lin2DTria6()
Destroys this lin2DTria6 object.
Eigen::MatrixXd ComputeInitialStiffnessMatrix() const
Compute the initial stiffness matrix of the element using gauss-integration.
Virtual class for creating an element in a mesh.
Definition: Element.hpp:51
Eigen::MatrixXd GetStrain() const
Gets the material/section (generalised) strain.
Eigen::MatrixXd GetStress() const
Gets the material/section (generalised) stress.
This file contains the "Node object" declarations, which is stores the coordinates, state variables, degrees-of-freedom, and total-degree of freedom lists of a node in a finite element mesh.
Eigen::VectorXd ComputeDomainReductionForces(const std::shared_ptr< Load > &drm, unsigned int k)
Compute the domain reduction forces acting on the element.
Eigen::MatrixXd ComputePMLMatrix()
Compute the PML history matrix using gauss-integration.
This file contains the "Damping object" declarations to update damping matrix to account for Caughey-...
Eigen::MatrixXd GetStrainAt(double x3, double x2) const
Gets the material strain in section at coordinate (x3,x2).
void SetDomain(std::map< unsigned int, std::shared_ptr< Node > > &nodes)
Sets the finite element dependance among objects.
Class for creating a 2D linearized six-node triangular element in a mesh.
Definition: lin2DTria6.hpp:55
std::vector< unsigned int > GetTotalDegreeOfFreedom() const
Gets the list of total-degree of freedom of this Element.
Eigen::MatrixXd ComputeStiffnessMatrix()
Compute the stiffness matrix of the element using gauss-integration.
void UpdateState()
Update the material states in the element.
Eigen::VectorXd ComputeStrain(const Eigen::MatrixXd &Bij) const
Update strain in the element.
Eigen::VectorXd ComputeInternalDynamicForces()
Compute the elastic, inertial, and viscous forces acting on the element.
This file contains the abstract "Quadrature Rule" declarations to integrate quantities in the iso-par...
lin2DTria6(const std::vector< unsigned int > nodes, std::unique_ptr< Material > &material, const double th, const std::string quadrature="GAUSS", const unsigned int nGauss=4)
Creates a lin2DTria6 in a finite element Mesh.
Eigen::VectorXd ComputeBodyForces(const std::shared_ptr< Load > &body, unsigned int k=0)
Compute the body forces acting on the element.
Eigen::VectorXd GetVTKResponse(std::string response) const
Gets the element internal response in VTK format for Paraview display.
Eigen::MatrixXd ComputeMassMatrix()
Compute the lumped/consistent mass matrix of the element.
This file contains the "Element" object declarations, which defines an element in a finite element me...