[ccwmpfa] Update stencils and related boundary volVars

dumux/discretization/cellcentered/wmpfa/elementvolumevariables.hh

@@ -112,6 +112,17 @@ namespace CCWMpfa {
                     volVars.emplace_back(std::move(dirichletVolVars));
                     volVarIndices.push_back(scvf.outsideScvIdx());
                 }
+                else
+                {
+                    VolumeVariables dirichletVolVars;
+                    dirichletVolVars.update(elementSolution<FVElemGeom>(problem.dirichlet(element, scvf)),
+                                            problem,
+                                            element,
+                                            scvI);
+
+                    volVars.emplace_back(std::move(dirichletVolVars));
+                    volVarIndices.push_back(scvf.outsideScvIdx());
+                }
             }
         }
         return std::make_pair(volVars, volVarIndices);
@@ -175,30 +186,30 @@ public:
               const FVElementGeometry& fvGeometry,
               const SolutionVector& sol)
     {
-        if (!fvGeometry.hasBoundaryScvf())
-            return;
-
         clear_();
 
-        auto&& [volVars, indices] = CCWMpfa::boundaryVolVars<VolumeVariables, std::size_t>(gridVolVars().problem(), element, fvGeometry);
-        std::move(volVars.begin(), volVars.end(), std::back_inserter(boundaryVolumeVariables_));
-        std::move(indices.begin(), indices.end(), std::back_inserter(boundaryVolVarIndices_));
+        if (fvGeometry.hasBoundaryScvf())
+        {
+            auto&& [volVars, indices] = CCWMpfa::boundaryVolVars<VolumeVariables, std::size_t>(gridVolVars().problem(), element, fvGeometry);
+            std::move(volVars.begin(), volVars.end(), std::back_inserter(boundaryVolumeVariables_));
+            std::move(indices.begin(), indices.end(), std::back_inserter(boundaryVolVarIndices_));
+        }
+
+        const auto& gridGeometry = fvGeometry.gridGeometry();
+        const auto globalI = gridGeometry.elementMapper().index(element);
+        const auto& assemblyMapI = gridGeometry.connectivityMap()[globalI];
 
         // add boundary volVars of neighbors
-        for (const auto& scvf : scvfs(fvGeometry))
+        for (const auto& dataJ : assemblyMapI)
         {
-            if (scvf.boundary())
-                continue;
-
-            const auto outsideScvIdx = scvf.outsideScvIdx();
-            const auto outsideElement = fvGeometry.gridGeometry().element(outsideScvIdx);
-            auto fvGeometryJ = localView(fvGeometry.gridGeometry());
-            fvGeometryJ.bind(outsideElement);
+            const auto& elementJ = gridGeometry.element(dataJ.globalJ);
+            auto fvGeometryJ = localView(gridGeometry);
+            fvGeometryJ.bind(elementJ);
 
             if (!fvGeometryJ.hasBoundaryScvf())
                 continue;
 
-            auto&& [volVars, indices]  = CCWMpfa::boundaryVolVars<VolumeVariables, std::size_t>(gridVolVars().problem(), outsideElement, fvGeometryJ);
+            auto&& [volVars, indices]  = CCWMpfa::boundaryVolVars<VolumeVariables, std::size_t>(gridVolVars().problem(), elementJ, fvGeometryJ);
             std::move(volVars.begin(), volVars.end(), std::back_inserter(boundaryVolumeVariables_));
             std::move(indices.begin(), indices.end(), std::back_inserter(boundaryVolVarIndices_));
         }
@@ -229,6 +240,10 @@ private:
     int getLocalIdx_(const int volVarIdx) const
     {
         auto it = std::find(boundaryVolVarIndices_.begin(), boundaryVolVarIndices_.end(), volVarIdx);
+        if(it == boundaryVolVarIndices_.end())
+        {
+            int shit = 1;
+        }
         assert(it != boundaryVolVarIndices_.end() && "Could not find the current volume variables for volVarIdx!");
         return std::distance(boundaryVolVarIndices_.begin(), it);
     }

dumux/discretization/fluxstencil.hh

@@ -147,12 +147,12 @@ public:
 
     //! Returns the indices of the elements required for flux calculation on an scvf.
     static Stencil stencil(const Element& element,
-                                  const FVElementGeometry& fvGeometry,
-                                  const SubControlVolumeFace& scvf)
+                           const FVElementGeometry& fvGeometry,
+                           const SubControlVolumeFace& scvf)
     {
         const auto& gridGeometry = fvGeometry.gridGeometry();
 
-        //ToDo correct stencil! At the moment simply all neighors of neighbors are added
+        //ToDo correct stencil! At the moment simply all neighors of the neighbor are added
         Stencil stencil({scvf.insideScvIdx()});
         for (const auto& scvf : scvfs(fvGeometry))
         {
@@ -160,7 +160,10 @@ public:
                 continue;
 
             stencil.push_back(scvf.outsideScvIdx());
+        }
 
+        if (!scvf.boundary())
+        {
             const auto outsideScvIdx = scvf.outsideScvIdx();
             const auto outsideElement = gridGeometry.element(outsideScvIdx);
             auto fvGeometryJ = localView(gridGeometry);

dumux/porousmediumflow/fluxvariablescachefiller.hh

@@ -732,10 +732,9 @@ class PorousMediumFluxVariablesCacheFillerImplementation<TypeTag, Discretization
     static constexpr bool heatConductionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentHeatConduction>();
 
 public:
-    //! This cache filler is always solution-dependent, as it updates the
-    //! vectors of cell unknowns with which the transmissibilities have to be
-    //! multiplied in order to obtain the fluxes.
-    static constexpr bool isSolDependent = true;
+    static constexpr bool isSolDependent = (advectionEnabled && advectionIsSolDependent) ||
+                                           (diffusionEnabled && diffusionIsSolDependent) ||
+                                           (heatConductionEnabled && heatConductionIsSolDependent);
 
     //! The constructor. Sets problem pointer.
     PorousMediumFluxVariablesCacheFillerImplementation(const Problem& problem)
@@ -817,6 +816,7 @@ private:
     template< class Handle, class IntOP, class TensorFunction >
     void prepareHandle_(Handle& handle, const IntOP& intOp, const TensorFunction& tensor, const FVElementGeometry& fvGeometry, const SubControlVolumeFace& scvf)
     {
+        handle.prepare();
         handle.decompose(intOp, tensor, fvGeometry, scvf);
     }
 
@@ -830,11 +830,11 @@ private:
                         const ElementVolumeVariables& elemVolVars,
                         const SubControlVolumeFace& scvf)
     {
-        // using AdvectionType = GetPropType<TypeTag, Properties::AdvectionType>;
-        // using AdvectionFiller = typename AdvectionType::Cache::Filler;
+        using AdvectionType = GetPropType<TypeTag, Properties::AdvectionType>;
+        using AdvectionFiller = typename AdvectionType::Cache::Filler;
 
-        // // forward to the filler for the advective quantities
-        // AdvectionFiller::fill(scvfFluxVarsCache, problem(), element, fvGeometry, elemVolVars, scvf, *this);
+        // forward to the filler for the advective quantities
+        AdvectionFiller::fill(scvfFluxVarsCache, problem(), element, fvGeometry, elemVolVars, scvf, *this);
     }
 
     // //! method to fill the diffusive quantities