- Home
- About DUNE
- Download
- Documentation
- Community
- Development
00001 // -*- tab-width: 4; indent-tabs-mode: nil -*- 00002 #ifndef DUNE_Q1_LOCALBASIS_HH 00003 #define DUNE_Q1_LOCALBASIS_HH 00004 00005 #include <dune/common/fmatrix.hh> 00006 00007 #include <dune/localfunctions/common/localbasis.hh> 00008 00009 namespace Dune 00010 { 00022 template<class D, class R, int dim> 00023 class Q1LocalBasis 00024 { 00025 public: 00026 typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>, 00027 Dune::FieldMatrix<R,1,dim> > Traits; 00028 00030 unsigned int size () const 00031 { 00032 return 1<<dim; 00033 } 00034 00036 inline void evaluateFunction (const typename Traits::DomainType& in, 00037 std::vector<typename Traits::RangeType>& out) const 00038 { 00039 out.resize(size()); 00040 00041 for (size_t i=0; i<size(); i++) { 00042 00043 out[i] = 1; 00044 00045 for (int j=0; j<dim; j++) 00046 // if j-th bit of i is set multiply with in[j], else with 1-in[j] 00047 out[i] *= (i & (1<<j)) ? in[j] : 1-in[j]; 00048 00049 } 00050 00051 } 00052 00054 inline void 00055 evaluateJacobian (const typename Traits::DomainType& in, // position 00056 std::vector<typename Traits::JacobianType>& out) const // return value 00057 { 00058 out.resize(size()); 00059 00060 // Loop over all shape functions 00061 for (size_t i=0; i<size(); i++) { 00062 00063 // Loop over all coordinate directions 00064 for (int j=0; j<dim; j++) { 00065 00066 // Initialize: the overall expression is a product 00067 // if j-th bit of i is set to -1, else 1 00068 out[i][0][j] = (i & (1<<j)) ? 1 : -1; 00069 00070 for (int k=0; k<dim; k++) { 00071 00072 if (j!=k) 00073 // if k-th bit of i is set multiply with in[j], else with 1-in[j] 00074 out[i][0][j] *= (i & (1<<k)) ? in[k] : 1-in[k]; 00075 00076 } 00077 00078 } 00079 00080 } 00081 00082 } 00083 00085 unsigned int order () const 00086 { 00087 return 1; 00088 } 00089 }; 00090 } 00091 #endif
Generated on Fri Apr 29 2011 with Doxygen (ver 1.7.1) [doxygen-log,error-log].