Feellgood
Public Member Functions | Public Attributes | Private Member Functions | Private Attributes | List of all members
Mesh::mesh Class Reference

#include <mesh.h>

Public Member Functions

 mesh (Settings const &mySets)
 
int getNbNodes (void) const
 
int getNbFacs (void) const
 
int getNbTets (void) const
 
const Eigen::Vector3d getNode_p (const int i) const
 
const Eigen::Vector3d getNode_u (const int i) const
 
const Eigen::Vector3d getNode_v (const int i) const
 
double getProj_ep (const int i) const
 
double getProj_eq (const int i) const
 
void set_node_u0 (const int i, Eigen::Vector3d const &val)
 
void set_node_zero_v (const int i)
 
void infos (void) const
 
void setBasis (const double r)
 
void updateNodes (Eigen::Ref< Eigen::VectorXd > X, const double dt)
 
void evolution (void)
 
double readSol (bool VERBOSE, const std::string fileName)
 
void init_distrib (Settings const &mySets)
 
double avg (std::function< double(Nodes::Node, Nodes::index)> getter, Nodes::index d) const
 
void savesol (const int precision, const std::string fileName, std::string const &metadata) const
 
bool savesol (const int precision, const std::string fileName, std::string const &metadata, std::vector< double > const &val) const
 
void set (const int i, std::function< void(Nodes::Node &, const double)> what_to_set, const double val)
 

Public Attributes

Eigen::Vector3d c
 
Eigen::Vector3d l
 
double vol
 
std::vector< Facette::Facfac
 
std::vector< Tetra::Tettet
 

Private Member Functions

void checkMeshFile (Settings const &mySets)
 
void readNodes (Settings const &mySets)
 
void readTetraedrons (Settings const &mySets)
 
void readTriangles (Settings const &mySets)
 
void readMesh (Settings const &mySets)
 
double doOnNodes (const double init_val, const Nodes::index coord, std::function< bool(double, double)> whatToDo) const
 
double minNodes (const Nodes::index coord) const
 
double maxNodes (const Nodes::index coord) const
 
void indexReorder (std::vector< Tetra::prm > const &prmTetra)
 
void sortNodes (Nodes::index long_axis)
 

Private Attributes

std::vector< Nodes::Nodenode
 
std::vector< int > node_index
 

Detailed Description

class for storing the mesh, including mesh geometry values, containers for the nodes, triangular faces and tetrahedrons. nodes data are private. They are accessible only through getter and setter.

Constructor & Destructor Documentation

◆ mesh()

Mesh::mesh::mesh ( Settings const &  mySets)
inline

constructor : read mesh file, reorder indices and computes some values related to the mesh : center and length along coordinates,full volume

Parameters
[in]mySets

Member Function Documentation

◆ avg()

double mesh::avg ( std::function< double(Nodes::Node, Nodes::index)>  getter,
Nodes::index  d 
) const

average component of either u or v through getter on the whole set of tetetrahedron

Parameters
[in]getter
[in]d

◆ checkMeshFile()

void Mesh::mesh::checkMeshFile ( Settings const &  mySets)
private

test if mesh file contains surfaces and regions mentionned in yaml settings and their dimensions

Parameters
[in]mySets

◆ doOnNodes()

double mesh::doOnNodes ( const double  init_val,
const Nodes::index  coord,
std::function< bool(double, double)>  whatToDo 
) const
private

loop on nodes to apply predicate 'whatTodo'

Parameters
[in]init_val
[in]coord
[in]whatToDo

◆ evolution()

void Mesh::mesh::evolution ( void  )
inline

call evolution for all the nodes

◆ getNbFacs()

int Mesh::mesh::getNbFacs ( void  ) const
inline

return number of triangular fac

◆ getNbNodes()

int Mesh::mesh::getNbNodes ( void  ) const
inline

return number of nodes

◆ getNbTets()

int Mesh::mesh::getNbTets ( void  ) const
inline

return number of tetrahedrons

◆ getNode_p()

const Eigen::Vector3d Mesh::mesh::getNode_p ( const int  i) const
inline

getter : return node.p

◆ getNode_u()

const Eigen::Vector3d Mesh::mesh::getNode_u ( const int  i) const
inline

getter : return node.u

◆ getNode_v()

const Eigen::Vector3d Mesh::mesh::getNode_v ( const int  i) const
inline

getter : return node.v

◆ getProj_ep()

double Mesh::mesh::getProj_ep ( const int  i) const
inline

return projection of speed at node i along ep

◆ getProj_eq()

double Mesh::mesh::getProj_eq ( const int  i) const
inline

return projection of speed at node i along eq

◆ indexReorder()

void mesh::indexReorder ( std::vector< Tetra::prm > const &  prmTetra)
private

