Commit 6209fad8 authored by Vishal Jambhekar's avatar Vishal Jambhekar
Browse files

Updated necessory files for 2pncmin model (Reviewed by Fetzer).

- Units in the fluidsystem/brinevariablesalinity and components/nacl checked. Additional commented units removed.
- Commented cout statements removed from computefromfugacities2pncmin, computefromreferencephase2pnc and miscible2pnccomposition
- The ccommented out removed from problem file and input file corrected.
- The test case is checked for both debug and optim model with clang compiler.

   

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@15203 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent d27214ec
...@@ -76,25 +76,25 @@ public: ...@@ -76,25 +76,25 @@ public:
* \brief Returns the critical temperature \f$\mathrm{[K]}\f$ of brine. Here, it is assumed to be equal to that of pure water. * \brief Returns the critical temperature \f$\mathrm{[K]}\f$ of brine. Here, it is assumed to be equal to that of pure water.
*/ */
static Scalar criticalTemperature() static Scalar criticalTemperature()
{ return H2O::criticalTemperature(); /* [K] */ } { return H2O::criticalTemperature(); }
/*! /*!
* \brief Returns the critical pressure \f$\mathrm{[Pa]}\f$ of brine. Here, it is assumed to be equal to that of pure water. * \brief Returns the critical pressure \f$\mathrm{[Pa]}\f$ of brine. Here, it is assumed to be equal to that of pure water.
*/ */
static Scalar criticalPressure() static Scalar criticalPressure()
{ return H2O::criticalPressure(); /* [N/m^2] */ } { return H2O::criticalPressure(); }
/*! /*!
* \brief Returns the temperature \f$\mathrm{[K]}\f$ at brine's triple point. Here, it is assumed to be equal to that of pure water. * \brief Returns the temperature \f$\mathrm{[K]}\f$ at brine's triple point. Here, it is assumed to be equal to that of pure water.
*/ */
static Scalar tripleTemperature() static Scalar tripleTemperature()
{ return H2O::tripleTemperature(); /* [K] */ } { return H2O::tripleTemperature(); }
/*! /*!
* \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at brine's triple point. Here, it is assumed to be equal to that of pure water. * \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at brine's triple point. Here, it is assumed to be equal to that of pure water.
*/ */
static Scalar triplePressure() static Scalar triplePressure()
{ return H2O::triplePressure(); /* [N/m^2] */ } { return H2O::triplePressure(); }
/*! /*!
* \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of pure brine * \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of pure brine
...@@ -103,7 +103,7 @@ public: ...@@ -103,7 +103,7 @@ public:
* \param T temperature of component in \f$\mathrm{[K]}\f$ * \param T temperature of component in \f$\mathrm{[K]}\f$
*/ */
static Scalar vaporPressure(Scalar T) static Scalar vaporPressure(Scalar T)
{ return H2O::vaporPressure(T); /* [N/m^2] */ } { return H2O::vaporPressure(T); }
/*! /*!
* \brief Specific enthalpy of gaseous brine \f$\mathrm{[J/kg]}\f$. * \brief Specific enthalpy of gaseous brine \f$\mathrm{[J/kg]}\f$.
...@@ -113,7 +113,7 @@ public: ...@@ -113,7 +113,7 @@ public:
*/ */
static const Scalar gasEnthalpy(Scalar temperature, static const Scalar gasEnthalpy(Scalar temperature,
Scalar pressure) Scalar pressure)
{ return H2O::gasEnthalpy(temperature, pressure); /* [J/kg] */ } { return H2O::gasEnthalpy(temperature, pressure); }
/*! /*!
* \brief Specific enthalpy of liquid brine \f$\mathrm{[J/kg]}\f$. * \brief Specific enthalpy of liquid brine \f$\mathrm{[J/kg]}\f$.
...@@ -124,6 +124,7 @@ public: ...@@ -124,6 +124,7 @@ public:
* Equations given in: - Palliser & McKibbin 1997 * Equations given in: - Palliser & McKibbin 1997
* - Michaelides 1981 * - Michaelides 1981
* - Daubert & Danner 1989 * - Daubert & Danner 1989
*
*/ */
static const Scalar liquidEnthalpy(Scalar T, static const Scalar liquidEnthalpy(Scalar T,
Scalar p, Scalar salinity) Scalar p, Scalar salinity)
...@@ -320,11 +321,11 @@ public: ...@@ -320,11 +321,11 @@ public:
salinity = std::abs(salinity); salinity = std::abs(salinity);
Scalar T_C = temperature - 273.15; Scalar T_C = temperature - 273.15;
if(salinity < 0.0) if(salinity < 0.0)
{salinity = 0.0;} {salinity = 0.0; }
Scalar A = (0.42*pow((pow(salinity, 0.8)-0.17), 2) + 0.045)*pow(T_C, 0.8); Scalar A = (0.42*pow((pow(salinity, 0.8)-0.17), 2) + 0.045)*pow(T_C, 0.8);
Scalar mu_brine = 0.1 + 0.333*salinity + (1.65+91.9*salinity*salinity*salinity)*exp(-A); Scalar mu_brine = 0.1 + 0.333*salinity + (1.65+91.9*salinity*salinity*salinity)*exp(-A);
assert(mu_brine > 0.0); assert(mu_brine > 0.0);
return mu_brine/1000.0; /* unit: Pa s */ return mu_brine/1000.0;
} }
}; };
} // end namespace } // end namespace
......
#ifndef NACL_HH // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
#define NACL_HH // vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*! /*!
* \file * \file
* *
* \brief A class for the NaCl properties * \ingroup Components
*
* \brief Material properties of pure salt \f$NaCl\f$.
*/ */
#ifndef DUMUX_NACL_HH #ifndef DUMUX_NACL_HH
#define DUMUX_NACL_HH #define DUMUX_NACL_HH
...@@ -33,15 +51,15 @@ public: ...@@ -33,15 +51,15 @@ public:
} }
/*! /*!
* \brief The mass in [kg] for one mole of NaCl. * \brief The molar mass of NaCl in \f$\mathrm{[kg/mol]}\f$.
*/ */
static Scalar molarMass() static Scalar molarMass()
{ {
return 58.4428e-3 ; return 58.4428e-3 ;
} // kg/mol }
/*! /*!
* \brief The diffusion Coefficient of NaCl in water. * \brief The diffusion Coefficient \f$\mathrm{[m^2/s]}\f$ of NaCl in water.
*/ */
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure) static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
{ {
...@@ -49,18 +67,15 @@ public: ...@@ -49,18 +67,15 @@ public:
} }
/*! /*!
* \brief The mass density of NaCl. * \brief The mass density \f$\mathrm{[kg/m^3]}\f$ of NaCl.
*/ */
static Scalar Density() static Scalar Density()
{ {
return (2165); /* 2165 kg/m³*/ return 2165.0;
} }
}; };
} // end namespace } // end namespace
#endif #endif
#endif // NACL_HH
...@@ -70,17 +70,19 @@ public: ...@@ -70,17 +70,19 @@ public:
return; return;
// Pure component fugacities // Pure component fugacities
for (int i = 0; i < numComponents; ++ i) { for (int i = 0; i < numComponents; ++ i)
//std::cout << f << " -> " << mutParams.fugacity(phaseIdx, i)/f << "\n"; {
fluidState.setMoleFraction(phaseIdx, fluidState.setMoleFraction(phaseIdx,i, 1.0/numComponents);
i,
1.0/numComponents);
} }
} }
/*! /*!
* \brief Calculates the chemical equilibrium from the component * \brief Calculates the chemical equilibrium from the component
* fugacities in a phase. This constraint solver is developed for drying scenarios where salt component is restricted to liquid phase and still for the sake for equilibrium calculation some residual salt must be considered in the gas phase. In such cases for existence of gas phase only, in the theoretical liquid phase, we set the mole fraction of salt to 1e-10. * fugacities in a phase. This constraint solver is developed for drying scenarios where
* salt component is restricted to liquid phase and still for the sake for equilibrium
* calculation some residual salt must be considered in the gas phase. In such cases for
* existence of gas phase only, in the theoretical liquid phase, we set the mole fraction
* of salt to 1e-10.
* \param fluidState Thermodynamic state of the fluids * \param fluidState Thermodynamic state of the fluids
* \param paramCache Container for cache parameters * \param paramCache Container for cache parameters
* \param phaseIdx The phase index * \param phaseIdx The phase index
...@@ -105,8 +107,7 @@ public: ...@@ -105,8 +107,7 @@ public:
else else
DUNE_THROW(NumericalProblem, "This constraint solver is not tested for non-ideal mixtures: Please refer computefromfugacities.hh for details" ); DUNE_THROW(NumericalProblem, "This constraint solver is not tested for non-ideal mixtures: Please refer computefromfugacities.hh for details" );
} }
protected: protected:
// update the phase composition in case the phase is an ideal // update the phase composition in case the phase is an ideal
// mixture, i.e. the component's fugacity coefficients are // mixture, i.e. the component's fugacity coefficients are
...@@ -118,7 +119,8 @@ protected: ...@@ -118,7 +119,8 @@ protected:
int phasePresence, int phasePresence,
const ComponentVector &fugacities) const ComponentVector &fugacities)
{ {
for (int i = 0; i < numComponents; ++ i) { for (int i = 0; i < numComponents; ++ i)
{
Scalar phi = FluidSystem::fugacityCoefficient(fluidState, Scalar phi = FluidSystem::fugacityCoefficient(fluidState,
paramCache, paramCache,
phaseIdx, phaseIdx,
...@@ -140,4 +142,4 @@ protected: ...@@ -140,4 +142,4 @@ protected:
}; };
} // end namespace Dumux } // end namespace Dumux
#endif #endif
\ No newline at end of file
...@@ -137,7 +137,8 @@ public: ...@@ -137,7 +137,8 @@ public:
refPhaseIdx)); refPhaseIdx));
// compute the fugacities of all components in the reference phase // compute the fugacities of all components in the reference phase
for (int compIdx = 0; compIdx < numComponents; ++compIdx) { for (int compIdx = 0; compIdx < numComponents; ++compIdx)
{
fluidState.setFugacityCoefficient(refPhaseIdx, fluidState.setFugacityCoefficient(refPhaseIdx,
compIdx, compIdx,
FluidSystem::fugacityCoefficient(fluidState, FluidSystem::fugacityCoefficient(fluidState,
...@@ -148,7 +149,8 @@ public: ...@@ -148,7 +149,8 @@ public:
} }
// compute all quantities for the non-reference phases // compute all quantities for the non-reference phases
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) { for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
{
if (phaseIdx == refPhaseIdx) if (phaseIdx == refPhaseIdx)
continue; // reference phase is already calculated continue; // reference phase is already calculated
......
...@@ -145,8 +145,7 @@ public: ...@@ -145,8 +145,7 @@ public:
} }
} }
//set the additional equations for the numComponents-numMajorComponents //set the additional equations for the numComponents-numMajorComponents
//Components, of which the molefractions are known, //Components, of which the molefractions are known, set to molefraction(knownCompIdx)=xKnown
//to molefraction(knownCompIdx)=xKnown
for(int knownCompIdx = 0; knownCompIdx < numComponents-numMajorComponents; ++knownCompIdx) for(int knownCompIdx = 0; knownCompIdx < numComponents-numMajorComponents; ++knownCompIdx)
{ {
int rowIdx = numComponents + numPhases + knownCompIdx; int rowIdx = numComponents + numPhases + knownCompIdx;
...@@ -195,24 +194,10 @@ public: ...@@ -195,24 +194,10 @@ public:
} }
} }
} }
// std::cout << "M_: "<< M << std::endl;
// std::cout << "b_: "<< b << std::endl;
// solve for all mole fractions // solve for all mole fractions
M.solve(x, b); M.solve(x, b);
// std::cout << "x_: "<< x << std::endl;
// In case known mole fraction is set to zero in problem the system determines unrealistic mole fractions
//resetting the mole fractions to zero again
// if(fluidState.moleFraction(knownPhaseIdx, numMajorComponents)== 0.0);
// for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
// int rowIdx = phaseIdx*numComponents + numMajorComponents;
// x[rowIdx]= 0.0;
// }
// std::cout << "x_1: "<< x << std::endl;
// set all mole fractions and the the additional quantities in // set all mole fractions and the the additional quantities in
// the fluid state // the fluid state
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) { for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
......
...@@ -147,46 +147,35 @@ public: ...@@ -147,46 +147,35 @@ public:
const GridView &gridView) const GridView &gridView)
: ParentType(timeManager, GridCreator::grid().leafGridView()) : ParentType(timeManager, GridCreator::grid().leafGridView())
{ {
try
{ outerSalinity_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterSalinity);
outerSalinity_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterSalinity); temperature_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, Temperature);
saltPorosity_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, SaltPorosity); reservoirPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, ReservoirPressure);
temperature_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, Temperature); initLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, LiquidSaturation);
reservoirPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, ReservoirPressure); initPrecipitatedSalt1_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InitPrecipitatedSalt1);
initLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, LiquidSaturation); initPrecipitatedSalt2_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InitPrecipitatedSalt2);
initPrecipitatedSalt1_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InitPrecipitatedSalt1);
initPrecipitatedSalt2_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InitPrecipitatedSalt2); outerLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterLiqSaturation);
innerLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerLiqSaturation);
outerLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterLiqSaturation); innerSalinity_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerSalinity);
innerLiqSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerLiqSaturation); innerPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerPressure);
innerSalinity_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerSalinity); outerPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterPressure);
innerPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerPressure); reservoirSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, reservoirSaturation);
outerPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterPressure);
reservoirSaturation_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, reservoirSaturation); nTemperature_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NTemperature);
nPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NPressure);
nTemperature_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NTemperature); pressureLow_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, PressureLow);
nPressure_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NPressure); pressureHigh_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, PressureHigh);
pressureLow_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, PressureLow); temperatureLow_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureLow);
pressureHigh_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, PressureHigh); temperatureHigh_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureHigh);
temperatureLow_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureLow); name_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, std::string, Problem, OutputName);
temperatureHigh_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureHigh); freqMassOutput_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, Output, FreqMassOutput);
name_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, std::string, Problem, OutputName); storageLastTimestep_ = Scalar(0);
freqMassOutput_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, Output, FreqMassOutput); lastMassOutputTime_ = Scalar(0);
storageLastTimestep_ = Scalar(0);
lastMassOutputTime_ = Scalar(0); outfile.open("evaporation.out");
outfile << "time; evaporationRate" << std::endl;
outfile.open("evaporation.out");
outfile << "time; evaporationRate" << std::endl;
}
catch (Dumux::ParameterException &e) {
std::cerr << e << ". Abort!\n";
exit(1) ;
}
catch (...) {
std::cerr << "Unknown exception thrown!\n";
exit(1);
}
FluidSystem::init(/*Tmin=*/temperatureLow_, FluidSystem::init(/*Tmin=*/temperatureLow_,
/*Tmax=*/temperatureHigh_, /*Tmax=*/temperatureHigh_,
/*nT=*/nTemperature_, /*nT=*/nTemperature_,
...@@ -233,7 +222,7 @@ public: ...@@ -233,7 +222,7 @@ public:
* This problem assumes a temperature of 10 degrees Celsius. * This problem assumes a temperature of 10 degrees Celsius.
*/ */
Scalar temperature() const Scalar temperature() const
{ return temperature_; }; { return temperature_; }
/*! /*!
* \name Boundary conditions * \name Boundary conditions
...@@ -354,7 +343,6 @@ public: ...@@ -354,7 +343,6 @@ public:
const FVElementGeometry &fvGeometry, const FVElementGeometry &fvGeometry,
int scvIdx, const ElementVolumeVariables &elemVolVars) const int scvIdx, const ElementVolumeVariables &elemVolVars) const
{ {
// std::cout<<"Hallo there :)"<<std::endl;
source = 0; source = 0;
const VolumeVariables &volVars = elemVolVars[scvIdx]; const VolumeVariables &volVars = elemVolVars[scvIdx];
Scalar moleFracNaCl_lPhase = volVars.fluidState().moleFraction(wPhaseIdx, NaClIdx); Scalar moleFracNaCl_lPhase = volVars.fluidState().moleFraction(wPhaseIdx, NaClIdx);
...@@ -362,6 +350,7 @@ public: ...@@ -362,6 +350,7 @@ public:
Scalar massFracNaCl_Max_lPhase = this->spatialParams().SolubilityLimit(); Scalar massFracNaCl_Max_lPhase = this->spatialParams().SolubilityLimit();
Scalar moleFracNaCl_Max_lPhase = massTomoleFrac_(massFracNaCl_Max_lPhase); Scalar moleFracNaCl_Max_lPhase = massTomoleFrac_(massFracNaCl_Max_lPhase);
Scalar moleFracNaCl_Max_gPhase = moleFracNaCl_Max_lPhase / volVars.fluidState().pressure(nPhaseIdx); Scalar moleFracNaCl_Max_gPhase = moleFracNaCl_Max_lPhase / volVars.fluidState().pressure(nPhaseIdx);
Scalar saltPorosity = this->spatialParams().porosityMin(element, fvGeometry, scvIdx);
// liquid phase // liquid phase
Scalar precipSalt = volVars.porosity() * volVars.molarDensity(wPhaseIdx) Scalar precipSalt = volVars.porosity() * volVars.molarDensity(wPhaseIdx)
...@@ -381,7 +370,7 @@ public: ...@@ -381,7 +370,7 @@ public:
if (precipSalt*this->timeManager().timeStepSize() + volVars.precipitateVolumeFraction(sPhaseIdx)* volVars.molarDensity(sPhaseIdx)< 0) if (precipSalt*this->timeManager().timeStepSize() + volVars.precipitateVolumeFraction(sPhaseIdx)* volVars.molarDensity(sPhaseIdx)< 0)
precipSalt = - volVars.precipitateVolumeFraction(sPhaseIdx)* volVars.molarDensity(sPhaseIdx)/this->timeManager().timeStepSize(); precipSalt = - volVars.precipitateVolumeFraction(sPhaseIdx)* volVars.molarDensity(sPhaseIdx)/this->timeManager().timeStepSize();
if (volVars.precipitateVolumeFraction(sPhaseIdx) >= volVars.InitialPorosity() - saltPorosity_ && precipSalt > 0) if (volVars.precipitateVolumeFraction(sPhaseIdx) >= volVars.InitialPorosity() - saltPorosity && precipSalt > 0)
precipSalt = 0; precipSalt = 0;
source[conti0EqIdx + NaClIdx] += -precipSalt; source[conti0EqIdx + NaClIdx] += -precipSalt;
...@@ -437,7 +426,6 @@ private: ...@@ -437,7 +426,6 @@ private:
Scalar pressureLow_, pressureHigh_; Scalar pressureLow_, pressureHigh_;
Scalar temperatureLow_, temperatureHigh_; Scalar temperatureLow_, temperatureHigh_;
Scalar outerSalinity_; Scalar outerSalinity_;
Scalar saltPorosity_;
Scalar reservoirPressure_; Scalar reservoirPressure_;
Scalar innerPressure_; Scalar innerPressure_;
Scalar outerPressure_; Scalar outerPressure_;
......
...@@ -105,15 +105,13 @@ public: ...@@ -105,15 +105,13 @@ public:
//set main diagonal entries of the permeability tensor to a value //set main diagonal entries of the permeability tensor to a value
//setting to one value means: isotropic, homogeneous //setting to one value means: isotropic, homogeneous
for (int i = 0; i < dim; i++) for (int i = 0; i < dim; i++)
// K_[i][i] = 1e-7; K_[i][i] = 2.23e-14;
K_[i][i] = 2.23e-14;//2.23e-14;
// residual saturations // residual saturations
materialParams_.setSwr(0.2); materialParams_.setSwr(0.2);
materialParams_.setSnr(1E-3); materialParams_.setSnr(1E-3);
//parameters of Brooks & Corey Law //parameters of Brooks & Corey Law
//materialParams_.setPe(500.0);
materialParams_.setPe(500); materialParams_.setPe(500);
materialParams_.setLambda(2); materialParams_.setLambda(2);
} }
...@@ -127,8 +125,7 @@ public: ...@@ -127,8 +125,7 @@ public:
* \param globalSolution The global solution vector * \param globalSolution The global solution vector
*/ */
void update(const SolutionVector &globalSolution) void update(const SolutionVector &globalSolution)
{ { };
};
/*! Intrinsic permeability tensor K \f$[m^2]\f$ depending /*! Intrinsic permeability tensor K \f$[m^2]\f$ depending
* on the position in the domain * on the position in the domain
...@@ -141,13 +138,15 @@ public: ...@@ -141,13 +138,15 @@ public:
* could be defined, where globalPos is the vector including the global coordinates * could be defined, where globalPos is the vector including the global coordinates
* of the finite volume. * of the finite volume.
*/ */
const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element, /*@\label{tutorial-coupled:permeability}@*/ const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element,
const FVElementGeometry &fvGeometry, const FVElementGeometry &fvGeometry,
const int scvIdx) const const int scvIdx) const
{ return K_; } {
return K_;
}
/*! /*!
* \brief Define the porosity \f$[-]\f$ of the spatial parameters * \brief Define the minimum porosity \f$[-]\f$ after salt precipitation
* *
* \param elemVolVars The data defined on the sub-control volume * \param elemVolVars The data defined on the sub-control volume
* \param element The finite element * \param element The finite element
...@@ -159,10 +158,10 @@ public: ...@@ -159,10 +158,10 @@ public:
const FVElementGeometry &fvGeometry, const FVElementGeometry &fvGeometry,
int scvIdx) const int scvIdx) const
{ {
return 1e-5; return 1e-5;
} }
/*! /*!
* \brief Define the minimum porosity \f$[-]\f$ after clogging caused by mineralization * \brief Define the minimum porosity \f$[-]\f$ after clogging caused by mineralization
* *
* \param elemVolVars The data defined on the sub-control volume * \param elemVolVars The data defined on the sub-control volume
...@@ -175,7 +174,7 @@ public: ...@@ -175,7 +174,7 @@ public:
const FVElementGeometry &fvGeometry, const FVElementGeometry &fvGeometry,
int scvIdx) const int scvIdx) const
{ {
return 0.11; return 0.11;
} }
...@@ -184,20 +183,20 @@ public: ...@@ -184,20 +183,20 @@ public:
int scvIdx) const int scvIdx) const
{ {
return 1 - 0.11; return 1 - 0.11;
} }