dune-localfunctions
2.2.0
|
00001 #ifndef DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH 00002 #define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH 00003 00004 #include <vector> 00005 00006 #include <dune/common/fmatrix.hh> 00007 00008 #include "../../common/localbasis.hh" 00009 00010 namespace Dune 00011 { 00021 template<class D, class R> 00022 class RT12DLocalBasis 00023 { 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; 00029 RT12DLocalBasis () 00030 { 00031 sign0 = sign1 = sign2 = 1.0; 00032 } 00033 00039 RT12DLocalBasis (unsigned int s) 00040 { 00041 sign0 = sign1 = sign2 = 1.0; 00042 if (s & 1) 00043 { 00044 sign0 = -1.0; 00045 } 00046 if (s & 2) 00047 { 00048 sign1 = -1.0; 00049 } 00050 if (s & 4) 00051 { 00052 sign2 = -1.0; 00053 } 00054 } 00055 00057 unsigned int size () const 00058 { 00059 return 8; 00060 } 00061 00068 inline void evaluateFunction (const typename Traits::DomainType& in, 00069 std::vector<typename Traits::RangeType>& out) const 00070 { 00071 out.resize(8); 00072 out[0][0] = sign0*(in[0] - 4.0*in[0]*in[1]); 00073 out[0][1] = sign0*(-1.0 + 5.0*in[1] - 4.0*in[1]*in[1]); 00074 out[1][0] = sign1*(-1.0 + 5.0*in[0] - 4.0*in[0]*in[0]); 00075 out[1][1] = sign1*(in[1] - 4.0*in[0]*in[1]); 00076 out[2][0] = sign2*(-3.0*in[0] + 4.0*in[0]*in[0] + 4.0*in[1]*in[0]); 00077 out[2][1] = sign2*(-3.0*in[1] + 4.0*in[0]*in[1] + 4.0*in[1]*in[1]); 00078 out[3][0] = -5.0*in[0] + 8.0*in[0]*in[0] + 4.0*in[1]*in[0]; 00079 out[3][1] = 3.0 - 6.0*in[0] - 7.0*in[1] + 8.0*in[0]*in[1] + 4.0*in[1]*in[1]; 00080 out[4][0] = -3.0 + 7.0*in[0] + 6.0*in[1] - 4.0*in[0]*in[0] - 8.0*in[1]*in[0]; 00081 out[4][1] = 5.0*in[1] - 4.0*in[0]*in[1] - 8.0*in[1]*in[1]; 00082 out[5][0] = in[0] - 4.0*in[0]*in[0] + 4.0*in[1]*in[0]; 00083 out[5][1] = -1.0*in[1] - 4.0*in[0]*in[1] + 4.0*in[1]*in[1]; 00084 out[6][0] = 16.0*in[0] - 16.0*in[0]*in[0] - 8.0*in[1]*in[0]; 00085 out[6][1] = 8.0*in[1] - 16.0*in[0]*in[1] - 8.0*in[1]*in[1]; 00086 out[7][0] = 8.0*in[0] - 8.0*in[0]*in[0] - 16.0*in[1]*in[0]; 00087 out[7][1] = 16.0*in[1] - 8.0*in[0]*in[1] - 16.0*in[1]*in[1]; 00088 } 00089 00096 inline void evaluateJacobian (const typename Traits::DomainType& in, 00097 std::vector<typename Traits::JacobianType>& out) const 00098 { 00099 out.resize(8); 00100 00101 out[0][0][0] = sign0*(1.0 - 4.0*in[1]); 00102 out[0][0][1] = sign0*(-4.0*in[0]); 00103 out[0][1][0] = 0.0; 00104 out[0][1][1] = sign0*(5.0 - 8.0*in[1]); 00105 00106 out[1][0][0] = sign1*(5.0 - 8.0*in[0]); 00107 out[1][0][1] = 0.0; 00108 out[1][1][0] = sign1*(-4.0*in[1]); 00109 out[1][1][1] = sign1*(1.0 - 4.0*in[0]); 00110 00111 out[2][0][0] = sign2*(-3.0 + 8.0*in[0] + 4.0*in[1]); 00112 out[2][0][1] = sign2*(4.0*in[0]); 00113 out[2][1][0] = sign2*(4.0*in[1]); 00114 out[2][1][1] = sign2*(-3.0 + 4.0*in[0] + 8.0*in[1]); 00115 00116 out[3][0][0] = -5.0 + 16.0*in[0] + 4.0*in[1]; 00117 out[3][0][1] = 4.0*in[0]; 00118 out[3][1][0] = -6.0 + 8.0*in[1]; 00119 out[3][1][1] = -7.0 + 8.0*in[0] + 8.0*in[1]; 00120 00121 out[4][0][0] = 7.0 - 8.0*in[0] - 8.0*in[1]; 00122 out[4][0][1] = 6.0 - 8.0*in[0]; 00123 out[4][1][0] = -4.0*in[1]; 00124 out[4][1][1] = 5.0 - 4.0*in[0] - 16.0*in[1]; 00125 00126 out[5][0][0] = 1.0 - 8.0*in[0] + 4*in[1]; 00127 out[5][0][1] = 4.0*in[0]; 00128 out[5][1][0] = -4.0*in[1]; 00129 out[5][1][1] = -1.0 - 4.0*in[0] + 8.0*in[1]; 00130 00131 out[6][0][0] = 16.0 - 32.0*in[0] - 8.0*in[1]; 00132 out[6][0][1] = -8.0*in[0]; 00133 out[6][1][0] = -16.0*in[1]; 00134 out[6][1][1] = 8.0 - 16.0*in[0] - 16.0*in[1]; 00135 00136 out[7][0][0] = 8.0 - 16.0*in[0] - 16.0*in[1]; 00137 out[7][0][1] = -16.0*in[0]; 00138 out[7][1][0] = -8.0*in[1]; 00139 out[7][1][1] = 16.0 - 8.0*in[0] - 32.0*in[1]; 00140 } 00141 00143 unsigned int order () const 00144 { 00145 return 2; 00146 } 00147 00148 private: 00149 R sign0, sign1, sign2; 00150 }; 00151 } 00152 #endif // DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH