Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
piola_util.h
Go to the documentation of this file.
1#ifndef _RHEOLEF_PIOLA_UTIL_H
2#define _RHEOLEF_PIOLA_UTIL_H
23//
24// piola transformation:
25//
26// F : hat_K ---> K
27// hat_x +--> x = F(hat_x)
28//
29// i.e. map any geo_element K from a reference_element hat_K
30//
31#include "rheolef/geo.h"
32#include "rheolef/basis_on_pointset.h"
33namespace rheolef {
34
35// ------------------------------------------
36// piola transformation:
37// x[q] = F(hat_x[q])
38// on all nodes of a quadrature formulae
39// ------------------------------------------
40template<class T, class M>
41void
43 const geo_basic<T,M>& omega,
44 const basis_on_pointset<T>& piola_on_pointset,
45 reference_element hat_K,
46 const std::vector<size_t>& dis_inod,
47 Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& x);
48
49// ------------------------------------------
50// jacobian of the piola transformation
51// DF(hat_x[q])
52// at an aritrarily point hat_x
53// ------------------------------------------
54template<class T, class M>
55void
57 const geo_basic<T,M>& omega,
58 const basis_basic<T>& piola_basis,
59 reference_element hat_K,
60 const std::vector<size_t>& dis_inod,
61 const point_basic<T>& hat_x,
62 tensor_basic<T>& DF);
63
64// ------------------------------------------
65// jacobian of the piola transformation
66// DF(hat_x[q])
67// on all nodes of a quadrature formulae
68// ------------------------------------------
69template<class T, class M>
70void
72 const geo_basic<T,M>& omega,
73 const basis_on_pointset<T>& piola_on_pointset,
74 reference_element hat_K,
75 const std::vector<size_t>& dis_inod,
76 Eigen::Matrix<tensor_basic<T>,Eigen::Dynamic,1>& DF);
77
78
79template <class T>
80T
81det_jacobian_piola_transformation (const tensor_basic<T>& DF, size_t d , size_t map_d);
82
83template<class T, class M>
84point_basic<T>
85normal_from_piola_transformation (const geo_basic<T,M>& omega, const geo_element& S, const tensor_basic<T>& DF, size_t d);
86
87// The pseudo inverse extend inv(DF) for face in 3d or edge in 2d
88// i.e. useful for Laplacian-Beltrami and others surfacic forms.
89//
90// pinvDF (hat_xq) = inv DF, if tetra in 3d, tri in 2d, etc
91// = pseudo-invese, when tri in 3d, edge in 2 or 3d
92// e.g. on 3d face : pinvDF*DF = [1, 0, 0; 0, 1, 0; 0, 0, 0]
93//
94// let DF = [u, v, w], where u, v, w are the column vectors of DF
95// then det(DF) = mixt(u,v,w)
96// and det(DF)*inv(DF)^T = [v^w, w^u, u^v] where u^v = vect(u,v)
97//
98// application:
99// if K=triangle(a,b,c) then u=ab=b-a, v=ac=c-a and w = n = u^v/|u^v|.
100// Thus DF = [ab,ac,n] and det(DF)=|ab^ac|
101// and inv(DF)^T = [ac^n/|ab^ac|, -ab^n/|ab^ac|, n]
102// The pseudo-inverse is obtained by remplacing the last column n by zero.
103//
104template<class T>
105tensor_basic<T>
107 const tensor_basic<T>& DF,
108 size_t d,
109 size_t map_d);
110
111// F^{-1}: x --> hat_x on K
112template<class T, class M>
113point_basic<T>
115 const geo_basic<T,M>& omega,
116 const reference_element& hat_K,
117 const std::vector<size_t>& dis_inod,
118 const point_basic<T>& x);
119
120// compute: P = I - nxn, the tangential projector on a map with unit normal n
121template <class T>
122void map_projector (const tensor_basic<T>& DF, size_t d, size_t map_d, tensor_basic<T>& P);
123
124// axisymetric weight ?
125// point_basic<T> xq = rheolef::piola_transformation (_omega, _piola_table, K, dis_inod, q);
126template<class T>
127T
128weight_coordinate_system (space_constant::coordinate_type sys_coord, const point_basic<T>& xq);
129
130// -------------------------------------------
131// weight integration: w = det_DF*wq
132// with optional axisymmetric r*dr factor
133// -------------------------------------------
134template<class T, class M>
135void
137 const geo_basic<T,M>& omega,
138 const basis_on_pointset<T>& piola_on_pointset,
139 reference_element hat_K,
140 const std::vector<size_t>& dis_inod,
141 bool ignore_sys_coord,
142 Eigen::Matrix<tensor_basic<T>,Eigen::Dynamic,1>& DF,
143 Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& x,
144 Eigen::Matrix<T,Eigen::Dynamic,1>& w);
145
146}// namespace rheolef
147#endif // _RHEOLEF_PIOLA_UTIL_H
Expr1::float_type T
Definition field_expr.h:230
This file is part of Rheolef.
point_basic< T > normal_from_piola_transformation(const geo_basic< T, M > &omega, const geo_element &S, const tensor_basic< T > &DF, size_t d)
void piola_transformation(const geo_basic< T, M > &omega, const basis_on_pointset< T > &piola_on_pointset, reference_element hat_K, const std::vector< size_t > &dis_inod, Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &x)
Definition piola_util.cc:43
tensor_basic< T > pseudo_inverse_jacobian_piola_transformation(const tensor_basic< T > &DF, size_t d, size_t map_d)
void map_projector(const tensor_basic< T > &DF, size_t d, size_t map_d, tensor_basic< T > &P)
T weight_coordinate_system(space_constant::coordinate_type sys_coord, const point_basic< T > &xq)
void piola_transformation_and_weight_integration(const geo_basic< T, M > &omega, const basis_on_pointset< T > &piola_on_quad, reference_element hat_K, const std::vector< size_t > &dis_inod, bool ignore_sys_coord, Eigen::Matrix< tensor_basic< T >, Eigen::Dynamic, 1 > &DF, Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &x, Eigen::Matrix< T, Eigen::Dynamic, 1 > &w)
void jacobian_piola_transformation(const geo_basic< T, M > &omega, const basis_on_pointset< T > &piola_on_pointset, reference_element hat_K, const std::vector< size_t > &dis_inod, Eigen::Matrix< tensor_basic< T >, Eigen::Dynamic, 1 > &DF)
Definition piola_util.cc:80
T det_jacobian_piola_transformation(const tensor_basic< T > &DF, size_t d, size_t map_d)
point_basic< T > inverse_piola_transformation(const geo_basic< T, M > &omega, const reference_element &hat_K, const std::vector< size_t > &dis_inod, const point_basic< T > &x)