Seismo-VLAB  1.3
An Open-Source Finite Element Software for Meso-Scale Simulations
Fib3DLineSection Class Reference

Class for creating a generic fiber section for 3D analysis with linear/nonlinear material for a frame elements. More...

#include <Fib3DLineSection.hpp>

Inheritance diagram for Fib3DLineSection:
Collaboration diagram for Fib3DLineSection:

Public Member Functions

 Fib3DLineSection (double h, double b, const std::vector< std::unique_ptr< Material > > &fibers, const std::vector< double > &zi, const std::vector< double > &yi, const std::vector< double > &Ai, double k2=1.0, double k3=1.0, double theta=0.0, unsigned int ip=10)
 Creates a Section to be specified at a Gauss-point in an Element. More...
 
 ~Fib3DLineSection ()
 Destroys this Fib3DLineSection object. More...
 
std::unique_ptr< SectionCopySection ()
 Clone the 'Fib3DLineSection' section. More...
 
Eigen::VectorXd GetStrain ()
 Returns the resultant (generalised) strain vector over the section. More...
 
Eigen::VectorXd GetStress ()
 Returns the resultant (generalised) stress vector over the section. More...
 
Eigen::MatrixXd GetDensity ()
 Access the section density matrix. More...
 
Eigen::MatrixXd GetTangentStiffness ()
 Returns the section stiffness matrix. More...
 
Eigen::MatrixXd GetInitialTangentStiffness ()
 Returns the section initial stiffness matrix. More...
 
Eigen::VectorXd GetStrainAt (double x3=0, double x2=0)
 Returns the section strain at given position. More...
 
Eigen::VectorXd GetStressAt (double x3=0, double x2=0)
 Returns the section stress at given position. More...
 
void CommitState ()
 Perform converged section state update. More...
 
void ReverseState ()
 Reverse the section states to previous converged state. More...
 
void InitialState ()
 Brings the section states to its initial state. More...
 
void UpdateState (const Eigen::VectorXd strain, const unsigned int cond=0)
 Update the section state for this iteration. More...
 
- Public Member Functions inherited from Section
 Section (std::string name)
 Creates a Section to be specified at a Gauss-point in an Element. More...
 
virtual ~Section ()=0
 Destroys this Section object. More...
 
std::string GetName ()
 Gets Section name. More...
 

Protected Member Functions

void ComputeSectionCenter ()
 Gets the section centroid. More...
 
unsigned int FiberIndexAt (double x3=0, double x2=0)
 Return the fiber strain closer to (x3, x2) More...
 
- Protected Member Functions inherited from Section
Eigen::MatrixXd GetLineRotationMatrix (double theta)
 Gets Local Axis Rotation for Line Section according to provided angle. More...
 
Eigen::MatrixXd GetAreaRotationMatrix (double Theta)
 Gets Local Axis Rotation for Area Section according to provided angle. More...
 
Eigen::MatrixXd GetLineTranslationMatrix (double h, double b, double zc, double yc, unsigned int ip)
 Gets centroid translation axis for Line Section according to insertion point. More...
 
void InsertionPointCoordinates (double &x3, double &x2, double h, double b, double zc, double yc, unsigned int ip)
 Gets the coordinate according to insertion point. More...
 
Eigen::MatrixXd ComputeLineLocalAxes (double h, double b, double zcm, double ycm, double angle, unsigned int ip)
 Transforms generalised strain/stresses from Element to Section local coordinate. More...
 

Private Attributes

double h
 Height. More...
 
double b
 Width. More...
 
double zcm
 Rotation angle. More...
 
double ycm
 Rotation angle. More...
 
double kappa2
 Shear area factor (As2/A) along the x2-axis. More...
 
double kappa3
 Shear area factor (As3/A) along the x3-axis. More...
 
double Theta
 Rotation angle. More...
 
unsigned int InsertPoint
 Local axis coordinate. More...
 
unsigned int NumberOfFibers
 Local axis coordinate. More...
 
Eigen::VectorXd Strain
 Generalized Strain vector. More...
 
std::vector< double > zi
 The horizontal (x3) local coordinates of each fiber. More...
 
std::vector< double > yi
 The vertical (x2) local coordinates of each fiber. More...
 
std::vector< double > Ai
 The discretized area value of each fiber. More...
 
std::vector< std::unique_ptr< Material > > theFiber
 The section's fiber list. More...
 

Detailed Description

Class for creating a generic fiber section for 3D analysis with linear/nonlinear material for a frame elements.

See also
Section.hpp

Constructor & Destructor Documentation

◆ Fib3DLineSection()

Fib3DLineSection::Fib3DLineSection ( double  h,
double  b,
const std::vector< std::unique_ptr< Material > > &  fibers,
const std::vector< double > &  zi,
const std::vector< double > &  yi,
const std::vector< double > &  Ai,
double  k2 = 1.0,
double  k3 = 1.0,
double  theta = 0.0,
unsigned int  ip = 10 
)

Creates a Section to be specified at a Gauss-point in an Element.

Parameters
fibersA vector of fibers for each discretized area.
ziA vector with horizontal coordinate (x3) of each fiber.
yiA vector with vertical coordinate (x2) of each fiber.
AiA vector with the area values of each fiber.
k2Shear area factor along the x2-axis.
k3Shear area factor along the x2-axis.
thetaThe angle of rotation about x3-axis.
ipThe insertion point where the section axes are located.
Note
Details about insertion point and local axis rotation can be found at Insertion Point and Local Axes.
See also
Fib3DLineSection::h, Fib3DLineSection::b, Fib3DLineSection::theMaterial, Fib3DLineSection::Theta, Fib3DLineSection::InsertPoint.

