From 892e61d03de2ebb341c211ca42fcf30c4a0ecfc2 Mon Sep 17 00:00:00 2001
From: Benjamin Faigle <benjamin.faigle@posteo.de>
Date: Mon, 20 Feb 2012 10:21:49 +0000
Subject: [PATCH] remove volume-derivatives functions of the old model which
was totally accidentially comitted. This will also remove the strange and
unnecessary warning.
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@7823 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
.../decoupled/2p2c/fvpressurecompositional.hh | 121 ------------------
1 file changed, 121 deletions(-)
diff --git a/dumux/decoupled/2p2c/fvpressurecompositional.hh b/dumux/decoupled/2p2c/fvpressurecompositional.hh
index 4b3b43f5ea..d3d03aa7f6 100644
--- a/dumux/decoupled/2p2c/fvpressurecompositional.hh
+++ b/dumux/decoupled/2p2c/fvpressurecompositional.hh
@@ -605,8 +605,6 @@ void FVPressureCompositional<TypeTag>::initialMaterialLaws(bool compositional)
* \param globalPos The global position of the current element
* \param ep A pointer to the current element
*/
-#if 1
-#warning decide which one!!!
template<class TypeTag>
void FVPressureCompositional<TypeTag>::volumeDerivatives(const GlobalPosition& globalPos, const Element& element)
{
@@ -717,125 +715,6 @@ void FVPressureCompositional<TypeTag>::volumeDerivatives(const GlobalPosition& g
}
cellData.confirmVolumeDerivatives();
}
-#else
-template<class TypeTag>
-void FVPressureCompositional<TypeTag>::volumeDerivatives(const GlobalPosition& globalPos, const Element& element)
-{
- // cell index
- int globalIdx = problem_.variables().index(element);
-
- CellData& cellData = problem_.variables().cellData(globalIdx);
-
- // get cell temperature
- Scalar temperature_ = problem_.temperatureAtPos(globalPos);
-
- // initialize an Fluid state for the update
- FluidState updFluidState;
-
- /**********************************
- * a) get necessary variables
- **********************************/
- //determine phase pressures from primary pressure variable
- PhaseVector pressure(0.);
- switch (pressureType)
- {
- //TODO: use pressure from cellData here
- case pw:
- {
- pressure[wPhaseIdx] = this->pressure()[globalIdx];
- pressure[nPhaseIdx] = this->pressure()[globalIdx]
- + cellData.capillaryPressure();
- break;
- }
- case pn:
- {
- pressure[wPhaseIdx] = this->pressure()[globalIdx]
- - cellData.capillaryPressure();
- pressure[nPhaseIdx] = this->pressure()[globalIdx];
- break;
- }
- }
-
- // mass of components inside the cell
- ComponentVector mass(0.);
- mass[0] = cellData.massConcentration(wCompIdx);
- mass[1] = cellData.massConcentration(nCompIdx);
-
- // shortcuts for density
- Scalar densityW = cellData.density(wPhaseIdx);
- Scalar densityNW = cellData.density(nPhaseIdx);
-
- // actual fluid volume
- Scalar volalt = (mass[0]+mass[1])
- * (cellData.phaseMassFraction(wPhaseIdx) / densityW
- + cellData.phaseMassFraction(nPhaseIdx) / densityNW);
-
- /**********************************
- * b) define increments
- **********************************/
- // increments for numerical derivatives
- ComponentVector massIncrement(0.);
- massIncrement[0] = updateEstimate_[wCompIdx][globalIdx];
- massIncrement[1] = updateEstimate_[nCompIdx][globalIdx];
- if(fabs(massIncrement[0]) < 1e-8 * densityW)
- massIncrement[0] = 1e-8* densityW;
- if(fabs(massIncrement[1]) < 1e-8 * densityNW)
- massIncrement[1] = 1e-8 * densityNW;
- Scalar incp = 1e-2;
-
-
- /**********************************
- * c) Secant method for derivatives
- **********************************/
-
- // numerical derivative of fluid volume with respect to pressure
- PhaseVector p_(incp);
- p_ += pressure;
- Scalar Z1 = mass[0] / (mass[0] + mass[1]);
- updFluidState.update(Z1,
- p_, problem_.spatialParameters().porosity(globalPos, element), temperature_);
- cellData.dv_dp() = (((mass[0]+mass[1]) * (updFluidState.phaseMassFraction(wPhaseIdx) /updFluidState.density(wPhaseIdx)
- + updFluidState.phaseMassFraction(nPhaseIdx) /updFluidState.density(nPhaseIdx))) - volalt) /incp;
-
- if (cellData.dv_dp()>0)
- {
- // dV_dp > 0 is unphysical: Try inverse increment for secant
- Dune::dinfo << "dv_dp larger 0 at Idx " << globalIdx << " , try and invert secant"<< std::endl;
-
- p_ -= 2*incp;
- updFluidState.update(Z1,
- p_, problem_.spatialParameters().porosity(globalPos, element), temperature_);
- cellData.dv_dp() = (((mass[0]+mass[1]) * (updFluidState.phaseMassFraction(wPhaseIdx) /updFluidState.density(wPhaseIdx)
- + updFluidState.phaseMassFraction(nPhaseIdx) /updFluidState.density(nPhaseIdx))) - volalt) /incp;
- // dV_dp > 0 is unphysical: Try inverse increment for secant
- if (cellData.dv_dp()>0)
- {
- Dune::dinfo << "dv_dp still larger 0 after inverting secant"<< std::endl;
- }
- }
-
- // numerical derivative of fluid volume with respect to mass of components
- for (int comp = 0; comp<numComponents; comp++)
- {
- mass[comp] += massIncrement[comp];
- Z1 = mass[0] / (mass[0] + mass[1]);
- updFluidState.update(Z1, pressure, problem_.spatialParameters().porosity(globalPos, element), temperature_);
-
- cellData.dv(comp) = ((mass[0]+mass[1])
- * (updFluidState.phaseMassFraction(wPhaseIdx) / updFluidState.density(wPhaseIdx)
- + updFluidState.phaseMassFraction(nPhaseIdx) / updFluidState.density(nPhaseIdx)) - volalt)
- / massIncrement[comp];
- mass[comp] -= massIncrement[comp];
-
- //check routines if derivatives are meaningful
- if (isnan(cellData.dv(comp)) || isinf(cellData.dv(comp)) )
- {
- DUNE_THROW(Dune::MathError, "NAN/inf of dV_dm. If that happens in first timestep, try smaller firstDt!");
- }
- }
- cellData.confirmVolumeDerivatives();
-}
-#endif
}//end namespace Dumux
#endif
--
GitLab