Commit 246f4ff0 authored by Thomas Fetzer's avatar Thomas Fetzer
Browse files

[implicit] reduced excessive line lengths, like proposed in #FS213, target line

length = 150, thanks to kristopherg

reviewed by fetzer


git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@12428 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent 8e211cee
......@@ -45,7 +45,8 @@ namespace Dumux
*
* and solves the mass continuity equation:
* \f[
\phi \frac{\partial \varrho}{\partial t} + \text{div} \left\lbrace - \varrho \frac{\textbf K}{\mu} \left( \textbf{grad}\, p -\varrho {\textbf g} \right) \right\rbrace = q,
\phi \frac{\partial \varrho}{\partial t} + \text{div} \left\lbrace
- \varrho \frac{\textbf K}{\mu} \left( \textbf{grad}\, p -\varrho {\textbf g} \right) \right\rbrace = q,
* \f]
* All equations are discretized using a vertex-centered finite volume (box)
* or cell-centered finite volume scheme as spatial
......
......@@ -45,7 +45,8 @@ struct OnePTwoCIndices
//! Component indices
static const int phaseCompIdx = phaseIdx;//!< The index of the main component of the considered phase
static const int transportCompIdx = (unsigned int)(1-phaseIdx); //!< The index of the transported (minor) component; ASSUMES phase indices of 0 and 1
//!< The index of the transported (minor) component; ASSUMES phase indices of 0 and 1
static const int transportCompIdx = (unsigned int)(1-phaseIdx);
// Equation indices
static const int conti0EqIdx = PVOffset + 0; //!< continuity equation index
......
......@@ -72,16 +72,23 @@ public:
static const int bothPhases = 2; //!< Both phases are present
// Primary variable indices
static const int pressureIdx = PVOffset + 0; //!< Index for wetting/non-wetting phase pressure (depending on the formulation) in a solution vector
static const int switchIdx = PVOffset + 1; //!< Index of either the saturation or the mass fraction of the non-wetting/wetting phase
//!< Index for wetting/non-wetting phase pressure (depending on the formulation) in a solution vector
static const int pressureIdx = PVOffset + 0;
//!< Index of either the saturation or the mass fraction of the non-wetting/wetting phase
static const int switchIdx = PVOffset + 1;
static const int pwIdx = pressureIdx; //!< Index for wetting phase pressure in a solution vector
static const int snOrXIdx = switchIdx; //!< Index of either the saturation of the non-wetting phase or the mass fraction secondary component in the only phase
//!< Index for wetting phase pressure in a solution vector
static const int pwIdx = pressureIdx;
//!< Index of either the saturation of the non-wetting phase or the mass fraction secondary component in the only phase
static const int snOrXIdx = switchIdx;
// equation indices
static const int conti0EqIdx = PVOffset; //!< Index of the mass conservation equation for the first component
static const int contiWEqIdx = conti0EqIdx + wCompIdx; //!< Index of the mass conservation equation for the primary component of the wetting phase
static const int contiNEqIdx = conti0EqIdx + nCompIdx; //!< Index of the mass conservation equation for the primary component of the non-wetting phase
//!< Index of the mass conservation equation for the first component
static const int conti0EqIdx = PVOffset;
//!< Index of the mass conservation equation for the primary component of the wetting phase
static const int contiWEqIdx = conti0EqIdx + wCompIdx;
//!< Index of the mass conservation equation for the primary component of the non-wetting phase
static const int contiNEqIdx = conti0EqIdx + nCompIdx;
};
/*!
......@@ -112,16 +119,23 @@ public:
static const int bothPhases = 3; //!< Both phases are present
// Primary variable indices
static const int pressureIdx = PVOffset + 0; //!< Index for wetting/non-wetting phase pressure (depending on the formulation) in a solution vector
static const int switchIdx = PVOffset + 1; //!< Index of either the saturation or the mass fraction of the non-wetting/wetting phase
//!< Index for wetting/non-wetting phase pressure (depending on the formulation) in a solution vector
static const int pressureIdx = PVOffset + 0;
//!< Index of either the saturation or the mass fraction of the non-wetting/wetting phase
static const int switchIdx = PVOffset + 1;
static const int pnIdx = pressureIdx; //!< Index for non-wetting phase pressure in a solution vector
static const int swOrXIdx = switchIdx; //!< Index of either the saturation of the liquid phase or the mass fraction of the secondary component in the only phase
//!< Index for non-wetting phase pressure in a solution vector
static const int pnIdx = pressureIdx;
//!< Index of either the saturation of the liquid phase or the mass fraction of the secondary component in the only phase
static const int swOrXIdx = switchIdx;
// Equation indices
static const int conti0EqIdx = PVOffset; //!< Index of the mass conservation equation for the first component
static const int contiWEqIdx = conti0EqIdx + wCompIdx; //!< Index of the mass conservation equation for the primary component of the wetting phase
static const int contiNEqIdx = conti0EqIdx + nCompIdx; //!< Index of the mass conservation equation for the primary component of the non-wetting phase
//!< Index of the mass conservation equation for the first component
static const int conti0EqIdx = PVOffset;
//!< Index of the mass conservation equation for the primary component of the wetting phase
static const int contiWEqIdx = conti0EqIdx + wCompIdx;
//!< Index of the mass conservation equation for the primary component of the non-wetting phase
static const int contiNEqIdx = conti0EqIdx + nCompIdx;
};
// \}
......
......@@ -66,7 +66,8 @@ NEW_PROP_TAG(ProblemEnableGravity); //!< Returns whether gravity is considered i
NEW_PROP_TAG(UseMoles); //!Defines whether mole (true) or mass (false) fractions are used
NEW_PROP_TAG(ImplicitMassUpwindWeight); //!< The value of the upwind weight for the mass conservation equations
NEW_PROP_TAG(ImplicitMobilityUpwindWeight); //!< Weight for the upwind mobility in the velocity calculation
NEW_PROP_TAG(ReplaceCompEqIdx); //!< The index of the total mass balance equation, if one component balance is replaced (ReplaceCompEqIdx < NumComponents)
NEW_PROP_TAG(ReplaceCompEqIdx); //!< The index of the total mass balance equation, if one component
//!< balance is replaced (ReplaceCompEqIdx < NumComponents)
NEW_PROP_TAG(VtkAddVelocity); //!< Returns whether velocity vectors are written into the VTK output
NEW_PROP_TAG(BaseFluxVariables); //! The base flux variables
NEW_PROP_TAG(SpatialParamsForchCoeff); //!< Property for the forchheimer coefficient
......
......@@ -76,7 +76,8 @@ namespace Dumux {
* default, the model uses \f$p_w\f$ and \f$S_n\f$.
* In case that only one phase (nonwetting or wetting phase) is present the second primary
* variable represents a mass fraction. The correct assignment of the second
* primary variable is performed by a phase state dependent primary variable switch. The phase state is stored for all nodes of the system. The following cases can be distinguished:
* primary variable is performed by a phase state dependent primary variable switch.
* The phase state is stored for all nodes of the system. The following cases can be distinguished:
* <ul>
* <li>
* Both phases are present: The saturation is used (either\f$S_n\f$ or \f$S_w\f$, dependent on the chosen formulation).
......
......@@ -68,9 +68,12 @@ public:
static const int switch1Idx = PVOffset + 1; //!< Index 1 of saturation or mole fraction
static const int switch2Idx = PVOffset + 2; //!< Index 2 of saturation or mole fraction
static const int pgIdx = pressureIdx; //!< Index for gas phase pressure in a solution vector
static const int sOrX1Idx = switch1Idx; //!< Index of the either the saturation of the gas phase or the mass fraction secondary component if a phase is not present
static const int sOrX2Idx = switch2Idx; //!< Index of the either the saturation of the gas phase or the mass fraction secondary component if a phase is not present
//!< Index for gas phase pressure in a solution vector
static const int pgIdx = pressureIdx;
//!< Index of the either the saturation of the gas phase or the mass fraction secondary component if a phase is not present
static const int sOrX1Idx = switch1Idx;
//!< Index of the either the saturation of the gas phase or the mass fraction secondary component if a phase is not present
static const int sOrX2Idx = switch2Idx;
// equation indices
static const int conti0EqIdx = PVOffset + wCompIdx; //!< Index of the mass conservation equation for the water component
......
......@@ -79,8 +79,11 @@ namespace Dumux
* An adaptive primary variable switch is included. The phase state is stored for all nodes
* of the system. The following cases can be distinguished:
* <ul>
* <li> All three phases are present: Primary variables are two saturations \f$(S_w\f$ and \f$S_n)\f$, and a pressure, in this case \f$p_g\f$. </li>
* <li> Only the water phase is present: Primary variables are now the mole fractions of air and contaminant in the water phase \f$(x_w^a\f$ and \f$x_w^c)\f$, as well as the gas pressure, which is, of course, in a case where only the water phase is present, just the same as the water pressure. </li>
* <li> All three phases are present: Primary variables are two saturations \f$(S_w\f$ and \f$S_n)\f$,
* and a pressure, in this case \f$p_g\f$. </li>
* <li> Only the water phase is present: Primary variables are now the mole fractions of air and
* contaminant in the water phase \f$(x_w^a\f$ and \f$x_w^c)\f$, as well as the gas pressure, which is,
* of course, in a case where only the water phase is present, just the same as the water pressure. </li>
* <li> Gas and NAPL phases are present: Primary variables \f$(S_n\f$, \f$x_g^w\f$, \f$p_g)\f$. </li>
* <li> Water and NAPL phases are present: Primary variables \f$(S_n\f$, \f$x_w^a\f$, \f$p_g)\f$. </li>
* <li> Only gas phase is present: Primary variables \f$(x_g^w\f$, \f$x_g^c\f$, \f$p_g)\f$. </li>
......
......@@ -730,7 +730,8 @@ private:
};
template <class TypeTag>
const typename ThreePThreeCVolumeVariables<TypeTag>::Scalar ThreePThreeCVolumeVariables<TypeTag>::R = Constants<typename GET_PROP_TYPE(TypeTag, Scalar)>::R;
const typename ThreePThreeCVolumeVariables<TypeTag>::Scalar ThreePThreeCVolumeVariables<TypeTag>::R
= Constants<typename GET_PROP_TYPE(TypeTag, Scalar)>::R;
} // end namespace
......
......@@ -53,8 +53,8 @@ namespace Dumux {
(\textbf{grad}\; p_\alpha - \varrho_\alpha \mbox{\bf g}) \right\}
\nonumber \\
\nonumber \\
&& - \sum\limits_\alpha \text{div} \left\{ D_{\alpha,\text{pm}}^\kappa \varrho_\alpha \frac{M^\kappa}{M_\alpha}
\textbf{grad} x^\kappa_{\alpha} \right\}
&& - \sum\limits_\alpha \text{div} \left\{ D_{\alpha,\text{pm}}^\kappa
\varrho_\alpha \frac{M^\kappa}{M_\alpha} \textbf{grad} x^\kappa_{\alpha} \right\}
- q^\kappa = 0 \qquad \forall \kappa , \; \forall \alpha
\f}
*
......@@ -88,8 +88,11 @@ namespace Dumux {
* An adaptive primary variable switch is included. The phase state is stored for all nodes
* of the system. The following cases can be distinguished:
* <ul>
* <li> All three phases are present: Primary variables are two saturations \f$(S_w\f$ and \f$S_n)\f$, a pressure, in this case \f$p_g\f$, and the temperature \f$T\f$. </li>
* <li> Only the water phase is present: Primary variables are now the mole fractions of air and contaminant in the water phase \f$(x_w^a\f$ and \f$x_w^c)\f$, as well as temperature and the gas pressure, which is, of course, in a case where only the water phase is present, just the same as the water pressure. </li>
* <li> All three phases are present: Primary variables are two saturations \f$(S_w\f$ and \f$S_n)\f$,
* a pressure, in this case \f$p_g\f$, and the temperature \f$T\f$. </li>
* <li> Only the water phase is present: Primary variables are now the mole fractions of air and
* contaminant in the water phase \f$(x_w^a\f$ and \f$x_w^c)\f$, as well as temperature and the gas pressure,
* which is, of course, in a case where only the water phase is present, just the same as the water pressure. </li>
* <li> Gas and NAPL phases are present: Primary variables \f$(S_n\f$, \f$x_g^w\f$, \f$p_g\f$, \f$T)\f$. </li>
* <li> Water and NAPL phases are present: Primary variables \f$(S_n\f$, \f$x_w^a\f$, \f$p_g\f$, \f$T)\f$. </li>
* <li> Only gas phase is present: Primary variables \f$(x_g^w\f$, \f$x_g^c\f$, \f$p_g\f$, \f$T)\f$. </li>
......
......@@ -344,7 +344,8 @@ class BoxFVElementGeometry
typedef typename LocalFiniteElement::Traits::LocalBasisType::Traits LocalBasisTraits;
typedef typename LocalBasisTraits::JacobianType ShapeJacobian;
Scalar quadrilateralArea(const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2, const GlobalPosition& p3)
Scalar quadrilateralArea(const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2,
const GlobalPosition& p3)
{
return 0.5*fabs((p3[0] - p1[0])*(p2[1] - p0[1]) - (p3[1] - p1[1])*(p2[0] - p0[0]));
}
......@@ -356,7 +357,8 @@ class BoxFVElementGeometry
c[2] = a[0]*b[1] - a[1]*b[0];
}
Scalar pyramidVolume (const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2, const GlobalPosition& p3, const GlobalPosition& p4)
Scalar pyramidVolume (const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2,
const GlobalPosition& p3, const GlobalPosition& p4)
{
GlobalPosition a(p2); a -= p0;
GlobalPosition b(p3); b -= p1;
......@@ -369,7 +371,8 @@ class BoxFVElementGeometry
return 1.0/6.0*(n*a);
}
Scalar prismVolume (const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2, const GlobalPosition& p3, const GlobalPosition& p4, const GlobalPosition& p5)
Scalar prismVolume (const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2,
const GlobalPosition& p3, const GlobalPosition& p4, const GlobalPosition& p5)
{
GlobalPosition a(p4);
for (int k = 0; k < dimWorld; ++k)
......@@ -405,7 +408,8 @@ class BoxFVElementGeometry
+ prismVolume(p0,p2,p3,p4,p6,p7);
}
void normalOfQuadrilateral3D(GlobalPosition &normal, const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2, const GlobalPosition& p3)
void normalOfQuadrilateral3D(GlobalPosition &normal, const GlobalPosition& p0, const GlobalPosition& p1, const GlobalPosition& p2,
const GlobalPosition& p3)
{
GlobalPosition a(p2);
for (int k = 0; k < dimWorld; ++k)
......@@ -592,12 +596,12 @@ public:
int numScvf; //!< number of inner-domain subcontrolvolume faces
int numNeighbors; //!< needed for compatibility with cc models
std::vector<ElementPointer> neighbors; //!< needed for compatibility with cc models
const LocalFiniteElementCache feCache_;
void updateInner(const Element& element) //!< needed for compatibility with cc models
{}
void update(const GridView& gridView, const Element& element)
{
const Geometry& geometry = element.geometry();
......@@ -772,7 +776,8 @@ public:
+ referenceElement.position(rightEdge, dim-1);
boundaryFace[bfIdx].ipLocal *= 0.25;
boundaryFace[bfIdx].area = quadrilateralArea3D(subContVol[vertInElement].global,
edgeCoordinates[rightEdge], faceCoordinates[faceIdx], edgeCoordinates[leftEdge]);
edgeCoordinates[rightEdge], faceCoordinates[faceIdx],
edgeCoordinates[leftEdge]);
break;
default:
DUNE_THROW(Dune::NotImplemented, "BoxFVElementGeometry for dim = " << dim);
......
......@@ -256,7 +256,8 @@ public:
* potentially solution dependent and requires some quantities that
* are specific to the fully-implicit method.
*
* \param values The neumann values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param values The neumann values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param element The finite element
* \param fvGeometry The finite-volume geometry
* \param intersection The intersection between element and boundary
......@@ -288,7 +289,8 @@ public:
* \brief Evaluate the boundary conditions for a neumann
* boundary segment.
*
* \param values The neumann values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param values The neumann values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param element The finite element
* \param fvGeometry The finite-volume geometry
* \param intersection The intersection between element and boundary
......@@ -313,7 +315,8 @@ public:
* \brief Evaluate the boundary conditions for a neumann
* boundary segment.
*
* \param values The neumann values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param values The neumann values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^2 \cdot s )] \f$
* \param globalPos The position of the boundary face's integration point in global coordinates
*
* For this method, the \a values parameter stores the mass flux
......@@ -338,7 +341,8 @@ public:
* potentially solution dependent and requires some quantities that
* are specific to the fully-implicit method.
*
* \param values The source and sink values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param values The source and sink values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param element The finite element
* \param fvGeometry The finite-volume geometry
* \param scvIdx The local subcontrolvolume index
......@@ -362,7 +366,8 @@ public:
* \brief Evaluate the source term for all phases within a given
* sub-control-volume.
*
* \param values The source and sink values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param values The source and sink values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param element The finite element
* \param fvGeometry The finite-volume geometry
* \param scvIdx The local subcontrolvolume index
......@@ -384,7 +389,8 @@ public:
* \brief Evaluate the source term for all phases within a given
* sub-control-volume.
*
* \param values The source and sink values for the conservation equations in units of \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param values The source and sink values for the conservation equations in units of
* \f$ [ \textnormal{unit of conserved quantity} / (m^3 \cdot s )] \f$
* \param globalPos The position of the center of the finite volume
* for which the source term ought to be
* specified in global coordinates
......
......@@ -71,7 +71,8 @@ public:
gasFlux_(fluxes, fluxVars, molarDensity);
else if ( FluidSystem::isLiquid(phaseIdx) ){
#if MACROSCALE_DIFFUSION_ONLY_GAS
return ; // in the case that only the diffusion in the gas phase is considered, the liquidFlux should not be called
return ; // in the case that only the diffusion in the gas phase is considered,
// the liquidFlux should not be called
#endif
liquidFlux_(fluxes, fluxVars, molarDensity);
}
......@@ -155,7 +156,8 @@ protected:
fluxes[compIdx] =
- xGrad *
molarDensity
* fluxVars.porousDiffCoeffG(compIdx, nCompIdx) ; // this is == 0 for nComp==comp, i.e. no diffusion of the main component of the phase
* fluxVars.porousDiffCoeffG(compIdx, nCompIdx) ; // this is == 0 for nComp==comp,
// i.e. no diffusion of the main component of the phase
}
}
}
......
......@@ -83,9 +83,12 @@ public:
* Imagine a displacement process, with the initial state being a fully saturated porous medium.
* If the residing phase is now displaced by another phase, what happens at the front of the process?
* At the very front there is one cell in which the invading phase enters and no invading (but residing) phase leaves.
* With enthalpy in the flux term and internal energy in the storage term, the difference (pv) has to be converted into temperature in order to fulfill energy conservation.
* -> A temperature peak at the front arises (if spatial discretization is sufficiently fine). This peak has a maximum value and does not increase with further refinement.
* -> Further evidence for this explanation: in a simple setting (constant parameters, few cells) the temperature peak can be correctly predicted a priori.
* With enthalpy in the flux term and internal energy in the storage term,
* the difference (pv) has to be converted into temperature in order to fulfill energy conservation.
* -> A temperature peak at the front arises (if spatial discretization is sufficiently fine).
* This peak has a maximum value and does not increase with further refinement.
* -> Further evidence for this explanation: in a simple setting (constant parameters,
* few cells) the temperature peak can be correctly predicted a priori.
* -> -> For those situations with a distinct displacement process the same quantity has to be stored and transported
* This is equivalent to neglecting volume changing work.
*
......@@ -196,7 +199,9 @@ public:
const VolumeVariables & up = elemVolVars[upIdx];
/* todo
* CAUTION: this is not exactly correct: does diffusion carry the upstream phase enthalpy? To be more precise this should be the components enthalpy. In the same vein: Counter current diffusion is not accounted for here.
* CAUTION: this is not exactly correct: does diffusion carry the upstream phase enthalpy?
* To be more precise this should be the components enthalpy.
* In the same vein: Counter current diffusion is not accounted for here.
*/
const Scalar enthalpy = up.fluidState().enthalpy(phaseIdx) ;
flux[energyEq0Idx + phaseIdx] += enthalpy * massFlux ;
......@@ -348,7 +353,8 @@ public:
// -> Therefore, this contribution needs to be added.
// -> the particle always brings the energy of the originating phase.
// -> Energy advectivly transported into a phase = the moles of a component that go into a phase * molMass * enthalpy of the component in the *originating* phase
// -> Energy advectivly transported into a phase = the moles of a component that go into a
// phase * molMass * enthalpy of the component in the *originating* phase
// The fluidsystem likes to get a fluidstate. ...
const FluidState & fluidState = volVars.fluidState();
......@@ -357,19 +363,19 @@ public:
switch (phaseIdx){
case wPhaseIdx:
source[energyEq0Idx + phaseIdx] += (componentIntoPhaseMassTransfer[wPhaseIdx][nCompIdx]
* FluidSystem::molarMass(nCompIdx)
* FluidSystem::componentEnthalpy(fluidState, nPhaseIdx, nCompIdx) );
* FluidSystem::molarMass(nCompIdx)
* FluidSystem::componentEnthalpy(fluidState, nPhaseIdx, nCompIdx) );
source[energyEq0Idx + phaseIdx] += (componentIntoPhaseMassTransfer[wPhaseIdx][wCompIdx]
* FluidSystem::molarMass(wCompIdx)
* FluidSystem::componentEnthalpy(fluidState, nPhaseIdx, wCompIdx));
* FluidSystem::molarMass(wCompIdx)
* FluidSystem::componentEnthalpy(fluidState, nPhaseIdx, wCompIdx));
break;
case nPhaseIdx:
source[energyEq0Idx + phaseIdx] += (componentIntoPhaseMassTransfer[nPhaseIdx][nCompIdx]
* FluidSystem::molarMass(nCompIdx)
* FluidSystem::componentEnthalpy(fluidState, wPhaseIdx, nCompIdx));
* FluidSystem::molarMass(nCompIdx)
* FluidSystem::componentEnthalpy(fluidState, wPhaseIdx, nCompIdx));
source[energyEq0Idx + phaseIdx] += (componentIntoPhaseMassTransfer[nPhaseIdx][wCompIdx]
* FluidSystem::molarMass(wCompIdx)
* FluidSystem::componentEnthalpy(fluidState, wPhaseIdx, wCompIdx));
* FluidSystem::molarMass(wCompIdx)
* FluidSystem::componentEnthalpy(fluidState, wPhaseIdx, wCompIdx));
break;
case sPhaseIdx:
break; // no sorption
......
......@@ -104,8 +104,8 @@ public:
this->resizePhaseBuffer_(prandtlNumber_);
this->resizePhaseBuffer_(nusseltNumber_);
if (velocityAveragingInModel and not velocityOutput/*only one of the two output options, otherwise paraview segfaults due to two times the same field name*/) {
/*only one of the two output options, otherwise paraview segfaults due to two times the same field name*/
if (velocityAveragingInModel and not velocityOutput) {
Scalar numVertices = this->problem_.gridView().size(dim);
for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
velocity_[phaseIdx].resize(numVertices);
......@@ -158,8 +158,10 @@ public:
aws_[globalIdx] = volVars.interfacialArea(wPhaseIdx, sPhaseIdx);
ans_[globalIdx] = volVars.interfacialArea(nPhaseIdx, sPhaseIdx);
if (velocityAveragingInModel and not velocityOutput/*only one of the two output options, otherwise paraview segfaults due to two times the same field name*/){
int numVertices = this->problem_.gridView().size(dim); // numVertices for vertexCentereed, numVolumes for volume centered
/*only one of the two output options, otherwise paraview segfaults due to two times the same field name*/
if (velocityAveragingInModel and not velocityOutput){
// numVertices for vertexCentereed, numVolumes for volume centered
int numVertices = this->problem_.gridView().size(dim);
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
for (int I = 0; I < numVertices; ++I)
velocity_[phaseIdx][I] = this->problem_.model().volumeDarcyVelocity(phaseIdx, I);
......@@ -196,7 +198,8 @@ public:
this->commitPhaseBuffer_(writer, "prandtlNumber_%s", prandtlNumber_);
if (nusseltOutput)
this->commitPhaseBuffer_(writer, "nusseltNumber_%s", nusseltNumber_);
if (velocityAveragingInModel and not velocityOutput/*only one of the two output options, otherwise paraview segfaults due to two timies the same field name*/){
/*only one of the two output options, otherwise paraview segfaults due to two timies the same field name*/
if (velocityAveragingInModel and not velocityOutput){
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
// commit the phase velocity
std::ostringstream oss;
......
......@@ -113,7 +113,8 @@ NEW_PROP_TAG(EnableKinetic);
//! Enable kinetic resolution of energy transfer processes?
NEW_PROP_TAG(EnableKineticEnergy);
//! Enable Maxwell Diffusion? (If false: use Fickian Diffusion) Maxwell incorporated the mutual influences of multiple diffusing components. However, Fick seems to be more robust.
//! Enable Maxwell Diffusion? (If false: use Fickian Diffusion) Maxwell incorporated the mutual
//! influences of multiple diffusing components. However, Fick seems to be more robust.
NEW_PROP_TAG(UseMaxwellDiffusion);
//! The model for the effective thermal conductivity
......
......@@ -99,17 +99,19 @@ public:
const unsigned int scvIdx)
{
// obtain (standard) material parameters (needed for the residual saturations)
const MaterialLawParams & materialParams = problem.spatialParams().materialLawParams(element,fvGeometry,scvIdx) ;
const MaterialLawParams & materialParams = problem.spatialParams().materialLawParams(element,fvGeometry,scvIdx);
//obtain parameters for interfacial area constitutive relations
const AwnSurfaceParams & aWettingNonWettingSurfaceParams = problem.spatialParams().aWettingNonWettingSurfaceParams(element,fvGeometry,scvIdx) ;
const AnsSurfaceParams & aNonWettingSolidSurfaceParams = problem.spatialParams().aNonWettingSolidSurfaceParams(element,fvGeometry,scvIdx) ;
const AwnSurfaceParams & aWettingNonWettingSurfaceParams
= problem.spatialParams().aWettingNonWettingSurfaceParams(element,fvGeometry,scvIdx);
const AnsSurfaceParams & aNonWettingSolidSurfaceParams
= problem.spatialParams().aNonWettingSolidSurfaceParams(element,fvGeometry,scvIdx);
Valgrind::CheckDefined(aWettingNonWettingSurfaceParams);
Valgrind::CheckDefined(aNonWettingSolidSurfaceParams);
const Scalar pc = fluidState.pressure(nPhaseIdx) - fluidState.pressure(wPhaseIdx) ;
const Scalar Sw = fluidState.saturation(wPhaseIdx) ;
const Scalar pc = fluidState.pressure(nPhaseIdx) - fluidState.pressure(wPhaseIdx);
const Scalar Sw = fluidState.saturation(wPhaseIdx);
Valgrind::CheckDefined(Sw);
Valgrind::CheckDefined(pc);
Scalar awn;
......@@ -168,9 +170,10 @@ if (AwnSurface::interfacialArea(aWettingNonWettingSurfaceParams, materialParams,
// positionString << " x"<< (i+1) << "=" << globalPos[i] << " " ;
// positionString << "\n";
//
// std::cout<<"a_{ns} > a_s, set a_{ns}=" << ans <<" to a_{ns}=a_s="<<solidSurface_ << " with S_w="<< Sw << " p_c= "<< pc << positionString.str() ;
// std::cout<<"a_{ns} > a_s, set a_{ns}=" << ans <<" to a_{ns}=a_s="<<solidSurface_ << "
// with S_w="<< Sw << " p_c= "<< pc << positionString.str() ;
//
// ans = solidSurface_ ;
// ans = solidSurface_ ;
// }
#endif
......@@ -186,7 +189,8 @@ if (AwnSurface::interfacialArea(aWettingNonWettingSurfaceParams, materialParams,
interfacialArea_[sPhaseIdx][wPhaseIdx] = interfacialArea_[wPhaseIdx][sPhaseIdx];
interfacialArea_[sPhaseIdx][sPhaseIdx] = 0. ;
#else
const AwsSurfaceParams & aWettingSolidSurfaceParams = problem.spatialParams().aWettingSolidSurfaceParams();
const AwsSurfaceParams & aWettingSolidSurfaceParams
= problem.spatialParams().aWettingSolidSurfaceParams();
Valgrind::CheckDefined(aWettingSolidSurfaceParams);
const Scalar aws = AwsSurface::interfacialArea(aWettingSolidSurfaceParams,materialParams, Sw, pc ); // 10.; //
interfacialArea_[wPhaseIdx][sPhaseIdx] = aws ;
......@@ -216,12 +220,8 @@ if (AwnSurface::interfacialArea(aWettingNonWettingSurfaceParams, materialParams,
const Scalar dynamicViscosity = fluidState.viscosity(phaseIdx);
const Scalar density = fluidState.density(phaseIdx);
const Scalar kinematicViscosity = dynamicViscosity / density;
const Scalar heatCapacity = FluidSystem::heatCapacity(fluidState,
paramCache,
phaseIdx);
const Scalar thermalConductivity = FluidSystem::thermalConductivity(fluidState,
paramCache,
phaseIdx);
const Scalar heatCapacity = FluidSystem::heatCapacity(fluidState, paramCache, phaseIdx);
const Scalar thermalConductivity = FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
// diffusion coefficient of non-wetting component in wetting phase
const Scalar diffCoeff = volVars.diffCoeff(phaseIdx, wCompIdx, nCompIdx) ;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment