interactionvolumedatahandle.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:
/*****************************************************************************
 *   See the file COPYING for full copying permissions.                      *
 *                                                                           *
 *   This program is free software: you can redistribute it and/or modify    *
 *   it under the terms of the GNU General Public License as published by    *
 *   the Free Software Foundation, either version 3 of the License, or       *
 *   (at your option) any later version.                                     *
 *                                                                           *
 *   This program is distributed in the hope that it will be useful,         *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
 *   GNU General Public License for more details.                            *
 *                                                                           *
 *   You should have received a copy of the GNU General Public License       *
 *   along with this program.  If not, see <http://www.gnu.org/licenses/>.   *
 *****************************************************************************/
#ifndef DUMUX_DISCRETIZATION_CC_MPFA_INTERACTIONVOLUMEDATAHANDLE_HH
#define DUMUX_DISCRETIZATION_CC_MPFA_INTERACTIONVOLUMEDATAHANDLE_HH
 
#include <cassert>
#include <vector>
 
#include <dune/common/dynvector.hh>
 
#include <dumux/common/parameters.hh>
 
namespace Dumux {
namespace CCMpfaDataHandleBases {
 
template<class MVT, int size1, int size2>
class SystemMatricesHandle
{
    using AMatrix = typename MVT::AMatrix;
    using BMatrix = typename MVT::BMatrix;
    using CMatrix = typename MVT::CMatrix;
    using TMatrix = typename MVT::TMatrix;
    using CellVector = typename MVT::CellVector;
    using OutsideTij = std::vector< std::vector<CellVector> >;
    using OmegaStorage = typename MVT::OmegaStorage;
 
public:
    const CMatrix& CA() const { return CA_[contextIdx1_][contextIdx2_]; }
    CMatrix& CA() { return CA_[contextIdx1_][contextIdx2_]; }
 
    const AMatrix& A() const { return A_[contextIdx1_][contextIdx2_]; }
    AMatrix& A() { return A_[contextIdx1_][contextIdx2_]; }
 
    const BMatrix& AB() const { return AB_[contextIdx1_][contextIdx2_]; }
    BMatrix& AB() { return AB_[contextIdx1_][contextIdx2_]; }
 
    const TMatrix& T() const { return T_[contextIdx1_][contextIdx2_]; }
    TMatrix& T() { return T_[contextIdx1_][contextIdx2_]; }
 
    const OutsideTij& tijOutside() const { return tijOutside_[contextIdx1_][contextIdx2_]; }
    OutsideTij& tijOutside() { return tijOutside_[contextIdx1_][contextIdx2_]; }
 
    const OmegaStorage& omegas() const { return wijk_[contextIdx1_][contextIdx2_]; }
    OmegaStorage& omegas() { return wijk_[contextIdx1_][contextIdx2_]; }
 
    void setContextIndex1(unsigned int idx) const { assert(idx < size1); contextIdx1_ = idx; }
    void setContextIndex2(unsigned int idx) const { assert(idx < size2); contextIdx2_ = idx; }
 
protected:
    mutable unsigned int contextIdx1_{0};
    mutable unsigned int contextIdx2_{0};
 
    std::array< std::array<OmegaStorage, size2>, size1 > wijk_;     
 
    std::array< std::array<TMatrix, size2>, size1 > T_;             
    std::array< std::array<AMatrix, size2>, size1 > A_;             
    std::array< std::array<BMatrix, size2>, size1 > AB_;            
    std::array< std::array<CMatrix, size2>, size1 > CA_;            
    std::array< std::array<OutsideTij, size2>, size1 > tijOutside_; 
};
 
template<class MVT, int size1, int size2>
class SystemVectorsHandle
{
    using CellVector = typename MVT::CellVector;
 
public:
    const CellVector& uj() const { return u_[contextIdx1_][contextIdx2_]; }
    CellVector& uj() { return u_[contextIdx1_][contextIdx2_]; }
 
protected:
    void setContextIndex1(unsigned int idx) const { assert(idx < size1); contextIdx1_ = idx; }
    void setContextIndex2(unsigned int idx) const { assert(idx < size2); contextIdx2_ = idx; }
 
    mutable unsigned int contextIdx1_{0};
    mutable unsigned int contextIdx2_{0};
 
    std::array< std::array<CellVector, size2>, size1 > u_;
};
 
} // end namespace CCMpfaDataHandleBases
 
class EmptyDataHandle {};
 
template<class MatVecTraits, class PhysicsTraits, bool EnableAdvection>
class AdvectionDataHandle
: public CCMpfaDataHandleBases::SystemMatricesHandle<MatVecTraits, 1, 1>
, public CCMpfaDataHandleBases::SystemVectorsHandle<MatVecTraits, PhysicsTraits::numPhases, 1>
{
    // we only have one local system for all phases since we
    // solve them w.r.t. permeability tensor (unique for all phases)
    using Base1 = CCMpfaDataHandleBases::SystemMatricesHandle<MatVecTraits, 1, 1>;
 
    // we do have cell/Dirichlet values for all phases though!
    static constexpr int numPhases = PhysicsTraits::numPhases;
    using Base2 = CCMpfaDataHandleBases::SystemVectorsHandle<MatVecTraits, numPhases, 1>;
 
    using UnknownVector = typename MatVecTraits::AMatrix::row_type;
    using FaceVector = typename MatVecTraits::FaceVector;
    using FaceScalar = typename FaceVector::value_type;
    using OutsideGravityStorage = std::vector< std::vector<FaceScalar> >;
 
public:
    void setPhaseIndex(unsigned int phaseIdx) const { Base2::setContextIndex1(phaseIdx); }
 
    const FaceVector& g() const { return g_[Base2::contextIdx1_]; }
    FaceVector& g() { return g_[Base2::contextIdx1_]; }
 
    const UnknownVector& deltaG() const { return deltaG_[Base2::contextIdx1_]; }
    UnknownVector& deltaG() { return deltaG_[Base2::contextIdx1_]; }
 
    const OutsideGravityStorage& gOutside() const { return outsideG_[Base2::contextIdx1_]; }
    OutsideGravityStorage& gOutside() { return outsideG_[Base2::contextIdx1_]; }
 
private:
    std::array< FaceVector, numPhases > g_; 
    std::array< UnknownVector, numPhases > deltaG_; 
    std::array< OutsideGravityStorage, numPhases > outsideG_;  
};
 
template<class MatVecTraits, class PhysicsTraits, bool EnableDiffusion>
class DiffusionDataHandle
: public CCMpfaDataHandleBases::SystemMatricesHandle<MatVecTraits, PhysicsTraits::numPhases, PhysicsTraits::numComponents>
, public CCMpfaDataHandleBases::SystemVectorsHandle<MatVecTraits, PhysicsTraits::numPhases, PhysicsTraits::numComponents>
{
    static constexpr int numPhases = PhysicsTraits::numPhases;
    static constexpr int numComponents = PhysicsTraits::numComponents;
    using Base1 = CCMpfaDataHandleBases::SystemMatricesHandle<MatVecTraits, numPhases, numComponents>;
    using Base2 = CCMpfaDataHandleBases::SystemVectorsHandle<MatVecTraits, numPhases, numComponents>;
 
public:
    void setPhaseIndex(unsigned int phaseIdx) const
    { Base1::setContextIndex1(phaseIdx); Base2::setContextIndex1(phaseIdx); }
    void setComponentIndex(unsigned int compIdx) const
    { Base1::setContextIndex2(compIdx); Base2::setContextIndex2(compIdx); }
};
 
template<class MatVecTraits, class PhysicsTraits, bool enableHeatConduction>
class HeatConductionDataHandle
: public CCMpfaDataHandleBases::SystemMatricesHandle<MatVecTraits, 1, 1>
, public CCMpfaDataHandleBases::SystemVectorsHandle<MatVecTraits, 1, 1>
{};
 
template<class MatVecTraits, class PhysicsTraits>
class AdvectionDataHandle<MatVecTraits, PhysicsTraits, false> : public EmptyDataHandle {};
template<class MatVecTraits, class PhysicsTraits>
class DiffusionDataHandle<MatVecTraits, PhysicsTraits, false> : public EmptyDataHandle {};
template<class MatVecTraits, class PhysicsTraits>
class HeatConductionDataHandle<MatVecTraits, PhysicsTraits, false> : public EmptyDataHandle {};
 
template<class MVT, class PT>
class InteractionVolumeDataHandle
{
 
public:
    using AdvectionHandle = AdvectionDataHandle<MVT, PT, PT::enableAdvection>;
    using DiffusionHandle = DiffusionDataHandle<MVT, PT, PT::enableMolecularDiffusion>;
    using HeatConductionHandle = HeatConductionDataHandle<MVT, PT, PT::enableHeatConduction>;
 
    const AdvectionHandle& advectionHandle() const { return advectionHandle_; }
    AdvectionHandle& advectionHandle() { return advectionHandle_; }
 
    const DiffusionHandle& diffusionHandle() const { return diffusionHandle_; }
    DiffusionHandle& diffusionHandle() { return diffusionHandle_; }
 
    const HeatConductionHandle& heatConductionHandle() const { return heatConductionHandle_; }
    HeatConductionHandle& heatConductionHandle() { return heatConductionHandle_; }
 
private:
    AdvectionHandle advectionHandle_;
    DiffusionHandle diffusionHandle_;
    HeatConductionHandle heatConductionHandle_;
};
 
} // end namespace Dumux
 
#endif