dune-localfunctions
2.2.0
|
00001 #ifndef DUNE_RT0Q2DALL_HH 00002 #define DUNE_RT0Q2DALL_HH 00003 00004 #include <cstddef> 00005 #include <vector> 00006 00007 #include <dune/common/fmatrix.hh> 00008 00009 #include <dune/localfunctions/common/localbasis.hh> 00010 #include <dune/localfunctions/common/localkey.hh> 00011 00012 namespace Dune 00013 { 00022 template<class D, class R> 00023 class RT0Q2DLocalBasis 00024 { 00025 public: 00026 typedef LocalBasisTraits<D,2,Dune::FieldVector<D,2>,R,2,Dune::FieldVector<R,2>, 00027 Dune::FieldMatrix<R,2,2> > Traits; 00028 00030 RT0Q2DLocalBasis () 00031 { 00032 sign0 = sign1 = sign2 = sign3 = 1.0; 00033 } 00034 00036 RT0Q2DLocalBasis (unsigned int s) 00037 { 00038 sign0 = sign1 = sign2 = sign3 = 1.0; 00039 if (s&1) sign0 = -1.0; 00040 if (s&2) sign1 = -1.0; 00041 if (s&4) sign2 = -1.0; 00042 if (s&8) sign3 = -1.0; 00043 } 00044 00046 unsigned int size () const 00047 { 00048 return 4; 00049 } 00050 00052 inline void evaluateFunction (const typename Traits::DomainType& in, 00053 std::vector<typename Traits::RangeType>& out) const 00054 { 00055 out.resize(4); 00056 out[0][0] = sign0*(in[0]-1.0); out[0][1]=0.0; 00057 out[1][0] = sign1*(in[0]); out[1][1]=0.0; 00058 out[2][0] = 0.0; out[2][1]=sign2*(in[1]-1.0); 00059 out[3][0] = 0.0; out[3][1]=sign3*(in[1]); 00060 } 00061 00063 inline void 00064 evaluateJacobian (const typename Traits::DomainType& in, // position 00065 std::vector<typename Traits::JacobianType>& out) const // return value 00066 { 00067 out.resize(4); 00068 out[0][0][0] = sign0; out[0][0][1] = 0; 00069 out[0][1][0] = 0; out[0][1][1] = 0; 00070 00071 out[1][0][0] = sign1; out[1][0][1] = 0; 00072 out[1][1][0] = 0; out[1][1][1] = 0; 00073 00074 out[2][0][0] = 0; out[2][0][1] = 0; 00075 out[2][1][0] = 0; out[2][1][1] = sign2; 00076 00077 out[3][0][0] = 0; out[3][0][1] = 0; 00078 out[3][1][0] = 0; out[3][1][1] = sign3; 00079 } 00080 00082 unsigned int order () const 00083 { 00084 return 1; 00085 } 00086 00087 private: 00088 R sign0, sign1, sign2, sign3; 00089 }; 00090 00091 00099 template<class LB> 00100 class RT0Q2DLocalInterpolation 00101 { 00102 public: 00103 00105 RT0Q2DLocalInterpolation () 00106 { 00107 sign0 = sign1 = sign2 = sign3 = 1.0; 00108 } 00109 00111 RT0Q2DLocalInterpolation (unsigned int s) 00112 { 00113 sign0 = sign1 = sign2 = sign3 = 1.0; 00114 if (s&1) sign0 *= -1.0; 00115 if (s&2) sign1 *= -1.0; 00116 if (s&4) sign2 *= -1.0; 00117 if (s&8) sign3 *= -1.0; 00118 00119 m0[0] = 0.0; m0[1] = 0.5; 00120 m1[0] = 1.0; m1[1] = 0.5; 00121 m2[0] = 0.5; m2[1] = 0.0; 00122 m3[0] = 0.5; m3[1] = 1.0; 00123 00124 n0[0] = -1.0; n0[1] = 0.0; 00125 n1[0] = 1.0; n1[1] = 0.0; 00126 n2[0] = 0.0; n2[1] = -1.0; 00127 n3[0] = 0.0; n3[1] = 1.0; 00128 } 00129 00130 template<typename F, typename C> 00131 void interpolate (const F& f, std::vector<C>& out) const 00132 { 00133 // f gives v*outer normal at a point on the edge! 00134 typename F::Traits::RangeType y; 00135 00136 out.resize(4); 00137 00138 f.evaluate(m0,y); out[0] = (y[0]*n0[0]+y[1]*n0[1])*sign0; 00139 f.evaluate(m1,y); out[1] = (y[0]*n1[0]+y[1]*n1[1])*sign1; 00140 f.evaluate(m2,y); out[2] = (y[0]*n2[0]+y[1]*n2[1])*sign2; 00141 f.evaluate(m3,y); out[3] = (y[0]*n3[0]+y[1]*n3[1])*sign3; 00142 } 00143 00144 private: 00145 typename LB::Traits::RangeFieldType sign0,sign1,sign2,sign3; 00146 typename LB::Traits::DomainType m0,m1,m2,m3; 00147 typename LB::Traits::DomainType n0,n1,n2,n3; 00148 }; 00149 00156 class RT0Q2DLocalCoefficients 00157 { 00158 public: 00160 RT0Q2DLocalCoefficients () : li(4) 00161 { 00162 for (std::size_t i=0; i<4; i++) 00163 li[i] = LocalKey(i,1,0); 00164 } 00165 00167 std::size_t size () const 00168 { 00169 return 4; 00170 } 00171 00173 const LocalKey& localKey (std::size_t i) const 00174 { 00175 return li[i]; 00176 } 00177 00178 private: 00179 std::vector<LocalKey> li; 00180 }; 00181 00182 } 00183 #endif