From 6900dea7c8518cba2530b421db789dd0a258d48a Mon Sep 17 00:00:00 2001
From: Markus Wolff <markus.wolff@twt-gmbh.de>
Date: Mon, 20 Aug 2012 13:13:52 +0000
Subject: [PATCH] Included support of VtkOutputLevel in decoupled 2p models
- OutputLevel 0 (default) only outputs the primary variables, level > 0 outputs
everything
- Changed corresponding tests and references
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@8901 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
CMakeLists.txt | 3 +
.../decoupled/2p/diffusion/fv/fvpressure2p.hh | 106 ++++++++----
.../decoupled/2p/diffusion/fv/fvvelocity2p.hh | 137 ++++++++--------
.../lmethod/fvmpfal2pfaboundpressure2p.hh | 112 ++++++++-----
.../fvmpfal2pfaboundpressure2padaptive.hh | 118 +++++++++-----
.../lmethod/fvmpfal2pfaboundvelocity2p.hh | 130 ++++++++-------
.../fvmpfal2pfaboundvelocity2padaptive.hh | 153 ++++++++++--------
.../omethod/fvmpfao2pfaboundpressure2p.hh | 112 ++++++++-----
.../omethod/fvmpfao2pfaboundvelocity2p.hh | 132 ++++++++-------
.../fvmpfa/omethod/fvmpfaopressure2p.hh | 24 ++-
.../2p/transport/fv/fvsaturation2p.hh | 45 +++++-
test/decoupled/2p/CMakeLists.txt | 3 +
test/decoupled/2p/test_impes-reference.vtu | 147 +----------------
test/decoupled/2p/test_impes.input | 3 +
test/decoupled/2p/test_mpfa2p-reference.vtu | 126 +++++++++++++++
test/decoupled/2p/test_mpfa2p.input | 3 +
test/decoupled/2p/test_mpfa2pproblem.hh | 17 +-
.../decoupled/2p/test_transport-reference.vtu | 11 --
test/decoupled/2p/test_transport.input | 3 +
19 files changed, 811 insertions(+), 574 deletions(-)
create mode 100644 test/decoupled/2p/test_mpfa2p-reference.vtu
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 064ce11c04..d1f43c6041 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -400,6 +400,8 @@ file(COPY test/decoupled/2p/test_transport-reference.vtu DESTINATION references)
file(COPY test/decoupled/2p/test_transport.input DESTINATION test/decoupled/2p)
file(COPY test/decoupled/2p/test_impes-reference.vtu DESTINATION references)
file(COPY test/decoupled/2p/test_impes.input DESTINATION test/decoupled/2p)
+file(COPY test/decoupled/2p/test_mpfa2p-reference.vtu DESTINATION references)
+file(COPY test/decoupled/2p/test_mpfa2p.input DESTINATION test/decoupled/2p)
file(COPY test/decoupled/2padaptive/test_2padaptive-reference.vtu DESTINATION references)
file(COPY test/decoupled/2padaptive/test_impesadaptive.input DESTINATION test/decoupled/2padaptive)
file(COPY test/decoupled/2p2c/test_dec2p2c-reference.vtu DESTINATION references)
@@ -451,6 +453,7 @@ add_test(test_diffusion bin/runTest.sh references/diffusion-reference.vtu
add_test(test_dec1p bin/runTest.sh references/test_1p-reference.vtu test_1p-00001.vtu test/decoupled/1p/test_dec1p -ParameterFile test/decoupled/1p/test_1p.input)
add_test(test_transport bin/runTest.sh references/test_transport-reference.vtu test_transport-00006.vtu test/decoupled/2p/test_transport -Grid.File test/decoupled/2p/grids/test_transport.dgf)
add_test(test_impes bin/runTest.sh references/test_impes-reference.vtu test_impes-00013.vtu test/decoupled/2p/test_impes)
+add_test(test_mpfa2p bin/runTest.sh references/test_mpfa2p-reference.vtu test_mpfa2p-00007.vtu test/decoupled/2p/test_mpfa2p)
add_test(test_impesadaptive bin/runTest.sh references/test_2padaptive-reference.vtu test_2padaptive-00007.vtu test/decoupled/2padaptive/test_impesadaptive)
add_test(test_dec2p2c bin/runTest.sh references/test_dec2p2c-reference.vtu test_dec2p2c-00021.vtu test/decoupled/2p2c/test_dec2p2c)
add_test(test_multiphysics2p2c bin/runTest.sh references/test_multiphysics2p2c-reference.vtu test_multiphysics2p2c-00021.vtu test/decoupled/2p2c/test_multiphysics2p2c)
diff --git a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
index 01bdd17de5..b615aa5a63 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
@@ -339,55 +339,97 @@ public:
void addOutputVtkFields(MultiWriter &writer)
{
int size = problem_.gridView().size(0);
- ScalarSolutionType *pressureW = writer.allocateManagedBuffer(size);
- ScalarSolutionType *pressureN = writer.allocateManagedBuffer(size);
- ScalarSolutionType *pC = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pressure = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pressureSecond = 0;
+ ScalarSolutionType *pC = 0;
+
+ if (vtkOutputLevel_ > 0)
+ {
+ pressureSecond = writer.allocateManagedBuffer(size);
+ pC = writer.allocateManagedBuffer(size);
+ }
for (int i = 0; i < size; i++)
{
CellData& cellData = problem_.variables().cellData(i);
storePressureSolution(i, cellData);
- (*pressureW)[i] = cellData.pressure(wPhaseIdx);
- (*pressureN)[i] = cellData.pressure(nPhaseIdx);
- (*pC)[i] = cellData.capillaryPressure();
- if (pressureType_ == pglobal)
+
+ if (pressureType_ == pw)
+ {
+ (*pressure)[i] = cellData.pressure(wPhaseIdx);
+ if (vtkOutputLevel_ > 0)
+ {
+ (*pressureSecond)[i] = cellData.pressure(nPhaseIdx);
+ }
+ }
+ else if (pressureType_ == pn)
{
- (*pressureW)[i] = cellData.globalPressure();
+ (*pressure)[i] = cellData.pressure(nPhaseIdx);
+ if (vtkOutputLevel_ > 0)
+ {
+ (*pressureSecond)[i] = cellData.pressure(wPhaseIdx);
+ }
+ }
+ else if (pressureType_ == pglobal)
+ {
+ (*pressure)[i] = cellData.globalPressure();
+ }
+ if (vtkOutputLevel_ > 0)
+ {
+ (*pC)[i] = cellData.capillaryPressure();
+ }
+ }
+
+ if (pressureType_ == pw)
+ {
+ writer.attachCellData(*pressure, "wetting pressure");
+ if (vtkOutputLevel_ > 0)
+ {
+ writer.attachCellData(*pressureSecond, "nonwetting pressure");
+ }
+ }
+ else if (pressureType_ == pn)
+ {
+ writer.attachCellData(*pressure, "nonwetting pressure");
+ if (vtkOutputLevel_ > 0)
+ {
+ writer.attachCellData(*pressureSecond, "wetting pressure");
}
}
if (pressureType_ == pglobal)
{
- writer.attachCellData(*pressureW, "global pressure");
+ writer.attachCellData(*pressure, "global pressure");
}
- else
+
+ if (vtkOutputLevel_ > 0)
{
- writer.attachCellData(*pressureW, "wetting pressure");
- writer.attachCellData(*pressureN, "nonwetting pressure");
+ writer.attachCellData(*pC, "capillary pressure");
}
- writer.attachCellData(*pC, "capillary pressure");
-
if (compressibility_)
{
- ScalarSolutionType *densityWetting = writer.allocateManagedBuffer(size);
- ScalarSolutionType *densityNonwetting = writer.allocateManagedBuffer(size);
- ScalarSolutionType *viscosityWetting = writer.allocateManagedBuffer(size);
- ScalarSolutionType *viscosityNonwetting = writer.allocateManagedBuffer(size);
-
- for (int i = 0; i < size; i++)
+ if (vtkOutputLevel_ > 0)
{
- CellData& cellData = problem_.variables().cellData(i);
- (*densityWetting)[i] = cellData.density(wPhaseIdx);
- (*densityNonwetting)[i] = cellData.density(nPhaseIdx);
- (*viscosityWetting)[i] = cellData.viscosity(wPhaseIdx);
- (*viscosityNonwetting)[i] = cellData.viscosity(nPhaseIdx);
- }
+ ScalarSolutionType *densityWetting = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *densityNonwetting = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *viscosityWetting = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *viscosityNonwetting = writer.allocateManagedBuffer(size);
- writer.attachCellData(*densityWetting, "wetting density");
- writer.attachCellData(*densityNonwetting, "nonwetting density");
- writer.attachCellData(*viscosityWetting, "wetting viscosity");
- writer.attachCellData(*viscosityNonwetting, "nonwetting viscosity");
+ for (int i = 0; i < size; i++)
+ {
+ CellData& cellData = problem_.variables().cellData(i);
+ (*densityWetting)[i] = cellData.density(wPhaseIdx);
+ (*densityNonwetting)[i] = cellData.density(nPhaseIdx);
+ (*viscosityWetting)[i] = cellData.viscosity(wPhaseIdx);
+ (*viscosityNonwetting)[i] = cellData.viscosity(nPhaseIdx);
+ }
+
+ writer.attachCellData(*densityWetting, "wetting density");
+ writer.attachCellData(*densityNonwetting, "nonwetting density");
+ writer.attachCellData(*viscosityWetting, "wetting viscosity");
+ writer.attachCellData(*viscosityNonwetting, "nonwetting viscosity");
+ }
}
return;
@@ -421,6 +463,8 @@ public:
density_[nPhaseIdx] = 0.;
viscosity_[wPhaseIdx] = 0.;
viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
private:
@@ -436,6 +480,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const bool compressibility_ = GET_PROP_VALUE(TypeTag, EnableCompressibility);
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$p_w\f$, \f$p_n\f$, \f$p_{global}\f$)
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$S_w\f$, \f$S_n\f$)
diff --git a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
index 02b36c210e..f261836989 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
@@ -131,6 +131,8 @@ public:
density_[nPhaseIdx] = 0.;
viscosity_[wPhaseIdx] = 0.;
viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
//! For initialization
@@ -178,84 +180,87 @@ public:
template<class MultiWriter>
void addOutputVtkFields(MultiWriter &writer)
{
- Dune::BlockVector < Dune::FieldVector<Scalar, dim> > &velocity = *(writer.template allocateManagedBuffer<Scalar,
- dim>(problem_.gridView().size(0)));
- Dune::BlockVector < Dune::FieldVector<Scalar, dim> > &velocitySecondPhase =
- *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
-
- // compute update vector
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ if (vtkOutputLevel_ > 0)
{
- // cell index
- int globalIdx = problem_.variables().index(*eIt);
-
- CellData& cellData = problem_.variables().cellData(globalIdx);
-
- Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
- Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
- // run through all intersections with neighbors and boundary
- IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
- for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
+ Dune::BlockVector < Dune::FieldVector<Scalar, dim> > &velocity = *(writer.template allocateManagedBuffer<Scalar,
+ dim>(problem_.gridView().size(0)));
+ Dune::BlockVector < Dune::FieldVector<Scalar, dim> > &velocitySecondPhase =
+ *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
+
+ // compute update vector
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
{
- int isIndex = isIt->indexInInside();
-
- fluxW[isIndex] += isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
- fluxNW[isIndex] += isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(nPhaseIdx, isIndex));
- }
-
- Dune::FieldVector < Scalar, dim > refVelocity(0);
- for (int i = 0; i < dim; i++)
+ // cell index
+ int globalIdx = problem_.variables().index(*eIt);
+
+ CellData& cellData = problem_.variables().cellData(globalIdx);
+
+ Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
+ Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
+ // run through all intersections with neighbors and boundary
+ IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
+ for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
+ {
+ int isIndex = isIt->indexInInside();
+
+ fluxW[isIndex] += isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
+ fluxNW[isIndex] += isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(nPhaseIdx, isIndex));
+ }
+
+ Dune::FieldVector < Scalar, dim > refVelocity(0);
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxW[2*i + 1] - fluxW[2*i]);
- const Dune::FieldVector<Scalar, dim>& localPos =
- ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
+ const Dune::FieldVector<Scalar, dim>& localPos =
+ ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
- // get the transposed Jacobian of the element mapping
- const DimMatrix& jacobianInv = eIt->geometry().jacobianInverseTransposed(localPos);
- DimMatrix jacobianT(jacobianInv);
- jacobianT.invert();
+ // get the transposed Jacobian of the element mapping
+ const DimMatrix& jacobianInv = eIt->geometry().jacobianInverseTransposed(localPos);
+ DimMatrix jacobianT(jacobianInv);
+ jacobianT.invert();
- // calculate the element velocity by the Piola transformation
- Dune::FieldVector < Scalar, dim > elementVelocity(0);
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ Dune::FieldVector < Scalar, dim > elementVelocity(0);
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocity[globalIdx] = elementVelocity;
+ velocity[globalIdx] = elementVelocity;
- refVelocity = 0;
- for (int i = 0; i < dim; i++)
+ refVelocity = 0;
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxNW[2*i + 1] - fluxNW[2*i]);
- // calculate the element velocity by the Piola transformation
- elementVelocity = 0;
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
-
- velocitySecondPhase[globalIdx] = elementVelocity;
- }
+ // calculate the element velocity by the Piola transformation
+ elementVelocity = 0;
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- //switch velocities
- switch (velocityType_)
- {
- case vw:
- {
- writer.attachCellData(velocity, "wetting-velocity", dim);
- writer.attachCellData(velocitySecondPhase, "non-wetting-velocity", dim);
- break;
- }
- case vn:
- {
- writer.attachCellData(velocity, "non-wetting-velocity", dim);
- writer.attachCellData(velocitySecondPhase, "wetting-velocity", dim);
- break;
+ velocitySecondPhase[globalIdx] = elementVelocity;
}
- case vt:
+
+ //switch velocities
+ switch (velocityType_)
{
- writer.attachCellData(velocity, "total velocity", dim);
- break;
+ case vw:
+ {
+ writer.attachCellData(velocity, "wetting-velocity", dim);
+ writer.attachCellData(velocitySecondPhase, "non-wetting-velocity", dim);
+ break;
+ }
+ case vn:
+ {
+ writer.attachCellData(velocity, "non-wetting-velocity", dim);
+ writer.attachCellData(velocitySecondPhase, "wetting-velocity", dim);
+ break;
+ }
+ case vt:
+ {
+ writer.attachCellData(velocity, "total velocity", dim);
+ break;
+ }
}
}
@@ -268,6 +273,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const int velocityType_ = GET_PROP_VALUE(TypeTag, VelocityFormulation);//!< gives kind of velocity used (\f$ 0 = v_w\f$, \f$ 1 = v_n\f$, \f$ 2 = v_t\f$)
static const bool compressibility_ = GET_PROP_VALUE(TypeTag, EnableCompressibility);
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation);//!< gives kind of pressure used (\f$p_w\f$, \f$p_n\f$, \f$p_{global}\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
index 46a799b075..8ef79f66f4 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
@@ -186,6 +186,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
updateMaterialLaws();
interactionVolumes_.resize(problem_.gridView().size(dim));
@@ -297,6 +309,8 @@ public:
/*! \brief Adds pressure output to the output file
*
* Adds the pressure, the potential and the capillary pressure to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -305,45 +319,72 @@ public:
template<class MultiWriter>
void addOutputVtkFields(MultiWriter &writer)
{
- ScalarSolutionType &pressure = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &potential = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &pC = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
+ int size = problem_.gridView().size(0);
+ ScalarSolutionType *potential = writer.allocateManagedBuffer(size);
- potential = this->pressure();
+ (*potential) = this->pressure();
- ElementIterator eItBegin = problem_.gridView().template begin<0>();
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ if (pressureType_ == pw)
{
- int idx = problem_.variables().index(*eIt);
+ writer.attachCellData(*potential, "wetting potential");
+ }
+
+ if (pressureType_ == pn)
+ {
+ writer.attachCellData(*potential, "nonwetting potential");
+ }
+
+ if (vtkOutputLevel_ > 0)
+ {
+ ScalarSolutionType *pressure = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pressureSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *potentialSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pC = writer.allocateManagedBuffer(size);
+
+ ElementIterator eItBegin = problem_.gridView().template begin<0>();
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ {
+ int idx = problem_.variables().index(*eIt);
+ CellData& cellData = problem_.variables().cellData(idx);
+
+ storePressureSolution(*eIt);
+ (*pC)[idx] = cellData.capillaryPressure();
+
+ if (pressureType_ == pw)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(nPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] + cellData.capillaryPressure();
+ }
+
+ if (pressureType_ == pn)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(wPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] - cellData.capillaryPressure();
+ }
+ }
+
if (pressureType_ == pw)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ writer.attachCellData(*pressure, "wetting pressure");
+ writer.attachCellData(*pressureSecond, "nonwetting pressure");
+ writer.attachCellData(*potentialSecond, "nonwetting potential");
}
if (pressureType_ == pn)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ writer.attachCellData(*pressure, "nonwetting pressure");
+ writer.attachCellData(*pressureSecond, "wetting pressure");
+ writer.attachCellData(*potentialSecond, "wetting potential");
}
- pC[idx] = problem_.variables().cellData(idx).capillaryPressure();
- }
- if (pressureType_ == pw)
- {
- writer.attachCellData(pressure, "wetting pressure");
- writer.attachCellData(potential, "wetting potential");
- }
-
- if (pressureType_ == pn)
- {
- writer.attachCellData(pressure, "nonwetting pressure");
- writer.attachCellData(potential, "nonwetting potential");
+ writer.attachCellData(*pC, "capillary pressure");
}
- writer.attachCellData(pC, "capillary pressure");
-
return;
}
@@ -384,17 +425,12 @@ public:
ErrorTermLowerBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermLowerBound);
ErrorTermUpperBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermUpperBound);
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
private:
@@ -425,6 +461,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$ 0 = S_w\f$, \f$ 1 = S_n\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
index c12186771c..37c5b3856b 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
@@ -199,6 +199,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
updateMaterialLaws();
updateInteractionVolumeInfo();
@@ -318,6 +330,8 @@ public:
/*! \brief Adds pressure output to the output file
*
* Adds the pressure, the potential and the capillary pressure to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -326,51 +340,72 @@ public:
template<class MultiWriter>
void addOutputVtkFields(MultiWriter &writer)
{
- ScalarSolutionType &pressure = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &potentialW = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &potentialNW = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &pC = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
+ int size = problem_.gridView().size(0);
+ ScalarSolutionType *potential = writer.allocateManagedBuffer(size);
-// potential = this->pressure();
+ (*potential) = this->pressure();
- ElementIterator eItBegin = problem_.gridView().template begin<0>();
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ if (pressureType_ == pw)
{
- int idx = problem_.variables().index(*eIt);
+ writer.attachCellData(*potential, "wetting potential");
+ }
+
+ if (pressureType_ == pn)
+ {
+ writer.attachCellData(*potential, "nonwetting potential");
+ }
+
+ if (vtkOutputLevel_ > 0)
+ {
+ ScalarSolutionType *pressure = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pressureSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *potentialSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pC = writer.allocateManagedBuffer(size);
+
+ ElementIterator eItBegin = problem_.gridView().template begin<0>();
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ {
+ int idx = problem_.variables().index(*eIt);
+ CellData& cellData = problem_.variables().cellData(idx);
+
+ storePressureSolution(*eIt);
+ (*pC)[idx] = cellData.capillaryPressure();
+
+ if (pressureType_ == pw)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(nPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] + cellData.capillaryPressure();
+ }
+
+ if (pressureType_ == pn)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(wPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] - cellData.capillaryPressure();
+ }
+ }
+
if (pressureType_ == pw)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
- potentialW[idx] = this->pressure()[idx][0];
- potentialNW[idx] = problem_.variables().cellData(idx).pressure(nPhaseIdx);
+ writer.attachCellData(*pressure, "wetting pressure");
+ writer.attachCellData(*pressureSecond, "nonwetting pressure");
+ writer.attachCellData(*potentialSecond, "nonwetting potential");
}
if (pressureType_ == pn)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
- potentialNW[idx] = this->pressure()[idx][0];
- potentialW[idx] = problem_.variables().cellData(idx).pressure(wPhaseIdx);
+ writer.attachCellData(*pressure, "nonwetting pressure");
+ writer.attachCellData(*pressureSecond, "wetting pressure");
+ writer.attachCellData(*potentialSecond, "wetting potential");
}
- pC[idx] = problem_.variables().cellData(idx).capillaryPressure();
- }
- if (pressureType_ == pw)
- {
- writer.attachCellData(pressure, "wetting pressure");
- writer.attachCellData(potentialW, "wetting potential");
- writer.attachCellData(potentialNW, "non-wetting potential");
- }
-
- if (pressureType_ == pn)
- {
- writer.attachCellData(pressure, "nonwetting pressure");
- writer.attachCellData(potentialNW, "nonwetting potential");
+ writer.attachCellData(*pC, "capillary pressure");
}
- writer.attachCellData(pC, "capillary pressure");
-
return;
}
@@ -417,17 +452,12 @@ public:
ErrorTermLowerBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermLowerBound);
ErrorTermUpperBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermUpperBound);
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
private:
@@ -469,6 +499,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$ 0 = S_w\f$, \f$ 1 = S_n\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2p.hh
index 03a3e284d9..76e5627f81 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2p.hh
@@ -139,17 +139,12 @@ public:
FVMPFAL2PFABoundVelocity2P(Problem& problem) :
ParentType(problem), problem_(problem), gravity_(problem.gravity())
{
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
//Calculates the velocities at all cell-cell interfaces.
@@ -163,6 +158,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
calculateVelocity();
return;
@@ -185,6 +192,8 @@ public:
/*! \brief Adds velocity output to the output file
*
* Adds the phase velocities or a total velocity (depending on the formulation) to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -195,68 +204,71 @@ public:
{
ParentType::addOutputVtkFields(writer);
- Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
- Scalar, dim>(problem_.gridView().size(0)));
- Dune::BlockVector < DimVector > &velocityNonwetting =
- *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
-
- // compute update vector
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ if (vtkOutputLevel_ > 0)
{
- // cell index
- int globalIdx = problem_.variables().index(*eIt);
+ Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
+ Scalar, dim>(problem_.gridView().size(0)));
+ Dune::BlockVector < DimVector > &velocityNonwetting =
+ *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
+
+ // compute update vector
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ // cell index
+ int globalIdx = problem_.variables().index(*eIt);
- CellData& cellData = problem_.variables().cellData(globalIdx);
+ CellData& cellData = problem_.variables().cellData(globalIdx);
- Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
- Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
- // run through all intersections with neighbors and boundary
- IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
- for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
- {
- int isIndex = isIt->indexInInside();
-
- fluxW[isIndex] = isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
- fluxNW[isIndex] =
- isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal()
- * cellData.fluxData().velocity(nPhaseIdx, isIndex));
- }
+ Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
+ Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
+ // run through all intersections with neighbors and boundary
+ IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
+ for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
+ {
+ int isIndex = isIt->indexInInside();
+
+ fluxW[isIndex] = isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
+ fluxNW[isIndex] =
+ isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal()
+ * cellData.fluxData().velocity(nPhaseIdx, isIndex));
+ }
- DimVector refVelocity(0);
- for (int i = 0; i < dim; i++)
+ DimVector refVelocity(0);
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxW[2*i + 1] - fluxW[2*i]);
- const DimVector& localPos =
- ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
+ const DimVector& localPos =
+ ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
- // get the transposed Jacobian of the element mapping
- const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
+ // get the transposed Jacobian of the element mapping
+ const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
- // calculate the element velocity by the Piola transformation
- DimVector elementVelocity(0);
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ DimVector elementVelocity(0);
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocityWetting[globalIdx] = elementVelocity;
+ velocityWetting[globalIdx] = elementVelocity;
- refVelocity = 0;
- for (int i = 0; i < dim; i++)
+ refVelocity = 0;
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxNW[2*i + 1] - fluxNW[2*i]);
- // calculate the element velocity by the Piola transformation
- elementVelocity = 0;
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ elementVelocity = 0;
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocityNonwetting[globalIdx] = elementVelocity;
- }
+ velocityNonwetting[globalIdx] = elementVelocity;
+ }
- writer.attachCellData(velocityWetting, "wetting-velocity", dim);
- writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ writer.attachCellData(velocityWetting, "wetting-velocity", dim);
+ writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ }
return;
}
@@ -267,6 +279,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int velocityType_ = GET_PROP_VALUE(TypeTag, VelocityFormulation); //!< gives kind of velocity used (\f$ 0 = v_w\f$, \f$ 1 = v_n\f$, \f$ 2 = v_t\f$)
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2padaptive.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2padaptive.hh
index b88b21d7bf..7ce8226574 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2padaptive.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundvelocity2padaptive.hh
@@ -147,17 +147,12 @@ public:
FVMPFAL2PFABoundVelocity2PAdaptive(Problem& problem) :
ParentType(problem), problem_(problem), gravity_(problem.gravity())
{
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
//Calculates the velocities at all cell-cell interfaces.
@@ -171,6 +166,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
calculateVelocity();
return;
@@ -193,80 +200,84 @@ public:
/*! \brief Adds velocity output to the output file
*
* Adds the phase velocities or a total velocity (depending on the formulation) to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
*
*/
- template<class MultiWriter>
- void addOutputVtkFields(MultiWriter &writer)
- {
- ParentType::addOutputVtkFields(writer);
-
- Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
- Scalar, dim>(problem_.gridView().size(0)));
- Dune::BlockVector < DimVector > &velocityNonwetting =
- *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
-
- // compute update vector
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
- {
- // cell index
- int globalIdx = problem_.variables().index(*eIt);
-
- CellData& cellData = problem_.variables().cellData(globalIdx);
-
- Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
- Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
- // run through all intersections with neighbors and boundary
- IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
- for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
- {
- int isIndex = isIt->indexInInside();
-
- fluxW[isIndex] += isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
- fluxNW[isIndex] +=
- isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal()
- * cellData.fluxData().velocity(nPhaseIdx, isIndex));
- }
-
- DimVector refVelocity(0);
- for (int i = 0; i < dim; i++)
+ template<class MultiWriter>
+ void addOutputVtkFields(MultiWriter &writer)
+ {
+ ParentType::addOutputVtkFields(writer);
+
+ if (vtkOutputLevel_ > 0)
+ {
+ Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
+ Scalar, dim>(problem_.gridView().size(0)));
+ Dune::BlockVector < DimVector > &velocityNonwetting =
+ *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
+
+ // compute update vector
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ // cell index
+ int globalIdx = problem_.variables().index(*eIt);
+
+ CellData& cellData = problem_.variables().cellData(globalIdx);
+
+ Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
+ Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
+ // run through all intersections with neighbors and boundary
+ IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
+ for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
+ {
+ int isIndex = isIt->indexInInside();
+
+ fluxW[isIndex] += isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
+ fluxNW[isIndex] +=
+ isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal()
+ * cellData.fluxData().velocity(nPhaseIdx, isIndex));
+ }
+
+ DimVector refVelocity(0);
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxW[2*i + 1] - fluxW[2*i]);
- const DimVector& localPos =
- ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
+ const DimVector& localPos =
+ ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
- // get the transposed Jacobian of the element mapping
- const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
+ // get the transposed Jacobian of the element mapping
+ const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
- // calculate the element velocity by the Piola transformation
- DimVector elementVelocity(0);
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ DimVector elementVelocity(0);
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocityWetting[globalIdx] = elementVelocity;
+ velocityWetting[globalIdx] = elementVelocity;
- refVelocity = 0;
- for (int i = 0; i < dim; i++)
+ refVelocity = 0;
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxNW[2*i + 1] - fluxNW[2*i]);
- // calculate the element velocity by the Piola transformation
- elementVelocity = 0;
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
-
- velocityNonwetting[globalIdx] = elementVelocity;
- }
+ // calculate the element velocity by the Piola transformation
+ elementVelocity = 0;
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- writer.attachCellData(velocityWetting, "wetting-velocity", dim);
- writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ velocityNonwetting[globalIdx] = elementVelocity;
+ }
- return;
- }
+ writer.attachCellData(velocityWetting, "wetting-velocity", dim);
+ writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ }
+ return;
+ }
private:
Problem& problem_;
@@ -275,6 +286,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int velocityType_ = GET_PROP_VALUE(TypeTag, VelocityFormulation); //!< gives kind of velocity used (\f$ 0 = v_w\f$, \f$ 1 = v_n\f$, \f$ 2 = v_t\f$)
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
index a2734ced08..2ac546010c 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
@@ -176,6 +176,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
updateMaterialLaws();
interactionVolumes_.resize(problem_.gridView().size(dim));
@@ -286,6 +298,8 @@ public:
/*! \brief Adds pressure output to the output file
*
* Adds the pressure, the potential and the capillary pressure to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -294,45 +308,72 @@ public:
template<class MultiWriter>
void addOutputVtkFields(MultiWriter &writer)
{
- ScalarSolutionType &pressure = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &potential = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
- ScalarSolutionType &pC = *(writer.template allocateManagedBuffer<Scalar, 1>(problem_.gridView().size(0)));
+ int size = problem_.gridView().size(0);
+ ScalarSolutionType *potential = writer.allocateManagedBuffer(size);
- potential = this->pressure();
+ (*potential) = this->pressure();
- ElementIterator eItBegin = problem_.gridView().template begin<0>();
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ if (pressureType_ == pw)
{
- int idx = problem_.variables().index(*eIt);
+ writer.attachCellData(*potential, "wetting potential");
+ }
+
+ if (pressureType_ == pn)
+ {
+ writer.attachCellData(*potential, "nonwetting potential");
+ }
+
+ if (vtkOutputLevel_ > 0)
+ {
+ ScalarSolutionType *pressure = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pressureSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *potentialSecond = writer.allocateManagedBuffer(size);
+ ScalarSolutionType *pC = writer.allocateManagedBuffer(size);
+
+ ElementIterator eItBegin = problem_.gridView().template begin<0>();
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = eItBegin; eIt != eItEnd; ++eIt)
+ {
+ int idx = problem_.variables().index(*eIt);
+ CellData& cellData = problem_.variables().cellData(idx);
+
+ storePressureSolution(*eIt);
+ (*pC)[idx] = cellData.capillaryPressure();
+
+ if (pressureType_ == pw)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(nPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] + cellData.capillaryPressure();
+ }
+
+ if (pressureType_ == pn)
+ {
+ (*pressure)[idx] = this->pressure()[idx][0]
+ - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ (*potentialSecond)[idx] = cellData.pressure(wPhaseIdx);
+ (*pressureSecond)[idx] = (*pressure)[idx] - cellData.capillaryPressure();
+ }
+ }
+
if (pressureType_ == pw)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[wPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ writer.attachCellData(*pressure, "wetting pressure");
+ writer.attachCellData(*pressureSecond, "nonwetting pressure");
+ writer.attachCellData(*potentialSecond, "nonwetting potential");
}
if (pressureType_ == pn)
{
- pressure[idx] = this->pressure()[idx][0]
- - density_[nPhaseIdx] * (gravity_ * (problem_.bboxMax() - eIt->geometry().center()));
+ writer.attachCellData(*pressure, "nonwetting pressure");
+ writer.attachCellData(*pressureSecond, "wetting pressure");
+ writer.attachCellData(*potentialSecond, "wetting potential");
}
- pC[idx] = problem_.variables().cellData(idx).capillaryPressure();
- }
- if (pressureType_ == pw)
- {
- writer.attachCellData(pressure, "wetting pressure");
- writer.attachCellData(potential, "wetting potential");
- }
-
- if (pressureType_ == pn)
- {
- writer.attachCellData(pressure, "nonwetting pressure");
- writer.attachCellData(potential, "nonwetting potential");
+ writer.attachCellData(*pC, "capillary pressure");
}
- writer.attachCellData(pC, "capillary pressure");
-
return;
}
@@ -365,17 +406,12 @@ public:
ErrorTermLowerBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermLowerBound);
ErrorTermUpperBound_ = GET_PARAM(TypeTag, Scalar, ErrorTermUpperBound);
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
private:
@@ -441,6 +477,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$ 0 = S_w\f$, \f$ 1 = S_n\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundvelocity2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundvelocity2p.hh
index 944f51a21a..af0c7a4f24 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundvelocity2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundvelocity2p.hh
@@ -139,17 +139,13 @@ public:
FVMPFAO2PFABoundVelocity2P(Problem& problem) :
ParentType(problem), problem_(problem), gravity_(problem.gravity())
{
- ElementIterator element = problem_.gridView().template begin<0>();
- FluidState fluidState;
- fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
- fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
- fluidState.setTemperature(problem_.temperature(*element));
- fluidState.setSaturation(wPhaseIdx, 1.);
- fluidState.setSaturation(nPhaseIdx, 0.);
- density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
- density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
- viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
- viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
+ density_[wPhaseIdx] = 0.;
+ density_[nPhaseIdx] = 0.;
+ viscosity_[wPhaseIdx] = 0.;
+ viscosity_[nPhaseIdx] = 0.;
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
//Calculates the velocities at all cell-cell interfaces.
@@ -163,6 +159,18 @@ public:
{
ParentType::initialize();
+ ElementIterator element = problem_.gridView().template begin<0>();
+ FluidState fluidState;
+ fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+ fluidState.setTemperature(problem_.temperature(*element));
+ fluidState.setSaturation(wPhaseIdx, 1.);
+ fluidState.setSaturation(nPhaseIdx, 0.);
+ density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+ density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+ viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+ viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+
calculateVelocity();
return;
@@ -185,6 +193,8 @@ public:
/*! \brief Adds velocity output to the output file
*
* Adds the phase velocities or a total velocity (depending on the formulation) to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -195,68 +205,70 @@ public:
{
ParentType::addOutputVtkFields(writer);
- Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
- Scalar, dim>(problem_.gridView().size(0)));
- Dune::BlockVector < DimVector > &velocityNonwetting =
- *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
-
- // compute update vector
- ElementIterator eItEnd = problem_.gridView().template end<0>();
- for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ if (vtkOutputLevel_ > 0)
{
- // cell index
- int globalIdx = problem_.variables().index(*eIt);
+ Dune::BlockVector < DimVector > &velocityWetting = *(writer.template allocateManagedBuffer<
+ Scalar, dim>(problem_.gridView().size(0)));
+ Dune::BlockVector < DimVector > &velocityNonwetting =
+ *(writer.template allocateManagedBuffer<Scalar, dim>(problem_.gridView().size(0)));
+
+ // compute update vector
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ // cell index
+ int globalIdx = problem_.variables().index(*eIt);
- CellData& cellData = problem_.variables().cellData(globalIdx);
+ CellData& cellData = problem_.variables().cellData(globalIdx);
- Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
- Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
- // run through all intersections with neighbors and boundary
- IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
- for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
- {
- int isIndex = isIt->indexInInside();
-
- fluxW[isIndex] = isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
- fluxNW[isIndex] =
- isIt->geometry().volume()
- * (isIt->centerUnitOuterNormal()
- * cellData.fluxData().velocity(nPhaseIdx, isIndex));
- }
+ Dune::FieldVector < Scalar, 2 * dim > fluxW(0);
+ Dune::FieldVector < Scalar, 2 * dim > fluxNW(0);
+ // run through all intersections with neighbors and boundary
+ IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
+ for (IntersectionIterator isIt = problem_.gridView().ibegin(*eIt); isIt != isItEnd; ++isIt)
+ {
+ int isIndex = isIt->indexInInside();
+
+ fluxW[isIndex] = isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal() * cellData.fluxData().velocity(wPhaseIdx, isIndex));
+ fluxNW[isIndex] =
+ isIt->geometry().volume()
+ * (isIt->centerUnitOuterNormal()
+ * cellData.fluxData().velocity(nPhaseIdx, isIndex));
+ }
- DimVector refVelocity(0);
- for (int i = 0; i < dim; i++)
+ DimVector refVelocity(0);
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxW[2*i + 1] - fluxW[2*i]);
- const DimVector& localPos =
- ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
+ const DimVector& localPos =
+ ReferenceElementContainer::general(eIt->geometry().type()).position(0, 0);
- // get the transposed Jacobian of the element mapping
- const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
+ // get the transposed Jacobian of the element mapping
+ const DimMatrix& jacobianT = eIt->geometry().jacobianTransposed(localPos);
- // calculate the element velocity by the Piola transformation
- DimVector elementVelocity(0);
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ DimVector elementVelocity(0);
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocityWetting[globalIdx] = elementVelocity;
+ velocityWetting[globalIdx] = elementVelocity;
- refVelocity = 0;
- for (int i = 0; i < dim; i++)
+ refVelocity = 0;
+ for (int i = 0; i < dim; i++)
refVelocity[i] = 0.5 * (fluxNW[2*i + 1] - fluxNW[2*i]);
- // calculate the element velocity by the Piola transformation
- elementVelocity = 0;
- jacobianT.umtv(refVelocity, elementVelocity);
- elementVelocity /= eIt->geometry().integrationElement(localPos);
+ // calculate the element velocity by the Piola transformation
+ elementVelocity = 0;
+ jacobianT.umtv(refVelocity, elementVelocity);
+ elementVelocity /= eIt->geometry().integrationElement(localPos);
- velocityNonwetting[globalIdx] = elementVelocity;
- }
-
- writer.attachCellData(velocityWetting, "wetting-velocity", dim);
- writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ velocityNonwetting[globalIdx] = elementVelocity;
+ }
+ writer.attachCellData(velocityWetting, "wetting-velocity", dim);
+ writer.attachCellData(velocityNonwetting, "non-wetting-velocity", dim);
+ }
return;
}
@@ -267,6 +279,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
static const Scalar threshold_ = 1e-15;
static const int velocityType_ = GET_PROP_VALUE(TypeTag, VelocityFormulation); //!< gives kind of velocity used (\f$ 0 = v_w\f$, \f$ 1 = v_n\f$, \f$ 2 = v_t\f$)
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$ 0 = p_w\f$, \f$ 1 = p_n\f$, \f$ 2 = p_{global}\f$)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
index c473828d60..9e0b123aa9 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
@@ -218,6 +218,8 @@ public:
/*! \brief Adds pressure output to the output file
*
* Adds the global pressure as well as the capillary pressure to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -232,16 +234,18 @@ public:
writer.attachCellData(*pressure, "global pressure");
- // output phase-dependent stuff
- ScalarSolutionType *pC = writer.allocateManagedBuffer (problem_.gridView().size(0));
- int size = problem_.gridView().size(0);
- for (int i = 0; i < size; i++)
+ if (vtkOutputLevel_ > 0)
{
- CellData& cellData = problem_.variables().cellData(i);
- (*pC)[i] = cellData.capillaryPressure();
+ // output phase-dependent stuff
+ ScalarSolutionType *pC = writer.allocateManagedBuffer (problem_.gridView().size(0));
+ int size = problem_.gridView().size(0);
+ for (int i = 0; i < size; i++)
+ {
+ CellData& cellData = problem_.variables().cellData(i);
+ (*pC)[i] = cellData.capillaryPressure();
+ }
+ writer.attachCellData(*pC, "capillary pressure");
}
- writer.attachCellData(*pC, "capillary pressure");
-
return;
}
@@ -260,6 +264,8 @@ public:
{
DUNE_THROW(Dune::NotImplemented, "MPFA method only implemented for 2-d!");
}
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
private:
@@ -267,6 +273,8 @@ private:
Scalar density_[numPhases];
Scalar viscosity_[numPhases];
+ int vtkOutputLevel_;
+
protected:
static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation); //!< gives kind of pressure used (\f$p_w\f$, \f$p_n\f$, \f$p_{global}\f$)
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation); //!< gives kind of saturation used (\f$S_w\f$, \f$S_n\f$)
diff --git a/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh b/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
index a32297ecd7..dd297f5cd4 100644
--- a/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
+++ b/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
@@ -266,6 +266,9 @@ public:
/*! \brief Adds saturation output to the output file
*
* Adds the phase saturation to the output. If the velocity is calculated in the transport model it is also added to the output.
+ * If the VtkOutputLevel is equal to zero (default) only primary variables are written,
+ * if it is larger than zero also secondary variables are written.
+ *
*
* \tparam MultiWriter Class defining the output writer
* \param writer The output writer (usually a <tt>VTKMultiWriter</tt> object)
@@ -275,21 +278,45 @@ public:
void addOutputVtkFields(MultiWriter &writer)
{
int size = problem_.gridView().size(0);
- TransportSolutionType *saturationW = writer.allocateManagedBuffer(size);
- TransportSolutionType *saturationN = writer.allocateManagedBuffer(size);
+ TransportSolutionType *saturation = writer.allocateManagedBuffer(size);
+ TransportSolutionType *saturationSecond = 0;
+
+ if (vtkOutputLevel_ > 0)
+ {
+ saturationSecond = writer.allocateManagedBuffer(size);
+ }
for (int i = 0; i < size; i++)
{
CellData& cellData = problem_.variables().cellData(i);
- (*saturationW)[i] = cellData.saturation(wPhaseIdx);
- (*saturationN)[i] = cellData.saturation(nPhaseIdx);
+
+ if (saturationType_ == Sw)
+ {
+ (*saturation)[i] = cellData.saturation(wPhaseIdx);
+ if (vtkOutputLevel_ > 0)
+ {
+ (*saturationSecond)[i] = cellData.saturation(nPhaseIdx);
+ }
+ }
+ else if (saturationType_ == Sn)
+ {
+ (*saturation)[i] = cellData.saturation(nPhaseIdx);
+ if (vtkOutputLevel_ > 0)
+ {
+ (*saturationSecond)[i] = cellData.saturation(wPhaseIdx);
+ }
+ }
}
- writer.attachCellData(*saturationW, "wetting saturation");
- writer.attachCellData(*saturationN, "nonwetting saturation");
+ writer.attachCellData(*saturation, "wetting saturation");
+
+ if (vtkOutputLevel_ > 0)
+ {
+ writer.attachCellData(*saturationSecond, "nonwetting saturation");
+ }
if (velocity().calculateVelocityInTransport())
- velocity().addOutputVtkFields(writer);
+ velocity().addOutputVtkFields(writer);
return;
}
@@ -371,6 +398,8 @@ public:
capillaryFlux_ = new CapillaryFlux(problem);
gravityFlux_ = new GravityFlux(problem);
velocity_ = new Velocity(problem);
+
+ vtkOutputLevel_ = GET_PARAM_FROM_GROUP(TypeTag, int, Vtk, OutputLevel);
}
//! Destructor
@@ -387,6 +416,8 @@ private:
CapillaryFlux* capillaryFlux_;
GravityFlux* gravityFlux_;
+ int vtkOutputLevel_;
+
static const bool compressibility_ = GET_PROP_VALUE(TypeTag, EnableCompressibility);
static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation);
static const int velocityType_ = GET_PROP_VALUE(TypeTag, VelocityFormulation);
diff --git a/test/decoupled/2p/CMakeLists.txt b/test/decoupled/2p/CMakeLists.txt
index 32a0c25b84..73869c1042 100644
--- a/test/decoupled/2p/CMakeLists.txt
+++ b/test/decoupled/2p/CMakeLists.txt
@@ -5,6 +5,9 @@ target_link_libraries("test_impes" ${DumuxLinkLibraries})
add_executable("test_transport" test_transport.cc)
target_link_libraries("test_transport" ${DumuxLinkLibraries})
+add_executable("test_mpfa2p" test_mpfa2p.cc)
+target_link_libraries("test_mpfa2p" ${DumuxLinkLibraries})
+
# add required libraries and includes to the build flags
link_directories(${DumuxLinkDirectories})
include_directories(${DumuxIncludeDirectories})
diff --git a/test/decoupled/2p/test_impes-reference.vtu b/test/decoupled/2p/test_impes-reference.vtu
index ce4a490d76..e62d98fdd1 100644
--- a/test/decoupled/2p/test_impes-reference.vtu
+++ b/test/decoupled/2p/test_impes-reference.vtu
@@ -2,7 +2,7 @@
<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">
<UnstructuredGrid>
<Piece NumberOfCells="180" NumberOfPoints="217">
- <CellData Scalars="wetting pressure" Vectors="wetting-velocity">
+ <CellData Scalars="wetting pressure">
<DataArray type="Float32" Name="wetting pressure" NumberOfComponents="1" format="ascii">
200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 199999
199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
@@ -20,40 +20,6 @@
200000 200000 200000 200000 200000 199999 199999 199999 199999 199999 199999 199999
199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
</DataArray>
- <DataArray type="Float32" Name="nonwetting pressure" NumberOfComponents="1" format="ascii">
- 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 200000 200000 200000 200000 200000 200000
- 200000 200000 200000 200000 200000 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 200000 200000 200000 200000 200000 200000
- 200000 200000 200000 200000 200000 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 200000 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 200000 200000 200000 200000 200000 200000
- 200000 200000 200000 200000 200000 199999 199999 199999 199999 199999 199999 199999
- 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999 199999
- </DataArray>
- <DataArray type="Float32" Name="capillary pressure" NumberOfComponents="1" format="ascii">
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- </DataArray>
<DataArray type="Float32" Name="wetting saturation" NumberOfComponents="1" format="ascii">
0.743888 0.661224 0.593957 0.300918 0.200013 0.2 0.2 0.2 0.2 0.2 0.2 0.2
0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2
@@ -71,117 +37,6 @@
0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2
0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2
</DataArray>
- <DataArray type="Float32" Name="nonwetting saturation" NumberOfComponents="1" format="ascii">
- 0.256112 0.338776 0.406043 0.699082 0.799987 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.256112 0.338776 0.406043 0.699082 0.799987 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.256112 0.338776 0.406043 0.699082 0.799987 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.256112 0.338776 0.406043 0.699082 0.799987 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.256112 0.338776 0.406043 0.699082 0.799987 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.256112 0.338776 0.406043 0.699082 0.799987 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
- </DataArray>
- <DataArray type="Float32" Name="wetting-velocity" NumberOfComponents="3" format="ascii">
- 2.99652e-07 0 0 2.90182e-07 6.46435e-16 0 2.43301e-07 1.21655e-16 0 1.02814e-07 6.35629e-20 0
- 4.30793e-11 1.34952e-38 0 2.78216e-30 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 2.99652e-07 -1.96372e-16 0 2.90182e-07 7.45886e-16 0
- 2.43301e-07 1.21655e-16 0 1.02814e-07 1.27126e-19 0 4.30793e-11 1.88932e-38 0 2.78216e-30 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 2.99652e-07 9.81862e-17 0 2.90182e-07 -5.46983e-16 0 2.43301e-07 -7.29928e-17 0 1.02814e-07 -1.27126e-19 0
- 4.30793e-11 0 0 2.78215e-30 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 2.99652e-07 8.83676e-16 0 2.90182e-07 -1.14369e-15 0
- 2.43301e-07 -2.6764e-16 0 1.02814e-07 -2.54252e-19 0 4.30793e-11 -1.61942e-38 0 2.78216e-30 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 2.99652e-07 -1.96372e-16 0 2.90182e-07 2.98355e-16 0 2.43301e-07 4.86618e-17 0 1.02814e-07 4.13159e-19 0
- 4.30793e-11 -1.61942e-38 0 2.78216e-30 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 2.99652e-07 -7.85489e-16 0 2.90182e-07 7.95612e-16 0
- 2.43301e-07 2.43309e-16 0 1.02814e-07 4.76722e-19 0 4.30793e-11 0 0 2.78215e-30 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- 0 0 0 0 0 0 0 0 0 0 0 0
- </DataArray>
- <DataArray type="Float32" Name="non-wetting-velocity" NumberOfComponents="3" format="ascii">
- 3.47586e-10 0 0 9.81759e-09 4.35618e-17 0 5.66992e-08 5.59075e-17 0 1.97186e-07 2.2126e-16 0
- 2.99957e-07 7.27593e-16 0 3e-07 4.36557e-16 0 3e-07 2.91038e-16 0 3e-07 -4.36557e-16 0
- 3e-07 -1.45519e-16 0 3e-07 -5.82077e-16 0 3e-07 -5.82077e-16 0 3e-07 -1.45519e-16 0
- 3e-07 1.45519e-16 0 3e-07 5.82077e-16 0 3e-07 -4.36557e-16 0 3e-07 2.91038e-16 0
- 3e-07 -1.45519e-16 0 3e-07 -1.45519e-16 0 3e-07 4.36557e-16 0 3e-07 2.91038e-16 0
- 3e-07 1.45519e-16 0 3e-07 2.91038e-16 0 3e-07 0 0 3e-07 -1.45519e-16 0
- 3e-07 -1.45519e-16 0 3e-07 2.91038e-16 0 3e-07 0 0 3e-07 2.91038e-16 0
- 3e-07 4.36557e-16 0 3e-07 7.27596e-16 0 3.47586e-10 -4.56099e-19 0 9.81759e-09 5.02636e-17 0
- 5.66992e-08 5.59075e-17 0 1.97186e-07 4.42519e-16 0 2.99957e-07 1.01863e-15 0 3e-07 -4.36557e-16 0
- 3e-07 -2.91038e-16 0 3e-07 0 0 3e-07 -4.36557e-16 0 3e-07 -4.36557e-16 0
- 3e-07 2.91038e-16 0 3e-07 -2.91038e-16 0 3e-07 5.82077e-16 0 3e-07 7.27596e-16 0
- 3e-07 -1.45519e-16 0 3e-07 0 0 3e-07 2.91038e-16 0 3e-07 -4.36557e-16 0
- 3e-07 2.91038e-16 0 3e-07 2.91038e-16 0 3e-07 2.91038e-16 0 3e-07 0 0
- 3e-07 1.45519e-16 0 3e-07 2.91038e-16 0 3e-07 -4.36557e-16 0 3e-07 -4.36557e-16 0
- 3e-07 -1.45519e-16 0 3e-07 -1.45519e-16 0 3e-07 0 0 3e-07 2.91038e-16 0
- 3.47586e-10 2.28049e-19 0 9.81759e-09 -3.686e-17 0 5.66992e-08 -3.35445e-17 0 1.97186e-07 -4.42519e-16 0
- 2.99957e-07 -2.91037e-16 0 3e-07 -1.16415e-15 0 3e-07 2.91038e-16 0 3e-07 -1.45519e-16 0
- 3e-07 -1.45519e-16 0 3e-07 7.27596e-16 0 3e-07 4.36557e-16 0 3e-07 0 0
- 3e-07 7.27596e-16 0 3e-07 -1.45519e-16 0 3e-07 0 0 3e-07 -2.91038e-16 0
- 3e-07 4.36557e-16 0 3e-07 -2.91038e-16 0 3e-07 0 0 3e-07 -4.36557e-16 0
- 3e-07 1.45519e-16 0 3e-07 -2.91038e-16 0 3e-07 4.36557e-16 0 3e-07 4.36557e-16 0
- 3e-07 1.45519e-16 0 3e-07 -4.36557e-16 0 3e-07 -1.45519e-16 0 3e-07 1.45519e-16 0
- 3e-07 -2.91038e-16 0 3e-07 -1.30967e-15 0 3.47586e-10 2.05244e-18 0 9.81759e-09 -7.70709e-17 0
- 5.66992e-08 -1.22996e-16 0 1.97186e-07 -8.85038e-16 0 2.99957e-07 -8.73111e-16 0 3e-07 -4.36557e-16 0
- 3e-07 2.91038e-16 0 3e-07 -4.36557e-16 0 3e-07 1.45519e-16 0 3e-07 4.36557e-16 0
- 3e-07 -4.36557e-16 0 3e-07 2.91038e-16 0 3e-07 0 0 3e-07 -2.91038e-16 0
- 3e-07 -4.36557e-16 0 3e-07 2.91038e-16 0 3e-07 -4.36557e-16 0 3e-07 1.45519e-16 0
- 3e-07 1.45519e-16 0 3e-07 -2.91038e-16 0 3e-07 -5.82077e-16 0 3e-07 4.36557e-16 0
- 3e-07 -1.45519e-16 0 3e-07 -4.36557e-16 0 3e-07 1.45519e-16 0 3e-07 5.82077e-16 0
- 3e-07 -1.45519e-16 0 3e-07 4.36557e-16 0 3e-07 0 0 3e-07 -7.27596e-16 0
- 3.47586e-10 -4.56099e-19 0 9.81759e-09 2.01054e-17 0 5.66992e-08 2.2363e-17 0 1.97186e-07 1.43819e-15 0
- 2.99957e-07 -8.73111e-16 0 3e-07 0 0 3e-07 -8.73115e-16 0 3e-07 0 0
- 3e-07 2.91038e-16 0 3e-07 -1.45519e-16 0 3e-07 -1.45519e-16 0 3e-07 1.45519e-16 0
- 3e-07 -2.91038e-16 0 3e-07 -5.82077e-16 0 3e-07 -5.82077e-16 0 3e-07 2.91038e-16 0
- 3e-07 -1.45519e-16 0 3e-07 2.91038e-16 0 3e-07 0 0 3e-07 2.91038e-16 0
- 3e-07 -5.82077e-16 0 3e-07 5.82077e-16 0 3e-07 -7.27596e-16 0 3e-07 -2.91038e-16 0
- 3e-07 -5.82077e-16 0 3e-07 0 0 3e-07 -1.45519e-16 0 3e-07 -4.36557e-16 0
- 3e-07 -2.91038e-16 0 3e-07 5.82077e-16 0 3.47586e-10 -1.82439e-18 0 9.81759e-09 5.36145e-17 0
- 5.66992e-08 1.11815e-16 0 1.97186e-07 1.65945e-15 0 2.99957e-07 -5.82074e-16 0 3e-07 1.45519e-16 0
- 3e-07 -2.91038e-16 0 3e-07 -1.45519e-16 0 3e-07 2.91038e-16 0 3e-07 0 0
- 3e-07 -1.45519e-16 0 3e-07 0 0 3e-07 0 0 3e-07 -5.82077e-16 0
- 3e-07 -4.36557e-16 0 3e-07 0 0 3e-07 2.91038e-16 0 3e-07 1.45519e-16 0
- 3e-07 0 0 3e-07 1.45519e-16 0 3e-07 0 0 3e-07 1.45519e-16 0
- 3e-07 -2.91038e-16 0 3e-07 1.45519e-16 0 3e-07 -2.91038e-16 0 3e-07 -2.91038e-16 0
- 3e-07 0 0 3e-07 -2.91038e-16 0 3e-07 -1.45519e-16 0 3e-07 4.36557e-16 0
- </DataArray>
</CellData>
<Points>
<DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="ascii">
diff --git a/test/decoupled/2p/test_impes.input b/test/decoupled/2p/test_impes.input
index 7812cf2c5c..4ca8ddf7a8 100644
--- a/test/decoupled/2p/test_impes.input
+++ b/test/decoupled/2p/test_impes.input
@@ -19,6 +19,9 @@ NumberOfCellsY = 6 # [-] resolution in y-direction
UpperRightX = 300 # [m] length of the domain
UpperRightY = 60 # [m] height of the domain
+[Vtk]
+#OutputLevel = 1 # 0 -> only primary variables (default), 1 -> also secondary variables
+
###############################################################
# Simulation restart
#
diff --git a/test/decoupled/2p/test_mpfa2p-reference.vtu b/test/decoupled/2p/test_mpfa2p-reference.vtu
new file mode 100644
index 0000000000..edd6c63da9
--- /dev/null
+++ b/test/decoupled/2p/test_mpfa2p-reference.vtu
@@ -0,0 +1,126 @@
+<?xml version="1.0"?>
+<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">
+ <UnstructuredGrid>
+ <Piece NumberOfCells="100" NumberOfPoints="121">
+ <CellData Scalars="wetting potential">
+ <DataArray type="Float32" Name="wetting potential" NumberOfComponents="1" format="ascii">
+ 100181 100435 100634 100798 100882 100756 100608 100458 100305 100140 100209 100448
+ 100643 100819 100934 100761 100609 100458 100308 100169 100285 100472 100657 100854
+ 101059 101172 101363 101820 101030 100235 100355 100494 100668 100886 101236 101588
+ 101927 103026 101587 100212 100285 101551 102377 101942 101450 101595 101721 101386
+ 100798 100174 100189 102476 103891 102964 101577 101602 101648 101151 100563 100150
+ 100131 100475 101013 103773 104038 104017 103200 100811 100424 100127 100115 100354
+ 100612 104806 106112 106110 104814 100581 100341 100112 100110 100324 100528 100884
+ 102171 102163 100861 100502 100311 100107 100109 100320 100528 100810 101241 101229
+ 100784 100504 100307 100105
+ </DataArray>
+ <DataArray type="Float32" Name="wetting saturation" NumberOfComponents="1" format="ascii">
+ 1 1 1 1 1 1 1 1 1 1 0.999944 1
+ 1 1 1 1 1 1 1 0.993708 0.964786 1 1 1
+ 1 1 0.99981 0.998789 0.999723 0.878958 0.809359 0.999943 1 1 0.998122 0.915612
+ 0.429419 0.240657 0.287893 0.776332 0.736877 0.999399 0.99824 0.999625 0.923737 0.999991 0.963424 0.56438
+ 0.73134 0.988434 0.708671 0.214496 0.207329 0.514697 0.823319 0.999915 0.976573 0.626967 0.844354 0.997189
+ 0.961564 0.617044 0.460227 0.996759 0.99406 0.99407 0.99628 0.616475 0.869689 0.998433 0.998359 0.887138
+ 0.62054 0.155622 0.157506 0.157504 0.155841 0.614653 0.898365 0.999191 0.999999 0.988317 0.807678 0.354803
+ 0.262015 0.261169 0.34195 0.803853 0.988146 1 1 0.99998 0.973048 0.850869 0.593308 0.593036
+ 0.847216 0.969937 0.999975 1
+ </DataArray>
+ </CellData>
+ <Points>
+ <DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="ascii">
+ 0 0 0 2 0 0 0 1 0 2 1 0
+ 4 0 0 4 1 0 6 0 0 6 1 0
+ 8 0 0 8 1 0 10 0 0 10 1 0
+ 12 0 0 12 1 0 14 0 0 14 1 0
+ 16 0 0 16 1 0 18 0 0 18 1 0
+ 20 0 0 20 1 0 0 2 0 2 2 0
+ 4 2 0 6 2 0 8 2 0 10 2 0
+ 12 2 0 14 2 0 16 2 0 18 2 0
+ 20 2 0 0 3 0 2 3 0 4 3 0
+ 6 3 0 8 3 0 10 3 0 12 3 0
+ 14 3 0 16 3 0 18 3 0 20 3 0
+ 0 4 0 2 4 0 4 4 0 6 4 0
+ 8 4 0 10 4 0 12 4 0 14 4 0
+ 16 4 0 18 4 0 20 4 0 0 5 0
+ 2 5 0 4 5 0 6 5 0 8 5 0
+ 10 5 0 12 5 0 14 5 0 16 5 0
+ 18 5 0 20 5 0 0 6 0 2 6 0
+ 4 6 0 6 6 0 8 6 0 10 6 0
+ 12 6 0 14 6 0 16 6 0 18 6 0
+ 20 6 0 0 7 0 2 7 0 4 7 0
+ 6 7 0 8 7 0 10 7 0 12 7 0
+ 14 7 0 16 7 0 18 7 0 20 7 0
+ 0 8 0 2 8 0 4 8 0 6 8 0
+ 8 8 0 10 8 0 12 8 0 14 8 0
+ 16 8 0 18 8 0 20 8 0 0 9 0
+ 2 9 0 4 9 0 6 9 0 8 9 0
+ 10 9 0 12 9 0 14 9 0 16 9 0
+ 18 9 0 20 9 0 0 10 0 2 10 0
+ 4 10 0 6 10 0 8 10 0 10 10 0
+ 12 10 0 14 10 0 16 10 0 18 10 0
+ 20 10 0
+ </DataArray>
+ </Points>
+ <Cells>
+ <DataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="ascii">
+ 0 1 3 2 1 4 5 3 4 6 7 5
+ 6 8 9 7 8 10 11 9 10 12 13 11
+ 12 14 15 13 14 16 17 15 16 18 19 17
+ 18 20 21 19 2 3 23 22 3 5 24 23
+ 5 7 25 24 7 9 26 25 9 11 27 26
+ 11 13 28 27 13 15 29 28 15 17 30 29
+ 17 19 31 30 19 21 32 31 22 23 34 33
+ 23 24 35 34 24 25 36 35 25 26 37 36
+ 26 27 38 37 27 28 39 38 28 29 40 39
+ 29 30 41 40 30 31 42 41 31 32 43 42
+ 33 34 45 44 34 35 46 45 35 36 47 46
+ 36 37 48 47 37 38 49 48 38 39 50 49
+ 39 40 51 50 40 41 52 51 41 42 53 52
+ 42 43 54 53 44 45 56 55 45 46 57 56
+ 46 47 58 57 47 48 59 58 48 49 60 59
+ 49 50 61 60 50 51 62 61 51 52 63 62
+ 52 53 64 63 53 54 65 64 55 56 67 66
+ 56 57 68 67 57 58 69 68 58 59 70 69
+ 59 60 71 70 60 61 72 71 61 62 73 72
+ 62 63 74 73 63 64 75 74 64 65 76 75
+ 66 67 78 77 67 68 79 78 68 69 80 79
+ 69 70 81 80 70 71 82 81 71 72 83 82
+ 72 73 84 83 73 74 85 84 74 75 86 85
+ 75 76 87 86 77 78 89 88 78 79 90 89
+ 79 80 91 90 80 81 92 91 81 82 93 92
+ 82 83 94 93 83 84 95 94 84 85 96 95
+ 85 86 97 96 86 87 98 97 88 89 100 99
+ 89 90 101 100 90 91 102 101 91 92 103 102
+ 92 93 104 103 93 94 105 104 94 95 106 105
+ 95 96 107 106 96 97 108 107 97 98 109 108
+ 99 100 111 110 100 101 112 111 101 102 113 112
+ 102 103 114 113 103 104 115 114 104 105 116 115
+ 105 106 117 116 106 107 118 117 107 108 119 118
+ 108 109 120 119
+ </DataArray>
+ <DataArray type="Int32" Name="offsets" NumberOfComponents="1" format="ascii">
+ 4 8 12 16 20 24 28 32 36 40 44 48
+ 52 56 60 64 68 72 76 80 84 88 92 96
+ 100 104 108 112 116 120 124 128 132 136 140 144
+ 148 152 156 160 164 168 172 176 180 184 188 192
+ 196 200 204 208 212 216 220 224 228 232 236 240
+ 244 248 252 256 260 264 268 272 276 280 284 288
+ 292 296 300 304 308 312 316 320 324 328 332 336
+ 340 344 348 352 356 360 364 368 372 376 380 384
+ 388 392 396 400
+ </DataArray>
+ <DataArray type="UInt8" Name="types" NumberOfComponents="1" format="ascii">
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9 9 9 9 9 9 9 9 9
+ 9 9 9 9
+ </DataArray>
+ </Cells>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
diff --git a/test/decoupled/2p/test_mpfa2p.input b/test/decoupled/2p/test_mpfa2p.input
index 2603d92ec5..ab99ef9f05 100644
--- a/test/decoupled/2p/test_mpfa2p.input
+++ b/test/decoupled/2p/test_mpfa2p.input
@@ -19,6 +19,9 @@ NumberOfCellsY = 10# [-] level 0 resolution in y-direction
UpperRightX = 20# [m] dimension of the grid
UpperRightY = 10# [m] dimension of the grid # 2D
+[Vtk]
+#OutputLevel = 1 # 0 -> only primary variables (default), 1 -> also secondary variables
+
############################################################
# Optional arguments
############################################################
diff --git a/test/decoupled/2p/test_mpfa2pproblem.hh b/test/decoupled/2p/test_mpfa2pproblem.hh
index 638636654a..20c10a54b7 100644
--- a/test/decoupled/2p/test_mpfa2pproblem.hh
+++ b/test/decoupled/2p/test_mpfa2pproblem.hh
@@ -123,8 +123,19 @@ NEW_TYPE_TAG(MPFALAdaptiveTwoPTestProblem, INHERITS_FROM(FVMPFALPressureTwoPAdap
* and free outflow for saturation are set. The domain is heterogeneous with a backround material and three lenses.
*
* To run the simulation execute the following line in shell:
- * <tt>./test_mpfa2p -parameterFile ./test_mpfa2p.input</tt>,
- * where the arguments define the parameter file..
+ *
+ * <tt>./test_mpfa2p</tt>,
+ *
+ * Additioinally, the numerical model can be switched by executing with the parameter "ModelType":
+ *
+ * <tt>./test_mpfa2p --ModelType=...</tt>,
+ *
+ * where ModelType can be:
+ * - FV (standard finite volume)
+ * - FVAdaptive (adaptive finite volume)
+ * - MPFAO (MPFA o-method)
+ * - MPFAL (MPFA l-method)
+ * - MPFALAdaptive (adaptive MPFA l-method)
*/
template<class TypeTag>
class MPFATwoPTestProblem: public IMPESProblem2P<TypeTag>
@@ -230,7 +241,7 @@ const char *name() const
}
else
{
- return "testmpfa2p";
+ return "test_mpfa2p";
}
}
diff --git a/test/decoupled/2p/test_transport-reference.vtu b/test/decoupled/2p/test_transport-reference.vtu
index d6547d6a9d..aa7a772856 100644
--- a/test/decoupled/2p/test_transport-reference.vtu
+++ b/test/decoupled/2p/test_transport-reference.vtu
@@ -14,17 +14,6 @@
1 0 0 0 0 0 1 1 1 1 1 0
0 0 0 0
</DataArray>
- <DataArray type="Float32" Name="nonwetting saturation" NumberOfComponents="1" format="ascii">
- -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16
- -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16
- 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1
- 1 1 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1
- 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1 1 1
- -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16
- -2.22045e-16 -2.22045e-16 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16
- 8.88178e-16 1 1 1 1 1 -2.22045e-16 -2.22045e-16 -2.22045e-16 -2.22045e-16 8.88178e-16 1
- 1 1 1 1
- </DataArray>
</CellData>
<Points>
<DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="ascii">
diff --git a/test/decoupled/2p/test_transport.input b/test/decoupled/2p/test_transport.input
index fd4d33c8fd..a94946cf55 100644
--- a/test/decoupled/2p/test_transport.input
+++ b/test/decoupled/2p/test_transport.input
@@ -15,6 +15,9 @@ DtInitial = 0 # [s]
[Grid]
File = ./grids/test_transport.dgf
+[Vtk]
+#OutputLevel = 1 # 0 -> only primary variables (default), 1 -> also secondary variables
+
###############################################################
# Simulation restart
#
--
GitLab