Feellgood
Public Member Functions | Public Attributes | Private Member Functions | Private Attributes | List of all members
electrostatSolver Class Reference

#include <electrostatSolver.h>

Inheritance diagram for electrostatSolver:
solver< DIM_PB_ELEC >

Public Member Functions

 electrostatSolver (Mesh::mesh &_msh, std::vector< Tetra::prm > &_pTetra, std::vector< Facette::prm > &_pFac, const double _tol, const bool v, const int max_iter)
 
bool checkBoundaryConditions (void) const
 
void infos (void)
 
void integrales (Tetra::Tet const &tet, Eigen::Ref< Eigen::Matrix< double, Tetra::N, Tetra::N > > AE)
 
void integrales (Facette::Fac const &fac, std::vector< double > &BE)
 
bool save (const std::string V_fileName, std::string const &metadata) const
 
double getSigma (Tetra::Tet const &tet) const
 
double getCurrentDensity (Facette::Fac const &fac) const
 
void compute (const bool verbose, const std::string V_fileName)
 
- Public Member Functions inherited from solver< DIM_PB_ELEC >
 solver (Mesh::mesh &_msh, std::vector< Tetra::prm > &_pTetra, std::vector< Facette::prm > &_pFac, const std::string name, const double _tol, const bool v, const int max_iter)
 

Public Attributes

std::vector< double > V
 

Private Member Functions

bool solve (void)
 

Private Attributes

const int precision = 8
 

Additional Inherited Members

- Protected Member Functions inherited from solver< DIM_PB_ELEC >
void buildMat (std::vector< int > &ind, Eigen::Matrix< double, DIM_PROBLEM *N, DIM_PROBLEM *N > &Ke, algebra::w_sparseMat &K)
 
void buildVect (std::vector< int > &ind, std::vector< double > &Le, std::vector< double > &L)
 
- Protected Attributes inherited from solver< DIM_PB_ELEC >
Mesh::meshmsh
 
const int NOD
 
const std::vector< Tetra::prm > & paramTet
 
const std::vector< Facette::prm > & paramFac
 
algebra::iteration< double > iter
 
- Static Protected Attributes inherited from solver< DIM_PB_ELEC >
static const int DIM_PB
 

Detailed Description

this class is containing both data and a solver to compute potential from dirichlet boundary conditions problem for the current density flowing in the sample.

Constructor & Destructor Documentation

◆ electrostatSolver()

electrostatSolver::electrostatSolver ( Mesh::mesh _msh,
std::vector< Tetra::prm > &  _pTetra,
std::vector< Facette::prm > &  _pFac,
const double  _tol,
const bool  v,
const int  max_iter 
)
inline

constructor

Parameters
[in]_mshreference to the mesh
[in]_pTetraref to vector of param tetra (volume region parameters)
[in]_pFacref to vector of param facette (surface region parameters)
[in]_tolcg_dir tolerance
[in]vverbose mode for iteration monitor
[in]max_itermaximum number of iterations

Member Function Documentation

◆ checkBoundaryConditions()

bool electrostatSolver::checkBoundaryConditions ( void  ) const

check boundary conditions: mesh and settings have to define a single surface with constant current density J (normal componant to the surface) and another single surface with constant potential V

◆ compute()

void electrostatSolver::compute ( const bool  verbose,
const std::string  V_fileName 
)

solves the potential and stores result in V, save to text file if needed if verbose set to true, some printing are sent to terminal

Parameters
[in]verbose
[in]V_fileNameoutput file name for V solution

◆ getCurrentDensity()

double electrostatSolver::getCurrentDensity ( Facette::Fac const &  fac) const

returns current density of the facette if it is defined in the boundary conditions, else zero

◆ getSigma()

double electrostatSolver::getSigma ( Tetra::Tet const &  tet) const

returns sigma of the tetraedron, (conductivity in (Ohm.m)^-1

◆ infos()

void electrostatSolver::infos ( void  )

basic informations on boundary conditions

◆ integrales() [1/2]

void electrostatSolver::integrales ( Facette::Fac const &  fac,
std::vector< double > &  BE 
)

compute integrales for vector coefficients, input from facette

◆ integrales() [2/2]

void electrostatSolver::integrales ( Tetra::Tet const &  tet,
Eigen::Ref< Eigen::Matrix< double, Tetra::N, Tetra::N > >  AE 
)

compute side problem (electrostatic potential on the nodes) integrales for matrix coefficients,inputs from tet

◆ save()

bool electrostatSolver::save ( const std::string  V_fileName,
std::string const &  metadata 
) const

text file (tsv) writing function for the solution V over all volume regions of the mesh, node indices are zero based

Parameters
[in]V_fileNameoutput file name
[in]metadata

◆ solve()

bool electrostatSolver::solve ( void  )
private

solver, using conjugate gradient with masking, with diagonal preconditionner and Dirichlet boundary conditions, returns true if has converged

Member Data Documentation

◆ precision

const int electrostatSolver::precision = 8
private

number of digits in the optional output file

◆ V

std::vector<double> electrostatSolver::V

electrostatic potential values for boundary conditions, V.size() is the size of the vector of nodes


The documentation for this class was generated from the following files: