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

#include <linear_algebra.h>

Public Member Functions

 LinAlgebra (Settings &s, Mesh::mesh &my_msh)
 
void prepareElements (Eigen::Vector3d const &Hext, timing const &t_prm)
 
void prepareElements (double const A_Hext, timing const &t_prm)
 
void buildInitGuess (Eigen::Ref< Eigen::VectorXd > G) const
 
int solver (timing const &t_prm)
 
void set_DW_vz (double vz)
 
double get_v_max (void)
 
void setExtSpaceField (Settings &s)
 
void base_projection ()
 

Private Attributes

Nodes::index idx_dir
 
const int NOD
 
const int MAXITER
 
const double TOL
 
double ILU_tol
 
double ILU_fill_factor
 
const int verbose
 
const std::vector< Tetra::prm > & prmTetra
 
const std::vector< Facette::prm > & prmFacette
 
Mesh::meshrefMsh
 
double DW_vz
 
double v_max
 
std::vector< Eigen::Matrix< double, Nodes::DIM, Tetra::NPI > > extSpaceField
 

Detailed Description

convenient class to grab altogether some part of the calculations involved using eigen BiCGSTAB solver at each timestep. The solver is handled by solver method, and is using Eigen::SparseMatrix, Row major matrix. This matrix is prepared in 'batch mode', using a vector of triplets (also called COO write sparse matrix).

Constructor & Destructor Documentation

◆ LinAlgebra()

LinAlgebra::LinAlgebra ( Settings s,
Mesh::mesh my_msh 
)
inline

constructor

Parameters
[in]s
[in]my_msh

Member Function Documentation

◆ base_projection()

void LinAlgebra::base_projection ( )

computes local vector basis {ep,eq} in the tangeant plane for projection on the elements

◆ buildInitGuess()

void LinAlgebra::buildInitGuess ( Eigen::Ref< Eigen::VectorXd >  G) const

build init guess for bicgstab solver

◆ get_v_max()

double LinAlgebra::get_v_max ( void  )
inline

getter for v_max

◆ prepareElements() [1/2]

void LinAlgebra::prepareElements ( double const  A_Hext,
timing const &  t_prm 
)

computes inner data structures of tetraedrons and triangular facettes (K matrices and L vectors)

Parameters
[in]A_Hextamplitude applied field
[in]t_prm

◆ prepareElements() [2/2]

void LinAlgebra::prepareElements ( Eigen::Vector3d const &  Hext,
timing const &  t_prm 
)

computes inner data structures of tetraedrons and triangular facettes (K matrices and L vectors)

Parameters
[in]Hextapplied field
[in]t_prm

◆ set_DW_vz()

void LinAlgebra::set_DW_vz ( double  vz)
inline

setter for DW_dz

Parameters
[in]vz

◆ setExtSpaceField()

void LinAlgebra::setExtSpaceField ( Settings s)

when external applied field is of field_type R4toR3 values of field_space are stored in spaceField

Parameters
[in]s

◆ solver()

int LinAlgebra::solver ( timing const &  t_prm)

solver, uses eigen stabilized biconjugate gradient solver (bicgstab) with ILU preconditionner, sparse matrix and vector are filled with multiThreading. Sparse matrix is row major.

Parameters
[in]t_prm

Member Data Documentation

◆ DW_vz

double LinAlgebra::DW_vz
private

speed of the domain wall

◆ extSpaceField

std::vector< Eigen::Matrix<double,Nodes::DIM,Tetra::NPI> > LinAlgebra::extSpaceField
private

external applied space field, values on gauss points, size is number of tetraedrons

◆ idx_dir

Nodes::index LinAlgebra::idx_dir
private

recentering index direction if any

◆ ILU_fill_factor

double LinAlgebra::ILU_fill_factor
private

ILU preconditionner filling factor

◆ ILU_tol

double LinAlgebra::ILU_tol
private

ILU preconditionner tolerance

◆ MAXITER

const int LinAlgebra::MAXITER
private

maximum number of iteration for bicgstab

◆ NOD

const int LinAlgebra::NOD
private

number of nodes, also an offset for filling sparseMatrix, initialized by constructor

◆ prmFacette

const std::vector<Facette::prm>& LinAlgebra::prmFacette
private

material parameters of the facettes

◆ prmTetra

const std::vector<Tetra::prm>& LinAlgebra::prmTetra
private

material parameters of the tetrahedrons

◆ refMsh

Mesh::mesh* LinAlgebra::refMsh
private

direct access to the mesh

◆ TOL

const double LinAlgebra::TOL
private

solver tolerance

◆ v_max

double LinAlgebra::v_max
private

maximum speed of the magnetization in the whole physical object

◆ verbose

const int LinAlgebra::verbose
private

verbosity


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