Skip to content
GitLab
Commit f72e5874 authored by Kilian Weishaupt's avatar Kilian Weishaupt
Browse files

[staggered] First attempt to introduce localFaceVars

parent 9a7ef39c
Hide whitespace changes
Inline Side-by-side
dumux/assembly/staggeredlocalassembler.hh
View file @ f72e5874
......@@ -31,6 +31,7 @@
#include <dumux/assembly/diffmethod.hh>
#include <dumux/implicit/staggered/primaryvariables.hh>
#include <dumux/discretization/staggered/facesolution.hh>
namespace Dumux {
......@@ -77,9 +78,6 @@ class StaggeredLocalAssembler<TypeTag,
using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
using FacePrimaryVariables = typename GET_PROP_TYPE(TypeTag, FacePrimaryVariables);
using CellCenterPrimaryVariables = typename GET_PROP_TYPE(TypeTag, CellCenterPrimaryVariables);
// using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
// typename DofTypeIndices::CellCenterIdx cellCenterIdx;
// typename DofTypeIndices::FaceIdx faceIdx;
static constexpr bool enableGlobalFluxVarsCache = GET_PROP_VALUE(TypeTag, EnableGlobalFluxVariablesCache);
......@@ -144,7 +142,6 @@ public:
for(auto&& scvf : scvfs(fvGeometry))
{
// res[faceIdx][scvf.dofIndex()] = 0.0;
faceResidualCache[scvf.localFaceIdx()] = localResidual.evalFace(problem,
element,
fvGeometry,
......@@ -371,25 +368,30 @@ static void dCCdFace_(Assembler& assembler,
const auto& problem = assembler.problem();
auto& localResidual = assembler.localResidual();
// auto& gridVariables = assembler.gridVariables();
// build derivatives with for cell center dofs w.r.t. cell center dofs
const auto cellCenterGlobalI = assembler.fvGridGeometry().elementMapper().index(element);
const auto& connectivityMap = assembler.fvGridGeometry().connectivityMap();
// build derivatives with for cell center dofs w.r.t. cell center dofs
const auto cellCenterGlobalI = assembler.fvGridGeometry().elementMapper().index(element);
for(const auto& globalJ : connectivityMap(cellCenterIdx, faceIdx, cellCenterGlobalI))
for(const auto& scvfJ : scvfs(fvGeometry))
{
const auto globalJ = scvfJ.dofIndex();
// get the faceVars of the face with respect to which we are going to build the derivative
auto origFaceVars = curGlobalFaceVars.faceVars(globalJ);
auto& curFaceVars = curGlobalFaceVars.faceVars(globalJ);
auto origFaceVars = curGlobalFaceVars.faceVars(scvfJ.index());
auto& curFaceVars = curGlobalFaceVars.faceVars(scvfJ.index());
for(auto pvIdx : PriVarIndices(faceIdx))
{
PrimaryVariables priVars(CellCenterPrimaryVariables(0.0), FacePrimaryVariables(curSol[faceIdx][globalJ]));
auto faceSolution = StaggeredFaceSolution<TypeTag>(scvfJ, curSol[faceIdx], assembler.fvGridGeometry());
const Scalar eps = numericEpsilon(priVars[pvIdx], cellCenterIdx, faceIdx);
priVars[pvIdx] += eps;
curFaceVars.update(priVars[faceIdx]);
faceSolution[globalJ][pvIdx] += eps;
curFaceVars.update(scvfJ, faceSolution);
const auto deflectedResidual = localResidual.evalCellCenter(problem, element, fvGeometry,
prevElemVolVars, curElemVolVars,
......@@ -507,16 +509,17 @@ static void dFacedFace_(Assembler& assembler,
for(const auto& globalJ : connectivityMap(faceIdx, faceIdx, scvf.index()))
{
// get the faceVars of the face with respect to which we are going to build the derivative
auto origFaceVars = curGlobalFaceVars.faceVars(globalJ);
auto& curFaceVars = curGlobalFaceVars.faceVars(globalJ);
auto origFaceVars = curGlobalFaceVars.faceVars(scvf.index());
auto& curFaceVars = curGlobalFaceVars.faceVars(scvf.index());
for(auto pvIdx : PriVarIndices(faceIdx))
{
PrimaryVariables priVars(CellCenterPrimaryVariables(0.0), FacePrimaryVariables(curSol[faceIdx][globalJ]));
auto faceSolution = StaggeredFaceSolution<TypeTag>(scvf, curSol[faceIdx], assembler.fvGridGeometry());
const Scalar eps = numericEpsilon(priVars[pvIdx], faceIdx, faceIdx);
priVars[pvIdx] += eps;
curFaceVars.update(priVars[faceIdx]);
const Scalar eps = numericEpsilon(faceSolution[globalJ][pvIdx], faceIdx, faceIdx);
faceSolution[globalJ][pvIdx] += eps;
curFaceVars.update(scvf, faceSolution);
const auto deflectedResidual = localResidual.evalFace(problem, element, fvGeometry, scvf,
prevElemVolVars, curElemVolVars,
......
dumux/discretization/staggered/facesolution.hh 0 → 100644
View file @ f72e5874
// -*- 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 2 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/>. *
*****************************************************************************/
/*!
* \file
* \brief The global volume variables class for cell centered models
*/
#ifndef DUMUX_DISCRETIZATION_STAGGERED_FACE_SOLUTION_HH
#define DUMUX_DISCRETIZATION_STAGGERED_FACE_SOLUTION_HH
#include <dumux/common/basicproperties.hh>
#include <dumux/discretization/staggered/elementvolumevariables.hh>
namespace Dumux
{
template<class TypeTag>
class StaggeredFaceSolution
{
using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
using Element = typename GridView::template Codim<0>::Entity;
using FaceSolutionVector = typename GET_PROP_TYPE(TypeTag, FaceSolutionVector);
using FacePrimaryVariables = typename GET_PROP_TYPE(TypeTag, FacePrimaryVariables);
using FVGridGeometry = typename GET_PROP_TYPE(TypeTag, FVGridGeometry);
using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
typename DofTypeIndices::CellCenterIdx cellCenterIdx;
typename DofTypeIndices::FaceIdx faceIdx;
public:
StaggeredFaceSolution(const SubControlVolumeFace& scvf, const FaceSolutionVector& sol,
const FVGridGeometry& fvGridGeometry)
{
const auto& connectivityMap = fvGridGeometry.connectivityMap();
const auto& stencil = connectivityMap(faceIdx, faceIdx, scvf.index());
facePriVars_.clear();
map_.clear();
for(const auto dofJ : stencil)
{
map_.push_back(dofJ);
facePriVars_.push_back(sol[dofJ]);
}
}
//! bracket operator const access
template<typename IndexType>
const FacePrimaryVariables& operator [](IndexType globalFaceDofIdx) const
{
const auto pos = std::find(map_.begin(), map_.end(), globalFaceDofIdx);
assert (pos != map_.end());
return facePriVars_[pos - map_.begin()];
}
//! bracket operator
template<typename IndexType>
FacePrimaryVariables& operator [](IndexType globalFaceDofIdx)
{
const auto pos = std::find(map_.begin(), map_.end(), globalFaceDofIdx);
assert (pos != map_.end());
return facePriVars_[pos - map_.begin()];
}
private:
std::vector<FacePrimaryVariables> facePriVars_;
std::vector<unsigned int> map_;
};
} // end namespace
#endif
dumux/discretization/staggered/freeflow/facevariables.hh
View file @ f72e5874
......@@ -27,25 +27,105 @@
namespace Dumux
{
namespace Properties
{
NEW_PROP_TAG(StaggeredFaceSolution);
}
template<class TypeTag>
class StaggeredFaceVariables
{
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using FaceSolutionVector = typename GET_PROP_TYPE(TypeTag, FaceSolutionVector);
using FacePrimaryVariables = typename GET_PROP_TYPE(TypeTag, FacePrimaryVariables);
using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
using FaceSolution = typename GET_PROP_TYPE(TypeTag, StaggeredFaceSolution);
struct SubFaceData
{
Scalar velocityNormalInside;
Scalar velocityNormalOutside;
Scalar velocityParallelInside;
Scalar velocityParallelOutside;
};
public:
void update(const FacePrimaryVariables &facePrivars)
{
velocity_ = facePrivars[0];
velocitySelf_ = facePrivars;
}
void update(const SubControlVolumeFace& scvf, const FaceSolutionVector& sol)
{
velocitySelf_ = sol[scvf.dofIndex()];
velocityOpposite_ = sol[scvf.dofIndexOpposingFace()];
subFaceVelocities_.resize(scvf.pairData().size());
for(int i = 0; i < scvf.pairData().size(); ++i)
{
subFaceVelocities_[i].velocityNormalInside = sol[scvf.pairData(i).normalPair.first];
if(scvf.pairData(i).normalPair.second >= 0)
subFaceVelocities_[i].velocityNormalOutside = sol[scvf.pairData(i).normalPair.second];
subFaceVelocities_[i].velocityParallelInside = sol[scvf.dofIndex()];
if(scvf.pairData(i).outerParallelFaceDofIdx >= 0)
subFaceVelocities_[i].velocityParallelOutside = sol[scvf.pairData(i).outerParallelFaceDofIdx];
}
}
void update(const SubControlVolumeFace& scvf, const FaceSolution& faceSol)
{
velocitySelf_ = faceSol[scvf.dofIndex()];
velocityOpposite_ = faceSol[scvf.dofIndexOpposingFace()];
subFaceVelocities_.resize(scvf.pairData().size());
for(int i = 0; i < scvf.pairData().size(); ++i)
{
subFaceVelocities_[i].velocityNormalInside = faceSol[scvf.pairData(i).normalPair.first];
if(scvf.pairData(i).normalPair.second >= 0)
subFaceVelocities_[i].velocityNormalOutside = faceSol[scvf.pairData(i).normalPair.second];
subFaceVelocities_[i].velocityParallelInside = faceSol[scvf.dofIndex()];
if(scvf.pairData(i).outerParallelFaceDofIdx >= 0)
subFaceVelocities_[i].velocityParallelOutside = faceSol[scvf.pairData(i).outerParallelFaceDofIdx];
}
}
Scalar velocitySelf() const
{
return velocitySelf_;
}
Scalar velocityOpposite() const
{
return velocityOpposite_;
}
auto& subFaceData() const
{
return subFaceVelocities_;
}
Scalar velocity() const
auto& subFaceData(const int localSubFaceIdx) const
{
return velocity_;
return subFaceVelocities_[localSubFaceIdx];
}
private:
Scalar velocity_;
Scalar velocitySelf_;
Scalar velocityOpposite_;
std::vector<SubFaceData> subFaceVelocities_;
};
} // end namespace
......
dumux/discretization/staggered/globalfacevariables.hh
View file @ f72e5874
......@@ -24,6 +24,7 @@
#define DUMUX_DISCRETIZATION_STAGGERED_GLOBAL_FACEVARIABLES_HH
#include <dumux/common/basicproperties.hh>
#include <dumux/discretization/staggered/facesolution.hh>
namespace Dumux
{
......@@ -51,15 +52,36 @@ public:
const auto& faceSol = sol[faceIdx];
const auto numFaceDofs = fvGridGeometry.gridView().size(1);
faceVariables_.resize(numFaceDofs);
assert(faceVariables_.size() == faceSol.size());
for(int i = 0; i < numFaceDofs; ++i)
// const auto numFaceDofs = fvGridGeometry.gridView().size(1);
// faceVariables_.resize(numFaceDofs);
faceVariables_.resize(fvGridGeometry.numScvf());
for(auto&& element : elements(fvGridGeometry.gridView()))
{
faceVariables_[i].update(faceSol[i]);
auto fvGeometry = localView(fvGridGeometry);
fvGeometry.bindElement(element);
for(auto&& scvf : scvfs(fvGeometry))
{
faceVariables_[scvf.index()].update(scvf, faceSol);
auto test = StaggeredFaceSolution<TypeTag>(scvf, faceSol, fvGridGeometry);
// std::cout << "test" << test[scvf.dofIndex()] << std::endl;
}
}
// for(auto& faceVariable : faceVariables_)
// {
// faceVariable.update()
// }
// // assert(faceVariables_.size() == faceSol.size());
//
// for(int i = 0; i < numFaceDofs; ++i)
// {
// faceVariables_[i].update(faceSol[i]);
// }
}
const FaceVariables& faceVars(const IndexType facetIdx) const
......
dumux/discretization/staggered/properties.hh
View file @ f72e5874
......@@ -44,6 +44,7 @@
#include <dumux/discretization/staggered/fvgridgeometry.hh>
#include <dumux/discretization/staggered/fvelementgeometry.hh>
#include <dumux/discretization/staggered/globalfacevariables.hh>
#include <dumux/discretization/staggered/facesolution.hh>
#include <dumux/common/intersectionmapper.hh>
#include <dune/istl/multitypeblockvector.hh>
......@@ -60,6 +61,7 @@ namespace Properties
NEW_PROP_TAG(CellCenterSolutionVector);
NEW_PROP_TAG(FaceSolutionVector);
NEW_PROP_TAG(StaggeredFaceSolution);
//! Type tag for the box scheme.
NEW_TYPE_TAG(StaggeredModel, INHERITS_FROM(FiniteVolumeModel, NumericModel));
......@@ -109,6 +111,8 @@ SET_TYPE_PROP(StaggeredModel, BaseLocalResidual, StaggeredLocalResidual<TypeTag>
SET_TYPE_PROP(StaggeredModel, IntersectionMapper, Dumux::ConformingGridIntersectionMapper<TypeTag>);
SET_TYPE_PROP(StaggeredModel, StaggeredFaceSolution, StaggeredFaceSolution<TypeTag>);
//! Definition of the indices for cell center and face dofs in the global solution vector
SET_PROP(StaggeredModel, DofTypeIndices)
{
......
dumux/freeflow/staggered/fluxvariables.hh
View file @ f72e5874
......@@ -111,7 +111,8 @@ public:
const auto& outsideVolVars = scvf.boundary() ? insideVolVars : elemVolVars[scvf.outsideScvIdx()];
CellCenterPrimaryVariables flux(0.0);
const Scalar velocity = globalFaceVars.faceVars(scvf.dofIndex()).velocity();
// const Scalar velocity = globalFaceVars.faceVars(scvf.dofIndex()).velocity();
const Scalar velocity = globalFaceVars.faceVars(scvf.index()).velocitySelf();
const bool insideIsUpstream = sign(scvf.outerNormalScalar()) == sign(velocity) ? true : false;
const auto& upstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
......@@ -201,8 +202,8 @@ public:
{
const auto insideScvIdx = scvf.insideScvIdx();
const auto& insideVolVars = elemVolVars[insideScvIdx];
const Scalar velocitySelf = globalFaceVars.faceVars(scvf.dofIndex()).velocity() ;
const Scalar velocityOpposite = globalFaceVars.faceVars(scvf.dofIndexOpposingFace()).velocity();
const Scalar velocitySelf = globalFaceVars.faceVars(scvf.index()).velocitySelf() ;
const Scalar velocityOpposite = globalFaceVars.faceVars(scvf.index()).velocityOpposite();
FacePrimaryVariables normalFlux(0.0);
if(navierStokes)
......@@ -253,20 +254,24 @@ public:
{
FacePrimaryVariables tangentialFlux(0.0);
// convenience function to get the velocity on a face
auto velocity = [&globalFaceVars](const int dofIdx)
{
return globalFaceVars.faceVars(dofIdx).velocity();
};
// // convenience function to get the velocity on a face
// auto velocity = [&globalFaceVars](const int dofIdx)
// {
// return globalFaceVars.faceVars(scvf.dofIdx()).velocity();
// };
auto& faceVars = globalFaceVars.faceVars(scvf.index());
const int numSubFaces = scvf.pairData().size();
// account for all sub-faces
for(const auto& subFaceData : scvf.pairData())
// for(const auto& subFaceData : scvf.pairData())
for(int localSubFaceIdx = 0; localSubFaceIdx < numSubFaces; ++localSubFaceIdx)
{
const auto eIdx = scvf.insideScvIdx();
const auto& normalFace = fvGeometry.scvf(eIdx, subFaceData.localNormalFaceIdx);
const auto& normalFace = fvGeometry.scvf(eIdx, scvf.pairData()[localSubFaceIdx].localNormalFaceIdx);
// Check if we have a symmetry boundary condition. If yes, the tangental part of the momentum flux can be neglected.
if(subFaceData.outerParallelFaceDofIdx < 0)
if(scvf.pairData()[localSubFaceIdx].outerParallelFaceDofIdx < 0)
{
// lambda to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest
auto makeGhostFace = [eIdx] (const GlobalPosition& pos)
......@@ -275,32 +280,33 @@ public:
};
// use the ghost face to check if there is a symmetry boundary condition and skip any further steps if yes
const auto bcTypes = problem.boundaryTypes(element, makeGhostFace(subFaceData.virtualOuterParallelFaceDofPos));
const auto bcTypes = problem.boundaryTypes(element, makeGhostFace(scvf.pairData()[localSubFaceIdx].virtualOuterParallelFaceDofPos));
if(bcTypes.isSymmetry())
continue;
}
// if there is no symmetry boundary condition, proceed to calculate the tangential momentum flux
if(navierStokes)
tangentialFlux += computeAdvectivePartOfTangentialMomentumFlux_(problem, element, scvf, normalFace, subFaceData, elemVolVars, velocity);
tangentialFlux += computeAdvectivePartOfTangentialMomentumFlux_(problem, element, scvf, normalFace, elemVolVars, faceVars, localSubFaceIdx);
tangentialFlux += computeDiffusivePartOfTangentialMomentumFlux_(problem, element, scvf, normalFace, subFaceData, elemVolVars, velocity);
tangentialFlux += computeDiffusivePartOfTangentialMomentumFlux_(problem, element, scvf, normalFace, elemVolVars, faceVars, localSubFaceIdx);
}
return tangentialFlux;
}
private:
template<class SubFaceData, class VelocityHelper>
template<class FaceVars>
FacePrimaryVariables computeAdvectivePartOfTangentialMomentumFlux_(const Problem& problem,
const Element& element,
const SubControlVolumeFace& scvf,
const SubControlVolumeFace& normalFace,
const SubFaceData& subFaceData,
const ElementVolumeVariables& elemVolVars,
VelocityHelper velocity)
const FaceVars& faceVars,
const int localSubFaceIdx)
{
const Scalar transportingVelocity = velocity(subFaceData.normalPair.first);
// const Scalar transportingVelocity = velocity(subFaceData.normalPair.first);
const Scalar transportingVelocity = faceVars.subFaceData(localSubFaceIdx).velocityNormalInside;
const auto insideScvIdx = normalFace.insideScvIdx();
const auto outsideScvIdx = normalFace.outsideScvIdx();
......@@ -317,12 +323,14 @@ private:
Scalar transportedVelocity(0.0);
if(innerElementIsUpstream)
transportedVelocity = velocity(scvf.dofIndex());
transportedVelocity = faceVars.subFaceData(localSubFaceIdx).velocityParallelInside;
// transportedVelocity = velocity(scvf.dofIndex());
else
{
const int outerDofIdx = subFaceData.outerParallelFaceDofIdx;
const int outerDofIdx = scvf.pairData(localSubFaceIdx).outerParallelFaceDofIdx;
if(outerDofIdx >= 0)
transportedVelocity = velocity(outerDofIdx);
// transportedVelocity = velocity(outerDofIdx);
transportedVelocity = faceVars.subFaceData(localSubFaceIdx).velocityParallelOutside;
else // this is the case when the outer parallal dof would lie outside the domain TODO: discuss which one is better
// transportedVelocity = problem.dirichlet(makeGhostFace(subFaceData.virtualOuterParallelFaceDofPos))[faceIdx][scvf.directionIndex()];
transportedVelocity = problem.dirichlet(element, scvf)[faceIdx][scvf.directionIndex()];
......@@ -334,14 +342,14 @@ private:
return transportingVelocity * momentum * sgn * normalFace.area() * 0.5;
}
template<class SubFaceData, class VelocityHelper>
template<class FaceVars>
FacePrimaryVariables computeDiffusivePartOfTangentialMomentumFlux_(const Problem& problem,
const Element& element,
const SubControlVolumeFace& scvf,
const SubControlVolumeFace& normalFace,
const SubFaceData& subFaceData,
const ElementVolumeVariables& elemVolVars,
VelocityHelper velocity)
const FaceVars& faceVars,
const int localSubFaceIdx)
{
FacePrimaryVariables tangentialDiffusiveFlux(0.0);
......@@ -362,35 +370,38 @@ private:
const Scalar muAvg = (insideVolVars.viscosity() + outsideVolVars.viscosity()) * 0.5;
// the normal derivative
const int innerNormalVelocityIdx = subFaceData.normalPair.first;
const int outerNormalVelocityIdx = subFaceData.normalPair.second;
// const int innerNormalVelocityIdx = subFaceData.normalPair.first;
const int outerNormalVelocityIdx = scvf.pairData(localSubFaceIdx).normalPair.second;
const Scalar innerNormalVelocity = velocity(innerNormalVelocityIdx);
// const Scalar innerNormalVelocity = velocity(innerNormalVelocityIdx);
const Scalar innerNormalVelocity = faceVars.subFaceData(localSubFaceIdx).velocityNormalInside;
const Scalar outerNormalVelocity = outerNormalVelocityIdx >= 0 ?
velocity(outerNormalVelocityIdx) :
problem.dirichlet(element, makeGhostFace(subFaceData.virtualOuterNormalFaceDofPos))[faceIdx][normalDirIdx];
faceVars.subFaceData(localSubFaceIdx).velocityNormalOutside :
problem.dirichlet(element, makeGhostFace(scvf.pairData(localSubFaceIdx).virtualOuterNormalFaceDofPos))[faceIdx][normalDirIdx];
const Scalar normalDeltaV = scvf.normalInPosCoordDir() ?
(outerNormalVelocity - innerNormalVelocity) :
(innerNormalVelocity - outerNormalVelocity);
const Scalar normalDerivative = normalDeltaV / subFaceData.normalDistance;
const Scalar normalDerivative = normalDeltaV / scvf.pairData(localSubFaceIdx).normalDistance;
tangentialDiffusiveFlux -= muAvg * normalDerivative;
// the parallel derivative
const Scalar innerParallelVelocity = velocity(scvf.dofIndex());
// const Scalar innerParallelVelocity = velocity(scvf.dofIndex());
const Scalar innerParallelVelocity = faceVars.subFaceData(localSubFaceIdx).velocityParallelInside;
const int outerParallelFaceDofIdx = subFaceData.outerParallelFaceDofIdx;
const int outerParallelFaceDofIdx = scvf.pairData(localSubFaceIdx).outerParallelFaceDofIdx;
const Scalar outerParallelVelocity = outerParallelFaceDofIdx >= 0 ?
velocity(outerParallelFaceDofIdx) :
problem.dirichlet(element, makeGhostFace(subFaceData.virtualOuterParallelFaceDofPos))[faceIdx][scvf.directionIndex()];
faceVars.subFaceData(localSubFaceIdx).velocityParallelOutside :
// velocity(outerParallelFaceDofIdx) :
problem.dirichlet(element, makeGhostFace(scvf.pairData(localSubFaceIdx).virtualOuterParallelFaceDofPos))[faceIdx][scvf.directionIndex()];
const Scalar parallelDeltaV = normalFace.normalInPosCoordDir() ?
(outerParallelVelocity - innerParallelVelocity) :
(innerParallelVelocity - outerParallelVelocity);
const Scalar parallelDerivative = parallelDeltaV / subFaceData.parallelDistance;
const Scalar parallelDerivative = parallelDeltaV / scvf.pairData(localSubFaceIdx).parallelDistance;
tangentialDiffusiveFlux -= muAvg * parallelDerivative;
const Scalar sgn = sign(normalFace.outerNormalScalar());
......
dumux/freeflow/staggered/localresidual.hh
View file @ f72e5874
......@@ -182,7 +182,7 @@ public:
const GlobalFaceVars& globalFaceVars)
{
FacePrimaryVariables storage(0.0);
const Scalar velocity = globalFaceVars.faceVars(scvf.dofIndex()).velocity();
const Scalar velocity = globalFaceVars.faceVars(scvf.index()).velocity();
storage[0] = volVars.density() * velocity;
return storage;
}
......@@ -327,7 +327,8 @@ protected:
// set a fixed value for the velocity for Dirichlet boundary conditions
if(bcTypes.isDirichlet(momentumBalanceIdx))
{
const Scalar velocity = faceVars.faceVars(scvf.dofIndex()).velocity();
// const Scalar velocity = faceVars.faceVars(scvf.dofIndex()).velocity();
const Scalar velocity = faceVars.faceVars(scvf.index()).velocitySelf();
const Scalar dirichletValue = problem.dirichlet(element, scvf)[faceIdx][scvf.directionIndex()];
residual = velocity - dirichletValue;
}
......@@ -336,7 +337,8 @@ protected:
// we therefore treat it like a Dirichlet boundary conditions with zero velocity
if(bcTypes.isSymmetry())
{
const Scalar velocity = faceVars.faceVars(scvf.dofIndex()).velocity();
// const Scalar velocity = faceVars.faceVars(scvf.dofIndex()).velocity();
const Scalar velocity = faceVars.faceVars(scvf.index()).velocitySelf();
const Scalar fixedValue = 0.0;
residual = velocity - fixedValue;
}
......
dumux/freeflow/staggered/velocityoutput.hh
View file @ f72e5874
......@@ -113,9 +113,9 @@ public:
for (auto&& scvf : scvfs(fvGeometry))
{
auto& origFaceVars = gridVariables_.curGridFaceVars().faceVars(scvf.dofIndex());
auto& origFaceVars = gridVariables_.curGridFaceVars().faceVars(scvf.index());
auto dirIdx = scvf.directionIndex();
velocity[dofIdxGlobal][dirIdx] += 0.5*origFaceVars.velocity();
velocity[dofIdxGlobal][dirIdx] += 0.5*origFaceVars.velocitySelf();
}
}