3 #ifndef DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE3D_ALL_HH 4 #define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE3D_ALL_HH 9 #include <dune/common/fmatrix.hh> 24 template<
class D,
class R>
34 sign0 = sign1 = sign2 = sign3 = sign4 = sign5 = 1.0;
40 sign0 = sign1 = sign2 = sign3 = sign4 = sign5 = 1.0;
41 if (s&1) sign0 = -1.0;
42 if (s&2) sign1 = -1.0;
43 if (s&4) sign2 = -1.0;
44 if (s&8) sign3 = -1.0;
45 if (s&16) sign4 = -1.0;
46 if (s&32) sign5 = -1.0;
57 std::vector<typename Traits::RangeType>& out)
const 60 out[0][0] = sign0*(in[0]-1.0); out[0][1]=0.0; out[0][2]=0.0;
61 out[1][0] = sign1*(in[0]); out[1][1]=0.0; out[1][2]=0.0;
62 out[2][0] = 0.0; out[2][1]=sign2*(in[1]-1.0); out[2][2]=0.0;
63 out[3][0] = 0.0; out[3][1]=sign3*(in[1]); out[3][2]=0.0;
64 out[4][0] = 0.0; out[4][1]=0.0; out[4][2]=sign4*(in[2]-1.0);
65 out[5][0] = 0.0; out[5][1]=0.0; out[5][2]=sign5*(in[2]);
71 std::vector<typename Traits::JacobianType>& out)
const 74 out[0][0][0] = sign0; out[0][0][1] = 0; out[0][0][2] = 0;
75 out[0][1][0] = 0; out[0][1][1] = 0; out[0][1][2] = 0;
76 out[0][2][0] = 0; out[0][2][1] = 0; out[0][2][2] = 0;
78 out[1][0][0] = sign1; out[1][0][1] = 0; out[1][0][2] = 0;
79 out[1][1][0] = 0; out[1][1][1] = 0; out[1][1][2] = 0;
80 out[1][2][0] = 0; out[1][2][1] = 0; out[1][2][2] = 0;
82 out[2][0][0] = 0; out[2][0][1] = 0; out[2][0][2] = 0;
83 out[2][1][0] = 0; out[2][1][1] = sign2; out[2][1][2] = 0;
84 out[2][2][0] = 0; out[2][2][1] = 0; out[2][2][2] = 0;
86 out[3][0][0] = 0; out[3][0][1] = 0; out[3][0][2] = 0;
87 out[3][1][0] = 0; out[3][1][1] = sign3; out[3][1][2] = 0;
88 out[3][2][0] = 0; out[3][2][1] = 0; out[3][2][2] = 0;
90 out[4][0][0] = 0; out[4][0][1] = 0; out[4][0][2] = 0;
91 out[4][1][0] = 0; out[4][1][1] = 0; out[4][1][2] = 0;
92 out[4][2][0] = 0; out[4][2][1] = 0; out[4][2][2] = sign4;
94 out[5][0][0] = 0; out[5][0][1] = 0; out[5][0][2] = 0;
95 out[5][1][0] = 0; out[5][1][1] = 0; out[5][1][2] = 0;
96 out[5][2][0] = 0; out[5][2][1] = 0; out[5][2][2] = sign5;
106 R sign0, sign1, sign2, sign3, sign4, sign5;
125 sign0 = sign1 = sign2 = sign3 = sign4 = sign5 = 1.0;
131 sign0 = sign1 = sign2 = sign3 = sign4 = sign5 = 1.0;
132 if (s&1) sign0 *= -1.0;
133 if (s&2) sign1 *= -1.0;
134 if (s&4) sign2 *= -1.0;
135 if (s&8) sign3 *= -1.0;
136 if (s&16) sign4 *= -1.0;
137 if (s&32) sign5 *= -1.0;
139 m0[0] = 0.0; m0[1] = 0.5; m0[2] = 0.5;
140 m1[0] = 1.0; m1[1] = 0.5; m1[2] = 0.5;
141 m2[0] = 0.5; m2[1] = 0.0; m2[2] = 0.5;
142 m3[0] = 0.5; m3[1] = 1.0; m3[2] = 0.5;
143 m4[0] = 0.5; m4[1] = 0.5; m4[2] = 0.0;
144 m5[0] = 0.5; m5[1] = 0.5; m5[2] = 1.0;
146 n0[0] = -1.0; n0[1] = 0.0; n0[2] = 0.0;
147 n1[0] = 1.0; n1[1] = 0.0; n1[2] = 0.0;
148 n2[0] = 0.0; n2[1] = -1.0; n2[2] = 0.0;
149 n3[0] = 0.0; n3[1] = 1.0; n3[2] = 0.0;
150 n4[0] = 0.0; n4[1] = 0.0; n4[2] =-1.0;
151 n5[0] = 0.0; n5[1] = 0.0; n5[2] = 1.0;
154 template<
typename F,
typename C>
158 typename F::Traits::RangeType y;
162 f.evaluate(m0,y); out[0] = (y[0]*n0[0]+y[1]*n0[1]+y[2]*n0[2])*sign0;
163 f.evaluate(m1,y); out[1] = (y[0]*n1[0]+y[1]*n1[1]+y[2]*n1[2])*sign1;
164 f.evaluate(m2,y); out[2] = (y[0]*n2[0]+y[1]*n2[1]+y[2]*n2[2])*sign2;
165 f.evaluate(m3,y); out[3] = (y[0]*n3[0]+y[1]*n3[1]+y[2]*n3[2])*sign3;
166 f.evaluate(m4,y); out[4] = (y[0]*n4[0]+y[1]*n4[1]+y[2]*n4[2])*sign4;
167 f.evaluate(m5,y); out[5] = (y[0]*n5[0]+y[1]*n5[1]+y[2]*n5[2])*sign5;
171 typename LB::Traits::RangeFieldType sign0,sign1,sign2,sign3,sign4,sign5;
172 typename LB::Traits::DomainType m0,m1,m2,m3,m4,m5;
173 typename LB::Traits::DomainType n0,n1,n2,n3,n4,n5;
188 for (std::size_t i=0; i<6; i++)
205 std::vector<LocalKey> li;
209 #endif // DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_CUBE3D_ALL_HH std::size_t size() const
number of coefficients
Definition: raviartthomas0cube3dall.hh:193
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: raviartthomas0cube3dall.hh:56
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: raviartthomas0cube3dall.hh:70
RT0Cube3DLocalInterpolation()
Standard constructor.
Definition: raviartthomas0cube3dall.hh:123
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: raviartthomas0cube3dall.hh:199
RT0Cube3DLocalBasis()
Standard constructor.
Definition: raviartthomas0cube3dall.hh:32
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 3, Dune::FieldVector< R, 3 >, Dune::FieldMatrix< R, 3, 3 > > Traits
Definition: raviartthomas0cube3dall.hh:29
RT0Cube3DLocalCoefficients()
Standard constructor.
Definition: raviartthomas0cube3dall.hh:186
Describe position of one degree of freedom.
Definition: localkey.hh:21
void interpolate(const F &f, std::vector< C > &out) const
Definition: raviartthomas0cube3dall.hh:155
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
unsigned int order() const
Polynomial order of the shape functions.
Definition: raviartthomas0cube3dall.hh:100
RT0Cube3DLocalBasis(unsigned int s)
Make set numer s, where 0<=s<64.
Definition: raviartthomas0cube3dall.hh:38
Layout map for RT0 elements on quadrilaterals.
Definition: raviartthomas0cube3dall.hh:182
unsigned int size() const
number of shape functions
Definition: raviartthomas0cube3dall.hh:50
RT0Cube3DLocalInterpolation(unsigned int s)
Make set numer s, where 0<=s<64.
Definition: raviartthomas0cube3dall.hh:129
D DomainType
domain type
Definition: localbasis.hh:49
Lowest order Raviart-Thomas shape functions on the reference hexahedron.
Definition: raviartthomas0cube3dall.hh:118
Lowest order Raviart-Thomas shape functions on the reference hexahedron.
Definition: raviartthomas0cube3dall.hh:25