dune-localfunctions  2.3.1-rc1
refinedsimplexlocalbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_REFINED_SIMPLEX_LOCALBASIS_HH
4 #define DUNE_REFINED_SIMPLEX_LOCALBASIS_HH
5 
10 #include <dune/common/fvector.hh>
11 #include <dune/common/exceptions.hh>
13 
14 namespace Dune
15 {
16  template<class D, int dim>
18  {
19  protected:
21  {
22  DUNE_THROW(Dune::NotImplemented,"RefinedSimplexLocalBasis not implemented for dim > 3.");
23  }
24  };
25 
33  template<class D>
35  {
36  protected:
37 
40 
51  static int getSubElement(const FieldVector<D,1>& global)
52  {
53  if (global[0] <= 0.5)
54  return 0;
55  else if (global[0] <= 1.0)
56  return 1;
57 
58  DUNE_THROW(InvalidStateException, "no subelement defined");
59  }
60 
67  static void getSubElement(const FieldVector<D,1>& global,
68  int& subElement,
69  FieldVector<D,1>& local)
70  {
71  if (global[0] <= 0.5) {
72  subElement = 0;
73  local[0] = 2.0 * global[0];
74  return;
75  }
76 
77  subElement = 1;
78  local[0] = 2.0 * global[0] - 1.0;
79  }
80 
81  };
82 
83 
94  template<class D>
96  {
97  protected:
98 
101 
116  static int getSubElement(const FieldVector<D,2>& global)
117  {
118  if (global[0] + global[1] <= 0.5)
119  return 0;
120  else if (global[0] >= 0.5)
121  return 1;
122  else if (global[1] >= 0.5)
123  return 2;
124 
125  return 3;
126  }
127 
134  static void getSubElement(const FieldVector<D,2>& global,
135  int& subElement,
136  FieldVector<D,2>& local)
137  {
138  if (global[0] + global[1] <= 0.5) {
139  subElement = 0;
140  local[0] = 2*global[0];
141  local[1] = 2*global[1];
142  return;
143  } else if (global[0] >= 0.5) {
144  subElement = 1;
145  local[0] = 2*global[0]-1;
146  local[1] = 2*global[1];
147  return;
148  } else if (global[1] >= 0.5) {
149  subElement = 2;
150  local[0] = 2*global[0];
151  local[1] = 2*global[1]-1;
152  return;
153  }
154 
155  subElement = 3;
156  local[0] = -2 * global[0] + 1;
157  local[1] = -2 * global[1] + 1;
158 
159  }
160 
161 
162  };
163 
174  template<class D>
176  {
177  protected:
178 
181 
212  static int getSubElement(const FieldVector<D,3>& global)
213  {
214  if (global[0] + global[1] + global[2] <= 0.5)
215  return 0;
216  else if (global[0] >= 0.5)
217  return 1;
218  else if (global[1] >= 0.5)
219  return 2;
220  else if (global[2] >= 0.5)
221  return 3;
222  else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] <= 0.5))
223  return 4;
224  else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] <= 0.5))
225  return 5;
226  else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] >= 0.5))
227  return 6;
228  else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] >= 0.5))
229  return 7;
230 
231  DUNE_THROW(InvalidStateException, "no subelement defined");
232 
233  }
240  static void getSubElement(const FieldVector<D,3>& global,
241  int& subElement,
242  FieldVector<D,3>& local)
243  {
244  if (global[0] + global[1] + global[2] <= 0.5) {
245  subElement = 0;
246  local = global;
247  local *= 2.0;
248  return;
249  } else if (global[0] >= 0.5) {
250  subElement = 1;
251  local = global;
252  local[0] -= 0.5;
253  local *= 2.0;
254  return;
255  } else if (global[1] >= 0.5) {
256  subElement = 2;
257  local = global;
258  local[1] -= 0.5;
259  local *= 2.0;
260  return;
261  } else if (global[2] >= 0.5) {
262  subElement = 3;
263  local = global;
264  local[2] -= 0.5;
265  local *= 2.0;
266  return;
267  } else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] <= 0.5)) {
268  subElement = 4;
269  local[0] = 2.0 * global[1];
270  local[1] = 2.0 * (0.5 - global[0] - global[1]);
271  local[2] = 2.0 * (-0.5 + global[0] + global[1] + global[2]);
272  // Dune::FieldMatrix<double,3,3> A(0.0);
273  // A[0][1] = 2.0;
274  // A[1][0] = -2.0;
275  // A[1][1] = -2.0;
276  // A[2][0] = 2.0;
277  // A[2][1] = 2.0;
278  // A[2][2] = 2.0;
279  // A.mv(global,local);
280  // local[1] += 1.0;
281  // local[2] -= 1.0;
282  return;
283  } else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] <= 0.5)) {
284  subElement = 5;
285  local[0] = 2.0 * (0.5 - global[0]);
286  local[1] = 2.0 * (0.5 - global[1] - global[2]);
287  local[2] = 2.0 * global[2];
288  // Dune::FieldMatrix<double,3,3> A(0.0);
289  // A[0][0] = -2.0;
290  // A[1][1] = -2.0;
291  // A[1][2] = -2.0;
292  // A[2][2] = 2.0;
293  // A.mv(global,local);
294  // local[0] += 1.0;
295  // local[1] += 1.0;
296  return;
297  } else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] >= 0.5)) {
298  subElement = 6;
299  local[0] = 2.0 * (0.5 - global[0] - global[1]);
300  local[1] = 2.0 * global[0];
301  local[2] = 2.0 * (-0.5 + global[1] + global[2]);
302  // Dune::FieldMatrix<double,3,3> A(0.0);
303  // A[0][0] = -2.0;
304  // A[0][1] = -2.0;
305  // A[1][0] = 2.0;
306  // A[2][1] = 2.0;
307  // A[2][2] = 2.0;
308  // A.mv(global,local);
309  // local[0] += 1.0;
310  // local[2] -= 1.0;
311  return;
312  } else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] >= 0.5)) {
313  subElement = 7;
314  local[0] = 2.0 * (-0.5 + global[1] + global[2]);
315  local[1] = 2.0 * (0.5 - global[1]);
316  local[2] = 2.0 * (-0.5 + global[0] + global[1]);
317  // Dune::FieldMatrix<double,3,3> A(0.0);
318  // A[0][1] = 2.0;
319  // A[0][2] = 2.0;
320  // A[1][1] = -2.0;
321  // A[2][0] = 2.0;
322  // A[2][1] = 2.0;
323  // A.mv(global,local);
324  // local[0] -= 1.0;
325  // local[1] += 1.0;
326  // local[2] -= 1.0;
327  return;
328  }
329 
330  DUNE_THROW(InvalidStateException, "no subelement defined");
331 
332  }
333 
334  };
335 
336 
337 }
338 
339 #endif
RefinedSimplexLocalBasis()
Definition: refinedsimplexlocalbasis.hh:20
static void getSubElement(const FieldVector< D, 3 > &global, int &subElement, FieldVector< D, 3 > &local)
Get local coordinates in the subsimplex.
Definition: refinedsimplexlocalbasis.hh:240
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:100
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:39
static void getSubElement(const FieldVector< D, 2 > &global, int &subElement, FieldVector< D, 2 > &local)
Get local coordinates in the subtriangle.
Definition: refinedsimplexlocalbasis.hh:134
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:180
static int getSubElement(const FieldVector< D, 1 > &global)
Get the number of the subelement containing a given point.
Definition: refinedsimplexlocalbasis.hh:51
Definition: refinedsimplexlocalbasis.hh:17
static void getSubElement(const FieldVector< D, 1 > &global, int &subElement, FieldVector< D, 1 > &local)
Get local coordinates in the subelement.
Definition: refinedsimplexlocalbasis.hh:67
static int getSubElement(const FieldVector< D, 2 > &global)
Get the number of the subtriangle containing a given point.
Definition: refinedsimplexlocalbasis.hh:116
static int getSubElement(const FieldVector< D, 3 > &global)
Get the number of the subsimplex containing a given point in the reference element.
Definition: refinedsimplexlocalbasis.hh:212