/usr/include/dune/localfunctions/dualmortarbasis/dualp1/dualp1localbasis.hh is in libdune-localfunctions-dev 2.4.1-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 | // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_DUAL_P1_LOCALBASIS_HH
#define DUNE_DUAL_P1_LOCALBASIS_HH
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/localfunctions/common/localbasis.hh>
namespace Dune
{
/**@ingroup LocalBasisImplementation
\brief Dual Lagrange shape functions on the simplex.
Defines the linear dual shape functions on the simplex.
Note that if the dual functions are chosen to be dual on the faces,
the integrated product of a Lagrange \f$\lambda_p\f$ and dual
function \f$\theta_q\f$ over faces not containing \f$q\f$ does in
general not vanish.
\tparam D Type to represent the field in the domain.
\tparam R Type to represent the field in the range.
\tparam dim The dimension of the simplex
\tparam faceDual If set, the basis functions are bi-orthogonal only on faces containing the corresponding vertex.
\nosubgrouping
*/
template<class D, class R, int dim, bool faceDualT=false>
class DualP1LocalBasis
{
public:
//! Determines if the basis is only biorthogonal on adjacent faces
static const bool faceDual = faceDualT;
//! \brief export type traits for function signature
typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
Dune::FieldMatrix<R,1,dim> > Traits;
//! \brief number of shape functions
unsigned int size () const
{
return dim+1;
}
//! \brief Evaluate all shape functions
inline void evaluateFunction (const typename Traits::DomainType& in,
std::vector<typename Traits::RangeType>& out) const
{
// evaluate P1 basis functions
std::vector<typename Traits::RangeType> p1Values(size());
p1Values[0] = 1.0;
for (int i=0; i<dim; i++) {
p1Values[0] -= in[i];
p1Values[i+1] = in[i];
}
// compute dual basis function values as a linear combination of the Lagrange values
out.resize(size());
for (int i=0; i<=dim; i++) {
out[i] = (dim+!faceDual)*p1Values[i];
for (int j=0; j<i; j++)
out[i] -= p1Values[j];
for (int j=i+1; j<=dim; j++)
out[i] -= p1Values[j];
}
}
//! \brief Evaluate Jacobian of all shape functions
inline void
evaluateJacobian (const typename Traits::DomainType& in,
std::vector<typename Traits::JacobianType>& out) const
{
// evaluate P1 jacobians
std::vector<typename Traits::JacobianType> p1Jacs(size());
for (int i=0; i<dim; i++)
p1Jacs[0][0][i] = -1;
for (int i=0; i<dim; i++)
for (int j=0; j<dim; j++)
p1Jacs[i+1][0][j] = (i==j);
// compute dual basis jacobians as linear combination of the Lagrange jacobians
out.resize(size());
for (size_t i=0; i<=dim; i++) {
out[i][0] = 0;
out[i][0].axpy(dim+!faceDual,p1Jacs[i][0]);
for (size_t j=0; j<i; j++)
out[i][0] -= p1Jacs[j][0];
for (int j=i+1; j<=dim; j++)
out[i][0] -= p1Jacs[j][0];
}
}
//! \brief Polynomial order of the shape functions
unsigned int order () const
{
return 1;
}
};
}
#endif
|