redefine orientation of triangular faces in accordance with the tetrahedron reorientation of the tetrahedrons if needed; definition of Ms on facette elements Indices and orientation convention :

                v
              .
            ,/
           /
        2(ic)                                 2
      ,/|`\                                 ,/|`\
    ,/  |  `\                             ,/  |  `\
  ,/    '.   `\                         ,6    '.   `5
,/       |     `\                     ,/       8     `\

,/ | \ ,/ |\ 0(ia)----—'.-----—1(ib) --> u 0-----—4–'.-----—1 \. | ,/. | ,/ \. | ,/. | ,9 ‘. ’. ,/ ‘7. ’. ,/ \. |/. |/ 3(id)3 \. w

◆ infos()

void mesh::infos ( void  ) const

basic informations on the mesh

◆ init_distrib()

void Mesh::mesh::init_distrib ( Settings const &  mySets)
inline

computes an analytical initial magnetization distribution as a starting point for the simulation

Parameters
[in]mySets

◆ maxNodes()

double Mesh::mesh::maxNodes ( const Nodes::index  coord) const
inlineprivate

return the maximum of all nodes coordinate along coord axis

Parameters
[in]coord

◆ minNodes()

double Mesh::mesh::minNodes ( const Nodes::index  coord) const
inlineprivate

return the minimum of all nodes coordinate along coord axis

Parameters
[in]coord

◆ readMesh()

void Mesh::mesh::readMesh ( Settings const &  mySets)
private

reading mesh format 2.2 text file function

Parameters
[in]mySets

◆ readNodes()

void Mesh::mesh::readNodes ( Settings const &  mySets)
private

read Nodes from mesh file

Parameters
[in]mySets

◆ readSol()

double Mesh::mesh::readSol ( bool  VERBOSE,
const std::string  fileName 
)

read a solution from a file (tsv formated) and initialize fem struct to restart computation from that distribution, return time

Parameters
[in]VERBOSE
[in]fileNameinput .sol text file

◆ readTetraedrons()

void Mesh::mesh::readTetraedrons ( Settings const &  mySets)
private

read tetraedrons of the settings volume regions

Parameters
[in]mySets

◆ readTriangles()

void Mesh::mesh::readTriangles ( Settings const &  mySets)
private

read facettes of the settings surface regions

Parameters
[in]mySets

◆ savesol() [1/2]

void Mesh::mesh::savesol ( const int  precision,
const std::string  fileName,
std::string const &  metadata 
) const

text file (tsv) writing function for a solution

Parameters
[in]precisionnumeric precision in .sol output text file
[in]fileName
[in]metadata

◆ savesol() [2/2]

bool Mesh::mesh::savesol ( const int  precision,
const std::string  fileName,
std::string const &  metadata,
std::vector< double > const &  val 
) const

text file (tsv) writing function for a solution of a side problem, used by electrostatSolver

Parameters
[in]precision
[in]fileName
[in]metadata
[in]val

◆ set()

void Mesh::mesh::set ( const int  i,
std::function< void(Nodes::Node &, const double)>  what_to_set,
const double  val 
)
inline

setter for node[i]; what_to_set will fix what is the part of the node struct to set (usefull for fmm_demag.h)

Parameters
[in]i
[in]what_to_set
[in]val

◆ set_node_u0()

void Mesh::mesh::set_node_u0 ( const int  i,
Eigen::Vector3d const &  val 
)
inline

setter for u0

◆ set_node_zero_v()

void Mesh::mesh::set_node_zero_v ( const int  i)
inline

fix to zero node[i].v

◆ setBasis()

void Mesh::mesh::setBasis ( const double  r)
inline

call setBasis for all nodes, and update P matrix for all elements

◆ sortNodes()

void mesh::sortNodes ( Nodes::index  long_axis)
private

Sort the nodes along the longest axis of the sample. This should reduce the bandwidth of the matrix we will have to solve for.

Parameters
[in]long_axis

◆ updateNodes()

void mesh::updateNodes ( Eigen::Ref< Eigen::VectorXd >  X,
const double  dt 
)

make_evol on all nodes

Member Data Documentation

◆ c

Eigen::Vector3d Mesh::mesh::c

isobarycenter

◆ fac

std::vector<Facette::Fac> Mesh::mesh::fac

face container

◆ l

Eigen::Vector3d Mesh::mesh::l

lengths along x,y,z axis

◆ node

std::vector<Nodes::Node> Mesh::mesh::node
private

node container: not initialized by constructor, but later while reading the mesh by member function init_node

◆ node_index

std::vector<int> Mesh::mesh::node_index
private

Index of a node in the node vector. This vector is itself indexed by the node position in the *.msh and *.sol files. In other words, the node found at file_idx in a file is stored as node[node_index[file_idx]].

This is the inverse of the permutation we applied when sorting the nodes.

◆ tet

std::vector<Tetra::Tet> Mesh::mesh::tet

tetrahedron container

◆ vol

double Mesh::mesh::vol

total volume of the mesh


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