From 33fd3a35737b8dedb317eab87c68fe96de2b1387 Mon Sep 17 00:00:00 2001 From: Martin Schneider <martin.schneider@iws.uni-stuttgart.de> Date: Thu, 22 Nov 2018 09:56:44 +0100 Subject: [PATCH] [velocityoutput] Divide velocity by extrusion factor --- .../discretization/box/subcontrolvolumeface.hh | 2 +- dumux/io/vtkoutputmodule.hh | 4 ++-- dumux/porousmediumflow/velocityoutput.hh | 17 +++++++++++------ 3 files changed, 14 insertions(+), 9 deletions(-) diff --git a/dumux/discretization/box/subcontrolvolumeface.hh b/dumux/discretization/box/subcontrolvolumeface.hh index 17a00578bf..9eca2a0652 100644 --- a/dumux/discretization/box/subcontrolvolumeface.hh +++ b/dumux/discretization/box/subcontrolvolumeface.hh @@ -195,7 +195,7 @@ public: return scvIndices_[1]; } - //! The global index of this sub control volume face + //! The local index of this sub control volume face GridIndexType index() const { return scvfIndex_; diff --git a/dumux/io/vtkoutputmodule.hh b/dumux/io/vtkoutputmodule.hh index bbe20095b9..94b74580ba 100644 --- a/dumux/io/vtkoutputmodule.hh +++ b/dumux/io/vtkoutputmodule.hh @@ -258,7 +258,7 @@ private: //! (1) Assemble all variable fields and add to writer ////////////////////////////////////////////////////////////// - // instatiate the velocity output + // instantiate the velocity output using VelocityVector = typename VelocityOutput::VelocityVector; std::vector<VelocityVector> velocity(velocityOutput_->numPhases()); @@ -419,7 +419,7 @@ private: writer_->clear(); } - //! Assembles the fields and adds them to the writer (conforming output) + //! Assembles the fields and adds them to the writer (nonconforming output) void writeNonConforming_(double time, Dune::VTK::OutputType type) { if(!isBox) diff --git a/dumux/porousmediumflow/velocityoutput.hh b/dumux/porousmediumflow/velocityoutput.hh index 320cf6ff54..3cf52dfe19 100644 --- a/dumux/porousmediumflow/velocityoutput.hh +++ b/dumux/porousmediumflow/velocityoutput.hh @@ -182,15 +182,17 @@ public: jacobianT1.mv(globalNormal, localNormal); localNormal /= localNormal.two_norm(); - // insantiate the flux variables + // instantiate the flux variables FluxVariables fluxVars; fluxVars.init(problem_, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache); // get the volume flux divided by the area of the // subcontrolvolume face in the reference element - // TODO: Divide by extrusion factor!!? Scalar localArea = scvfReferenceArea_(geomType, scvf.index()); Scalar flux = fluxVars.advectiveFlux(phaseIdx, upwindTerm) / localArea; + flux /= problem_.extrusionFactor(element, + fvGeometry.scv(scvf.insideScvIdx()), + elementSolution(element, elemVolVars, fvGeometry)); // transform the volume flux into a velocity vector Velocity tmpVelocity = localNormal; @@ -230,17 +232,19 @@ public: // first we extract the corner indices for each scv for the CIV method // for network grids there might be multiple intersection with the same geometryInInside - // we indentify those by the indexInInside for now (assumes conforming grids at branching facets) + // we identify those by the indexInInside for now (assumes conforming grids at branching facets) // here we keep track of them std::vector<bool> handledScvf; - if (dim < dimWorld) handledScvf.resize(element.subEntities(1), false); + if (dim < dimWorld) + handledScvf.resize(element.subEntities(1), false); // find the local face indices of the scvfs (for conforming meshes) std::vector<unsigned int> scvfIndexInInside(fvGeometry.numScvf()); int localScvfIdx = 0; for (const auto& intersection : intersections(fvGridGeometry_.gridView(), element)) { - if (dim < dimWorld) if (handledScvf[intersection.indexInInside()]) continue; + if (dim < dimWorld && handledScvf[intersection.indexInInside()]) + continue; if (intersection.neighbor() || intersection.boundary()) { @@ -248,7 +252,8 @@ public: scvfIndexInInside[localScvfIdx++] = intersection.indexInInside(); // for surface and network grids mark that we handled this face - if (dim < dimWorld) handledScvf[intersection.indexInInside()] = true; + if (dim < dimWorld) + handledScvf[intersection.indexInInside()] = true; } } -- GitLab