Toast++  1.0.2 (r.539)
Forward and inverse modelling in optical tomography
tri3old.h
1 // -*-C++-*-
2 // ==========================================================================
3 // Module libfe
4 // File tri3.h
5 // Declaration of unstructured 2-D element class Triangle3old
6 // (3-noded triangle, 1st order shape functions)
7 //
8 // Inheritance:
9 // ------------
10 // Element
11 // ---> Element_Unstructured
12 // ---> Element_Unstructured_2D
13 // ---> Triangle3old
14 // ==========================================================================
15 
16 #ifndef __TRI3OLD_H
17 #define __TRI3OLD_H
18 
19 #include "toastdef.h"
20 
21 class Surface;
22 
24 public:
25 
26  Triangle3old();
27  Triangle3old (const Triangle3old& el);
28  ~Triangle3old();
29  // constructors, destructor
30 
34  Element *Copy();
35 
36  void Initialise (const NodeList& nlist);
37 
38  BYTE Type () const { return ELID_TRI3OLD; }
39  // returns element type id
40 
41  unsigned long GetCaps () const { return ELCAPS_SUBSAMPLING; }
42  // Return element capability flags
43 
44  int nNode() const { return 3; }
45  // returns number of nodes per element
46 
47  int nSide () const { return 3; };
48  // number of sides per element
49 
50  int nSideNode (int /*side*/) const { return 2; };
51  // returns number of nodes attached to side 'side'
52 
53  int SideNode (int side, int node) const;
54  // returns number of 'node'-th node of side 'side'
55 
56  double SideSize (int sd, const NodeList &nlist) const;
57  // returns length of side 'sd'
58 
59  Point Local (const NodeList& nlist, const Point& glob) const;
60  // returns the local coordinate corresponding to global coordinate 'glob'
61 
62  Point NodeLocal (int node) const;
63  // returns local coordinates of node 'node'
64 
65  Point SurfToLocal (int side, const Point &p) const;
66  // Map 1-D surface point p (0-1) along side 'side' into 2-D local coords
67 
68  RVector DirectionCosine (int side, RDenseMatrix &jacin);
69  // returns a vector of direction cosines in global coordinates of the
70  // normal to side 'side'.
71 
72  const RVector &LNormal (int side) const;
73 
74  bool LContains (const Point& loc, bool pad = true) const;
75  // returns TRUE if point 'loc' (in local coordinates of the element) is
76  // within the standard triangle
77 
78  bool GContains (const Point& glob, const NodeList& nlist) const;
79  // returns TRUE if point 'glob' (in global coordinates) is within the
80  // element. nlist is the mesh node coordinate list
81 
82  RVector LocalShapeF (const Point &loc) const;
83  // vector of shape functions u_i(loc) at point 'loc' given in local
84  // element coordinates
85 
86  RDenseMatrix LocalShapeD (const Point &loc) const;
87  // 2x3 matrix of shape function derivatives (d u_j)/(d x_i) at point 'loc'
88  // given in local element coordinates
89 
90  RVector GlobalShapeF (const NodeList& nlist, const Point& glob) const;
91  // returns the shape functions for each node at a point, given in global
92  // coordinates
93 
94  RDenseMatrix GlobalShapeD (const NodeList& nlist, const Point& glob) const;
95  // returns derivatives of shape functions in global coordinates at point
96  // `glob', given in global coordinates
97 
98  double IntF (int i) const;
99 
100  RSymMatrix IntFF () const;
101  // Return integral over element of product of shape functions:
102  // FF = Int_el { F_i(r) F_j(r) } dr
103 
104  int QuadRule (int order, const double **wght, const Point **absc) const;
105 
106  double IntFF (int i, int j) const;
107  // Return a single element of IntFF
108 
109  double IntFFF (int i, int j, int k) const;
110  // returns a single element of integral over element of product of three
111  // shape functions:
112  // IntFFF = Int_el { F_i(r) * F_j(r) * F_k(r) } dr
113 
114 #ifdef UNDEF
115  // MS 29.6.99 Removed because incompatible with quadratic shape functions
116 
117  void IntFFF (double &iii, double &iij, double &ijk) const;
118  // returns the values of the FFF tensor for:
119  // all indices equal (iii)
120  // two indices equal (iij)
121  // all indices different (ijk)
122 #endif
123 
124  RSymMatrix IntPFF (const RVector& P) const;
125  // Returns integral over element of product of two shape functions and a
126  // function P defined on nodes:
127  // PFF = Int_el { P(r) F_i(r) F_J(r) } dr
128  // Vector P contains the complete nodal solution for the mesh
129 
130  double IntPFF (int i, int j, const RVector& P) const;
131  // Returns a single element of IntPFF
132 
133  RVector IntFD (int i, int j) const;
134  // Returns single (vector) element of FD integral over element:
135  // IntFD = Int_el { F_i(r) D_j(r) } dr
136 
137  double IntFDD (int i, int j, int k) const
138  { return intdd.Get(j,k) * 0.3333333333; }
139  // returns a single element of integral over element:
140  // IntFDD = Int_el { F_i(r) * D_j(r) * D_k(r) } dr
141 
142  RSymMatrix IntPDD (const RVector& P) const
143  { return intdd * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
144  // Returns integral over element of product of two shape derivatives and a
145  // function P defined in nodal basis:
146  // PDD = Int_el { P(r) D_i(r) D_j(r) } dr
147  // Vector P contains the complete nodal solution for the mesh
148 
149  double IntPDD (int i, int j, const RVector &P) const
150  { return intdd.Get(i,j) * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
151  // Returns a single element of IntPDD
152 
153  double BndIntFFSide (int i, int j, int sd);
154 
155  RSymMatrix BndIntPFF (const RVector &P) const
156  { return intbff * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
157  // This is a hack! Really P would have to be taken into the integral,
158  // rather than just averaged.
159 
160  double BndIntPFF (int i, int j, const RVector &P) const
161  { return intbff.Get(i,j) * ((P[Node[0]]+P[Node[1]]+P[Node[2]])/3.0); }
162  // Returns a single element of BndIntPFF
163 
164  RSymMatrix Intdd() const;
165  // returns matrix of mixed derivatives
166 
167  RVector BndIntFX (int side, double (*func)(const Point&),
168  const NodeList &nlist) const;
169  // See element.h for description
170 
171  RVector BndIntFCos (int side, const Surface *surf, const RVector &cntcos,
172  double sigma, double sup, const NodeList &nlist) const;
173  // See element.h for description
174 
175  RVector BndIntFCos (int side, const RVector &cntcos, double a,
176  const NodeList &nlist) const;
177  // See element.h for description
178 
179  RVector BndIntFDelta (int side, const Surface *surf, const RVector &pos,
180  const NodeList &nlist) const;
181  // See element.h for description
182 
183  int GetLocalSubsampleAbsc (const Point *&absc) const;
184  // returns abscissae in local coords for numerical integration by uniform
185  // subsampling. Return value is the number of samples returned in 'absc'
186 
187  int GetBndSubsampleAbsc (int side, const Point *&absc) const;
188  // abscissae for numerical integration over boundary side in local
189  // coordinates of boundary element (dim-1). Use SurfToLocal to convert
190  // into local element coordinates
191 
192  int GlobalIntersection (const NodeList &nlist, const Point &p1,
193  const Point &p2, Point **list);
194  // Same as 'Intersection' but p1 and p2 given in global coords
195  // The return list however is still in local coords
196 
197  int Intersection (const Point &p1, const Point &p2, Point** pi);
198  // creates a list of points where the line defined by p1 and p2 intersects
199  // the element (in local coordinates) or starts/ends within the element.
200  // Returns the length of the list
201 
202 protected:
203 
204  double ComputeSize (const NodeList &nlist) const;
205  // area of triangle in global coordinates
206 
207  RSymMatrix ComputeIntDD (const NodeList &nlist) const;
208  // Returns integral over element of product of shape derivatives:
209  // DD = Int_el { D_i(r) D_j(r) } dr
210 
211  void ComputeIntFD (const NodeList &nlist);
212  // Calculates the FD integrals over the element:
213  // IntFD = Int_el { F_i(r) D_j(r) } dr
214 
215  RSymMatrix ComputeBndIntFF (const NodeList& nlist) const;
216  // Returns line integral of product of two shape functions along sides of
217  // the element which belong to the mesh boundary
218  // If the element does not contain any boundary sides then the returned
219  // matrix is empty.
220 
221 private:
222 
223  RDenseMatrix jac;
224  // 2x2 matrix of derivatives of global coordinates with respect to local
225  // ones
226 
227  double djac;
228  // determinant of jac
229 
230  double *intfd_0, *intfd_1;
231  // IntFD matrix storage
232 
233 #ifdef TRI3_STORE_INTFF
234  RSymMatrix intff;
235  // stores integral over element of product of shape functions
236  // Int_el { F_i(r) F_j(r) } dr
237  // set by Initialise
238 #endif
239 
240 #ifdef TRI3_STORE_COORDS
241  double n0x, n0y, n1x, n1y, n2x, n2y;
242 #endif
243 };
244 
245 #endif
Definition: ndlist.h:14
double IntFDD(int i, int j, int k) const
Integral of a product of a shape function and two shape function derivatives over the element...
Definition: tri3old.h:137
RSymMatrix BndIntPFF(const RVector &P) const
Surface integrals of all products of a nodal function and two shape functions over all boundary sides...
Definition: tri3old.h:155
Templated vector class.
Definition: vector.h:39
virtual bool GContains(const Point &glob, const NodeList &nlist) const
Checks if a global point coordinate is inside the element.
#define ELID_TRI3OLD
old-style 3-noded triangle
Definition: element.h:31
virtual RSymMatrix IntFF() const =0
Integrals of all products of two shape functions over the element.
unsigned long GetCaps() const
Returns element capability flags.
Definition: tri3old.h:41
Definition: node.h:39
Definition: surface.h:21
virtual RVector DirectionCosine(int side, RDenseMatrix &jacin)=0
Returns the direction cosines of a side normal.
virtual RVector GlobalShapeF(const NodeList &nlist, const Point &glob) const
Returns the values of the shape functions at a global point.
Definition: element.h:490
virtual double IntFFF(int i, int j, int k) const =0
Integral of a product of three shape functions over the element.
int nSide() const
Returns the number of element sides.
Definition: tri3old.h:47
Definition: point.h:18
virtual double BndIntFFSide(int i, int j, int sd)=0
Surface integral of a product of two shape functions over one of the sides of the element...
virtual RDenseMatrix LocalShapeD(const Point &loc) const =0
Returns the values of the shape function derivatives at a local point.
virtual const RVector & LNormal(int side) const =0
Returns a side normal in local coordinates.
virtual RSymMatrix IntPFF(const RVector &P) const =0
Integrals of all products of two shape functions and a nodal function over the element.
virtual void Initialise(const NodeList &nlist)
Element initialisation.
virtual RDenseMatrix GlobalShapeD(const NodeList &nlist, const Point &glob) const
Returns the values of the shape function derivatives at a global point.
Definition: element.h:502
RSymMatrix IntPDD(const RVector &P) const
All integrals of products of a nodal function and two shape function derivatives over the element...
Definition: tri3old.h:142
virtual RVector IntFD(int i, int j) const
Integral of a product of a shape function and a shape function derivative over the element...
Definition: element.h:637
Base class for finite element types.
Definition: element.h:84
double IntPDD(int i, int j, const RVector &P) const
Integrals of a product of a nodal function and two shape function derivatives over the element...
Definition: tri3old.h:149
Base class for all 2-D unstructured element types.
Definition: element.h:1234
virtual double IntF(int i) const =0
Integral of a shape function over the element.
double BndIntPFF(int i, int j, const RVector &P) const
Surface integrals of a product of a nodal function and two shape functions over all boundary sides of...
Definition: tri3old.h:160
virtual double SideSize(int side, const NodeList &nlist) const
Returns the size of an element side.
Definition: element.h:397
virtual RVector LocalShapeF(const Point &loc) const =0
Returns the values of the shape functions at a local point.
virtual Point SurfToLocal(int side, const Point &p) const
Maps a point from surface coordinates to local element coordinates.
Definition: element.h:302
virtual RSymMatrix Intdd() const
Integral of the product of two partial shape function derivatives over the element.
Definition: element.h:836
Definition: tri3old.h:23
virtual Element * Copy()=0
Create a copy of the element and return a pointer to it.
Dense matrix class.
Definition: crmatrix.h:38
int nNode() const
Returns the number of nodes associated with the element.
Definition: tri3old.h:44
virtual Point Local(const NodeList &nlist, const Point &glob) const =0
Maps a point from global to local element coordinates.
int nSideNode(int) const
Returns the number of vertices associated with a side.
Definition: tri3old.h:50
virtual Point NodeLocal(int node) const =0
Returns the local coordinates of an element node.
#define ELCAPS_SUBSAMPLING
element implements IntFFF and IntFDD by subsampling
Definition: element.h:61
BYTE Type() const
Returns an element type identifier.
Definition: tri3old.h:38
virtual int SideNode(int side, int node) const =0
Returns relative node index for a side vertex.
virtual int QuadRule(int order, const double **wght, const Point **absc) const
Returns the weights and abscissae of quadrature rules over the element.
Definition: element.h:521
virtual bool LContains(const Point &loc, bool pad=true) const =0
Checks if a local point coordinate is inside the element.