fem/interpolationpointdata.hh

Go to the documentation of this file.

// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
//
// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
// SPDX-License-Identifier: GPL-3.0-or-later
//
#ifndef DUMUX_DISCRETIZATION_FE_INTEGRATION_POINT_DATA_HH
#define DUMUX_DISCRETIZATION_FE_INTEGRATION_POINT_DATA_HH
 
#include <vector>
 
namespace Dumux {
 
template<class GlobalPosition, class LocalBasis>
class FEInterpolationPointData
{
    using LocalPosition = typename LocalBasis::Traits::DomainType;
    using RangeType = typename LocalBasis::Traits::RangeType;
    using JacobianType = typename LocalBasis::Traits::JacobianType;
 
    using ShapeValues = std::vector<RangeType>;
    using ShapeJacobians = std::vector<JacobianType>;
    using Gradients = std::vector<GlobalPosition>;
 
public:
    // The default constructor
    FEInterpolationPointData() = delete;
 
    // The constructor
    template<class Geometry>
    FEInterpolationPointData(const Geometry& geometry,
                             const LocalPosition& local,
                             const GlobalPosition& global,
                             const LocalBasis& localBasis)
    : local_(local),
      global_(global)
    {
        auto numLocalDofs = localBasis.size();
 
        // set the shape values
        shapeValues_.resize(numLocalDofs);
        localBasis.evaluateFunction(local, shapeValues_);
 
        // the local shape function gradients
        thread_local ShapeJacobians shapeJacobians;
        localBasis.evaluateJacobian(local, shapeJacobians);
 
        // the global shape function gradients
        const auto jacInvT = geometry.jacobianInverseTransposed(local);
        gradients_.resize(numLocalDofs, GlobalPosition(0.0));
        for (unsigned int i = 0; i < numLocalDofs; ++i)
            jacInvT.umv(shapeJacobians[i][0], gradients_[i]);
    }
 
    // The constructor
    template<class Geometry>
    FEInterpolationPointData(const Geometry& geometry,
                             const LocalPosition& local,
                             const LocalBasis& localBasis)
    : FEInterpolationPointData(geometry, local, geometry.global(local), localBasis)
    { }

    const ShapeValues& shapeValues() const
    { return shapeValues_; }

    const RangeType& shapeValue(int i) const
    { return shapeValues_[i]; }

    const GlobalPosition& gradN(int i) const
    { return gradients_[i]; }

    const LocalPosition& local() const
    { return local_; }

    const GlobalPosition& global() const
    { return global_; }
 
private:
    LocalPosition local_;
    GlobalPosition global_;
 
    ShapeValues shapeValues_;
    Gradients gradients_;
};

template<class BaseClass, class BoundaryFlag, class LocalIndex>
class FEFaceInterpolationPointData : public BaseClass
{
public:
    using GlobalPosition = std::remove_cvref_t<decltype(std::declval<BaseClass>().global())>;
    using LocalPosition = std::remove_cvref_t<decltype(std::declval<BaseClass>().local())>;
 
    template<class... Args>
    FEFaceInterpolationPointData(GlobalPosition&& n, BoundaryFlag&& bFlag, LocalIndex index, Args&&... args)
    : BaseClass(std::forward<Args>(args)...), normal_(std::move(n)), boundaryFlag_(std::move(bFlag)), index_(index) {}
 
    template<class... Args>
    FEFaceInterpolationPointData(const GlobalPosition& n, const BoundaryFlag& bFlag, LocalIndex index, Args&&... args)
    : BaseClass(std::forward<Args>(args)...), normal_(n), boundaryFlag_(bFlag), index_(index) {}

    const GlobalPosition& unitOuterNormal() const
    { return normal_; }

    typename BoundaryFlag::value_type boundaryFlag() const
    { return boundaryFlag_.get(); }

    LocalIndex intersectionIndex() const
    { return index_; }
 
private:
    GlobalPosition normal_;
    BoundaryFlag boundaryFlag_;
    LocalIndex index_;
};

template<class BaseClass, class LocalIndex>
class FEBoundaryFaceInterpolationPointData : public BaseClass
{
public:
    using GlobalPosition = std::remove_cvref_t<decltype(std::declval<BaseClass>().global())>;
    using LocalPosition = std::remove_cvref_t<decltype(std::declval<BaseClass>().local())>;
 
    template<class... Args>
    FEBoundaryFaceInterpolationPointData(GlobalPosition&& n, LocalIndex index, Args&&... args)
    : BaseClass(std::forward<Args>(args)...), normal_(std::move(n)), boundaryFaceIndex_(index) {}
 
    template<class... Args>
    FEBoundaryFaceInterpolationPointData(const GlobalPosition& n, LocalIndex index, Args&&... args)
    : BaseClass(std::forward<Args>(args)...), normal_(n), boundaryFaceIndex_(index) {}

    const GlobalPosition& unitOuterNormal() const
    { return normal_; }

    LocalIndex boundaryFaceIndex() const
    { return boundaryFaceIndex_; }
 
private:
    GlobalPosition normal_;
    LocalIndex boundaryFaceIndex_;
};
 
} // end namespace Dumux
 
#endif