◆ ~Fib3DLineSection()

Fib3DLineSection::~Fib3DLineSection ( )

Destroys this Fib3DLineSection object.

Member Function Documentation

◆ CommitState()

void Fib3DLineSection::CommitState ( )
virtual

Perform converged section state update.

Note
This function sets the trail stress and strain as converged.

Implements Section.

◆ ComputeSectionCenter()

void Fib3DLineSection::ComputeSectionCenter ( )
protected

Gets the section centroid.

See also
The geometric properties can be revisited at Fib3DLineSection.

◆ CopySection()

std::unique_ptr<Section> Fib3DLineSection::CopySection ( )
virtual

Clone the 'Fib3DLineSection' section.

Returns
A Section pointer to the derived class.
See also
Section.

Implements Section.

◆ FiberIndexAt()

unsigned int Fib3DLineSection::FiberIndexAt ( double  x3 = 0,
double  x2 = 0 
)
protected

Return the fiber strain closer to (x3, x2)

Parameters
x3Position on the x3-axis where the stress is evaluated.
x2Position on the x2-axis where the stress is evaluated.
Returns
Index of the fiber closer to coordinate (x3,x2).

◆ GetDensity()

Eigen::MatrixXd Fib3DLineSection::GetDensity ( )
virtual

Access the section density matrix.

Returns
The section density matrix.

Implements Section.

◆ GetInitialTangentStiffness()

Eigen::MatrixXd Fib3DLineSection::GetInitialTangentStiffness ( )
virtual

Returns the section initial stiffness matrix.

Returns
Matrix with the initial section tangent stiffness matrix.
Note
The initial tangent stiffness matrix is computed when the generalized strain vector is zero.

Implements Section.

◆ GetStrain()

Eigen::VectorXd Fib3DLineSection::GetStrain ( )
virtual

Returns the resultant (generalised) strain vector over the section.

Returns
Vector with the resultant strain components.

Implements Section.

◆ GetStrainAt()

Eigen::VectorXd Fib3DLineSection::GetStrainAt ( double  x3 = 0,
double  x2 = 0 
)
virtual

Returns the section strain at given position.

Parameters
x3Position on the x3-axis where the strain is evaluated.
x2Position on the x2-axis where the strain is evaluated.
Returns
Vector with the strain components at coordinate (x3,x2).

Implements Section.

◆ GetStress()

Eigen::VectorXd Fib3DLineSection::GetStress ( )
virtual

Returns the resultant (generalised) stress vector over the section.

Returns
Vector with the resultant stress components.

Implements Section.

◆ GetStressAt()

Eigen::VectorXd Fib3DLineSection::GetStressAt ( double  x3 = 0,
double  x2 = 0 
)
virtual

Returns the section stress at given position.

Parameters
x3Position on the x3-axis where the stress is evaluated.
x2Position on the x2-axis where the stress is evaluated.
Returns
Vector with the stress components at coordinate (x3,x2).

Implements Section.

◆ GetTangentStiffness()

Eigen::MatrixXd Fib3DLineSection::GetTangentStiffness ( )
virtual

Returns the section stiffness matrix.

Returns
Matrix with the section consistent tangent stiffness matrix.

Implements Section.

◆ InitialState()

void Fib3DLineSection::InitialState ( )
virtual

Brings the section states to its initial state.

Note
This function returns the fiber states to the beginning.

Implements Section.

◆ ReverseState()

void Fib3DLineSection::ReverseState ( )
virtual

Reverse the section states to previous converged state.

Note
This function returns the fiber states to previous converged states.

Implements Section.

◆ UpdateState()

void Fib3DLineSection::UpdateState ( const Eigen::VectorXd  strain,
const unsigned int  cond = 0 
)
virtual

Update the section state for this iteration.

Parameters
strainVector with the strain components at this Gauss-point.
condIf the the elastic/plastic fiber components will be updated.
Note
This function computes the strain and tanget stiffness matrix once the trial strain converged.

Implements Section.

Member Data Documentation

◆ Ai

std::vector<double> Fib3DLineSection::Ai
private

The discretized area value of each fiber.

◆ b

double Fib3DLineSection::b
private

Width.

◆ h

double Fib3DLineSection::h
private

Height.

◆ InsertPoint

unsigned int Fib3DLineSection::InsertPoint
private

Local axis coordinate.

◆ kappa2

double Fib3DLineSection::kappa2
private

Shear area factor (As2/A) along the x2-axis.

◆ kappa3

double Fib3DLineSection::kappa3
private

Shear area factor (As3/A) along the x3-axis.

◆ NumberOfFibers

unsigned int Fib3DLineSection::NumberOfFibers
private

Local axis coordinate.

◆ Strain

Eigen::VectorXd Fib3DLineSection::Strain
private

Generalized Strain vector.

◆ theFiber

std::vector<std::unique_ptr<Material> > Fib3DLineSection::theFiber
private

The section's fiber list.

◆ Theta

double Fib3DLineSection::Theta
private

Rotation angle.

◆ ycm

double Fib3DLineSection::ycm
private

Rotation angle.

◆ yi

std::vector<double> Fib3DLineSection::yi
private

The vertical (x2) local coordinates of each fiber.

◆ zcm

double Fib3DLineSection::zcm
private

Rotation angle.

◆ zi

std::vector<double> Fib3DLineSection::zi
private

The horizontal (x3) local coordinates of each fiber.