62 using FVElementGeometry =
typename GridGeometry::LocalView;
63 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
68 static const int dimWorld = GridView::dimensionworld;
69 using DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
76 template <
class ElementFluxVariablesCache>
78 const SubControlVolumeFace& scvf,
79 const FVElementGeometry& fvGeometry,
80 const ElementVolumeVariables& elemVolVars,
81 const ElementFluxVariablesCache& elemFluxVarsCache,
86 DUNE_THROW(Dune::NotImplemented,
"Scheidegger dispersion tensors are only implemented for single phase flows.");
89 auto velocity = dispersionVelocity_(problem, scvf, fvGeometry, elemVolVars, elemFluxVarsCache);
92 std::array<Scalar,2> dispersivity = problem.spatialParams().dispersionAlphas(scvf.center(), phaseIdx, compIdx);
94 return scheideggerDispersionTensor_(dispersivity, velocity);
97 template <
class ElementFluxVariablesCache>
99 const SubControlVolumeFace& scvf,
100 const FVElementGeometry& fvGeometry,
101 const ElementVolumeVariables& elemVolVars,
102 const ElementFluxVariablesCache& elemFluxVarsCache,
106 DUNE_THROW(Dune::NotImplemented,
"Scheidegger dispersion tensors are only implemented for single phase flows.");
109 auto velocity = dispersionVelocity_(problem, scvf, fvGeometry, elemVolVars, elemFluxVarsCache);
112 std::array<Scalar,2> dispersivity = problem.spatialParams().dispersionAlphas(scvf.center(), phaseIdx);
114 return scheideggerDispersionTensor_(dispersivity, velocity);
119 template <
class ElementFluxVariablesCache>
120 static Dune::FieldVector<Scalar, dimWorld> dispersionVelocity_(
const Problem& problem,
121 const SubControlVolumeFace& scvf,
122 [[maybe_unused]]
const FVElementGeometry& fvGeometry,
123 [[maybe_unused]]
const ElementVolumeVariables& elemVolVars,
124 [[maybe_unused]]
const ElementFluxVariablesCache& elemFluxVarsCache)
127 Dune::FieldVector<Scalar, dimWorld> velocity(0.0);
128 if constexpr (stationaryVelocityField)
131 DUNE_THROW(Dune::NotImplemented,
"\n Please provide the stationary velocity field in the spatialparams via a velocity function.");
133 velocity = problem.spatialParams().velocity(scvf);
139 const auto& fluxVarsCache = elemFluxVarsCache[scvf];
140 const auto& shapeValues = fluxVarsCache.shapeValues();
143 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
144 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
146 outsideVolVars.permeability(),
147 scvf.unitOuterNormal());
150 Dune::FieldVector<Scalar, dimWorld> gradP(0.0);
153 for (
auto&& scv : scvs(fvGeometry))
155 const auto& volVars = elemVolVars[scv];
158 rho += volVars.density(0)*shapeValues[scv.indexInElement()][0];
160 gradP.axpy(volVars.pressure(0), fluxVarsCache.gradN(scv.indexInElement()));
164 gradP.axpy(-rho, problem.spatialParams().gravity(scvf.center()));
169 velocity /= -0.5 * (insideVolVars.viscosity() + outsideVolVars.viscosity());
172 DUNE_THROW(Dune::NotImplemented,
"\n Scheidegger Dispersion for compositional models without given constant velocity field is only implemented using the Box method.");
178 static DimWorldMatrix scheideggerDispersionTensor_(
const std::array<Scalar,2>& dispersivity,
179 const Dune::FieldVector<Scalar, dimWorld>& velocity)
181 DimWorldMatrix dispersionTensor(0.0);
184 for (
int i=0; i < dimWorld; i++)
185 for (
int j = 0; j < dimWorld; j++)
186 dispersionTensor[i][j] = velocity[i]*velocity[j];
189 Scalar vNorm = velocity.two_norm();
191 dispersionTensor /= vNorm;
193 dispersionTensor = 0;
196 dispersionTensor *= (dispersivity[0] - dispersivity[1]);
199 for (
int i = 0; i < dimWorld; i++)
200 dispersionTensor[i][i] += vNorm*dispersivity[1];
202 return dispersionTensor;
static DimWorldMatrix compositionalDispersionTensor(const Problem &problem, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const int phaseIdx, const int compIdx)
Definition scheidegger.hh:77
static DimWorldMatrix thermalDispersionTensor(const Problem &problem, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const int phaseIdx)
Definition scheidegger.hh:98
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Scalar faceTensorAverage(const Scalar T1, const Scalar T2, const Dune::FieldVector< Scalar, dim > &normal)
Average of a discontinuous scalar field at discontinuity interface (for compatibility reasons with th...
Definition facetensoraverage.hh:29