Builds finite element arrays for thermal diffusion and convection in 2-D domains using 3-Node triangles. More...

#include <DC2DT3.h>

Inheritance diagram for DC2DT3:
Equa_Therm< 3, 3, 2, 2 > Equation< NEN_, NEE_, NSN_, NSE_ > Equa

Public Member Functions

 DC2DT3 ()
 Default Constructor. Constructs an empty equation.
 
 DC2DT3 (Mesh &ms)
 Constructor using Mesh data.
 
 DC2DT3 (Mesh &ms, Vect< real_t > &u)
 Constructor using Mesh and initial condition.
 
 ~DC2DT3 ()
 Destructor.
 
void LCapacity (real_t coef=1)
 Add lumped capacity matrix to element matrix after multiplying it by coefficient coef
 
void Capacity (real_t coef=1)
 Add Consistent capacity matrix to element matrix after multiplying it by coefficient coef
 
void Diffusion (real_t coef=1)
 Add diffusion matrix to element matrix after multiplying it by coefficient coef
 
void Diffusion (const LocalMatrix< real_t, 2, 2 > &diff, real_t coef=1)
 Add diffusion matrix to element matrix after multiplying it by coefficient coef
 
void Convection (const Point< real_t > &v, real_t coef=1)
 Add convection matrix to element matrix after multiplying it by coefficient coef
 
void Convection (const Vect< real_t > &v, real_t coef=1)
 Add convection matrix to element matrix after multiplying it by coefficient coef
 
void Convection (real_t coef=1)
 Add convection matrix to element matrix after multiplying it by coefficient coef
 
void LinearExchange (real_t coef, real_t T)
 Add an edge linear exchange term to left and right-hand sides.
 
void BodyRHS (const Vect< real_t > &f)
 Add body right-hand side term to right hand side.
 
void BodyRHS (real_t f)
 Add body right-hand side term to right hand side.
 
void BoundaryRHS (real_t flux)
 Add boundary right-hand side flux to right hand side.
 
void BoundaryRHS (const Vect< real_t > &f)
 Add boundary right-hand side term to right hand side after multiplying it by coefficient coef
 
void Periodic (real_t coef=1.e20)
 Add contribution of periodic boundary condition (by a penalty technique).
 
Point< real_t > & Flux () const
 Return (constant) heat flux in element.
 
void Grad (Vect< Point< real_t > > &g)
 Compute gradient of solution.
 
Point< real_t > & Grad (const Vect< real_t > &u) const
 Return gradient of a vector in element.
 
void setInput (EType opt, Vect< real_t > &u)
 Set equation input data.
 
void JouleHeating (const Vect< real_t > &sigma, const Vect< real_t > &psi)
 Set Joule heating term as source.
 
- Public Member Functions inherited from Equa_Therm< 3, 3, 2, 2 >
 Equa_Therm ()
 Default constructor.
 
virtual ~Equa_Therm ()
 Destructor.
 
virtual void setStab ()
 Set stabilized formulation.
 
void build ()
 Build the linear system of equations.
 
void build (TimeStepping &s)
 Build the linear system of equations.
 
void build (EigenProblemSolver &e)
 Build the linear system for an eigenvalue problem.
 
void setRho (const real_t &rho)
 Set Density (constant)
 
void setCp (const real_t &cp)
 Set Specific heat (constant)
 
void setConductivity (const real_t &diff)
 Set (constant) thermal conductivity.
 
- Public Member Functions inherited from Equation< NEN_, NEE_, NSN_, NSE_ >
 Equation ()
 
 Equation (Mesh &mesh)
 Constructor with mesh instance.
 
 Equation (Mesh &mesh, Vect< real_t > &u)
 Constructor with mesh instance and solution vector.
 
 Equation (Mesh &mesh, Vect< real_t > &u, real_t &init_time, real_t &final_time, real_t &time_step)
 Constructor with mesh instance, matrix and right-hand side.
 
 ~Equation ()
 Destructor.
 
void updateBC (const Element &el, const Vect< real_t > &bc)
 Update Right-Hand side by taking into account essential boundary conditions.
 
void DiagBC (DOFSupport dof_type=NODE_DOF, int dof=0)
 Update element matrix to impose bc by diagonalization technique.
 
void LocalNodeVector (Vect< real_t > &b)
 Localize element vector from a Vect instance.
 
void ElementNodeVector (const Vect< real_t > &b, LocalVect< real_t, NEE_ > &be)
 Localize element vector from a Vect instance.
 
void SideNodeVector (const Vect< real_t > &b, LocalVect< real_t, NSE_ > &bs)
 Localize side vector from a Vect instance.
 
void SideSideVector (const Vect< real_t > &b, vector< real_t > &bs)
 Localize side vector from a Vect instance.
 
void ElementNodeVectorSingleDOF (const Vect< real_t > &b, LocalVect< real_t, NEN_ > &be)
 Localize Element Vector from a Vect instance.
 
