From 1edee53600be0629c932f20e892b89bb85c004c4 Mon Sep 17 00:00:00 2001
From: Markus Wolff <markus.wolff@twt-gmbh.de>
Date: Wed, 18 Jan 2012 15:08:58 +0000
Subject: [PATCH] same clean-up and small modifications
- clean-up in 2p models
- added possibility to fix a pressure somewhere in the domain
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@7430 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
.../decoupled/2p/diffusion/fv/fvpressure2p.hh | 5 ++-
.../decoupled/2p/diffusion/fv/fvvelocity2p.hh | 6 +++-
.../2p/transport/fv/fvsaturation2p.hh | 2 +-
.../decoupled/2p/transport/fv/gravitypart.hh | 2 +-
dumux/decoupled/common/fv/fvpressure.hh | 32 ++++++++++++++++++-
5 files changed, 42 insertions(+), 5 deletions(-)
diff --git a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
index 073f9a196c..c5d04d53ee 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
@@ -163,6 +163,9 @@ public:
}
// std::cout<<"Pressure defect = "<<pressureNorm<<"; "<<numIter<<" Iterations needed for initial pressure field"<<std::endl;
}
+
+ storePressureSolution();
+
return;
}
@@ -334,7 +337,7 @@ public:
if (!compressibility_)
{
- const Element& element = *(problem_.gridView().template end<0> ());
+ const Element& element = *(problem_.gridView().template begin<0> ());
FluidState fluidState;
fluidState.setPressure(wPhaseIdx, problem_.referencePressure(element));
fluidState.setPressure(nPhaseIdx, problem_.referencePressure(element));
diff --git a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
index deebb8c043..7de45cacf0 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh
@@ -125,7 +125,7 @@ public:
if (!compressibility_)
{
- const Element& element = *(problem_.gridView().template end<0> ());
+ const Element& element = *(problem_.gridView().template begin<0> ());
FluidState fluidState;
fluidState.setPressure(wPhaseIdx, problem_.referencePressure(element));
fluidState.setPressure(nPhaseIdx, problem_.referencePressure(element));
@@ -654,6 +654,10 @@ void FVVelocity2P<TypeTag>::calculateVelocityOnBoundary(const Intersection& inte
velocityNW /= density_[nPhaseIdx];
}
+ //store potential gradients for further calculations
+ cellDataI.fluxData().setPotential(wPhaseIdx, isIndex, boundValues[wPhaseIdx]);
+ cellDataI.fluxData().setPotential(nPhaseIdx, isIndex, boundValues[nPhaseIdx]);
+
cellDataI.fluxData().setVelocity(wPhaseIdx, isIndex, velocityW);
cellDataI.fluxData().setVelocity(nPhaseIdx, isIndex, velocityNW);
cellDataI.fluxData().setVelocityMarker(isIndex);
diff --git a/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh b/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
index 710cdaf306..dc5120e7e3 100644
--- a/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
+++ b/dumux/decoupled/2p/transport/fv/fvsaturation2p.hh
@@ -340,7 +340,7 @@ public:
if (!compressibility_)
{
- const Element& element = *(problem_.gridView().template end<0> ());
+ const Element& element = *(problem_.gridView().template begin<0> ());
FluidState fluidState;
fluidState.setPressure(wPhaseIdx, problem_.referencePressure(element));
fluidState.setPressure(nPhaseIdx, problem_.referencePressure(element));
diff --git a/dumux/decoupled/2p/transport/fv/gravitypart.hh b/dumux/decoupled/2p/transport/fv/gravitypart.hh
index ce6b65daad..6af0d9b17e 100644
--- a/dumux/decoupled/2p/transport/fv/gravitypart.hh
+++ b/dumux/decoupled/2p/transport/fv/gravitypart.hh
@@ -179,7 +179,7 @@ public:
GravityPart (Problem& problem)
: ConvectivePart<TypeTag>(problem), problem_(problem), preComput_(GET_PROP_VALUE(TypeTag, PrecomputedConstRels))
{
- const Element& element = *(problem_.gridView().template end<0> ());
+ const Element& element = *(problem_.gridView().template begin<0> ());
FluidState fluidState;
fluidState.setPressure(wPhaseIdx, problem_.referencePressure(element));
fluidState.setPressure(nPhaseIdx, problem_.referencePressure(element));
diff --git a/dumux/decoupled/common/fv/fvpressure.hh b/dumux/decoupled/common/fv/fvpressure.hh
index 8082259b15..46d244df53 100644
--- a/dumux/decoupled/common/fv/fvpressure.hh
+++ b/dumux/decoupled/common/fv/fvpressure.hh
@@ -142,13 +142,32 @@ public:
}
//@}
+ /*! set a pressure to be fixed at a certain cell */
+ void setFixPressureAtIndex(Scalar pressure, int globalIdx)
+ {
+ hasFixPressureAtIndex_ = true;
+ fixPressureAtIndex_ = pressure;
+ idxFixPressureAtIndex_ = globalIdx;
+ }
+
+ /*! unset a fixed pressure at a certain cell */
+ void unsetFixPressureAtIndex(int globalIdx)
+ {
+ hasFixPressureAtIndex_ = false;
+ fixPressureAtIndex_ = 0.0;
+ idxFixPressureAtIndex_ = 0.0;
+ }
+
//! Constructs a FVPressure object
/**
* \param problem a problem class object
*/
FVPressure(Problem& problem) :
problem_(problem), size_(problem.gridView().size(0)),
- pressure_(size_), A_(size_, size_, (2 * dim + 1) * size_, Matrix::random), f_(size_)
+ pressure_(size_), A_(size_, size_, (2 * dim + 1) * size_, Matrix::random), f_(size_),
+ fixPressureAtIndex_(0),
+ idxFixPressureAtIndex_(0),
+ hasFixPressureAtIndex_(false)
{}
private:
@@ -170,6 +189,9 @@ private:
protected:
Matrix A_;
RHSVector f_;
+ Scalar fixPressureAtIndex_;
+ Scalar idxFixPressureAtIndex_;
+ bool hasFixPressureAtIndex_;
};
//! initializes the matrix to store the system of equations
@@ -325,6 +347,14 @@ void FVPressure<TypeTag>::solve()
if (verboseLevelSolver)
std::cout << __FILE__ <<": solve for pressure" << std::endl;
+ //set a fixed pressure for a certain cell
+ if (hasFixPressureAtIndex_)
+ {
+ A_[idxFixPressureAtIndex_] = 0;
+ A_[idxFixPressureAtIndex_][idxFixPressureAtIndex_] = 1;
+ f_[idxFixPressureAtIndex_] = fixPressureAtIndex_;
+ }
+
Solver solver(problem_);
solver.solve(A_, pressure_, f_);
// printmatrix(std::cout, A_, "global stiffness matrix", "row", 11, 3);
--
GitLab