diff --git a/dumux/implicit/1p/1pmodel.hh b/dumux/implicit/1p/1pmodel.hh
index 132ab642345683abd97b2fa6df79b575a8a05630..8e5ad1e622b78f95579454b64b4a5ac17ab291bc 100644
--- a/dumux/implicit/1p/1pmodel.hh
+++ b/dumux/implicit/1p/1pmodel.hh
@@ -134,8 +134,6 @@ public:
velocityOutput.calculateVelocity(*velocity, elemVolVars, fvGeometry, *elemIt, /*phaseIdx=*/0);
}
- velocityOutput.completeVelocityCalculation(*velocity);
-
writer.attachDofData(*p, "p", isBox);
writer.attachDofData(*K, "K", isBox);
if (velocityOutput.enableOutput())
diff --git a/dumux/implicit/1p2c/1p2cmodel.hh b/dumux/implicit/1p2c/1p2cmodel.hh
index 110070e5401a7b8da2f7bea835b5d4e31049e980..be916cee7c3523a50c5086feedc6095b19270607 100644
--- a/dumux/implicit/1p2c/1p2cmodel.hh
+++ b/dumux/implicit/1p2c/1p2cmodel.hh
@@ -175,8 +175,6 @@ public:
velocityOutput.calculateVelocity(*velocity, elemVolVars, fvGeometry, *elemIt, phaseIdx);
}
- velocityOutput.completeVelocityCalculation(*velocity);
-
writer.attachDofData(pressure, "P", isBox);
writer.attachDofData(delp, "delp", isBox);
if (velocityOutput.enableOutput())
diff --git a/dumux/implicit/2p/2pmodel.hh b/dumux/implicit/2p/2pmodel.hh
index 3f39e856880ea0fc06226097b8ca503807e553f0..3d64da3563e44488cdd0006103425fa6d32c7c80 100644
--- a/dumux/implicit/2p/2pmodel.hh
+++ b/dumux/implicit/2p/2pmodel.hh
@@ -187,9 +187,6 @@ public:
velocityOutput.calculateVelocity(*velocityN, elemVolVars, fvGeometry, *elemIt, nPhaseIdx);
}
- velocityOutput.completeVelocityCalculation(*velocityW);
- velocityOutput.completeVelocityCalculation(*velocityN);
-
writer.attachDofData(*Sn, "Sn", isBox);
writer.attachDofData(*Sw, "Sw", isBox);
writer.attachDofData(*pN, "pn", isBox);
diff --git a/dumux/implicit/2p2c/2p2cmodel.hh b/dumux/implicit/2p2c/2p2cmodel.hh
index 2af0d88b85a9cf15f0c8f39ac1547e6bc97da490..97ab736e738d2cb8f39398f3c5b9507035e1c9b9 100644
--- a/dumux/implicit/2p2c/2p2cmodel.hh
+++ b/dumux/implicit/2p2c/2p2cmodel.hh
@@ -375,9 +375,6 @@ public:
} // loop over elements
- velocityOutput.completeVelocityCalculation(*velocityW);
- velocityOutput.completeVelocityCalculation(*velocityN);
-
writer.attachDofData(*sN, "Sn", isBox);
writer.attachDofData(*sW, "Sw", isBox);
writer.attachDofData(*pN, "pN", isBox);
diff --git a/dumux/implicit/3p3c/3p3cmodel.hh b/dumux/implicit/3p3c/3p3cmodel.hh
index 2add58906cd0b481fcf0cbc8b233c6060ce15808..b171e3283cf7e64aca7f9228a84f0f850098169d 100644
--- a/dumux/implicit/3p3c/3p3cmodel.hh
+++ b/dumux/implicit/3p3c/3p3cmodel.hh
@@ -379,10 +379,6 @@ public:
}
- velocityOutput.completeVelocityCalculation(*velocityW);
- velocityOutput.completeVelocityCalculation(*velocityN);
- velocityOutput.completeVelocityCalculation(*velocityG);
-
writer.attachDofData(*saturation[wPhaseIdx], "Sw", isBox);
writer.attachDofData(*saturation[nPhaseIdx], "Sn", isBox);
writer.attachDofData(*saturation[gPhaseIdx], "Sg", isBox);
diff --git a/dumux/implicit/common/implicitvelocityoutput.hh b/dumux/implicit/common/implicitvelocityoutput.hh
index a8b752bdb8113fa403023d79659bc749114767e6..455aadcebf7cedc56411a70dae3470f9554b97de 100644
--- a/dumux/implicit/common/implicitvelocityoutput.hh
+++ b/dumux/implicit/common/implicitvelocityoutput.hh
@@ -70,26 +70,41 @@ public:
* \param problem The problem to be solved
*/
ImplicitVelocityOutput(const Problem& problem)
- : problem_(problem)
+ : problem_(problem)
{
// check, if velocity output can be used (works only for cubes so far)
velocityOutput_ = GET_PARAM_FROM_GROUP(TypeTag, bool, Vtk, AddVelocity);
+ if (velocityOutput_ && isBox)
+ {
+ cellNum_.assign(problem_.gridView().size(dofCodim), 0);
+ }
+
ElementIterator elemIt = problem_.gridView().template begin<0>();
ElementIterator elemEndIt = problem_.gridView().template end<0>();
for (; elemIt != elemEndIt; ++elemIt)
{
- if (elemIt->geometry().type().isCube() == false){
+ if (elemIt->geometry().type().isCube() == false)
+ {
velocityOutput_ = false;
}
+
+ if (isBox)
+ {
+ FVElementGeometry fvGeometry;
+ fvGeometry.update(problem_.gridView(), *elemIt);
+
+ // transform vertex velocities from local to global coordinates
+ for (int scvIdx = 0; scvIdx < fvGeometry.numScv; ++scvIdx)
+ {
+ int globalIdx = problem_.vertexMapper().map(*elemIt, scvIdx, dofCodim);
+
+ cellNum_[globalIdx] += 1;
+ }
+ }
}
if (velocityOutput_ != GET_PARAM_FROM_GROUP(TypeTag, bool, Vtk, AddVelocity))
std::cout << "ATTENTION: Velocity output only works for cubes and is set to false for simplices\n";
-
- if (velocityOutput_ && isBox)
- {
- cellNum_.resize(problem_.gridView().size(dofCodim), 0);
- }
}
bool enableOutput()
@@ -137,7 +152,7 @@ public:
fvGeometry,
faceIdx,
elemVolVars);
-
+
const GlobalPosition globalNormal = fluxVars.face().normal;
GlobalPosition localNormal(0);
@@ -149,8 +164,7 @@ public:
// Get the Darcy velocities. The Darcy velocities are divided by the area of the subcontrolvolume
// face in the reference element.
- Scalar flux;
- flux = fluxVars.volumeFlux(phaseIdx) / localArea;
+ Scalar flux = fluxVars.volumeFlux(phaseIdx) / localArea;
// transform the normal Darcy velocity into a vector
tmpVelocity = localNormal;
@@ -168,11 +182,9 @@ public:
Dune::FieldVector<CoordScalar, dim> scvVelocity(0);
jacobianT2.mtv(scvVelocities[scvIdx], scvVelocity);
- scvVelocity /= element.geometry().integrationElement(localPos);
+ scvVelocity /= element.geometry().integrationElement(localPos)*cellNum_[globalIdx];
// add up the wetting phase subcontrolvolume velocities for each vertex
velocity[globalIdx] += scvVelocity;
-
- cellNum_[globalIdx] +=1;
}
}
else
@@ -206,12 +218,12 @@ public:
fvGeometry,
faceIdx,
elemVolVars,true);
-
+
Scalar flux = fluxVars.volumeFlux(phaseIdx);
scvVelocities[faceIdx] = flux;
}
}
-
+
Dune::FieldVector < CoordScalar, dim > refVelocity(0);
for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (scvVelocities[2*i + 1] - scvVelocities[2*i]);
@@ -228,32 +240,10 @@ public:
} // velocity output
}
- template<class VelocityVector>
- void completeVelocityCalculation(VelocityVector& velocity)
- {
- if (velocityOutput_ && isBox)
- {
- unsigned int size = velocity.size();
-
- assert(size != cellNum_.size());
-
-
- // divide the vertex velocities by the number of adjacent scvs i.e. cells
- for(unsigned int globalIdx = 0; globalIdx < size; ++globalIdx)
- {
- velocity[globalIdx] /= cellNum_[globalIdx];
- }
- fluxVariables_.clear();
- }
- }
protected:
const Problem& problem_;
-
- // parameters given in constructor
bool velocityOutput_;
std::vector<int> cellNum_;
-
- std::vector<FluxVariables> fluxVariables_;
};
}
diff --git a/dumux/implicit/mpnc/mpncvtkwritercommon.hh b/dumux/implicit/mpnc/mpncvtkwritercommon.hh
index c386545c1caff6e54750360f697e067de51c2efa..df346904bc112e0446509c37874222efee245cdb 100644
--- a/dumux/implicit/mpnc/mpncvtkwritercommon.hh
+++ b/dumux/implicit/mpnc/mpncvtkwritercommon.hh
@@ -65,7 +65,7 @@ class MPNCVtkWriterCommon : public MPNCVtkWriterModule<TypeTag>
typedef Dune::array<DimField, numPhases> PhaseDimField;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
enum { dofCodim = isBox ? dim : 0 };
-
+
public:
MPNCVtkWriterCommon(const Problem &problem)
: ParentType(problem), velocityOutput_(problem)
@@ -99,7 +99,7 @@ public:
if (velocityOutput_.enableOutput()) {
Scalar nDofs = this->problem_.model().numDofs();
- for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
+ for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
velocity_[phaseIdx].resize(nDofs);
velocity_[phaseIdx] = 0;
}
@@ -127,7 +127,7 @@ public:
for (int scvIdx = 0; scvIdx < fvGeometry.numScv; ++scvIdx)
{
int globalIdx = this->problem_.model().dofMapper().map(element, scvIdx, dofCodim);
-
+
const VolumeVariables &volVars = elemVolVars[scvIdx];
if (porosityOutput_) porosity_[globalIdx] = volVars.porosity();
@@ -203,7 +203,6 @@ public:
if (velocityOutput_.enableOutput()) {
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
- velocityOutput_.completeVelocityCalculation(velocity_[phaseIdx]);
std::ostringstream oss;
oss << "velocity_" << FluidSystem::phaseName(phaseIdx);
writer.attachDofData(velocity_[phaseIdx],
diff --git a/dumux/implicit/richards/richardsmodel.hh b/dumux/implicit/richards/richardsmodel.hh
index 59ac1d62ef74fd9d605b8d4e4dd9ef752de27f75..bc75f3279cb13965a2b5f8044348c22fefef1abe 100644
--- a/dumux/implicit/richards/richardsmodel.hh
+++ b/dumux/implicit/richards/richardsmodel.hh
@@ -206,8 +206,6 @@ public:
velocityOutput.calculateVelocity(*velocity, elemVolVars, fvGeometry, *elemIt, /*phaseIdx=*/0);
}
- velocityOutput.completeVelocityCalculation(*velocity);
-
writer.attachDofData(*Sn, "Sn", isBox);
writer.attachDofData(*Sw, "Sw", isBox);
writer.attachDofData(*pN, "pn", isBox);