void ElementNodeVector (const Vect< real_t > &b, LocalVect< real_t, NEN_ > &be, int dof)
 Localize Element Vector from a Vect instance.
 
void ElementSideVector (const Vect< real_t > &b, LocalVect< real_t, NSE_ > &be)
 Localize Element Vector from a Vect instance.
 
void ElementVector (const Vect< real_t > &b, DOFSupport dof_type=NODE_DOF, int flag=0)
 Localize element vector.
 
void SideVector (const Vect< real_t > &b, vector< real_t > &sb)
 Localize side vector.
 
void ElementNodeCoordinates ()
 Localize coordinates of element nodes.
 
void SideNodeCoordinates ()
 Localize coordinates of side nodes.
 
void ElementAssembly (Matrix< real_t > *A)
 Assemble element matrix into global one.
 
void ElementAssembly (BMatrix< real_t > &A)
 Assemble element matrix into global one.
 
void ElementAssembly (SkSMatrix< real_t > &A)
 Assemble element matrix into global one.
 
void ElementAssembly (SkMatrix< real_t > &A)
 Assemble element matrix into global one.
 
void ElementAssembly (SpMatrix< real_t > &A)
 Assemble element matrix into global one.
 
void ElementAssembly (TrMatrix< real_t > &A)
 Assemble element matrix into global one.
 
void DGElementAssembly (Matrix< real_t > *A)
 Assemble element matrix into global one for the Discontinuous Galerkin approximation.
 
void DGElementAssembly (SkSMatrix< real_t > &A)
 Assemble element matrix into global one for the Discontinuous Galerkin approximation.
 
void DGElementAssembly (SkMatrix< real_t > &A)
 Assemble element matrix into global one for the Discontinuous Galerkin approximation.
 
void DGElementAssembly (SpMatrix< real_t > &A)
 Assemble element matrix into global one for the Discontinuous Galerkin approximation.
 
void DGElementAssembly (TrMatrix< real_t > &A)
 Assemble element matrix into global one for the Discontinuous Galerkin approximation.
 
void SideAssembly (Matrix< real_t > *A)
 Assemble side (edge or face) matrix into global one.
 
void SideAssembly (SkSMatrix< real_t > &A)
 Assemble side (edge or face) matrix into global one.
 
void SideAssembly (SkMatrix< real_t > &A)
 Assemble side (edge or face) matrix into global one.
 
void SideAssembly (SpMatrix< real_t > &A)
 Assemble side (edge or face) matrix into global one.
 
void ElementAssembly (Vect< real_t > &v)
 Assemble element vector into global one.
 
void SideAssembly (Vect< real_t > &v)
 Assemble side (edge or face) vector into global one.
 
void AxbAssembly (const Element &el, const Vect< real_t > &x, Vect< real_t > &b)
 Assemble product of element matrix by element vector into global vector.
 
void AxbAssembly (const Side &sd, const Vect< real_t > &x, Vect< real_t > &b)
 Assemble product of side matrix by side vector into global vector.
 
size_t getNbNodes () const
 Return number of element nodes.
 
size_t getNbEq () const
 Return number of element equations.
 
real_t setMaterialProperty (const string &exp, const string &prop)
 Define a material property by an algebraic expression.
 
- Public Member Functions inherited from Equa
 Equa ()
 Default constructor.
 
virtual ~Equa ()
 Destructor.
 
void setMesh (Mesh &m)
 Define mesh and renumber DOFs after removing imposed ones.
 
MeshgetMesh () const
 Return reference to Mesh instance.
 
LinearSolvergetLinearSolver ()
 Return reference to linear solver instance.
 
Matrix< real_t > * getMatrix () const
 Return pointer to matrix.
 
void setSolver (Iteration ls, Preconditioner pc=IDENT_PREC)
 Choose solver for the linear system.
 
void setMatrixType (int t)
 Choose type of matrix.
 
int solveLinearSystem (Matrix< real_t > *A, Vect< real_t > &b, Vect< real_t > &x)
 Solve the linear system with given matrix and right-hand side.
 
int solveLinearSystem (Vect< real_t > &b, Vect< real_t > &x)
 Solve the linear system with given right-hand side.
 
void LinearSystemInfo ()
 Print info on linear system solver.
 

Additional Inherited Members

- Protected Member Functions inherited from Equa_Therm< 3, 3, 2, 2 >
void setMaterial ()
 Set material properties.
 

Detailed Description

Builds finite element arrays for thermal diffusion and convection in 2-D domains using 3-Node triangles.

Note that members calculating element arrays have as an argument a real coef that will be multiplied by the contribution of the current element. This makes possible testing different algorithms.

Constructor & Destructor Documentation

◆ DC2DT3() [1/2]

DC2DT3 ( Mesh ms)

Constructor using Mesh data.

Parameters
[in]msMesh instance

◆ DC2DT3() [2/2]

DC2DT3 ( Mesh ms,
Vect< real_t > &  u 
)

