diff --git a/dumux/implicit/1p/1pmodel.hh b/dumux/implicit/1p/1pmodel.hh
index 40b419d33f8f7507715b649a72fe68ef61965928..01603e78c22407c15dd68fed3c78d0b11caf29ed 100644
--- a/dumux/implicit/1p/1pmodel.hh
+++ b/dumux/implicit/1p/1pmodel.hh
@@ -60,6 +60,7 @@ class OnePBoxModel : public GET_PROP_TYPE(TypeTag, BaseModel)
enum { dim = GridView::dimension };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -98,11 +99,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox)
- globalIdx = this->vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*elemIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
@@ -113,23 +110,14 @@ public:
const SpatialParams &spatialParams = this->problem_().spatialParams();
(*p)[globalIdx] = volVars.pressure();
- (*K)[globalIdx]= spatialParams.intrinsicPermeability(*elemIt,
+ (*K)[globalIdx] = spatialParams.intrinsicPermeability(*elemIt,
fvGeometry,
scvIdx);
}
}
- if (isBox)
- {
- writer.attachVertexData(*p, "p");
- writer.attachVertexData(*K, "K");
- }
- else
- {
- writer.attachCellData(*p, "p");
- writer.attachCellData(*K, "K");
- }
-
+ writer.attachDofData(*p, "p", isBox);
+ writer.attachDofData(*K, "K", isBox);
writer.attachCellData(*rank, "process rank");
}
};
diff --git a/dumux/implicit/1p2c/1p2cmodel.hh b/dumux/implicit/1p2c/1p2cmodel.hh
index d565a1fa1a9fd1d29125bbb52be510b02bf11591..75d726463875f205a2a7a9fc33cdc529fd028c82 100644
--- a/dumux/implicit/1p2c/1p2cmodel.hh
+++ b/dumux/implicit/1p2c/1p2cmodel.hh
@@ -91,6 +91,7 @@ class OnePTwoCBoxModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef Dune::FieldVector<Scalar, dim> DimVector;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -175,11 +176,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox)
- globalIdx = this->vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*elemIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
@@ -245,7 +242,7 @@ public:
//use vertical faces (horizontal noraml vector) to calculate vx
//in case of heterogeneities it seams to be better to define intrinisc permeability elementwise
- if(xDir > yDir)//(fluxVars.face().normal[0] > 1e-10 || fluxVars.face().normal[0] < -1e-10)// (xDir > yDir)
+ if (xDir > yDir)//(fluxVars.face().normal[0] > 1e-10 || fluxVars.face().normal[0] < -1e-10)// (xDir > yDir)
{
// get darcy velocity
//calculate (v n) n/A
@@ -309,50 +306,28 @@ public:
}
}
- if (isBox)
- {
- writer.attachVertexData(pressure, "P");
- writer.attachVertexData(delp, "delp");
- if (velocityOutput)
- {
- writer.attachVertexData(velocityX, "Vx");
- writer.attachVertexData(velocityY, "Vy");
- if (dim > 2)
- writer.attachVertexData(velocityZ, "Vz");
- }
- char nameMoleFraction0[42], nameMoleFraction1[42];
- snprintf(nameMoleFraction0, 42, "x_%s", FluidSystem::componentName(0));
- snprintf(nameMoleFraction1, 42, "x_%s", FluidSystem::componentName(1));
- writer.attachVertexData(moleFraction0, nameMoleFraction0);
- writer.attachVertexData(moleFraction1, nameMoleFraction1);
-
- char nameMassFraction0[42], nameMassFraction1[42];
- snprintf(nameMassFraction0, 42, "X_%s", FluidSystem::componentName(0));
- snprintf(nameMassFraction1, 42, "X_%s", FluidSystem::componentName(1));
- writer.attachVertexData(massFraction0, nameMassFraction0);
- writer.attachVertexData(massFraction1, nameMassFraction1);
- writer.attachVertexData(rho, "rho");
- writer.attachVertexData(mu, "mu");
- }
- else
+ writer.attachDofData(pressure, "P", isBox);
+ writer.attachDofData(delp, "delp", isBox);
+ if (velocityOutput)
{
- writer.attachCellData(pressure, "P");
- writer.attachCellData(delp, "delp");
- char nameMoleFraction0[42], nameMoleFraction1[42];
- snprintf(nameMoleFraction0, 42, "x_%s", FluidSystem::componentName(0));
- snprintf(nameMoleFraction1, 42, "x_%s", FluidSystem::componentName(1));
- writer.attachCellData(moleFraction0, nameMoleFraction0);
- writer.attachCellData(moleFraction1, nameMoleFraction1);
-
- char nameMassFraction0[42], nameMassFraction1[42];
- snprintf(nameMassFraction0, 42, "X_%s", FluidSystem::componentName(0));
- snprintf(nameMassFraction1, 42, "X_%s", FluidSystem::componentName(1));
- writer.attachCellData(massFraction0, nameMassFraction0);
- writer.attachCellData(massFraction1, nameMassFraction1);
- writer.attachCellData(rho, "rho");
- writer.attachCellData(mu, "mu");
+ writer.attachDofData(velocityX, "Vx", isBox);
+ writer.attachDofData(velocityY, "Vy", isBox);
+ if (dim > 2)
+ writer.attachDofData(velocityZ, "Vz", isBox);
}
-
+ char nameMoleFraction0[42], nameMoleFraction1[42];
+ snprintf(nameMoleFraction0, 42, "x_%s", FluidSystem::componentName(0));
+ snprintf(nameMoleFraction1, 42, "x_%s", FluidSystem::componentName(1));
+ writer.attachDofData(moleFraction0, nameMoleFraction0, isBox);
+ writer.attachDofData(moleFraction1, nameMoleFraction1, isBox);
+
+ char nameMassFraction0[42], nameMassFraction1[42];
+ snprintf(nameMassFraction0, 42, "X_%s", FluidSystem::componentName(0));
+ snprintf(nameMassFraction1, 42, "X_%s", FluidSystem::componentName(1));
+ writer.attachDofData(massFraction0, nameMassFraction0, isBox);
+ writer.attachDofData(massFraction1, nameMassFraction1, isBox);
+ writer.attachDofData(rho, "rho", isBox);
+ writer.attachDofData(mu, "mu", isBox);
writer.attachCellData(rank, "process rank");
}
diff --git a/dumux/implicit/2p/2pmodel.hh b/dumux/implicit/2p/2pmodel.hh
index 40922141833b1ec0cc674b392341c20f404fe972..a045ae796e34d9244f3cf5327a5716c94509511a 100644
--- a/dumux/implicit/2p/2pmodel.hh
+++ b/dumux/implicit/2p/2pmodel.hh
@@ -99,6 +99,7 @@ class TwoPModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -140,7 +141,7 @@ public:
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
- if(velocityOutput) // check if velocity output is demanded
+ if (velocityOutput) // check if velocity output is demanded
{
// initialize velocity fields
for (unsigned int i = 0; i < numDofs; ++i)
@@ -161,7 +162,7 @@ public:
ElementIterator elemEndIt = this->gridView_().template end<0>();
for (; elemIt != elemEndIt; ++elemIt)
{
- if(velocityOutput && !elemIt->geometry().type().isCube()){
+ if (velocityOutput && !elemIt->geometry().type().isCube()){
DUNE_THROW(Dune::InvalidStateException,
"Currently, velocity output only works for cubes. "
"Please set the EnableVelocityOutput property to false!");
@@ -173,11 +174,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox)
- globalIdx = this->vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*elemIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
@@ -197,13 +194,13 @@ public:
(*mobN)[globalIdx] = volVars.mobility(nPhaseIdx);
(*poro)[globalIdx] = volVars.porosity();
(*Te)[globalIdx] = volVars.temperature();
- if(velocityOutput)
+ if (velocityOutput)
{
(*cellNum)[globalIdx] += 1;
}
}
- if(velocityOutput)
+ if (velocityOutput)
{
// calculate vertex velocities
GlobalPosition tmpVelocity[numPhases];
@@ -261,11 +258,11 @@ public:
tmpVelocity[phaseIdx] = localNormal;
tmpVelocity[phaseIdx] *= q[phaseIdx];
- if(phaseIdx == wPhaseIdx){
+ if (phaseIdx == wPhaseIdx){
scvVelocityW[fluxVars.face().i] += tmpVelocity[phaseIdx];
scvVelocityW[fluxVars.face().j] += tmpVelocity[phaseIdx];
}
- else if(phaseIdx == nPhaseIdx){
+ else if (phaseIdx == nPhaseIdx){
scvVelocityN[fluxVars.face().i] += tmpVelocity[phaseIdx];
scvVelocityN[fluxVars.face().j] += tmpVelocity[phaseIdx];
}
@@ -280,18 +277,18 @@ public:
= elemIt->geometry().jacobianTransposed(localPos);
// transform vertex velocities from local to global coordinates
- for (int i = 0; i < fvGeometry.numSCV; ++i)
+ for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx = this->vertexMapper().map(*elemIt, i, dim);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
// calculate the subcontrolvolume velocity by the Piola transformation
Dune::FieldVector<CoordScalar, dim> scvVelocity(0);
- jacobianT2.mtv(scvVelocityW[i], scvVelocity);
+ jacobianT2.mtv(scvVelocityW[scvIdx], scvVelocity);
scvVelocity /= elemIt->geometry().integrationElement(localPos);
// add up the wetting phase subcontrolvolume velocities for each vertex
(*velocityW)[globalIdx] += scvVelocity;
- jacobianT2.mtv(scvVelocityN[i], scvVelocity);
+ jacobianT2.mtv(scvVelocityN[scvIdx], scvVelocity);
scvVelocity /= elemIt->geometry().integrationElement(localPos);
// add up the nonwetting phase subcontrolvolume velocities for each vertex
(*velocityN)[globalIdx] += scvVelocity;
@@ -299,43 +296,26 @@ public:
}
}
- if (isBox) // vertex data
+ writer.attachDofData(*Sn, "Sn", isBox);
+ writer.attachDofData(*Sw, "Sw", isBox);
+ writer.attachDofData(*pN, "pn", isBox);
+ writer.attachDofData(*pW, "pw", isBox);
+ writer.attachDofData(*pC, "pc", isBox);
+ writer.attachDofData(*rhoW, "rhoW", isBox);
+ writer.attachDofData(*rhoN, "rhoN", isBox);
+ writer.attachDofData(*mobW, "mobW", isBox);
+ writer.attachDofData(*mobN, "mobN", isBox);
+ writer.attachDofData(*poro, "porosity", isBox);
+ writer.attachDofData(*Te, "temperature", isBox);
+ if (velocityOutput) // check if velocity output is demanded
{
- writer.attachVertexData(*Sn, "Sn");
- writer.attachVertexData(*Sw, "Sw");
- writer.attachVertexData(*pN, "pn");
- writer.attachVertexData(*pW, "pw");
- writer.attachVertexData(*pC, "pc");
- writer.attachVertexData(*rhoW, "rhoW");
- writer.attachVertexData(*rhoN, "rhoN");
- writer.attachVertexData(*mobW, "mobW");
- writer.attachVertexData(*mobN, "mobN");
- writer.attachVertexData(*poro, "porosity");
- writer.attachVertexData(*Te, "temperature");
- if(velocityOutput) // check if velocity output is demanded
- {
- // divide the vertex velocities by the number of adjacent scvs i.e. cells
- for(unsigned int globalIdx = 0; globalIdx < velocityW->size(); ++globalIdx){
- (*velocityW)[globalIdx] /= (*cellNum)[globalIdx];
- (*velocityN)[globalIdx] /= (*cellNum)[globalIdx];
- }
- writer.attachVertexData(*velocityW, "velocityW", dim);
- writer.attachVertexData(*velocityN, "velocityN", dim);
+ // divide the vertex velocities by the number of adjacent scvs i.e. cells
+ for(unsigned int globalIdx = 0; globalIdx < velocityW->size(); ++globalIdx){
+ (*velocityW)[globalIdx] /= (*cellNum)[globalIdx];
+ (*velocityN)[globalIdx] /= (*cellNum)[globalIdx];
}
- }
- else // cell data
- {
- writer.attachCellData(*Sn, "Sn");
- writer.attachCellData(*Sw, "Sw");
- writer.attachCellData(*pN, "pn");
- writer.attachCellData(*pW, "pw");
- writer.attachCellData(*pC, "pc");
- writer.attachCellData(*rhoW, "rhoW");
- writer.attachCellData(*rhoN, "rhoN");
- writer.attachCellData(*mobW, "mobW");
- writer.attachCellData(*mobN, "mobN");
- writer.attachCellData(*poro, "porosity");
- writer.attachCellData(*Te, "temperature");
+ writer.attachDofData(*velocityW, "velocityW", isBox, dim);
+ writer.attachDofData(*velocityN, "velocityN", isBox, dim);
}
writer.attachCellData(*rank, "process rank");
}
diff --git a/dumux/implicit/2p2c/2p2cmodel.hh b/dumux/implicit/2p2c/2p2cmodel.hh
index 19e9475954ab17c3fcfda0fcce3a7217ba7b6c02..d7c295649d7c243265c7ff2df0f98f0af1d31fe7 100644
--- a/dumux/implicit/2p2c/2p2cmodel.hh
+++ b/dumux/implicit/2p2c/2p2cmodel.hh
@@ -135,6 +135,7 @@ class TwoPTwoCModel: public GET_PROP_TYPE(TypeTag, BaseModel)
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -154,11 +155,11 @@ public:
// check, if velocity output can be used (works only for cubes so far)
velocityOutput_ = GET_PARAM_FROM_GROUP(TypeTag, bool, Vtk, AddVelocity);
- ElementIterator eIt = this->gridView_().template begin<0>();
- ElementIterator eEndIt = this->gridView_().template end<0>();
- for (; eIt != eEndIt; ++eIt)
+ ElementIterator elemIt = this->gridView_().template begin<0>();
+ ElementIterator elemEndIt = this->gridView_().template end<0>();
+ for (; elemIt != elemEndIt; ++elemIt)
{
- if(eIt->geometry().type().isCube() == false){
+ if (elemIt->geometry().type().isCube() == false){
velocityOutput_ = false;
}
@@ -166,8 +167,8 @@ public:
{
velocityOutput_ = false;
- int globalIdx = this->dofMapper().map(*eIt);
- const GlobalPosition &globalPos = eIt->geometry().center();
+ int globalIdx = this->dofMapper().map(*elemIt);
+ const GlobalPosition &globalPos = elemIt->geometry().center();
// initialize phase presence
staticDat_[globalIdx].phasePresence
@@ -215,10 +216,10 @@ public:
{
storage = 0;
- ElementIterator eIt = this->gridView_().template begin<0>();
- const ElementIterator eEndIt = this->gridView_().template end<0>();
- for (; eIt != eEndIt; ++eIt) {
- this->localResidual().evalPhaseStorage(*eIt, phaseIdx);
+ ElementIterator elemIt = this->gridView_().template begin<0>();
+ const ElementIterator elemEndIt = this->gridView_().template end<0>();
+ for (; elemIt != elemEndIt; ++elemIt) {
+ this->localResidual().evalPhaseStorage(*elemIt, phaseIdx);
for (unsigned int i = 0; i < this->localResidual().storageTerm().size(); ++i)
storage += this->localResidual().storageTerm()[i];
@@ -320,7 +321,7 @@ public:
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
- if(velocityOutput_) // check if velocity output is demanded
+ if (velocityOutput_) // check if velocity output is demanded
{
// initialize velocity fields
for (unsigned int i = 0; i < numDofs; ++i)
@@ -337,25 +338,21 @@ public:
FVElementGeometry fvGeometry;
VolumeVariables volVars;
- ElementIterator eIt = this->gridView_().template begin<0>();
- ElementIterator eEndIt = this->gridView_().template end<0>();
- for (; eIt != eEndIt; ++eIt)
+ ElementIterator elemIt = this->gridView_().template begin<0>();
+ ElementIterator elemEndIt = this->gridView_().template end<0>();
+ for (; elemIt != elemEndIt; ++elemIt)
{
- int idx = this->elementMapper().map(*eIt);
+ int idx = this->elementMapper().map(*elemIt);
(*rank)[idx] = this->gridView_().comm().rank();
- fvGeometry.update(this->gridView_(), *eIt);
+ fvGeometry.update(this->gridView_(), *elemIt);
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (numDofs == numElements) // element data
- globalIdx = idx;
- else
- globalIdx = this->vertexMapper().map(*eIt, scvIdx, dim);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
scvIdx,
false);
@@ -380,13 +377,13 @@ public:
(*temperature)[globalIdx] = volVars.temperature();
(*phasePresence)[globalIdx]
= staticDat_[globalIdx].phasePresence;
- if(velocityOutput_)
+ if (velocityOutput_)
{
(*cellNum)[globalIdx] += 1;
}
}
- if(velocityOutput_)
+ if (velocityOutput_)
{
// calculate vertex velocities
GlobalPosition tmpVelocity[numPhases];
@@ -405,7 +402,7 @@ public:
ElementVolumeVariables elemVolVars;
elemVolVars.update(this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
false /* oldSol? */);
@@ -413,7 +410,7 @@ public:
{
FluxVariables fluxVars(this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
faceIdx,
elemVolVars);
@@ -425,7 +422,7 @@ public:
// Transformation of the global normal vector to normal vector in the reference element
const typename Element::Geometry::JacobianTransposed jacobianT1 =
- eIt->geometry().jacobianTransposed(localPosIP);
+ elemIt->geometry().jacobianTransposed(localPosIP);
const GlobalPosition globalNormal = fluxVars.face().normal;
GlobalPosition localNormal(0);
@@ -444,11 +441,11 @@ public:
tmpVelocity[phaseIdx] = localNormal;
tmpVelocity[phaseIdx] *= flux[phaseIdx];
- if(phaseIdx == wPhaseIdx){
+ if (phaseIdx == wPhaseIdx){
scvVelocityW[fluxVars.face().i] += tmpVelocity[phaseIdx];
scvVelocityW[fluxVars.face().j] += tmpVelocity[phaseIdx];
}
- else if(phaseIdx == nPhaseIdx){
+ else if (phaseIdx == nPhaseIdx){
scvVelocityN[fluxVars.face().i] += tmpVelocity[phaseIdx];
scvVelocityN[fluxVars.face().j] += tmpVelocity[phaseIdx];
}
@@ -457,33 +454,33 @@ public:
typedef Dune::GenericReferenceElements<Scalar, dim> ReferenceElements;
const Dune::FieldVector<Scalar, dim>& localPos =
- ReferenceElements::general(eIt->geometry().type()).position(0, 0);
+ ReferenceElements::general(elemIt->geometry().type()).position(0, 0);
// get the transposed Jacobian of the element mapping
const typename Element::Geometry::JacobianTransposed jacobianT2 =
- eIt->geometry().jacobianTransposed(localPos);
+ elemIt->geometry().jacobianTransposed(localPos);
// transform vertex velocities from local to global coordinates
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx = this->vertexMapper().map(*eIt, scvIdx, dim);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
// calculate the subcontrolvolume velocity by the Piola transformation
Dune::FieldVector<CoordScalar, dim> scvVelocity(0);
jacobianT2.mtv(scvVelocityW[scvIdx], scvVelocity);
- scvVelocity /= eIt->geometry().integrationElement(localPos);
+ scvVelocity /= elemIt->geometry().integrationElement(localPos);
// add up the wetting phase subcontrolvolume velocities for each vertex
(*velocityW)[globalIdx] += scvVelocity;
jacobianT2.mtv(scvVelocityN[scvIdx], scvVelocity);
- scvVelocity /= eIt->geometry().integrationElement(localPos);
+ scvVelocity /= elemIt->geometry().integrationElement(localPos);
// add up the nonwetting phase subcontrolvolume velocities for each vertex
(*velocityN)[globalIdx] += scvVelocity;
}
} // velocity output
} // loop over elements
- if(velocityOutput_)
+ if (velocityOutput_)
{
// divide the vertex velocities by the number of adjacent scvs i.e. cells
for(unsigned int globalIdx = 0; globalIdx < numDofs; ++globalIdx){
@@ -492,59 +489,32 @@ public:
}
}
- if (isBox) // vertex data
+ writer.attachDofData(*sN, "Sn", isBox);
+ writer.attachDofData(*sW, "Sw", isBox);
+ writer.attachDofData(*pN, "pN", isBox);
+ writer.attachDofData(*pW, "pW", isBox);
+ writer.attachDofData(*pC, "pC", isBox);
+ writer.attachDofData(*rhoW, "rhoW", isBox);
+ writer.attachDofData(*rhoN, "rhoN", isBox);
+ writer.attachDofData(*mobW, "mobW", isBox);
+ writer.attachDofData(*mobN, "mobN", isBox);
+ for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
{
- writer.attachVertexData(*sN, "Sn");
- writer.attachVertexData(*sW, "Sw");
- writer.attachVertexData(*pN, "pN");
- writer.attachVertexData(*pW, "pW");
- writer.attachVertexData(*pC, "pC");
- writer.attachVertexData(*rhoW, "rhoW");
- writer.attachVertexData(*rhoN, "rhoN");
- writer.attachVertexData(*mobW, "mobW");
- writer.attachVertexData(*mobN, "mobN");
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- std::ostringstream oss;
- oss << "X_" << FluidSystem::phaseName(phaseIdx) << "^" << FluidSystem::componentName(compIdx);
- writer.attachVertexData(*massFrac[phaseIdx][compIdx], oss.str());
- }
- }
- writer.attachVertexData(*poro, "porosity");
- writer.attachVertexData(*temperature, "temperature");
- writer.attachVertexData(*phasePresence, "phase presence");
-
- if(velocityOutput_) // check if velocity output is demanded
+ for (int compIdx = 0; compIdx < numComponents; ++compIdx)
{
- writer.attachVertexData(*velocityW, "velocityW", dim);
- writer.attachVertexData(*velocityN, "velocityN", dim);
+ std::ostringstream oss;
+ oss << "X_" << FluidSystem::phaseName(phaseIdx) << "^" << FluidSystem::componentName(compIdx);
+ writer.attachDofData(*massFrac[phaseIdx][compIdx], oss.str(), isBox);
}
}
- else // cell data
+ writer.attachDofData(*poro, "porosity", isBox);
+ writer.attachDofData(*temperature, "temperature", isBox);
+ writer.attachDofData(*phasePresence, "phase presence", isBox);
+
+ if (velocityOutput_) // check if velocity output is demanded
{
- writer.attachCellData(*sN, "Sn");
- writer.attachCellData(*sW, "Sw");
- writer.attachCellData(*pN, "pN");
- writer.attachCellData(*pW, "pW");
- writer.attachCellData(*pC, "pC");
- writer.attachCellData(*rhoW, "rhoW");
- writer.attachCellData(*rhoN, "rhoN");
- writer.attachCellData(*mobW, "mobW");
- writer.attachCellData(*mobN, "mobN");
- for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
- {
- for (int compIdx = 0; compIdx < numComponents; ++compIdx)
- {
- std::ostringstream oss;
- oss << "X_" << FluidSystem::phaseName(phaseIdx) << "^" << FluidSystem::componentName(compIdx);
- writer.attachCellData(*massFrac[phaseIdx][compIdx], oss.str());
- }
- }
- writer.attachCellData(*poro, "porosity");
- writer.attachCellData(*temperature, "temperature");
- writer.attachCellData(*phasePresence, "phase presence");
+ writer.attachDofData(*velocityW, "velocityW", isBox, dim);
+ writer.attachDofData(*velocityN, "velocityN", isBox, dim);
}
writer.attachCellData(*rank, "process rank");
@@ -609,19 +579,14 @@ public:
FVElementGeometry fvGeometry;
static VolumeVariables volVars;
- ElementIterator eIt = this->gridView_().template begin<0> ();
- const ElementIterator &eEndIt = this->gridView_().template end<0> ();
- for (; eIt != eEndIt; ++eIt)
+ ElementIterator elemIt = this->gridView_().template begin<0> ();
+ const ElementIterator &elemEndIt = this->gridView_().template end<0> ();
+ for (; elemIt != elemEndIt; ++elemIt)
{
- fvGeometry.update(this->gridView_(), *eIt);
+ fvGeometry.update(this->gridView_(), *elemIt);
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
-
- if (isBox)
- globalIdx = this->vertexMapper().map(*eIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*eIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
if (staticDat_[globalIdx].visited)
continue;
@@ -629,7 +594,7 @@ public:
staticDat_[globalIdx].visited = true;
volVars.update(curGlobalSol[globalIdx],
this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
scvIdx,
false);
diff --git a/dumux/implicit/2p2c/2p2cvolumevariables.hh b/dumux/implicit/2p2c/2p2cvolumevariables.hh
index 9591cae3f738fc3650f0394aa321e0e5278667c0..cc1d8bfd661556c405a470f189b599eca5bde395 100644
--- a/dumux/implicit/2p2c/2p2cvolumevariables.hh
+++ b/dumux/implicit/2p2c/2p2cvolumevariables.hh
@@ -102,6 +102,7 @@ class TwoPTwoCVolumeVariables : public ImplicitVolumeVariables<TypeTag>
typedef Dumux::ComputeFromReferencePhase<Scalar, FluidSystem> ComputeFromReferencePhase;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
//! The type of the object returned by the fluidState() method
@@ -184,14 +185,7 @@ public:
fvGeometry, scvIdx);
fluidState.setTemperature(t);
- unsigned numDofs = problem.model().numDofs();
- unsigned numVertices = problem.gridView().size(dim);
-
- int globalIdx;
- if (isBox) // vertex data
- globalIdx = problem.model().dofMapper().map(element, scvIdx, dim);
- else
- globalIdx = problem.model().dofMapper().map(element);
+ int globalIdx = problem.model().dofMapper().map(element, scvIdx, dofCodim);
int phasePresence = problem.model().phasePresence(globalIdx, isOldSol);
diff --git a/dumux/implicit/3p3c/3p3cmodel.hh b/dumux/implicit/3p3c/3p3cmodel.hh
index 9d8d48cd64ae30b8ec8b88f2b1c8e9ea29fada30..b973451b1c9e6e80ebb1a2988586dccc102fe5d3 100644
--- a/dumux/implicit/3p3c/3p3cmodel.hh
+++ b/dumux/implicit/3p3c/3p3cmodel.hh
@@ -140,6 +140,7 @@ class ThreePThreeCModel: public GET_PROP_TYPE(TypeTag, BaseModel)
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -206,10 +207,10 @@ public:
{
storage = 0;
- ElementIterator eIt = this->gridView_().template begin<0>();
- const ElementIterator eEndIt = this->gridView_().template end<0>();
- for (; eIt != eEndIt; ++eIt) {
- this->localResidual().evalPhaseStorage(*eIt, phaseIdx);
+ ElementIterator elemIt = this->gridView_().template begin<0>();
+ const ElementIterator elemEndIt = this->gridView_().template end<0>();
+ for (; elemIt != elemEndIt; ++elemIt) {
+ this->localResidual().evalPhaseStorage(*elemIt, phaseIdx);
for (unsigned int i = 0; i < this->localResidual().storageTerm().size(); ++i)
storage += this->localResidual().storageTerm()[i];
@@ -325,11 +326,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox) // vertex data
- globalIdx = this->vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = idx;
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
@@ -362,64 +359,32 @@ public:
}
- if (isBox) // vertex data
- {
- writer.attachVertexData(*saturation[wPhaseIdx], "Sw");
- writer.attachVertexData(*saturation[nPhaseIdx], "Sn");
- writer.attachVertexData(*saturation[gPhaseIdx], "Sg");
- writer.attachVertexData(*pressure[wPhaseIdx], "pw");
- writer.attachVertexData(*pressure[nPhaseIdx], "pn");
- writer.attachVertexData(*pressure[gPhaseIdx], "pg");
- writer.attachVertexData(*density[wPhaseIdx], "rhow");
- writer.attachVertexData(*density[nPhaseIdx], "rhon");
- writer.attachVertexData(*density[gPhaseIdx], "rhog");
-
- for (int i = 0; i < numPhases; ++i)
- {
- for (int j = 0; j < numComponents; ++j)
- {
- std::ostringstream oss;
- oss << "x^"
- << FluidSystem::phaseName(i)
- << "_"
- << FluidSystem::componentName(j);
- writer.attachVertexData(*moleFraction[i][j], oss.str().c_str());
- }
- }
- writer.attachVertexData(*poro, "porosity");
- writer.attachVertexData(*perm, "permeability");
- writer.attachVertexData(*temperature, "temperature");
- writer.attachVertexData(*phasePresence, "phase presence");
- }
- else // cell data
+ writer.attachDofData(*saturation[wPhaseIdx], "Sw", isBox);
+ writer.attachDofData(*saturation[nPhaseIdx], "Sn", isBox);
+ writer.attachDofData(*saturation[gPhaseIdx], "Sg", isBox);
+ writer.attachDofData(*pressure[wPhaseIdx], "pw", isBox);
+ writer.attachDofData(*pressure[nPhaseIdx], "pn", isBox);
+ writer.attachDofData(*pressure[gPhaseIdx], "pg", isBox);
+ writer.attachDofData(*density[wPhaseIdx], "rhow", isBox);
+ writer.attachDofData(*density[nPhaseIdx], "rhon", isBox);
+ writer.attachDofData(*density[gPhaseIdx], "rhog", isBox);
+
+ for (int i = 0; i < numPhases; ++i)
{
- writer.attachCellData(*saturation[wPhaseIdx], "Sw");
- writer.attachCellData(*saturation[nPhaseIdx], "Sn");
- writer.attachCellData(*saturation[gPhaseIdx], "Sg");
- writer.attachCellData(*pressure[wPhaseIdx], "pw");
- writer.attachCellData(*pressure[nPhaseIdx], "pn");
- writer.attachCellData(*pressure[gPhaseIdx], "pg");
- writer.attachCellData(*density[wPhaseIdx], "rhow");
- writer.attachCellData(*density[nPhaseIdx], "rhon");
- writer.attachCellData(*density[gPhaseIdx], "rhog");
-
- for (int i = 0; i < numPhases; ++i)
- {
- for (int j = 0; j < numComponents; ++j)
- {
- std::ostringstream oss;
- oss << "x^"
- << FluidSystem::phaseName(i)
- << "_"
- << FluidSystem::componentName(j);
- writer.attachCellData(*moleFraction[i][j], oss.str().c_str());
- }
+ for (int j = 0; j < numComponents; ++j)
+ {
+ std::ostringstream oss;
+ oss << "x^"
+ << FluidSystem::phaseName(i)
+ << "_"
+ << FluidSystem::componentName(j);
+ writer.attachDofData(*moleFraction[i][j], oss.str().c_str(), isBox);
}
- writer.attachCellData(*poro, "porosity");
- writer.attachCellData(*perm, "permeability");
- writer.attachCellData(*temperature, "temperature");
- writer.attachCellData(*phasePresence, "phase presence");
}
+ writer.attachDofData(*poro, "porosity", isBox);
+ writer.attachDofData(*perm, "permeability", isBox);
+ writer.attachDofData(*temperature, "temperature", isBox);
+ writer.attachDofData(*phasePresence, "phase presence", isBox);
writer.attachCellData(*rank, "process rank");
}
@@ -483,19 +448,14 @@ public:
FVElementGeometry fvGeometry;
static VolumeVariables volVars;
- ElementIterator eIt = this->gridView_().template begin<0> ();
- const ElementIterator &eEndIt = this->gridView_().template end<0> ();
- for (; eIt != eEndIt; ++eIt)
+ ElementIterator elemIt = this->gridView_().template begin<0> ();
+ const ElementIterator &elemEndIt = this->gridView_().template end<0> ();
+ for (; elemIt != elemEndIt; ++elemIt)
{
- fvGeometry.update(this->gridView_(), *eIt);
+ fvGeometry.update(this->gridView_(), *elemIt);
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
-
- if (isBox)
- globalIdx = this->vertexMapper().map(*eIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*eIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
if (staticDat_[globalIdx].visited)
continue;
@@ -503,7 +463,7 @@ public:
staticDat_[globalIdx].visited = true;
volVars.update(curGlobalSol[globalIdx],
this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
scvIdx,
false);
diff --git a/dumux/implicit/3p3c/3p3cvolumevariables.hh b/dumux/implicit/3p3c/3p3cvolumevariables.hh
index edbbf9ebc81228fc59e27e11112458c6fb53a570..52f0a55448a0931786d8f3903a2f988e64185aed 100644
--- a/dumux/implicit/3p3c/3p3cvolumevariables.hh
+++ b/dumux/implicit/3p3c/3p3cvolumevariables.hh
@@ -101,6 +101,7 @@ class ThreePThreeCVolumeVariables : public ImplicitVolumeVariables<TypeTag>
static const Scalar R; // universial gas constant
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
//! The type of the object returned by the fluidState() method
@@ -130,11 +131,7 @@ public:
const MaterialLawParams &materialParams =
problem.spatialParams().materialLawParams(element, fvGeometry, scvIdx);
- int globalIdx;
- if (isBox) // vertex data
- globalIdx = problem.model().dofMapper().map(element, scvIdx, dim);
- else
- globalIdx = problem.model().dofMapper().map(element);
+ int globalIdx = problem.model().dofMapper().map(element, scvIdx, dofCodim);
int phasePresence = problem.model().phasePresence(globalIdx, isOldSol);
diff --git a/dumux/implicit/co2/co2model.hh b/dumux/implicit/co2/co2model.hh
index c5f098d011eaaea27d9eaf9ad9ed5f4b255d5e8d..31b2fa649f897ac9de24ae12d44a2b3df779b616 100644
--- a/dumux/implicit/co2/co2model.hh
+++ b/dumux/implicit/co2/co2model.hh
@@ -100,8 +100,10 @@ class CO2Model: public TwoPTwoCModel<TypeTag>
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef Dune::FieldVector<Scalar, numPhases> PhasesVector;
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
+ enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
- public:
+public:
/*!
@@ -118,24 +120,16 @@ class CO2Model: public TwoPTwoCModel<TypeTag>
for (unsigned i = 0; i < ParentType::staticDat_.size(); ++i)
ParentType::staticDat_[i].visited = false;
- unsigned numDofs = this->numDofs();
- unsigned numVertices = this->problem_().gridView().size(dim);
-
FVElementGeometry fvGeometry;
static VolumeVariables volVars;
- ElementIterator eIt = this->gridView_().template begin<0> ();
- const ElementIterator &eEndIt = this->gridView_().template end<0> ();
- for (; eIt != eEndIt; ++eIt)
+ ElementIterator elemIt = this->gridView_().template begin<0> ();
+ const ElementIterator &elemEndIt = this->gridView_().template end<0> ();
+ for (; elemIt != elemEndIt; ++elemIt)
{
- fvGeometry.update(this->gridView_(), *eIt);
+ fvGeometry.update(this->gridView_(), *elemIt);
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
-
- if (numDofs != numVertices)
- globalIdx = this->elementMapper().map(*eIt);
- else
- globalIdx = this->vertexMapper().map(*eIt, scvIdx, dim);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
if (ParentType::staticDat_[globalIdx].visited)
continue;
@@ -143,17 +137,17 @@ class CO2Model: public TwoPTwoCModel<TypeTag>
ParentType::staticDat_[globalIdx].visited = true;
volVars.update(curGlobalSol[globalIdx],
this->problem_(),
- *eIt,
+ *elemIt,
fvGeometry,
scvIdx,
false);
- const GlobalPosition &globalPos = eIt->geometry().corner(scvIdx);
+ const GlobalPosition &globalPos = elemIt->geometry().corner(scvIdx);
if (primaryVarSwitch_(curGlobalSol,
volVars,
globalIdx,
globalPos))
{
- this->jacobianAssembler().markVertexRed(globalIdx);
+ this->jacobianAssembler().markDofRed(globalIdx);
wasSwitched = true;
}
}
diff --git a/dumux/implicit/co2/co2volumevariables.hh b/dumux/implicit/co2/co2volumevariables.hh
index 3d346443ec22f302830551b3f7d7ea366ba606c9..e494aa65108e20465aa860cddb8ec536f8f1a496 100644
--- a/dumux/implicit/co2/co2volumevariables.hh
+++ b/dumux/implicit/co2/co2volumevariables.hh
@@ -87,6 +87,8 @@ class CO2VolumeVariables: public TwoPTwoCVolumeVariables<TypeTag>
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GridView::template Codim<0>::Entity Element;
enum { dim = GridView::dimension};
+ enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
static const Scalar R; // universial nonwetting constant
@@ -112,7 +114,7 @@ public:
const int scvIdx,
const bool isOldSol)
{
- // Update BoxVolVars but not 2p2cvolvars
+ // Update BoxVolVars but not 2p2cvolvars
// ToDo: Is BaseClassType the right name?
BaseClassType::update(priVars,
problem,
@@ -121,14 +123,7 @@ public:
scvIdx,
isOldSol);
- unsigned numDofs = problem.model().numDofs();
- unsigned numVertices = problem.gridView().size(dim);
-
- int globalIdx;
- if (numDofs != numVertices) // element data
- globalIdx = problem.model().dofMapper().map(element);
- else
- globalIdx = problem.model().dofMapper().map(element, scvIdx, dim);
+ int globalIdx = problem.model().dofMapper().map(element, scvIdx, dofCodim);
int phasePresence = problem.model().phasePresence(globalIdx, isOldSol);
diff --git a/dumux/implicit/common/implicitmodel.hh b/dumux/implicit/common/implicitmodel.hh
index ceafafa3d16a2bb18e0324a060d9c2fb14c920c1..1b45cb14e940b3b6f87a0a2485bc4c74c25770e5 100644
--- a/dumux/implicit/common/implicitmodel.hh
+++ b/dumux/implicit/common/implicitmodel.hh
@@ -78,6 +78,7 @@ class ImplicitModel
typedef typename Dune::GenericReferenceElement<CoordScalar, dim> ReferenceElement;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
// copying a model is not a good idea
ImplicitModel(const ImplicitModel &);
@@ -884,18 +885,8 @@ protected:
// loop over all element vertices, i.e. sub control volumes
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; scvIdx++)
{
- int globalIdx;
-
- if (isBox)
- {
- // map the local vertex index to the global one
- globalIdx = vertexMapper().map(*eIt, scvIdx, dim);
- }
- else
- {
- // get the global index of the element
- globalIdx = elementMapper().map(*eIt);
- }
+ // get the global index of the degree of freedom
+ int globalIdx = dofMapper().map(*eIt, scvIdx, dofCodim);
// let the problem do the dirty work of nailing down
// the initial solution.
diff --git a/dumux/implicit/mpnc/energy/mpncvtkwriterenergy.hh b/dumux/implicit/mpnc/energy/mpncvtkwriterenergy.hh
index 79be576f16f333a08e9738ad8eba607da512955c..ab71e760174b7a5c5be9362c1d625d57d89e10c1 100644
--- a/dumux/implicit/mpnc/energy/mpncvtkwriterenergy.hh
+++ b/dumux/implicit/mpnc/energy/mpncvtkwriterenergy.hh
@@ -62,6 +62,7 @@ class MPNCVtkWriterEnergy : public MPNCVtkWriterModule<TypeTag>
typedef typename ParentType::ScalarVector ScalarVector;
typedef typename ParentType::PhaseVector PhaseVector;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
MPNCVtkWriterEnergy(const Problem &problem)
@@ -84,15 +85,14 @@ public:
* \brief Modify the internal buffers according to the volume
* variables seen on an element
*/
- void processElement(const Element &elem,
+ void processElement(const Element &element,
const FVElementGeometry &fvGeometry,
const ElementVolumeVariables &elemVolVars,
const ElementBoundaryTypes &elemBcTypes)
{
- int numLocalVertices = elem.geometry().corners();
- for (int localVertexIdx = 0; localVertexIdx < numLocalVertices; ++localVertexIdx) {
- const unsigned int globalIdx = this->problem_.vertexMapper().map(elem, localVertexIdx, dim);
- const VolumeVariables &volVars = elemVolVars[localVertexIdx];
+ for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx) {
+ const unsigned int globalIdx = this->problem_.model().dofMapper().map(element, scvIdx, dofCodim);
+ const VolumeVariables &volVars = elemVolVars[scvIdx];
if (temperatureOutput_)
temperature_[globalIdx] = volVars.fluidState().temperature(/*phaseIdx=*/0);
@@ -144,6 +144,7 @@ class MPNCVtkWriterEnergy<TypeTag, /* enableEnergy = */ true, /* enableKineticEn
typedef typename ParentType::ScalarVector ScalarVector;
typedef typename ParentType::PhaseVector PhaseVector;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
MPNCVtkWriterEnergy(const Problem &problem)
@@ -170,15 +171,14 @@ public:
* \brief Modify the internal buffers according to the volume
* variables seen on an element
*/
- void processElement(const Element &elem,
+ void processElement(const Element &element,
const FVElementGeometry &fvGeometry,
const ElementVolumeVariables &elemVolVars,
const ElementBoundaryTypes &elemBcTypes)
{
- const unsigned int numLocalVertices = elem.geometry().corners();
- for (int localVertexIdx = 0; localVertexIdx < numLocalVertices; ++localVertexIdx) {
- int gobalIdx = this->problem_.vertexMapper().map(elem, localVertexIdx, dim);
- const VolumeVariables &volVars = elemVolVars[localVertexIdx];
+ for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx) {
+ int gobalIdx = this->problem_.model().dofMapper().map(element, scvIdx, dofCodim);
+ const VolumeVariables &volVars = elemVolVars[scvIdx];
if (temperatureOutput_) temperature_[gobalIdx] = volVars.fluidState().temperature(/*phaseIdx=*/0);
for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
diff --git a/dumux/implicit/mpnc/mass/mpncvtkwritermass.hh b/dumux/implicit/mpnc/mass/mpncvtkwritermass.hh
index 3968f741df3804564d8c61bea6767292baefc308..e20077b6d30a32eca085eebbac38d70b0a9daa02 100644
--- a/dumux/implicit/mpnc/mass/mpncvtkwritermass.hh
+++ b/dumux/implicit/mpnc/mass/mpncvtkwritermass.hh
@@ -57,6 +57,7 @@ class MPNCVtkWriterMass : public MPNCVtkWriterModule<TypeTag>
enum { dim = GridView::dimension };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
typedef typename ParentType::ComponentVector ComponentVector;
bool fugacityOutput_;
@@ -82,15 +83,14 @@ public:
* \brief Modify the internal buffers according to the volume
* variables seen on an element
*/
- void processElement(const Element &elem,
- const FVElementGeometry &fvElemGeom,
+ void processElement(const Element &element,
+ const FVElementGeometry &fvGeometry,
const ElementVolumeVariables &elemVolVars,
const ElementBoundaryTypes &elemBcTypes)
{
- int numLocalVertices = elem.geometry().corners();
- for (int localVertexIdx = 0; localVertexIdx < numLocalVertices; ++localVertexIdx) {
- int globalIdx = this->problem_.vertexMapper().map(elem, localVertexIdx, dim);
- const VolumeVariables &volVars = elemVolVars[localVertexIdx];
+ for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx) {
+ const unsigned int globalIdx = this->problem_.model().dofMapper().map(element, scvIdx, dofCodim);
+ const VolumeVariables &volVars = elemVolVars[scvIdx];
if (fugacityOutput_) {
for (int compIdx = 0; compIdx < numComponents; ++compIdx) {
diff --git a/dumux/implicit/mpnc/mpncvtkwritercommon.hh b/dumux/implicit/mpnc/mpncvtkwritercommon.hh
index c66de24e9588afa8c997aae04b22887730238663..dcb98d07fed7c1e32e8ad3c93c4f61d8ff6cddf6 100644
--- a/dumux/implicit/mpnc/mpncvtkwritercommon.hh
+++ b/dumux/implicit/mpnc/mpncvtkwritercommon.hh
@@ -63,6 +63,7 @@ class MPNCVtkWriterCommon : public MPNCVtkWriterModule<TypeTag>
typedef Dune::BlockVector<DimVector> DimField;
typedef Dune::array<DimField, numPhases> PhaseDimField;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
MPNCVtkWriterCommon(const Problem &problem)
@@ -130,11 +131,7 @@ public:
{
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox) // vertex data
- globalIdx = this->problem_.vertexMapper().map(element, scvIdx, dim);
- else
- globalIdx = this->problem_.elementMapper().map(element);
+ int globalIdx = this->problem_.model().dofMapper().map(element, scvIdx, dofCodim);
const VolumeVariables &volVars = elemVolVars[scvIdx];
diff --git a/dumux/implicit/richards/richardsmodel.hh b/dumux/implicit/richards/richardsmodel.hh
index 8ed12767a9132b9f9e71c1f86bcc0519220557b2..7551403fc735da888b7106bf83d97c6b72a5bbf6 100644
--- a/dumux/implicit/richards/richardsmodel.hh
+++ b/dumux/implicit/richards/richardsmodel.hh
@@ -120,6 +120,7 @@ class RichardsModel : public GET_PROP_TYPE(TypeTag, BaseModel)
enum { dim = GridView::dimension };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -168,11 +169,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox)
- globalIdx = this->vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = this->elementMapper().map(*elemIt);
+ int globalIdx = this->dofMapper().map(*elemIt, scvIdx, dofCodim);
volVars.update(sol[globalIdx],
this->problem_(),
@@ -195,34 +192,17 @@ public:
}
}
- if (isBox) // vertex data
- {
- writer.attachVertexData(*Sn, "Sn");
- writer.attachVertexData(*Sw, "Sw");
- writer.attachVertexData(*pN, "pn");
- writer.attachVertexData(*pW, "pw");
- writer.attachVertexData(*pC, "pc");
- writer.attachVertexData(*rhoW, "rhoW");
- writer.attachVertexData(*rhoN, "rhoN");
- writer.attachVertexData(*mobW, "mobW");
- writer.attachVertexData(*mobN, "mobN");
- writer.attachVertexData(*poro, "porosity");
- writer.attachVertexData(*Te, "temperature");
- }
- else // cell data
- {
- writer.attachCellData(*Sn, "Sn");
- writer.attachCellData(*Sw, "Sw");
- writer.attachCellData(*pN, "pn");
- writer.attachCellData(*pW, "pw");
- writer.attachCellData(*pC, "pc");
- writer.attachCellData(*rhoW, "rhoW");
- writer.attachCellData(*rhoN, "rhoN");
- writer.attachCellData(*mobW, "mobW");
- writer.attachCellData(*mobN, "mobN");
- writer.attachCellData(*poro, "porosity");
- writer.attachCellData(*Te, "temperature");
- }
+ writer.attachDofData(*Sn, "Sn", isBox);
+ writer.attachDofData(*Sw, "Sw", isBox);
+ writer.attachDofData(*pN, "pn", isBox);
+ writer.attachDofData(*pW, "pw", isBox);
+ writer.attachDofData(*pC, "pc", isBox);
+ writer.attachDofData(*rhoW, "rhoW", isBox);
+ writer.attachDofData(*rhoN, "rhoN", isBox);
+ writer.attachDofData(*mobW, "mobW", isBox);
+ writer.attachDofData(*mobN, "mobN", isBox);
+ writer.attachDofData(*poro, "porosity", isBox);
+ writer.attachDofData(*Te, "temperature", isBox);
writer.attachCellData(*rank, "process rank");
}
};
diff --git a/dumux/implicit/richards/richardsnewtoncontroller.hh b/dumux/implicit/richards/richardsnewtoncontroller.hh
index 2f0b035d6d34e99a81ba72e0b0e24628428c4bd4..5513dcf5af5c7efedf351a58bbefb53474d5d5a1 100644
--- a/dumux/implicit/richards/richardsnewtoncontroller.hh
+++ b/dumux/implicit/richards/richardsnewtoncontroller.hh
@@ -55,6 +55,7 @@ class RichardsNewtonController : public NewtonController<TypeTag>
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
enum { dim = GridView::dimension };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+ enum { dofCodim = isBox ? dim : 0 };
public:
/*!
@@ -96,11 +97,7 @@ public:
fvGeometry.update(gridView, *elemIt);
for (int scvIdx = 0; scvIdx < fvGeometry.numSCV; ++scvIdx)
{
- int globalIdx;
- if (isBox)
- globalIdx = this->problem_().vertexMapper().map(*elemIt, scvIdx, dim);
- else
- globalIdx = this->problem_().elementMapper().map(*elemIt);
+ int globalIdx = this->model_().dofMapper().map(*elemIt, scvIdx, dofCodim);
// calculate the old wetting phase saturation
const SpatialParams &spatialParams = this->problem_().spatialParams();
diff --git a/dumux/io/vtkmultiwriter.hh b/dumux/io/vtkmultiwriter.hh
index c4c65b0c80e8036ff7e1e6378b7adb52aebbaffe..00c0cb0149b078eba618f57485ef7d10c1c32f49 100644
--- a/dumux/io/vtkmultiwriter.hh
+++ b/dumux/io/vtkmultiwriter.hh
@@ -219,6 +219,22 @@ public:
curWriter_->addCellData(fnPtr);
}
+ /*!
+ * \brief Add data associated with degrees of freedom to the output.
+ *
+ * This is a wrapper for the functions attachVertexData and attachCelldata.
+ * Depending on the value of \param isVertexData, the appropriate function
+ * is selected.
+ */
+ template <class DataBuffer>
+ void attachDofData(DataBuffer &buf, std::string name, bool isVertexData, int nComps = 1)
+ {
+ if (isVertexData)
+ attachVertexData(buf, name, nComps);
+ else
+ attachCellData(buf, name, nComps);
+ }
+
/*!
* \brief Finalizes the current writer.
*