From fb178a9f52805169bbd45a02ba581b9ac57ed867 Mon Sep 17 00:00:00 2001
From: Klaus Mosthaf <klmos@env.dtu.dk>
Date: Tue, 27 Aug 2013 14:31:01 +0000
Subject: [PATCH] Model for the calculation of the effective diffusivity in the
porous medium can be changed as property. Default diffusivity model is set to
DiffusivityMillingtonQuirk (as it was before).
Reviewed by MartinS
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@11277 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
dumux/implicit/2p2c/2p2cfluxvariables.hh | 44 +++++++++++----------
dumux/implicit/2p2c/2p2cproperties.hh | 1 +
dumux/implicit/2p2c/2p2cpropertydefaults.hh | 13 +++++-
3 files changed, 35 insertions(+), 23 deletions(-)
diff --git a/dumux/implicit/2p2c/2p2cfluxvariables.hh b/dumux/implicit/2p2c/2p2cfluxvariables.hh
index 64ea218dba..2239ae8633 100644
--- a/dumux/implicit/2p2c/2p2cfluxvariables.hh
+++ b/dumux/implicit/2p2c/2p2cfluxvariables.hh
@@ -51,6 +51,7 @@ class TwoPTwoCFluxVariables : public GET_PROP_TYPE(TypeTag, BaseFluxVariables)
typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
typedef typename GET_PROP_TYPE(TypeTag, SpatialParams) SpatialParams;
+ typedef typename GET_PROP_TYPE(TypeTag, EffectiveDiffusivityModel) EffectiveDiffusivityModel;
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
@@ -72,7 +73,7 @@ class TwoPTwoCFluxVariables : public GET_PROP_TYPE(TypeTag, BaseFluxVariables)
typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
typedef typename FVElementGeometry::SubControlVolumeFace SCVFace;
-public:
+ public:
/*!
* \brief The constructor
*
@@ -89,7 +90,7 @@ public:
const int faceIdx,
const ElementVolumeVariables &elemVolVars,
const bool onBoundary = false)
- : BaseFluxVariables(problem, element, fvGeometry, faceIdx, elemVolVars, onBoundary)
+ : BaseFluxVariables(problem, element, fvGeometry, faceIdx, elemVolVars, onBoundary)
{
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
density_[phaseIdx] = Scalar(0);
@@ -100,7 +101,7 @@ public:
calculateValues_(problem, element, elemVolVars);
}
-protected:
+ protected:
void calculateValues_(const Problem &problem,
const Element &element,
const ElementVolumeVariables &elemVolVars)
@@ -111,7 +112,7 @@ protected:
idx < this->face().numFap;
idx++) // loop over adjacent vertices
{
- // index for the element volume variables
+ // index for the element volume variables
int volVarsIdx = this->face().fapIndices[idx];
for (int phaseIdx = 0; phaseIdx < numPhases; phaseIdx++)
@@ -140,7 +141,7 @@ protected:
// FE gradient at vertex idx
const DimVector &feGrad = this->face().grad[idx];
- // index for the element volume variables
+ // index for the element volume variables
int volVarsIdx = this->face().fapIndices[idx];
// the mole fraction gradient of the wetting phase
@@ -177,36 +178,37 @@ protected:
const VolumeVariables &volVarsI = elemVolVars[this->face().i];
const VolumeVariables &volVarsJ = elemVolVars[this->face().j];
+ // the effective diffusion coefficients at vertex i and j
+ Scalar diffCoeffI;
+ Scalar diffCoeffJ;
+
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
{
// make sure to only calculate diffusion coefficients
// for phases which exist in both finite volumes
- if (volVarsI.saturation(phaseIdx) <= 0 ||
- volVarsJ.saturation(phaseIdx) <= 0)
+ if (volVarsI.saturation(phaseIdx) <= 0 || volVarsJ.saturation(phaseIdx) <= 0)
{
porousDiffCoeff_[phaseIdx] = 0.0;
continue;
}
- // calculate tortuosity at the nodes i and j needed
- // for porous media diffusion coefficient
- Scalar tauI =
- 1.0/(volVarsI.porosity() * volVarsI.porosity()) *
- pow(volVarsI.porosity() * volVarsI.saturation(phaseIdx), 7.0/3);
- Scalar tauJ =
- 1.0/(volVarsJ.porosity() * volVarsJ.porosity()) *
- pow(volVarsJ.porosity() * volVarsJ.saturation(phaseIdx), 7.0/3);
- // Diffusion coefficient in the porous medium
+ diffCoeffI = EffectiveDiffusivityModel::effectiveDiffusivity(volVarsI.porosity(),
+ volVarsI.saturation(phaseIdx),
+ volVarsI.diffCoeff(phaseIdx));
+
+ diffCoeffJ = EffectiveDiffusivityModel::effectiveDiffusivity(volVarsJ.porosity(),
+ volVarsJ.saturation(phaseIdx),
+ volVarsJ.diffCoeff(phaseIdx));
// -> harmonic mean
- porousDiffCoeff_[phaseIdx] = harmonicMean(volVarsI.porosity() * volVarsI.saturation(phaseIdx) * tauI * volVarsI.diffCoeff(phaseIdx),
- volVarsJ.porosity() * volVarsJ.saturation(phaseIdx) * tauJ * volVarsJ.diffCoeff(phaseIdx));
+ porousDiffCoeff_[phaseIdx] = harmonicMean(diffCoeffI, diffCoeffJ);
}
}
-public:
+ public:
/*!
- * \brief The binary diffusion coefficient for each fluid phase.
+ * \brief The effective diffusion coefficient \f$\mathrm{[m^2/s]}\f$
+ * for each fluid phase in the porous medium.
*/
Scalar porousDiffCoeff(int phaseIdx) const
{ return porousDiffCoeff_[phaseIdx]; };
@@ -229,7 +231,7 @@ public:
const DimVector &moleFractionGrad(int phaseIdx) const
{ return moleFractionGrad_[phaseIdx]; };
-protected:
+ protected:
// mole fraction gradients
DimVector moleFractionGrad_[numPhases];
diff --git a/dumux/implicit/2p2c/2p2cproperties.hh b/dumux/implicit/2p2c/2p2cproperties.hh
index 739639c2c7..1c87bd9ea8 100644
--- a/dumux/implicit/2p2c/2p2cproperties.hh
+++ b/dumux/implicit/2p2c/2p2cproperties.hh
@@ -59,6 +59,7 @@ NEW_PROP_TAG(FluidSystem); //!< Type of the multi-component relations
NEW_PROP_TAG(MaterialLaw); //!< The material law which ought to be used (extracted from the spatial parameters)
NEW_PROP_TAG(MaterialLawParams); //!< The parameters of the material law (extracted from the spatial parameters)
+NEW_PROP_TAG(EffectiveDiffusivityModel); //!< The employed model for the computation of the effective diffusivity
NEW_PROP_TAG(ProblemEnableGravity); //!< Returns whether gravity is considered in the problem
NEW_PROP_TAG(ImplicitMassUpwindWeight); //!< The value of the upwind weight for the mass conservation equations
diff --git a/dumux/implicit/2p2c/2p2cpropertydefaults.hh b/dumux/implicit/2p2c/2p2cpropertydefaults.hh
index 6aa279c55a..15814127ab 100644
--- a/dumux/implicit/2p2c/2p2cpropertydefaults.hh
+++ b/dumux/implicit/2p2c/2p2cpropertydefaults.hh
@@ -36,6 +36,7 @@
#include "2p2clocalresidual.hh"
#include "2p2cnewtoncontroller.hh"
+#include <dumux/material/fluidmatrixinteractions/2p/diffusivitymillingtonquirk.hh>
#include <dumux/implicit/common/implicitdarcyfluxvariables.hh>
#include <dumux/material/spatialparams/implicitspatialparams.hh>
@@ -132,17 +133,25 @@ SET_SCALAR_PROP(TwoPTwoC, ImplicitMassUpwindWeight, 1.0);
SET_SCALAR_PROP(TwoPTwoC, ImplicitMobilityUpwindWeight, 1.0);
//! The indices required by the isothermal 2p2c model
-SET_PROP(TwoPTwoC, Indices)
+SET_PROP(TwoPTwoC, Indices)
{ private:
enum { Formulation = GET_PROP_VALUE(TypeTag, Formulation) };
public:
typedef TwoPTwoCIndices<TypeTag, Formulation, 0> type;
};
-//! The spatial parameters to be employed.
+//! The spatial parameters to be employed.
//! Use ImplicitSpatialParams by default.
SET_TYPE_PROP(TwoPTwoC, SpatialParams, ImplicitSpatialParams<TypeTag>);
+//! The model after Millington (1961) is used for the effective diffusivity
+SET_PROP(TwoPTwoC, EffectiveDiffusivityModel)
+{ private :
+ typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+ public:
+ typedef DiffusivityMillingtonQuirk<Scalar> type;
+};
+
// disable velocity output by default
SET_BOOL_PROP(TwoPTwoC, VtkAddVelocity, false);
--
GitLab