Commit b49e1cc6 authored by Dennis Gläser's avatar Dennis Gläser Committed by Timo Koch
Browse files

[disc][scv] rename indexInElement() to localDofIndex()

parent 7b811b59
......@@ -226,8 +226,8 @@ public:
// Get bcTypes and maybe mark for refinement on Dirichlet boundaries
for (const auto& scv : scvs(fvGeometry))
{
bcTypes[scv.indexInElement()] = problem_->boundaryTypes(element, scv);
if (refineAtDirichletBC_ && bcTypes[scv.indexInElement()].hasDirichlet())
bcTypes[scv.localDofIndex()] = problem_->boundaryTypes(element, scv);
if (refineAtDirichletBC_ && bcTypes[scv.localDofIndex()].hasDirichlet())
{
indicatorVector_[eIdx] = true;
break; // element is marked, escape scv loop
......
......@@ -86,13 +86,13 @@ public:
{
const auto residual = this->asImp_().evalLocalResidual(); // forward to the internal implementation
for (const auto& scv : scvs(this->fvGeometry()))
res[scv.dofIndex()] += residual[scv.indexInElement()];
res[scv.dofIndex()] += residual[scv.localDofIndex()];
}
else
{
const auto residual = this->asImp_().assembleJacobianAndResidualImpl(jac, gridVariables, partialReassembler); // forward to the internal implementation
for (const auto& scv : scvs(this->fvGeometry()))
res[scv.dofIndex()] += residual[scv.indexInElement()];
res[scv.dofIndex()] += residual[scv.localDofIndex()];
}
auto applyDirichlet = [&] (const auto& scvI,
......@@ -104,7 +104,7 @@ public:
for (const auto& scvJ : scvs(this->fvGeometry()))
{
jac[scvI.dofIndex()][scvJ.dofIndex()][eqIdx] = 0.0;
if (scvI.indexInElement() == scvJ.indexInElement())
if (scvI.localDofIndex() == scvJ.localDofIndex())
jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
}
};
......@@ -129,7 +129,7 @@ public:
for (const auto& scvJ : scvs(this->fvGeometry()))
{
jac[scvI.dofIndex()][scvJ.dofIndex()][eqIdx] = 0.0;
if (scvI.indexInElement() == scvJ.indexInElement())
if (scvI.localDofIndex() == scvJ.localDofIndex())
jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
}
};
......@@ -146,7 +146,7 @@ public:
const auto residual = this->evalLocalResidual();
for (const auto& scv : scvs(this->fvGeometry()))
res[scv.dofIndex()] += residual[scv.indexInElement()];
res[scv.dofIndex()] += residual[scv.localDofIndex()];
auto applyDirichlet = [&] (const auto& scvI,
const auto& dirichletValues,
......@@ -171,7 +171,7 @@ public:
{
for (const auto& scvI : scvs(this->fvGeometry()))
{
const auto bcTypes = this->elemBcTypes()[scvI.indexInElement()];
const auto bcTypes = this->elemBcTypes()[scvI.localDofIndex()];
if (bcTypes.hasDirichlet())
{
const auto dirichletValues = this->problem().dirichlet(this->element(), scvI);
......@@ -281,7 +281,7 @@ public:
auto evalResiduals = [&](Scalar priVar)
{
// update the volume variables and compute element residual
elemSol[scv.indexInElement()][pvIdx] = priVar;
elemSol[scv.localDofIndex()][pvIdx] = priVar;
curVolVars.update(elemSol, this->problem(), element, scv);
return this->evalLocalResidual();
};
......@@ -289,8 +289,8 @@ public:
// derive the residuals numerically
static const NumericEpsilon<Scalar, numEq> eps_{GET_PROP_VALUE(TypeTag, ModelParameterGroup)};
static const int numDiffMethod = getParamFromGroup<int>(GET_PROP_VALUE(TypeTag, ModelParameterGroup), "Assembly.NumericDifferenceMethod");
NumericDifferentiation::partialDerivative(evalResiduals, elemSol[scv.indexInElement()][pvIdx], partialDerivs, origResiduals,
eps_(elemSol[scv.indexInElement()][pvIdx], pvIdx), numDiffMethod);
NumericDifferentiation::partialDerivative(evalResiduals, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origResiduals,
eps_(elemSol[scv.localDofIndex()][pvIdx], pvIdx), numDiffMethod);
// update the global stiffness matrix with the current partial derivatives
for (auto&& scvJ : scvs(fvGeometry))
......@@ -305,7 +305,7 @@ public:
// the residual of equation 'eqIdx' at dof 'i'
// depending on the primary variable 'pvIdx' at dof
// 'col'.
A[scvJ.dofIndex()][dofIdx][eqIdx][pvIdx] += partialDerivs[scvJ.indexInElement()][eqIdx];
A[scvJ.dofIndex()][dofIdx][eqIdx][pvIdx] += partialDerivs[scvJ.localDofIndex()][eqIdx];
}
}
}
......@@ -314,7 +314,7 @@ public:
curVolVars = origVolVars;
// restore the original element solution
elemSol[scv.indexInElement()][pvIdx] = curSol[scv.dofIndex()][pvIdx];
elemSol[scv.localDofIndex()][pvIdx] = curSol[scv.dofIndex()][pvIdx];
// TODO additional dof dependencies
}
}
......@@ -401,7 +401,7 @@ public:
auto evalStorage = [&](Scalar priVar)
{
// auto partialDerivsTmp = partialDerivs;
elemSol[scv.indexInElement()][pvIdx] = priVar;
elemSol[scv.localDofIndex()][pvIdx] = priVar;
curVolVars.update(elemSol, this->problem(), element, scv);
return this->evalLocalStorageResidual();
};
......@@ -409,8 +409,8 @@ public:
// derive the residuals numerically
static const NumericEpsilon<Scalar, numEq> eps_{GET_PROP_VALUE(TypeTag, ModelParameterGroup)};
static const int numDiffMethod = getParamFromGroup<int>(GET_PROP_VALUE(TypeTag, ModelParameterGroup), "Assembly.NumericDifferenceMethod");
NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.indexInElement()][pvIdx], partialDerivs, origResiduals,
eps_(elemSol[scv.indexInElement()][pvIdx], pvIdx), numDiffMethod);
NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origResiduals,
eps_(elemSol[scv.localDofIndex()][pvIdx], pvIdx), numDiffMethod);
// update the global stiffness matrix with the current partial derivatives
for (int eqIdx = 0; eqIdx < numEq; eqIdx++)
......@@ -419,14 +419,14 @@ public:
// the residual of equation 'eqIdx' at dof 'i'
// depending on the primary variable 'pvIdx' at dof
// 'col'.
A[dofIdx][dofIdx][eqIdx][pvIdx] += partialDerivs[scv.indexInElement()][eqIdx];
A[dofIdx][dofIdx][eqIdx][pvIdx] += partialDerivs[scv.localDofIndex()][eqIdx];
}
// restore the original state of the scv's volume variables
curVolVars = origVolVars;
// restore the original element solution
elemSol[scv.indexInElement()][pvIdx] = curSol[scv.dofIndex()][pvIdx];
elemSol[scv.localDofIndex()][pvIdx] = curSol[scv.dofIndex()][pvIdx];
// TODO additional dof dependencies
}
}
......@@ -534,7 +534,7 @@ public:
else
{
const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
if (this->elemBcTypes()[insideScv.indexInElement()].hasNeumann())
if (this->elemBcTypes()[insideScv.localDofIndex()].hasNeumann())
{
// add flux term derivatives
this->localResidual().addRobinFluxDerivatives(A[insideScv.dofIndex()],
......
......@@ -73,13 +73,13 @@ public:
{
const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
residual[insideScv.indexInElement()] += flux;
residual[outsideScv.indexInElement()] -= flux;
residual[insideScv.localDofIndex()] += flux;
residual[outsideScv.localDofIndex()] -= flux;
}
else
{
const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
residual[insideScv.indexInElement()] += flux;
residual[insideScv.localDofIndex()] += flux;
}
}
......@@ -104,7 +104,7 @@ public:
else
{
const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
const auto& bcTypes = elemBcTypes[scv.indexInElement()];
const auto& bcTypes = elemBcTypes[scv.localDofIndex()];
// Neumann and Robin ("solution dependent Neumann") boundary conditions
if (bcTypes.hasNeumann() && !bcTypes.hasDirichlet())
......
......@@ -64,7 +64,7 @@ public:
const SubControlVolumeFace& scvf) const
{
const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
const auto localScvIdx = scv.indexInElement();
const auto localScvIdx = scv.localDofIndex();
residual[localScvIdx] += this->asImp().evalFlux(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
}
......
......@@ -112,7 +112,7 @@ public:
// all sub control volumes
for (auto&& scv : scvs(fvGeometry))
{
auto localScvIdx = scv.indexInElement();
auto localScvIdx = scv.localDofIndex();
const auto& volVars = elemVolVars[scv];
storage[localScvIdx] = asImp().computeStorage(scv, volVars);
storage[localScvIdx] *= scv.volume() * volVars.extrusionFactor();
......@@ -382,7 +382,7 @@ public:
storage *= scv.volume();
storage /= timeLoop_->timeStepSize();
residual[scv.indexInElement()] += storage;
residual[scv.localDofIndex()] += storage;
}
/*!
......@@ -411,7 +411,7 @@ public:
source *= scv.volume()*curVolVars.extrusionFactor();
//! subtract source from local rate (sign convention in user interface)
residual[scv.indexInElement()] -= source;
residual[scv.localDofIndex()] -= source;
}
/*!
......
......@@ -407,7 +407,7 @@ public:
const SubControlVolume &scv) const
{
NumEqVector source(0);
auto scvIdx = scv.indexInElement();
auto scvIdx = scv.localDofIndex();
auto key = std::make_pair(fvGridGeometry_->elementMapper().index(element), scvIdx);
if (pointSourceMap_.count(key))
{
......
......@@ -303,7 +303,7 @@ public:
for (auto&& scv : scvs(fvGeometry))
{
if (intersectsPointGeometry(globalPos, scv.geometry()))
scvIndices.push_back(scv.indexInElement());
scvIndices.push_back(scv.localDofIndex());
}
// for all scvs that where tested positiv add the point sources
// to the element/scv to point source map
......
......@@ -98,10 +98,10 @@ public:
const auto& volVars = elemVolVars[scv];
if (enableGravity)
rho += volVars.density(phaseIdx)*shapeValues[scv.indexInElement()][0];
rho += volVars.density(phaseIdx)*shapeValues[scv.localDofIndex()][0];
// the global shape function gradient
gradP.axpy(volVars.pressure(phaseIdx), fluxVarCache.gradN(scv.indexInElement()));
gradP.axpy(volVars.pressure(phaseIdx), fluxVarCache.gradN(scv.localDofIndex()));
}
if (enableGravity)
......@@ -135,8 +135,8 @@ public:
std::vector<Scalar> ti(fvGeometry.numScv());
for (const auto& scv : scvs(fvGeometry))
ti[scv.indexInElement()] =
-1.0*scvf.area()*vtmv(scvf.unitOuterNormal(), K, fluxVarCache.gradN(scv.indexInElement()));
ti[scv.localDofIndex()] =
-1.0*scvf.area()*vtmv(scvf.unitOuterNormal(), K, fluxVarCache.gradN(scv.localDofIndex()));
return ti;
}
......
......@@ -62,7 +62,7 @@ public:
for (const auto& scv : scvs(fvGeometry))
{
int scvIdxLocal = scv.indexInElement();
int scvIdxLocal = scv.localDofIndex();
vertexBCTypes_[scvIdxLocal].reset();
if (fvGeometry.fvGridGeometry().dofOnBoundary(scv.dofIndex()))
......
......@@ -63,7 +63,7 @@ public:
const auto numVert = element.subEntities(GridView::dimension);
priVars_.resize(numVert);
for (const auto& scv : scvs(fvGeometry))
priVars_[scv.indexInElement()] = elemVolVars[scv].priVars();
priVars_[scv.localDofIndex()] = elemVolVars[scv].priVars();
}
//! extract the element solution from the solution vector using a mapper
......@@ -85,7 +85,7 @@ public:
const auto numVert = element.subEntities(GridView::dimension);
priVars_.resize(numVert);
for (const auto& scv : scvs(fvGeometry))
priVars_[scv.indexInElement()] = sol[scv.dofIndex()];
priVars_[scv.localDofIndex()] = sol[scv.dofIndex()];
}
//! bracket operator const access
......
......@@ -63,7 +63,7 @@ public:
template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, int> = 0>
const VolumeVariables& operator [](const SubControlVolume& scv) const
{ return gridVolVars().volVars(eIdx_, scv.indexInElement()); }
{ return gridVolVars().volVars(eIdx_, scv.localDofIndex()); }
// For compatibility reasons with the case of not storing the vol vars.
// function to be called before assembling an element, preparing the vol vars within the stencil
......@@ -133,7 +133,7 @@ public:
// resize volume variables to the required size
volumeVariables_.resize(fvGeometry.numScv());
for (auto&& scv : scvs(fvGeometry))
volumeVariables_[scv.indexInElement()].update(elemSol, gridVolVars().problem(), element, scv);
volumeVariables_[scv.localDofIndex()].update(elemSol, gridVolVars().problem(), element, scv);
}
const VolumeVariables& operator [](std::size_t scvIdx) const
......@@ -144,11 +144,11 @@ public:
template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, int> = 0>
const VolumeVariables& operator [](const SubControlVolume& scv) const
{ return volumeVariables_[scv.indexInElement()]; }
{ return volumeVariables_[scv.localDofIndex()]; }
template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, int> = 0>
VolumeVariables& operator [](const SubControlVolume& scv)
{ return volumeVariables_[scv.indexInElement()]; }
{ return volumeVariables_[scv.localDofIndex()]; }
//! The global volume variables object we are a restriction of
const GridVolumeVariables& gridVolVars() const
......
......@@ -93,7 +93,7 @@ public:
// density interpolation
Scalar rho(0.0);
for (auto&& scv : scvs(fvGeometry))
rho += elemVolVars[scv].molarDensity(phaseIdx)*shapeValues[scv.indexInElement()][0];
rho += elemVolVars[scv].molarDensity(phaseIdx)*shapeValues[scv.localDofIndex()][0];
for (int compIdx = 0; compIdx < numComponents; compIdx++)
{
......@@ -119,7 +119,7 @@ public:
// the mole/mass fraction gradient
GlobalPosition gradX(0.0);
for (auto&& scv : scvs(fvGeometry))
gradX.axpy(elemVolVars[scv].moleFraction(phaseIdx, compIdx), fluxVarsCache.gradN(scv.indexInElement()));
gradX.axpy(elemVolVars[scv].moleFraction(phaseIdx, compIdx), fluxVarsCache.gradN(scv.localDofIndex()));
// compute the diffusive flux
componentFlux[compIdx] = -1.0*rho*vtmv(scvf.unitOuterNormal(), D, gradX)*scvf.area();
......@@ -142,7 +142,7 @@ public:
Scalar rho(0.0);
const auto& shapeValues = fluxVarCache.shapeValues();
for (auto&& scv : scvs(fvGeometry))
rho += elemVolVars[scv].molarDensity(phaseIdx)*shapeValues[scv.indexInElement()][0];
rho += elemVolVars[scv].molarDensity(phaseIdx)*shapeValues[scv.localDofIndex()][0];
const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
......@@ -171,7 +171,7 @@ public:
ti[compIdx].resize(fvGeometry.numScv());
for (auto&& scv : scvs(fvGeometry))
ti[compIdx][scv.indexInElement()] = -rho*vtmv(scvf.unitOuterNormal(), D, fluxVarCache.gradN(scv.indexInElement()))*scvf.area();
ti[compIdx][scv.localDofIndex()] = -rho*vtmv(scvf.unitOuterNormal(), D, fluxVarCache.gradN(scv.localDofIndex()))*scvf.area();
}
return ti;
......
......@@ -96,7 +96,7 @@ public:
// compute the temperature gradient with the shape functions
GlobalPosition gradTemp(0.0);
for (auto&& scv : scvs(fvGeometry))
gradTemp.axpy(elemVolVars[scv].temperature(), fluxVarsCache.gradN(scv.indexInElement()));
gradTemp.axpy(elemVolVars[scv].temperature(), fluxVarsCache.gradN(scv.localDofIndex()));
// comute the heat conduction flux
return -1.0*vtmv(scvf.unitOuterNormal(), lambda, gradTemp)*scvf.area();
......
......@@ -122,9 +122,9 @@ public:
{
// compute the temperature gradient with the shape functions
if (phaseIdx < numEnergyEqFluid)
gradTemp.axpy(elemVolVars[scv].temperature(phaseIdx), fluxVarsCache.gradN(scv.indexInElement()));
gradTemp.axpy(elemVolVars[scv].temperature(phaseIdx), fluxVarsCache.gradN(scv.localDofIndex()));
else
gradTemp.axpy(elemVolVars[scv].temperatureSolid(), fluxVarsCache.gradN(scv.indexInElement()));
gradTemp.axpy(elemVolVars[scv].temperatureSolid(), fluxVarsCache.gradN(scv.localDofIndex()));
}
// comute the heat conduction flux
......
......@@ -88,17 +88,17 @@ public:
// update the volvars of the element
volumeVariables_[eIdx].resize(fvGeometry.numScv());
for (auto&& scv : scvs(fvGeometry))
volumeVariables_[eIdx][scv.indexInElement()].update(elemSol, problem(), element, scv);
volumeVariables_[eIdx][scv.localDofIndex()].update(elemSol, problem(), element, scv);
}
}
template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, int> = 0>
const VolumeVariables& volVars(const SubControlVolume& scv) const
{ return volumeVariables_[scv.elementIndex()][scv.indexInElement()]; }
{ return volumeVariables_[scv.elementIndex()][scv.localDofIndex()]; }
template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, int> = 0>
VolumeVariables& volVars(const SubControlVolume& scv)
{ return volumeVariables_[scv.elementIndex()][scv.indexInElement()]; }
{ return volumeVariables_[scv.elementIndex()][scv.localDofIndex()]; }
const VolumeVariables& volVars(const std::size_t eIdx, const std::size_t scvIdx) const
{ return volumeVariables_[eIdx][scvIdx]; }
......
......@@ -103,14 +103,14 @@ public:
const auto& volVars = elemVolVars[scv];
// density interpolation
rho += volVars.molarDensity(phaseIdx)*shapeValues[scv.indexInElement()][0];
rho += volVars.molarDensity(phaseIdx)*shapeValues[scv.localDofIndex()][0];
// the ansatz function gradient
jacInvT.mv(shapeJacobian[scv.indexInElement()][0], gradN[scv.indexInElement()]);
jacInvT.mv(shapeJacobian[scv.localDofIndex()][0], gradN[scv.localDofIndex()]);
//interpolate the mole fraction for the diffusion matrix
for (int compIdx = 0; compIdx < numComponents; compIdx++)
{
moleFrac[compIdx] += volVars.moleFraction(phaseIdx, compIdx)*shapeValues[scv.indexInElement()][0];
moleFrac[compIdx] += volVars.moleFraction(phaseIdx, compIdx)*shapeValues[scv.localDofIndex()][0];
}
}
......@@ -124,7 +124,7 @@ public:
const auto& volVars = elemVolVars[scv];
// the mole/mass fraction gradient
gradX.axpy(volVars.moleFraction(phaseIdx, compIdx), gradN[scv.indexInElement()]);
gradX.axpy(volVars.moleFraction(phaseIdx, compIdx), gradN[scv.localDofIndex()]);
}
normalX[compIdx] = gradX *scvf.unitOuterNormal();
......
......@@ -146,7 +146,7 @@ public:
}
//! The global index of this scv
LocalIndexType indexInElement() const
LocalIndexType localDofIndex() const
{
return scvIdx_;
}
......
......@@ -139,7 +139,7 @@ public:
}
//! The global index of this scv
LocalIndexType indexInElement() const
LocalIndexType localDofIndex() const
{
return 0;
}
......
......@@ -65,9 +65,9 @@ public:
}
//! The index of the dof this scv is embedded in (box)
LocalIndexType indexInElement() const
LocalIndexType localDofIndex() const
{
return asImp_().indexInElement();
return asImp_().localDofIndex();
}
// The position of the dof this scv is embedded in
......
......@@ -195,7 +195,7 @@ public:
/ curElemVolVars[scvf.insideScvIdx()].viscosity();
for (const auto& scv : scvs(fvGeometry))
{
auto d = up*ti[scv.indexInElement()];
auto d = up*ti[scv.localDofIndex()];
A[insideScv.dofIndex()][scv.dofIndex()][conti0EqIdx][pressureIdx] += d;
A[outsideScv.dofIndex()][scv.dofIndex()][conti0EqIdx][pressureIdx] -= d;
}
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment