Forked from
agnumpde / dune-fufem
1862 commits behind the upstream repository.
-
Oliver Sander authored
'ag-common' to dune/fufem. [[Imported from SVN: r3520]]
Oliver Sander authored'ag-common' to dune/fufem. [[Imported from SVN: r3520]]
Code owners
Assign users and groups as approvers for specific file changes. Learn more.
gradientassembler.hh 3.39 KiB
#ifndef GRADIENT_ASSEMBLER_HH
#define GRADIENT_ASSEMBLER_HH
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/grid/common/quadraturerules.hh>
#include <dune/grid/common/genericreferenceelements.hh>
#include <dune/istl/matrix.hh>
#include <dune/ag-common/functions/localbasisderivativefunction.hh>
#include <dune/ag-common/assemblers/localassembler.hh>
//! \todo Please doc me !
template <class GridType, class TrialLocalFE, class AnsatzLocalFE>
class GradientAssembler
: public LinearLocalAssembler < GridType, TrialLocalFE, AnsatzLocalFE, Dune::FieldMatrix<double,GridType::dimension,1> >
{
static const int dim = GridType::dimension;
public:
typedef typename Dune::FieldMatrix<double,dim,1> T;
typedef typename GridType::template Codim<0>::Entity::EntityPointer ElementPointer;
typedef typename LocalAssembler < GridType, TrialLocalFE, AnsatzLocalFE, T >::BoolMatrix BoolMatrix;
typedef typename LocalAssembler < GridType, TrialLocalFE, AnsatzLocalFE, T >::LocalMatrix LocalMatrix;
void indices(const ElementPointer& element, BoolMatrix& isNonZero, const TrialLocalFE& tFE, const AnsatzLocalFE& aFE) const
{
isNonZero = true;
return;
}
void assemble(const ElementPointer& element, LocalMatrix& localMatrix
, const TrialLocalFE& tFE, const AnsatzLocalFE& aFE) const
{
typedef typename Dune::template FieldVector<double,dim> FVdim;
typedef typename Dune::template FieldMatrix<double,dim,dim> FMdimdim;
typedef typename AnsatzLocalFE::Traits::LocalBasisType::Traits::JacobianType JacobianType;
typedef typename TrialLocalFE::Traits::LocalBasisType::Traits::DomainType DomainType;
typedef typename TrialLocalFE::Traits::LocalBasisType::Traits::RangeType RangeType;
localMatrix = 0.0;
// compute the element barycenter in local coordinates
const FVdim pos = Dune::GenericReferenceElements<double,dim>::general(element->type()).position(0,0);
// get transposed inverse of Jacobian of transformation
const FMdimdim& invJacobian = element->geometry().jacobianInverseTransposed(pos);
typedef typename Dune::LocalFiniteElementFunctionBase<TrialLocalFE>::type FunctionBaseClass;
LocalBasisDerivativeFunction<AnsatzLocalFE, FunctionBaseClass> derivative(aFE.localBasis());
std::vector<JacobianType> referenceGradients(tFE.localBasis().size());
std::vector<typename RangeType::field_type> partialDerivatives(tFE.localBasis().size());
for (int i = 0; i < aFE.localBasis().size(); ++i)
{
// interpolate all partial derivatives of ansatz function by test functions
derivative.setIndex(i);
for (int j = 0; j < dim; ++j)
{
derivative.setComponent(j);
tFE.localInterpolation().interpolate(derivative, partialDerivatives);
for (int k = 0; k < tFE.localBasis().size(); ++k)
referenceGradients[k][0][j] = partialDerivatives[k];
}
// transform gradients
for (int k = 0; k < tFE.localBasis().size(); ++k)
invJacobian.mv(referenceGradients[k][0], localMatrix[k][i]);
}
}
};
#endif