Toast++  1.0.2 (r.539)
Forward and inverse modelling in optical tomography
tri3D3.h
1 // -*-C++-*-
2 // ==========================================================================
3 // Module libfe
4 // File tri3D3.h
5 // Declaration of unstructured 2-D element class Triangle3D3
6 // (3-noded triangle, 1st order shape functions, nodes are three dimensional)
7 //
8 //
9 // Node arrangement: N2
10 // ----------------- + The local element has node coordinates
11 // |\ N0 = (0,0)
12 // | \ N1 = (1,0)
13 // | \ N2 = (0,1)
14 // | \
15 // | \ Sides: side 0 contains N0,N1 (y=0)
16 // | \ side 1 contains N1,N2 (x+y=1)
17 // | \ side 2 contains N2,N0 (x=0)
18 // N0+-------+N1
19 //
20 // Inheritance:
21 // ------------
22 // Element
23 // ---> Element_Unstructured
24 // ---> Element_Unstructured_2D
25 // ---> Triangle3
26 // ---> Triangle3D3
27 // ==========================================================================
28 
29 #ifndef __TRI3D3_H
30 #define __TRI3D3_H
31 
32 #include "toastdef.h"
33 
34 class Surface;
35 
36 class FELIB Triangle3D3: public Triangle3 {
37 public:
38 
39  Triangle3D3();
40  Triangle3D3 (const Triangle3D3& el);
41  ~Triangle3D3();
42  // constructors, destructor
43 
44  void Initialise (const NodeList& nlist);
45 
46  BYTE Type () const { return ELID_TRI3D3; }
47  // returns element type id
48 
49  int Dimension () const { return 3; }
50 
51  // unsigned long GetCaps () const { return ELCAPS_SUBSAMPLING; }
52  // Return element capability flags
53 
54  // int nNode() const { return 3; }
55  // returns number of nodes per element
56 
57  // int nSide () const { return 3; };
58  // number of sides per element
59 
60  // int nSideNode (int /*side*/) const { return 2; };
61  // returns number of nodes attached to side 'side'
62 
63  // int SideNode (int side, int node) const;
64  // returns number of 'node'-th node of side 'side'
65 
66  // double SideSize (int sd, const NodeList &nlist) const;
67  // returns length of side 'sd'
68 
69  Point Local (const NodeList& nlist, const Point& glob) const;
70  // returns the local coordinate corresponding to global coordinate 'glob'
71 
72  // Point NodeLocal (int node) const;
73  // returns local coordinates of node 'node'
74 
75  Point SurfToLocal (int side, const Point &p) const;
76  // Map 1-D surface point p (0-1) along side 'side' into 2-D local coords
77 
78  RVector DirectionCosine (int side, RDenseMatrix &jacin);
79  // returns a vector of direction cosines in global coordinates of the
80  // normal to side 'side'.
81 
82  // const RVector &LNormal (int side) const;
83 
84  // bool LContains (const Point& loc, bool pad = true) const;
85  // returns TRUE if point 'loc' (in local coordinates of the element) is
86  // within the standard triangle
87 
88  bool GContains (const Point& glob, const NodeList& nlist) const;
89  // returns TRUE if point 'glob' (in global coordinates) is within the
90  // element. nlist is the mesh node coordinate list
91 
92  // RVector LocalShapeF (const Point &loc) const;
93  // vector of shape functions u_i(loc) at point 'loc' given in local
94  // element coordinates
95 
96  // RDenseMatrix LocalShapeD (const Point &loc) const;
97  // 2x3 matrix of shape function derivatives (d u_j)/(d x_i) at point 'loc'
98  // given in local element coordinates
99 
100  RVector GlobalShapeF (const NodeList& nlist, const Point& glob) const;
101  // returns the shape functions for each node at a point, given in global
102  // coordinates
103 
104  RDenseMatrix GlobalShapeD (const NodeList& nlist, const Point& glob) const;
105  // returns derivatives of shape functions in global coordinates at point
106  // `glob', given in global coordinates
107 
108  // double IntF (int i) const;
109 
110  // RSymMatrix IntFF () const;
111  // Return integral over element of product of shape functions:
112  // FF = Int_el { F_i(r) F_j(r) } dr
113 
114  // int QuadRule (int order, const double **wght, const Point **absc) const;
115 
116  // double IntFF (int i, int j) const;
117  // Return a single element of IntFF
118 
119  // double IntFFF (int i, int j, int k) const;
120  // returns a single element of integral over element of product of three
121  // shape functions:
122  // IntFFF = Int_el { F_i(r) * F_j(r) * F_k(r) } dr
123 
124 #ifdef UNDEF
125  // MS 29.6.99 Removed because incompatible with quadratic shape functions
126 
127  // void IntFFF (double &iii, double &iij, double &ijk) const;
128  // returns the values of the FFF tensor for:
129  // all indices equal (iii)
130  // two indices equal (iij)
131  // all indices different (ijk)
132 #endif
133 
134  // RSymMatrix IntPFF (const RVector& P) const;
135  // Returns integral over element of product of two shape functions and a
136  // function P defined on nodes:
137  // PFF = Int_el { P(r) F_i(r) F_J(r) } dr
138  // Vector P contains the complete nodal solution for the mesh
139 
140  // double IntPFF (int i, int j, const RVector& P) const;
141  // Returns a single element of IntPFF
142 
143  // RVector IntFD (int i, int j) const;
144  // Returns single (vector) element of FD integral over element:
145  // IntFD = Int_el { F_i(r) D_j(r) } dr
146 
147  // double IntFDD (int i, int j, int k) const
148  // { return intdd.Get(j,k) * 0.3333333333; }
149  // returns a single element of integral over element:
150  // IntFDD = Int_el { F_i(r) * D_j(r) * D_k(r) } dr
151 
152  // RSymMatrix IntPDD (const RVector& P) const
153  // { return intdd * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
154  // Returns integral over element of product of two shape derivatives and a
155  // function P defined in nodal basis:
156  // PDD = Int_el { P(r) D_i(r) D_j(r) } dr
157  // Vector P contains the complete nodal solution for the mesh
158 
159  // double IntPDD (int i, int j, const RVector &P) const
160  // { return intdd.Get(i,j) * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
161  // Returns a single element of IntPDD
162 
163  // RSymMatrix BndIntPFF (const RVector &P) const
164  // { return intbff * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
165  // This is a hack! Really P would have to be taken into the integral,
166  // rather than just averaged.
167 
168  // double BndIntPFF (int i, int j, const RVector &P) const
169  // { return intbff.Get(i,j) * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
170  // Returns a single element of BndIntPFF
171 
172  // double IntFd (int i, int j, int k) const;
173  // Returns a single element of Int [u_i du_j/dx_k] dr
174 
175  // double IntPd (const RVector &P, int j, int k) const;
176  // Returns Int [p(r) du_j/dx_k] dr where p(r) is defined by its nodal
177  // basis coefficients P.
178 
179  // double IntFdd (int i, int j, int k, int l, int m) const;
180  // Int u_i du_j/dx_l du_k/dx_m dr
181 
182  // double IntPdd (const RVector &p, int j, int k, int l, int m) const;
183  // Int f(r) du_j/dx_l du_k/dx_m dr
184  // where f(r) is given as a nodal vector
185 
186  // RSymMatrix Intdd() const;
187  // returns matrix of mixed derivatives
188 
189  // double BndIntFFSide (int i, int j, int sd);
190  // Int [u_i u_j] dr along side sd
191 
192  // RVector BndIntFX (int side, double (*func)(const Point&),
193  // const NodeList &nlist) const;
194  // See element.h for description
195 
196 // RVector BndIntFCos (int side, const Surface *surf, const RVector &cntcos,
197 // double sigma, double sup, const NodeList &nlist) const;
198  // See element.h for description
199 
200 // RVector BndIntFCos (int side, const RVector &cntcos, double a,
201 // const NodeList &nlist) const;
202  // See element.h for description
203 
204 
205 // RVector BndIntFDelta (int side, const Surface *surf, const RVector &pos,
206 // const NodeList &nlist) const;
207  // See element.h for description
208 
209  // int GetLocalSubsampleAbsc (const Point *&absc) const;
210  // returns abscissae in local coords for numerical integration by uniform
211  // subsampling. Return value is the number of samples returned in 'absc'
212 
213  // int GetBndSubsampleAbsc (int side, const Point *&absc) const;
214  // abscissae for numerical integration over boundary side in local
215  // coordinates of boundary element (dim-1). Use SurfToLocal to convert
216  // into local element coordinates
217 
218  int GlobalIntersection (const NodeList &nlist, const Point &p1,
219  const Point &p2, Point **list);
220  // Same as 'Intersection' but p1 and p2 given in global coords
221  // The return list however is still in local coords
222 
223  // int Intersection (const Point &p1, const Point &p2, Point** pi);
224  // creates a list of points where the line defined by p1 and p2 intersects
225  // the element (in local coordinates) or starts/ends within the element.
226  // Returns the length of the list
227  RVector LocalShapeQF (const Point &loc) const;
228  RDenseMatrix LocalShapeQD (const Point &loc) const;
229  double IntUnitSphereP (const NodeList& nlist, const RVector& P) const;
230  double IntUnitSpherePFF(const NodeList& nlist, const int i, const int j, const RVector& P) const;
231  double IntUnitSphereFF (const NodeList& nlist, const int i, const int j) const;
232 
233  RDenseMatrix LocaltoGlobalMat () const
234  { return MapLocaltoGlobal; }
235 
236  RDenseMatrix GlobaltoLocalMat () const
237  { return inverse(MapLocaltoGlobal); }
238 
239  RDenseMatrix FTAMat () const ;
240 
241 
242 protected:
243 
244  // double ComputeSize (const NodeList &nlist) const;
245  // area of triangle in global coordinates
246 
247  // RSymMatrix ComputeIntDD (const NodeList &nlist) const;
248  // Returns integral over element of product of shape derivatives:
249  // DD = Int_el { D_i(r) D_j(r) } dr
250 
251  // void ComputeIntFD (const NodeList &nlist);
252  // Calculates the FD integrals over the element:
253  // IntFD = Int_el { F_i(r) D_j(r) } dr
254 
255  // RSymMatrix ComputeBndIntFF (const NodeList& nlist) const;
256  // Returns line integral of product of two shape functions along sides of
257  // the element which belong to the mesh boundary
258  // If the element does not contain any boundary sides then the returned
259  // matrix is empty.
260 
261 private:
262 
263  // double *intfd_0, *intfd_1;
264  // IntFD matrix storage
265 
266  // double a0, b0, c0, a1, b1, c1, a2, b2, c2;
267  // triangle geometry parameters
268 
269 #ifdef TRI3_STORE_INTFF
270  // RSymMatrix intff;
271  // stores integral over element of product of shape functions
272  // Int_el { F_i(r) F_j(r) } dr
273  // set by Initialise
274 #endif
275  RDenseMatrix Map2Dto3D;
276  RDenseMatrix Map3Dto2D;
277  RDenseMatrix MapLocaltoGlobal;
278  RVector *V2D;
279  RVector uspts[6]; // used by unit sphere integration
280 };
281 
282 #endif
Definition: ndlist.h:14
Definition: tri3D3.h:36
Templated vector class.
Definition: vector.h:39
void Initialise(const NodeList &nlist)
Initialises the triangle.
Point Local(const NodeList &nlist, const Point &glob) const
Maps a point from global to local element coordinates.
Definition: surface.h:21
RDenseMatrix GlobalShapeD(const NodeList &nlist, const Point &glob) const
Returns the values of the shape function derivatives at a global point.
Definition: point.h:18
RVector GlobalShapeF(const NodeList &nlist, const Point &glob) const
Returns the values of the shape functions at a global point.
BYTE Type() const
Returns an element type identifier.
Definition: tri3D3.h:46
A 3-noded (linear) 2-dimensional triangle element.
Definition: tri3.h:42
Point SurfToLocal(int side, const Point &p) const
Maps a point from surface coordinates to local element coordinates.
Dense matrix class.
Definition: crmatrix.h:38
#define ELID_TRI3D3
3-noded surface triangle
Definition: element.h:46
bool GContains(const Point &glob, const NodeList &nlist) const
Checks if a global point coordinate is inside the element.
RVector DirectionCosine(int side, RDenseMatrix &jacin)
Returns the direction cosines of a side normal.
int Dimension() const
Returns the spatial dimension of the element.
Definition: tri3D3.h:49