Constructor using Mesh and initial condition.

Parameters
[in]msMesh instance
[in]uVect instance containing initial solution

Member Function Documentation

◆ BodyRHS() [1/2]

void BodyRHS ( const Vect< real_t > &  f)
virtual

Add body right-hand side term to right hand side.

Parameters
[in]fVector containing source at nodes.

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ BodyRHS() [2/2]

void BodyRHS ( real_t  f)

Add body right-hand side term to right hand side.

Case where the body right-hand side is piecewise constant.

Parameters
[in]fValue of thermal source (Constant in element).

◆ BoundaryRHS() [1/2]

void BoundaryRHS ( const Vect< real_t > &  f)
virtual

Add boundary right-hand side term to right hand side after multiplying it by coefficient coef

Parameters
[in]fVector containing source at nodes

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ BoundaryRHS() [2/2]

void BoundaryRHS ( real_t  flux)

Add boundary right-hand side flux to right hand side.

Parameters
[in]fluxVector containing source at side nodes.

◆ Capacity()

void Capacity ( real_t  coef = 1)
virtual

Add Consistent capacity matrix to element matrix after multiplying it by coefficient coef

Parameters
[in]coefCoefficient to multiply by added term [Default: 1]

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ Convection() [1/3]

void Convection ( const Point< real_t > &  v,
real_t  coef = 1 
)

Add convection matrix to element matrix after multiplying it by coefficient coef

Parameters
[in]vConstant velocity vector
[in]coefCoefficient to multiply by added term [Default: 1]

◆ Convection() [2/3]

void Convection ( const Vect< real_t > &  v,
real_t  coef = 1 
)

Add convection matrix to element matrix after multiplying it by coefficient coef

Case where velocity field is given by a vector v

Parameters
[in]vVelocity vector
[in]coefCoefficient to multiply by added term (Default: 1]

◆ Convection() [3/3]

void Convection ( real_t  coef = 1)
virtual

Add convection matrix to element matrix after multiplying it by coefficient coef

Case where velocity field has been previouly defined

Parameters
[in]coefCoefficient to multiply by added term [Default: 1]

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ Diffusion() [1/2]

void Diffusion ( const LocalMatrix< real_t, 2, 2 > &  diff,
real_t  coef = 1 
)

Add diffusion matrix to element matrix after multiplying it by coefficient coef

Case where the diffusivity matrix is given as an argument.

Parameters
[in]diffDiffusion matrix (class LocalMatrix).
[in]coefCoefficient to multiply by added term [Default: 1]

◆ Diffusion() [2/2]

void Diffusion ( real_t  coef = 1)
virtual

Add diffusion matrix to element matrix after multiplying it by coefficient coef

Parameters
[in]coefCoefficient to multiply by added term [Default: 1]

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ Grad() [1/2]

Point< real_t > & Grad ( const Vect< real_t > &  u) const

Return gradient of a vector in element.

Parameters
[in]uGlobal vector for which gradient is computed. Vector u has as size the total number of nodes

◆ Grad() [2/2]

void Grad ( Vect< Point< real_t > > &  g)

Compute gradient of solution.

Parameters
[in]gElementwise vector containing gradient of solution.

◆ JouleHeating()

void JouleHeating ( const Vect< real_t > &  sigma,
const Vect< real_t > &  psi 
)

Set Joule heating term as source.

Parameters
[in]sigmaVect instance containing electric conductivity (elementwise)
[in]psiVect instance containing electric potential (elementwise)

◆ LCapacity()

void LCapacity ( real_t  coef = 1)
virtual

Add lumped capacity matrix to element matrix after multiplying it by coefficient coef

Parameters
[in]coefCoefficient to multiply by added term [Default: 1].

Reimplemented from Equa_Therm< 3, 3, 2, 2 >.

◆ LinearExchange()

void LinearExchange ( real_t  coef,
real_t  T 
)

Add an edge linear exchange term to left and right-hand sides.

Parameters
[in]coefCoefficient of exchange
[in]TExternal value for exchange
Remarks
This assumes a constant value of T

◆ Periodic()

void Periodic ( real_t  coef = 1.e20)

Add contribution of periodic boundary condition (by a penalty technique).

Boundary nodes where periodic boundary conditions are to be imposed must have codes equal to PERIODIC_A on one side and PERIODIC_B on the opposite side.

Parameters
[in]coefValue of penalty parameter [Default: 1.e20]

◆ setInput()

void setInput ( EType  opt,
Vect< real_t > &  u 
)

Set equation input data.

Parameters
[in]optParameter to select type of input (enumerated values)
  • INITIAL_FIELD: Initial temperature
  • BOUNDARY_CONDITION_DATA: Boundary condition (Dirichlet)
  • SOURCE_DATA: Heat source
  • FLUX_DATA: Heat flux (Neumann boundary condition)
  • VELOCITY_FIELD: Velocity vector (for the convection term)
[in]uVector containing input data