3 #ifndef DUNE_REFINED_SIMPLEX_LOCALBASIS_HH 4 #define DUNE_REFINED_SIMPLEX_LOCALBASIS_HH 10 #include <dune/common/fvector.hh> 11 #include <dune/common/exceptions.hh> 16 template<
class D,
int dim>
22 DUNE_THROW(Dune::NotImplemented,
"RefinedSimplexLocalBasis not implemented for dim > 3.");
55 else if (global[0] <= 1.0)
58 DUNE_THROW(InvalidStateException,
"no subelement defined");
69 FieldVector<D,1>& local)
71 if (global[0] <= 0.5) {
73 local[0] = 2.0 * global[0];
78 local[0] = 2.0 * global[0] - 1.0;
119 if (global[0] + global[1] <= 0.5)
121 else if (global[0] >= 0.5)
123 else if (global[1] >= 0.5)
137 FieldVector<D,2>& local)
139 if (global[0] + global[1] <= 0.5) {
141 local[0] = 2*global[0];
142 local[1] = 2*global[1];
144 }
else if (global[0] >= 0.5) {
146 local[0] = 2*global[0]-1;
147 local[1] = 2*global[1];
149 }
else if (global[1] >= 0.5) {
151 local[0] = 2*global[0];
152 local[1] = 2*global[1]-1;
157 local[0] = -2 * global[0] + 1;
158 local[1] = -2 * global[1] + 1;
215 if (global[0] + global[1] + global[2] <= 0.5)
217 else if (global[0] >= 0.5)
219 else if (global[1] >= 0.5)
221 else if (global[2] >= 0.5)
223 else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] <= 0.5))
225 else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] <= 0.5))
227 else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] >= 0.5))
229 else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] >= 0.5))
232 DUNE_THROW(InvalidStateException,
"no subelement defined");
243 FieldVector<D,3>& local)
245 if (global[0] + global[1] + global[2] <= 0.5) {
250 }
else if (global[0] >= 0.5) {
256 }
else if (global[1] >= 0.5) {
262 }
else if (global[2] >= 0.5) {
268 }
else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] <= 0.5)) {
270 local[0] = 2.0 * global[1];
271 local[1] = 2.0 * (0.5 - global[0] - global[1]);
272 local[2] = 2.0 * (-0.5 + global[0] + global[1] + global[2]);
284 }
else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] <= 0.5)) {
286 local[0] = 2.0 * (0.5 - global[0]);
287 local[1] = 2.0 * (0.5 - global[1] - global[2]);
288 local[2] = 2.0 * global[2];
298 }
else if ((global[0] + global[1] <= 0.5)and (global[1] + global[2] >= 0.5)) {
300 local[0] = 2.0 * (0.5 - global[0] - global[1]);
301 local[1] = 2.0 * global[0];
302 local[2] = 2.0 * (-0.5 + global[1] + global[2]);
313 }
else if ((global[0] + global[1] >= 0.5)and (global[1] + global[2] >= 0.5)) {
315 local[0] = 2.0 * (-0.5 + global[1] + global[2]);
316 local[1] = 2.0 * (0.5 - global[1]);
317 local[2] = 2.0 * (-0.5 + global[0] + global[1]);
331 DUNE_THROW(InvalidStateException,
"no subelement defined");
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:100
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 void getSubElement(const FieldVector< D, 3 > &global, int &subElement, FieldVector< D, 3 > &local)
Get local coordinates in the subsimplex.
Definition: refinedsimplexlocalbasis.hh:241
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:213
RefinedSimplexLocalBasis()
Definition: refinedsimplexlocalbasis.hh:20
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:39
static int getSubElement(const FieldVector< D, 1 > &global)
Get the number of the subelement containing a given point.
Definition: refinedsimplexlocalbasis.hh:51
RefinedSimplexLocalBasis()
Protected default constructor so this class can only be instantiated as a base class.
Definition: refinedsimplexlocalbasis.hh:181
static int getSubElement(const FieldVector< D, 2 > &global)
Get the number of the subtriangle containing a given point.
Definition: refinedsimplexlocalbasis.hh:117
Definition: refinedsimplexlocalbasis.hh:17
static void getSubElement(const FieldVector< D, 2 > &global, int &subElement, FieldVector< D, 2 > &local)
Get local coordinates in the subtriangle.
Definition: refinedsimplexlocalbasis.hh:135