dune-localfunctions
2.2.0
|
00001 #ifndef DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_HH 00002 #define DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_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 BDM1Q2DLocalBasis 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; 00028 00030 BDM1Q2DLocalBasis () 00031 { 00032 sign0 = sign1 = sign2 = sign3 = 1.0; 00033 } 00034 00040 BDM1Q2DLocalBasis (unsigned int s) 00041 { 00042 sign0 = sign1 = sign2 = sign3 = 1.0; 00043 if (s & 1) 00044 { 00045 sign0 = -1.0; 00046 } 00047 if (s & 2) 00048 { 00049 sign1 = -1.0; 00050 } 00051 if (s & 4) 00052 { 00053 sign2 = -1.0; 00054 } 00055 if (s & 8) 00056 { 00057 sign3 = -1.0; 00058 } 00059 } 00060 00062 unsigned int size () const 00063 { 00064 return 8; 00065 } 00066 00073 inline void evaluateFunction (const typename Traits::DomainType& in, 00074 std::vector<typename Traits::RangeType>& out) const 00075 { 00076 out.resize(8); 00077 00078 out[0][0] = sign0*(in[0] - 1.0); 00079 out[0][1] = 0.0; 00080 out[1][0] = 6.0*in[0]*in[1] - 3.0*in[0]-6*in[1] + 3.0; 00081 out[1][1] = -3.0*in[1]*in[1] + 3.0*in[1]; 00082 out[2][0] = sign1*(in[0]); 00083 out[2][1] = 0.0; 00084 out[3][0] = -6.0*in[0]*in[1] + 3.0*in[0]; 00085 out[3][1] = 3.0*in[1]*in[1] - 3.0*in[1]; 00086 out[4][0] = 0.0; 00087 out[4][1] = sign2*(in[1] - 1.0); 00088 out[5][0] = 3.0*in[0]*in[0] - 3.0*in[0]; 00089 out[5][1] = -6.0*in[0]*in[1] + 6.0*in[0] + 3.0*in[1] - 3.0; 00090 out[6][0] = 0.0; 00091 out[6][1] = sign3*(in[1]); 00092 out[7][0] = -3.0*in[0]*in[0] + 3.0*in[0]; 00093 out[7][1] = 6.0*in[0]*in[1] - 3.0*in[1]; 00094 } 00095 00102 inline void evaluateJacobian (const typename Traits::DomainType& in, 00103 std::vector<typename Traits::JacobianType>& out) const 00104 { 00105 out.resize(8); 00106 00107 out[0][0][0] = sign0; 00108 out[0][0][1] = 0.0; 00109 out[0][1][0] = 0.0; 00110 out[0][1][1] = 0.0; 00111 00112 out[1][0][0] = 6.0*in[1] - 3.0; 00113 out[1][0][1] = 6.0*in[0] - 6.0; 00114 out[1][1][0] = 0.0; 00115 out[1][1][1] = -6.0*in[1] + 3.0; 00116 00117 out[2][0][0] = sign1; 00118 out[2][0][1] = 0.0; 00119 out[2][1][0] = 0.0; 00120 out[2][1][1] = 0.0; 00121 00122 out[3][0][0] = -6.0*in[1] + 3.0; 00123 out[3][0][1] = -6.0*in[0]; 00124 out[3][1][0] = 0.0; 00125 out[3][1][1] = 6.0*in[1] - 3.0; 00126 00127 out[4][0][0] = 0.0; 00128 out[4][0][1] = 0.0; 00129 out[4][1][0] = 0.0; 00130 out[4][1][1] = sign2; 00131 00132 out[5][0][0] = 6.0*in[0] - 3.0; 00133 out[5][0][1] = 0.0; 00134 out[5][1][0] = -6.0*in[1] + 6.0; 00135 out[5][1][1] = -6.0*in[0] + 3.0; 00136 00137 out[6][0][0] = 0.0; 00138 out[6][0][1] = 0.0; 00139 out[6][1][0] = 0.0; 00140 out[6][1][1] = sign3; 00141 00142 out[7][0][0] = -6.0*in[0] + 3.0; 00143 out[7][0][1] = 0.0; 00144 out[7][1][0] = 6.0*in[1]; 00145 out[7][1][1] = 6.0*in[0] - 3.0; 00146 } 00147 00149 unsigned int order () const 00150 { 00151 return 2; 00152 } 00153 00154 private: 00155 R sign0, sign1, sign2, sign3; 00156 }; 00157 } 00158 #endif // DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_HH