diff --git a/dumux/discretization/staggered/freeflow/staggeredgeometryhelper.hh b/dumux/discretization/staggered/freeflow/staggeredgeometryhelper.hh index 31acc832cf5816ccc8ec35dd68fd5ed7f4d48fac..73ea53b26577dcea741ca4206763c9d4de302e13 100644 --- a/dumux/discretization/staggered/freeflow/staggeredgeometryhelper.hh +++ b/dumux/discretization/staggered/freeflow/staggeredgeometryhelper.hh @@ -57,7 +57,7 @@ struct PairData std::pair<GridIndexType, GridIndexType> lateralPair; SmallLocalIndexType localLateralFaceIdx; Scalar lateralDistance; - GlobalPosition virtualBoundaryFaceDofPos; + GlobalPosition lateralStaggeredFaceCenter; }; @@ -370,6 +370,10 @@ private: { const auto innerElementIntersectionIdx = innerElementIntersection.indexInInside(); setLateralPairFirstInfo_(innerElementIntersectionIdx, element_, numPairInnerLateralIdx); + + const auto distance = innerElementIntersection.geometry().center() - selfElementCenter; + pairData_[numPairInnerLateralIdx].lateralStaggeredFaceCenter = std::move(selfFacetCenter + distance); + numPairInnerLateralIdx++; } } @@ -399,8 +403,6 @@ private: if (facetIsNormal_(intersection.indexInInside(), intersection_.indexInInside())) { assert(!pairData_[numPairOuterLateralIdx].hasOuterLateral); - const auto boundaryDistanceOffset = intersection.geometry().center() - selfElementCenter; - pairData_[numPairOuterLateralIdx].virtualBoundaryFaceDofPos = std::move(boundaryDistanceOffset + selfElementCenter); const auto normalDistanceoffset = selfFacetCenter - selfElementCenter; pairData_[numPairOuterLateralIdx].lateralDistance = std::move(normalDistanceoffset.two_norm()); @@ -419,8 +421,6 @@ private: void fillParallelPairData_(std::false_type) { // set basic global positions and stencil size definitions - const auto& elementCenter = element_.geometry().center(); - // get the parallel Dofs const auto parallelLocalIdx = intersection_.indexInInside(); SmallLocalIndexType numPairParallelIdx = 0; @@ -438,15 +438,7 @@ private: pairData_[numPairParallelIdx].parallelDofs[0] = gridView_.indexSet().subIndex(outerElement, parallelLocalIdx, codimIntersection); pairData_[numPairParallelIdx].parallelCellWidths[0] = setParallelPairCellWidths_(outerElement, parallelAxisIdx); } - else // No parallel neighbor available - { - // If the intersection has no neighbor we have to deal with the virtual outer parallel dof - const auto& boundaryFacetCenter = intersection.geometry().center(); - const auto distance = boundaryFacetCenter - elementCenter; - const auto virtualBoundaryFaceDofPos = intersection_.geometry().center() + distance; - pairData_[numPairParallelIdx].virtualBoundaryFaceDofPos = std::move(virtualBoundaryFaceDofPos); - } numPairParallelIdx++; } } @@ -460,7 +452,6 @@ private: { // set basic global positions and stencil size definitions const auto numParallelFaces = pairData_[0].parallelCellWidths.size(); - const auto& elementCenter = element_.geometry().center(); // get the parallel Dofs const auto parallelLocalIdx = intersection_.indexInInside(); @@ -506,16 +497,7 @@ private: } } - else - { - // If the intersection has no neighbor we have to deal with the virtual outer parallel dof - const auto& boundaryFacetCenter = intersection.geometry().center(); - const auto distance = boundaryFacetCenter - elementCenter; - const auto virtualBoundaryFaceDofPos = intersection_.geometry().center() + distance; - - pairData_[numPairParallelIdx].virtualBoundaryFaceDofPos = std::move(virtualBoundaryFaceDofPos); - } numPairParallelIdx++; } } diff --git a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh index 0d955ec3f87b5f0af7ec7ba76e854c0ef000ac3f..2276366472f0922e72527b3f114df0aaa2b13362 100644 --- a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh +++ b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh @@ -287,7 +287,7 @@ public: { const auto eIdx = scvf.insideScvIdx(); // Get the face normal to the face the dof lives on. The staggered sub face conincides with half of this lateral face. - const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx); + const auto& lateralScvf = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx); // create a boundaryTypes object (will be empty if not at a boundary) Dune::Std::optional<BoundaryTypes> lateralFaceBoundaryTypes; @@ -298,7 +298,7 @@ public: { // Construct a temporary scvf which corresponds to the staggered sub face, featuring the location // the sub faces's center. - auto lateralBoundaryFaceCenter = scvf.pairData(localSubFaceIdx).virtualBoundaryFaceDofPos + lateralFace.center(); + auto lateralBoundaryFaceCenter = scvf.pairData(localSubFaceIdx).lateralStaggeredFaceCenter + lateralScvf.center(); lateralBoundaryFaceCenter *= 0.5; // ________________ @@ -310,7 +310,7 @@ public: // | || | o position at which the boundary conditions will be evaluated // ---------------- (lateralBoundaryFaceCenter) - const auto lateralBoundaryFace = lateralFace.makeBoundaryFace(lateralBoundaryFaceCenter); + const auto lateralBoundaryFace = lateralScvf.makeBoundaryFace(lateralBoundaryFaceCenter); // Retrieve the boundary types that correspond to the sub face. lateralFaceBoundaryTypes.emplace(problem.boundaryTypes(element, lateralBoundaryFace)); @@ -324,7 +324,7 @@ public: if (lateralFaceBoundaryTypes->isNeumann(Indices::velocity(scvf.directionIndex()))) { lateralFlux += problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, lateralBoundaryFace)[Indices::velocity(scvf.directionIndex())] - * elemVolVars[lateralFace.insideScvIdx()].extrusionFactor() * lateralFace.area() * 0.5; + * elemVolVars[lateralScvf.insideScvIdx()].extrusionFactor() * lateralScvf.area() * 0.5; continue; } @@ -457,7 +457,7 @@ private: { // Create a boundaryTypes object (will be empty if not at a boundary). Dune::Std::optional<BoundaryTypes> bcTypes; - const auto& boundaryFace = scvf.makeBoundaryFace(scvf.pairData(localSubFaceIdx).virtualBoundaryFaceDofPos); + const auto& boundaryFace = scvf.makeBoundaryFace(scvf.pairData(localSubFaceIdx).lateralStaggeredFaceCenter); // Get the boundary conditions if we are at a boundary. if (scvf.boundary()) @@ -612,7 +612,7 @@ private: // | || | o position at which the boundary conditions will be evaluated // ---------------- - return ownScvf.makeBoundaryFace(ownScvf.pairData(localSubFaceIdx).virtualBoundaryFaceDofPos); + return ownScvf.makeBoundaryFace(ownScvf.pairData(localSubFaceIdx).lateralStaggeredFaceCenter); }; //! helper function to get the averaged extrusion factor for a face diff --git a/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh b/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh index 7ec3833120e69572eec9268df480012e6d103e4c..29baf9ae1b31d94b02f13a176cdf2abe5094ad6b 100644 --- a/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh +++ b/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh @@ -426,7 +426,7 @@ private: momenta[2] = faceVars.velocityParallel(localSubFaceIdx, 1) * outsideVolVars.density(); else { - const Element& elementParallel = fvGeometry.fvGridGeometry().element(fvGeometry.scv(lateralFace.outsideScvIdx())); + const Element& elementParallel = fvGeometry.fvGridGeometry().element(lateralFace.outsideScvIdx()); const SubControlVolumeFace& firstParallelScvf = fvGeometry.scvf(lateralFace.outsideScvIdx(), ownScvf.localFaceIdx()); momenta[2] = getParallelVelocityFromOtherBoundary_(problem, fvGeometry, firstParallelScvf, localSubFaceIdx, elementParallel, @@ -626,8 +626,8 @@ private: const Scalar parallelVelocity) { // A ghost subface at the boundary is created, featuring the location of the sub face's center - const SubControlVolumeFace& lateralFace = fvGeometry.scvf(scvf.insideScvIdx(), scvf.pairData(localIdx).localLateralFaceIdx); - GlobalPosition lateralBoundaryFaceCenter = scvf.pairData(localIdx).virtualBoundaryFaceDofPos + lateralFace.center(); + const SubControlVolumeFace& lateralScvf = fvGeometry.scvf(scvf.insideScvIdx(), scvf.pairData(localIdx).localLateralFaceIdx); + GlobalPosition lateralBoundaryFaceCenter = scvf.pairData(localIdx).lateralStaggeredFaceCenter + lateralScvf.center(); lateralBoundaryFaceCenter *= 0.5; // ________________ @@ -639,7 +639,7 @@ private: // | || | o position at which the boundary conditions will be evaluated // ---------------- (lateralBoundaryFaceCenter) - const SubControlVolumeFace lateralBoundaryFace = lateralFace.makeBoundaryFace(lateralBoundaryFaceCenter); + const SubControlVolumeFace lateralBoundaryFace = lateralScvf.makeBoundaryFace(lateralBoundaryFaceCenter); // The boundary condition is checked, in case of symmetry or Dirichlet for the pressure // a gradient of zero is assumed in the direction normal to the bounadry, while if there is @@ -665,7 +665,7 @@ private: else if (bcTypes.isBJS(Indices::velocity(scvf.directionIndex()))) { const SubControlVolume& scv = fvGeometry.scv(scvf.insideScvIdx()); - return problem.bjsVelocity(boundaryElement, scv, lateralFace, parallelVelocity); + return problem.bjsVelocity(boundaryElement, scv, lateralScvf, parallelVelocity); } else { @@ -677,7 +677,7 @@ private: //! helper function to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest static SubControlVolumeFace makeParallelGhostFace_(const SubControlVolumeFace& ownScvf, const int localSubFaceIdx) { - return ownScvf.makeBoundaryFace(ownScvf.pairData(localSubFaceIdx).virtualBoundaryFaceDofPos); + return ownScvf.makeBoundaryFace(ownScvf.pairData(localSubFaceIdx).lateralStaggeredFaceCenter); }; };