diff --git a/dumux/decoupled/2p/2pfluidstate.hh b/dumux/decoupled/2p/2pfluidstate.hh
deleted file mode 100644
index 42d108c2a832458b2ef96ca0342a8e8f731c70e7..0000000000000000000000000000000000000000
--- a/dumux/decoupled/2p/2pfluidstate.hh
+++ /dev/null
@@ -1,184 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2009 by Markus Wolff *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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
- *
- * \brief Determines the phase state of the immiscible sequential 2p model.
- */
-#ifndef DUMUX_2P_FLUID_STATE_DECOUPLED_HH
-#define DUMUX_2P_FLUID_STATE_DECOUPLED_HH
-
-#include <dumux/material/old_fluidsystems/fluidstate.hh>
-#include <dumux/decoupled/2p/2pproperties_old.hh>
-
-namespace Dumux
-{
-/*!
- * \ingroup IMPES
- */
-/*!
- * \brief Calcultes the phase state from the primary variables in the sequential
- * 2p model.
- */
-template <class TypeTag>
-class DecoupledTwoPFluidState : public FluidState<typename GET_PROP_TYPE(TypeTag, Scalar),
-DecoupledTwoPFluidState<TypeTag> >
-{
- typedef DecoupledTwoPFluidState<TypeTag> ThisType;
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
-
- typedef typename GET_PROP_TYPE(TypeTag, TwoPIndices) Indices;
-
- enum {
- wPhaseIdx = Indices::wPhaseIdx,
- nPhaseIdx = Indices::nPhaseIdx,
-
- wCompIdx = Indices::wPhaseIdx,
- nCompIdx = Indices::nPhaseIdx
- };
-
-public:
- enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
-
-public:
- void update(Scalar Sw, Scalar pressW, Scalar pressN, Scalar temperature)
- {
- Sw_ = Sw;
- phasePressure_[wPhaseIdx] = pressW;
- phasePressure_[nPhaseIdx] = pressN;
- temperature_=temperature;
- density_[wPhaseIdx] = FluidSystem::phaseDensity(wPhaseIdx,
- temperature,
- pressW,
- *this);
- density_[nPhaseIdx] = FluidSystem::phaseDensity(nPhaseIdx,
- temperature,
- pressN,
- *this);
- }
-
- /*!
- * \brief Returns the saturation of a phase.
- */
- Scalar saturation(int phaseIdx) const
- {
- if (phaseIdx == wPhaseIdx)
- return Sw_;
- else
- return Scalar(1.0) - Sw_;
- };
-
- /*!
- * \brief Returns the mass fraction of a component in a phase.
- */
- Scalar massFrac(int phaseIdx, int compIdx) const
- {
- if (compIdx == phaseIdx)
- return 1.0;
- return 0;
- }
-
- /*!
- * \brief Returns the molar fraction of a component in a fluid phase.
- */
- Scalar moleFrac(int phaseIdx, int compIdx) const
- {
- return massFrac(phaseIdx, compIdx);
- }
-
- /*!
- * \brief Returns the total concentration of a phase \f$\mathrm{[mol/m^3]}\f$.
- *
- * This is equivalent to the sum of all component concentrations.
- */
- Scalar phaseConcentration(int phaseIdx) const
- {
- return density_[phaseIdx]/FluidSystem::molarMass(phaseIdx);
- };
-
- /*!
- * \brief Returns the concentration of a component in a phase \f$\mathrm{[mol/m^3]}\f$.
- */
- Scalar concentration(int phaseIdx, int compIdx) const
- {
- if (phaseIdx == compIdx)
- return phaseConcentration(phaseIdx);
- return 0;
- };
-
- /*!
- * \brief Returns the density of a phase \f$\mathrm{[kg/m^3]}\f$.
- */
- Scalar density(int phaseIdx) const
- { return density_[phaseIdx]; }
-
- /*!
- * \brief Returns mean molar mass of a phase \f$\mathrm{[kg/mol]}\f$.
- *
- * This is equivalent to the sum of all component molar masses
- * weighted by their respective mole fraction.
- */
- Scalar averageMolarMass(int phaseIdx) const
- { return FluidSystem::molarMass(phaseIdx); };
-
- /*!
- * \brief Returns the partial pressure of a component in the gas phase \f$\mathrm{[Pa]}\f$.
- */
- Scalar partialPressure(int compIdx) const
- {
- if (compIdx == wPhaseIdx)
- return 0;
- return phasePressure_[nPhaseIdx];
- }
-
- /*!
- * \brief Returns the pressure of a fluid phase \f$\mathrm{[Pa]}\f$.
- */
- Scalar phasePressure(int phaseIdx) const
- { return phasePressure_[phaseIdx]; }
-
- /*!
- * \brief Returns the capillary pressure \f$\mathrm{[Pa]}\f$
- */
- Scalar capillaryPressure() const
- { return phasePressure_[nPhaseIdx] - phasePressure_[wPhaseIdx]; }
-
- /*!
- * \brief Returns the temperature of the fluids \f$\mathrm{[K]}\f$.
- *
- * Note that we assume thermodynamic equilibrium, so all fluids
- * and the rock matrix exhibit the same temperature.
- */
- Scalar temperature() const
- { return temperature_; };
-
-public:
- Scalar density_[numPhases];
- Scalar phasePressure_[numPhases];
- Scalar temperature_;
- Scalar Sw_;
-};
-
-} // end namepace
-
-#endif
diff --git a/dumux/decoupled/2p/2pproperties_old.hh b/dumux/decoupled/2p/2pproperties_old.hh
deleted file mode 100644
index bc54f00aede8a906520c137bf6c4bdec24ec1709..0000000000000000000000000000000000000000
--- a/dumux/decoupled/2p/2pproperties_old.hh
+++ /dev/null
@@ -1,171 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2008-2009 by Markus Wolff *
- * Copyright (C) 2008-2009 by Andreas Lauser *
- * Copyright (C) 2008 by Bernd Flemisch *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-
-/*!
- * \ingroup IMPES
- * \ingroup Properties
- */
-/*!
- * \file
- *
- * \brief Defines the properties required for (immiscible) twophase sequential models.
- */
-
-#ifndef DUMUX_2PPROPERTIES_DECOUPLED_HH
-#define DUMUX_2PPROPERTIES_DECOUPLED_HH
-
-//Dumux-includes
-#include <dumux/decoupled/common/impetproperties.hh>
-#include <dumux/decoupled/common/transportproperties.hh>
-#include "2pindices.hh"
-
-namespace Dumux
-{
-
-////////////////////////////////
-// forward declarations
-////////////////////////////////
-
-template<class TypeTag>
-class VariableClass2P;
-
-template<class TypeTag>
-class FluidSystem2P;
-
-template<class TypeTag>
-class DecoupledTwoPFluidState;
-
-//template<class TypeTag>
-//struct DecoupledTwoPCommonIndices;
-
-////////////////////////////////
-// properties
-////////////////////////////////
-namespace Properties
-{
-//////////////////////////////////////////////////////////////////
-// Type tags
-//////////////////////////////////////////////////////////////////
-
-//! The type tag for the two-phase problems
-NEW_TYPE_TAG(DecoupledTwoP, INHERITS_FROM(IMPET, Transport))
-;
-
-//////////////////////////////////////////////////////////////////
-// Property tags
-//////////////////////////////////////////////////////////////////
-
-NEW_PROP_TAG ( TwoPIndices )
-;
-NEW_PROP_TAG( SpatialParameters )
-; //!< The type of the spatial parameters object
-NEW_PROP_TAG(MaterialLaw); //!< The material law which ought to be used (extracted from the spatial parameters)
-NEW_PROP_TAG(MaterialLawParams); //!< The context material law (extracted from the spatial parameters)
-NEW_PROP_TAG( EnableGravity)
-; //!< Returns whether gravity is considered in the problem
-NEW_PROP_TAG( Formulation);
-NEW_PROP_TAG( PressureFormulation)
-; //!< The formulation of the model
-NEW_PROP_TAG( SaturationFormulation)
-; //!< The formulation of the model
-NEW_PROP_TAG( VelocityFormulation)
-; //!< The formulation of the model
-NEW_PROP_TAG( EnableCompressibility)
-;// !< Returns whether compressibility is allowed
-NEW_PROP_TAG( WettingPhase)
-; //!< The wetting phase for two-phase models
-NEW_PROP_TAG( NonwettingPhase)
-; //!< The non-wetting phase for two-phase models
-NEW_PROP_TAG( FluidSystem )//!< Defines the fluid system
-;
-NEW_PROP_TAG( FluidState )//!< Defines the fluid state
-;
-
-NEW_PROP_TAG( ErrorTermFactor );
-NEW_PROP_TAG( ErrorTermLowerBound );
-NEW_PROP_TAG( ErrorTermUpperBound );
-
-//////////////////////////////////////////////////////////////////
-// Properties
-//////////////////////////////////////////////////////////////////
-
-SET_INT_PROP(DecoupledTwoP, NumEq, 2);
-
-SET_INT_PROP(DecoupledTwoP, NumPhases, 2);//!< The number of phases in the 2p model is 2
-
-SET_INT_PROP(DecoupledTwoP, NumComponents, 1); //!< Each phase consists of 1 pure component
-
-SET_INT_PROP(DecoupledTwoP,
- Formulation,
- DecoupledTwoPCommonIndices::pwSw);
-
-SET_PROP(DecoupledTwoP, TwoPIndices)
-{
-typedef DecoupledTwoPIndices<GET_PROP_VALUE(TypeTag, Formulation), 0> type;
-};
-
-//! Set the default formulation
-SET_INT_PROP(DecoupledTwoP,
- PressureFormulation,
- GET_PROP_TYPE(TypeTag, TwoPIndices)::pressureType);
-
-SET_INT_PROP(DecoupledTwoP,
- SaturationFormulation,
- GET_PROP_TYPE(TypeTag, TwoPIndices)::saturationType);
-
-SET_INT_PROP(DecoupledTwoP,
- VelocityFormulation,
- GET_PROP_TYPE(TypeTag, TwoPIndices)::velocityDefault);
-
-SET_BOOL_PROP(DecoupledTwoP, EnableCompressibility, false);
-
-SET_TYPE_PROP(DecoupledTwoP, Variables, VariableClass2P<TypeTag>);
-
-SET_TYPE_PROP(DecoupledTwoP, FluidSystem, FluidSystem2P<TypeTag>);
-
-SET_TYPE_PROP(DecoupledTwoP, FluidState, DecoupledTwoPFluidState<TypeTag>);
-
-/*!
- * \brief Set the property for the material parameters by extracting
- * it from the material law.
- */
-SET_PROP(DecoupledTwoP, MaterialLawParams)
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
-
-public:
- typedef typename MaterialLaw::Params type;
-};
-
-SET_SCALAR_PROP(DecoupledTwoP, ErrorTermFactor, 0.5);
-SET_SCALAR_PROP(DecoupledTwoP, ErrorTermLowerBound, 0.1);
-SET_SCALAR_PROP(DecoupledTwoP, ErrorTermUpperBound, 0.9);
-
-// \}
-}
-
-}
-
-#endif
diff --git a/dumux/decoupled/2p/variableclass2p_gridadapt_old.hh b/dumux/decoupled/2p/variableclass2p_gridadapt_old.hh
deleted file mode 100644
index ee17cb303b0743027475d8b18a356708b6bb7278..0000000000000000000000000000000000000000
--- a/dumux/decoupled/2p/variableclass2p_gridadapt_old.hh
+++ /dev/null
@@ -1,255 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2011 by Markus Wolff *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-#ifndef DUMUX_VariableClass2PGridAdapt_GRIDADAPT_HH
-#define DUMUX_VariableClass2PGridAdapt_GRIDADAPT_HH
-
-#include <dumux/decoupled/2p/variableclass2p.hh>
-#include <dune/grid/utility/persistentcontainer.hh>
-
-/**
- * @file
- * @brief Class including the variables and data of discretized data of the constitutive relations
- * @author Markus Wolff, Benjamin Faigle
- */
-
-namespace Dumux
-{
-/*!
- * \ingroup IMPES
- */
-//! Class including the variables and data of discretized data of the constitutive relations.
-/*! The variables of two-phase flow, which are one pressure and one saturation are stored in this class.
- *
- * \tparam TypeTag The Type Tag
- */
-template<class TypeTag>
-class VariableClass2PGridAdapt: public VariableClass2P<TypeTag>
-{
-private:
- typedef VariableClass2P<TypeTag> ParentType;
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
- typedef typename GridView::Grid Grid;
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
-
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
- typedef typename GET_PROP_TYPE(TypeTag, TransportSolutionType) TransportSolutionType;
-
- typedef typename GET_PROP_TYPE(TypeTag, TwoPIndices) Indices;
-
- typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
-
- enum
- {
- dim = GridView::dimension, dimWorld = GridView::dimensionworld
- };
- enum
- {
- pw = Indices::pressureW, pn = Indices::pressureNW, pglobal = Indices::pressureGlobal
- };
- enum
- {
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
- };
-
- static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation);
-
- typedef typename GridView::Traits::template Codim<0>::Entity Element;
- typedef typename Grid::template Codim<0>::EntityPointer ElementPointer;
-
- typedef Dune::FieldVector<Scalar, dim> LocalPosition;
- typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
- typedef typename GridView::IndexSet IndexSet;
- typedef typename Grid::LocalIdSet IdSet;
- typedef typename IdSet::IdType IdType;
-
- typedef typename Grid::LevelGridView LevelGridView;
- typedef typename LevelGridView::template Codim<0>::Iterator LevelIterator;
-
- typedef typename GridView::template Codim<0>::Iterator ElementIterator;
- typedef typename GridView::IntersectionIterator IntersectionIterator;
-
- //*******************************
- // TODO: Doc me! Besseren Namen finden!
- struct RestrictedValue
- {
- Scalar saturation;
- Scalar press;
- Scalar volCorr;
- int count;
- bool front;
- RestrictedValue()
- {
- saturation = 0.;
- press = 0.;
- count = 0;
- volCorr = 0;
- front = false;
- }
- };
-
- const Grid& grid_;
- Dune::PersistentContainer<Grid, RestrictedValue> restrictionmap_;
- std::vector<bool> front_;
-
-public:
-
- //! Constructs a VariableClass object
- /**
- * @param gridView a DUNE gridview object corresponding to diffusion and transport equation
- */
-
- VariableClass2PGridAdapt(const GridView& gridView) :
- ParentType(gridView), grid_(gridView.grid()), restrictionmap_(grid_, 0), front_(
- gridView.size(0), false) //codim 0 map
- {
- }
-
- bool isFront(int idx) const
- {
- return front_[idx];
- }
-
- void setFront(int idx)
- {
- front_[idx] = true;
- }
-
- void adaptVariableSize(int size)
- {
- front_.resize(size);
-
- ParentType::adaptVariableSize(size); //Also adapt pressure, velocity and potential
- }
-
- /*!
- * Store primary variables
- *
- * To reconstruct the solution in father elements, problem properties might
- * need to be accessed.
- *
- * @param problem The current problem
- */
- void storePrimVars(const Problem& problem)
- {
- // loop over all levels of the grid
- for (int level = grid_.maxLevel(); level >= 0; level--)
- {
- //get grid view on level grid
- LevelGridView levelView = grid_.levelView(level);
- for (LevelIterator it = levelView.template begin<0>(); it != levelView.template end<0>(); ++it)
- {
- //get your map entry
- RestrictedValue &rv = restrictionmap_[*it];
-
- // put your value in the map
- if (it->isLeaf())
- {
- // get index
- int indexI = this->index(*it);
-
- rv.saturation = this->saturation()[indexI][0];
- rv.press = this->pressure()[indexI][0];
- rv.volCorr = this->volumecorrection(indexI);
- rv.front = front_[indexI];
- rv.count = 1;
- }
- //Average in father
- if (it.level() > 0)
- {
- ElementPointer epFather = it->father();
- RestrictedValue& rvf = restrictionmap_[*epFather];
- rvf.saturation += rv.saturation / rv.count;
- rvf.press += rv.press / rv.count;
- rvf.volCorr += rv.volCorr / rv.count;
- rvf.count += 1;
- }
- }
- }
- }
-
- /*!
- * Reconstruct missing primary variables (where elements are created/deleted)
- *
- * To reconstruct the solution in father elements, problem properties might
- * need to be accessed.
- *
- * @param problem The current problem
- */
- void reconstructPrimVars(const Problem& problem)
- {
- restrictionmap_.reserve();
-
- for (int level = 0; level <= grid_.maxLevel(); level++)
- {
- LevelGridView levelView = grid_.levelView(level);
- for (LevelIterator it = levelView.template begin<0>(); it != levelView.template end<0>(); ++it)
- {
- if (!it->isNew())
- {
- //entry is in map, write in leaf
- if (it->isLeaf())
- {
- RestrictedValue &rv = restrictionmap_[*it];
- int newIdxI = this->index(*it);
- this->saturation()[newIdxI][0] = rv.saturation / rv.count;
- this->pressure()[newIdxI][0] = rv.press / rv.count;
- this->volumecorrection(newIdxI) = rv.volCorr / rv.count;
- front_[newIdxI] = rv.front;
- }
- }
- else
- {
- // value is not in map, interpolate from father element
- if (it.level() > 0)
- {
- ElementPointer ep = it->father();
- RestrictedValue& rvf = restrictionmap_[*ep];
- if (it->isLeaf())
- {
- int newIdxI = this->index(*it);
-
- this->saturation()[newIdxI][0] = rvf.saturation / rvf.count;
- this->pressure()[newIdxI][0] = rvf.press / rvf.count;
- this->volumecorrection(newIdxI) = rvf.volCorr / rvf.count;
- front_[newIdxI] = false;
- }
- else
- {
- //create new entry
- RestrictedValue& rv = restrictionmap_[*it];
- rv.saturation = rvf.saturation / rvf.count;
- rv.press = rvf.press / rvf.count;
- rv.volCorr = rvf.volCorr / rvf.count;
- rv.count = 1;
- }
- }
- }
- }
-
- }
- // reset entries in restrictionmap
- restrictionmap_.clear();
- }
-};
-}
-#endif
diff --git a/dumux/decoupled/2p/variableclass2p_old.hh b/dumux/decoupled/2p/variableclass2p_old.hh
deleted file mode 100644
index b8f07bdb26bc7ff8aaa5b629742a3531fa7b46a8..0000000000000000000000000000000000000000
--- a/dumux/decoupled/2p/variableclass2p_old.hh
+++ /dev/null
@@ -1,627 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2009 by Markus Wolff *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-#ifndef DUMUX_VARIABLECLASS2P_HH
-#define DUMUX_VARIABLECLASS2P_HH
-
-#include <dumux/decoupled/common/variableclass_old.hh>
-#include <dumux/decoupled/2p/2pfluidstate.hh>
-#include "2pproperties_old.hh"
-
-/**
- * @file
- * @brief Class including the variables and data of discretized data of the constitutive relations
- * @author Markus Wolff
- */
-
-namespace Dumux
-{
-/*!
- * \ingroup IMPES
- */
-//! Class including the variables and data of discretized data of the constitutive relations.
-/*! The variables of two-phase flow, which are one pressure and one saturation are stored in this class.
- * Additionally, a velocity needed in the transport part of the decoupled two-phase flow is stored, as well as discretized data of constitutive relationships like
- * mobilities, fractional flow functions and capillary pressure. Thus, they have to be callculated just once in every time step or every iteration step.
- *
- * @tparam TypeTag The Type Tag
-1*/
-template<class TypeTag>
-class VariableClass2P: public VariableClass<TypeTag>
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
-
- typedef typename GET_PROP_TYPE(TypeTag, TwoPIndices) Indices;
-
- typedef VariableClass<TypeTag> ParentType;
-
- enum
- {
- dim = GridView::dimension, dimWorld = GridView::dimensionworld
- };
- enum
- {
- pw = Indices::pressureW, pn = Indices::pressureNW, pglobal = Indices::pressureGlobal,
- Sw = Indices::saturationW,
- Sn = Indices::saturationNW
- };
- enum
- {
- wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
- };
-
- static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation);
- static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation);
-
- typedef typename SolutionTypes::ElementMapper ElementMapper;
-
- typedef typename GridView::Traits::template Codim<0>::Entity Element;
- typedef Dune::FieldVector<Scalar, dim> LocalPosition;
- typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
- typedef typename GridView::IndexSet IndexSet;
- typedef typename GridView::template Codim<0>::Iterator ElementIterator;
- typedef typename GridView::IntersectionIterator IntersectionIterator;
-
-public:
- typedef typename SolutionTypes::ScalarSolution ScalarSolutionType;//!<type for vector of scalars
- typedef typename SolutionTypes::PhaseProperty PhasePropertyType;//!<type for vector of phase properties
- typedef typename SolutionTypes::FluidProperty FluidPropertyType;//!<type for vector of fluid properties
- typedef typename SolutionTypes::PhasePropertyElemFace PhasePropertyElemFaceType;//!<type for vector of vectors (of size 2 x dimension) of scalars
- typedef typename SolutionTypes::DimVecElemFace DimVecElemFaceType;//!<type for vector of vectors (of size 2 x dimension) of vector (of size dimension) of scalars
-
-private:
- const int codim_;
-
- ScalarSolutionType saturation_;
- PhasePropertyType mobility_;//store lambda for efficiency reasons
- PhasePropertyType fracFlowFunc_;
- ScalarSolutionType capillaryPressure_;
- DimVecElemFaceType velocitySecondPhase_;
-
- FluidPropertyType density_;
- FluidPropertyType viscosity_;
-
- ScalarSolutionType volumecorrection_;
-
-public:
-
- //! Constructs a VariableClass object
- /**
- * @param gridView a DUNE gridview object corresponding to diffusion and transport equation
- * @param initialSat initial value for the saturation (only necessary if only diffusion part is solved)
- * @param initialVel initial value for the velocity (only necessary if only transport part is solved)
- */
-
- VariableClass2P(const GridView& gridView, Scalar& initialSat = *(new Scalar(1)),
- Dune::FieldVector<Scalar, dim>& initialVel = *(new Dune::FieldVector<Scalar, dim>(0))) :
- VariableClass<TypeTag> (gridView, initialVel), codim_(0)
- {
- initializeGlobalVariablesDiffPart(initialVel);
- initializeGlobalVariablesTransPart(initialSat);
- }
-
- //! Constructs a VariableClass object
- /**
- * @param gridView a DUNE gridview object corresponding to diffusion and transport equation
- * @param codim codimension of the entity of which data has to be strored
- * @param initialSat initial value for the saturation (only necessary if only diffusion part is solved)
- * @param initialVel initial value for the velocity (only necessary if only transport part is solved)
- */
- VariableClass2P(const GridView& gridView, int codim, Scalar& initialSat = *(new Scalar(1)), Dune::FieldVector<
- Scalar, dim>& initialVel = *(new Dune::FieldVector<Scalar, dim>(0))) :
- VariableClass<TypeTag> (gridView, codim, initialVel), codim_(codim)
- {
- initializeGlobalVariablesDiffPart(initialVel);
- initializeGlobalVariablesTransPart(initialSat);
- }
-
- // serialization methods
- //! Function needed for restart option.
- template<class Restarter>
- void serialize(Restarter &res)
- {
- res.template serializeEntities<0> (*this, this->gridView());
- }
- //! Function needed for restart option.
- template<class Restarter>
- void deserialize(Restarter &res)
- {
- res.template deserializeEntities<0> (*this, this->gridView());
- }
-
- //! Function needed for restart option.
- void serializeEntity(std::ostream &outstream, const Element &element)
- {
- int globalIdx = this->elementMapper().map(element);
- outstream << this->pressure()[globalIdx] << " " << saturation_[globalIdx];
- }
- //! Function needed for restart option.
- void deserializeEntity(std::istream &instream, const Element &element)
- {
- int globalIdx = this->elementMapper().map(element);
- instream >> this->pressure()[globalIdx] >> saturation_[globalIdx];
- }
-
-private:
- void initializeGlobalVariablesDiffPart(Dune::FieldVector<Scalar, dim>& initialVel)
- {
- //resize to grid size
- velocitySecondPhase_.resize(this->gridSize());//depends on pressure
- for (int i=0; i<2; i++) //for both phases
- {
- density_[i].resize(this->gridSize());//depends on pressure
- viscosity_[i].resize(this->gridSize());//depends on pressure
- }
- //initialise variables
- velocitySecondPhase_ = initialVel;
- }
- void initializeGlobalVariablesTransPart(Scalar initialSat)
- {
- //resize to grid size
- saturation_.resize(this->gridSize());
- for (int i=0; i<2; i++) //for both phases
- {
- mobility_[i].resize(this->gridSize());//lambda is dependent on saturation! ->choose same size
- fracFlowFunc_[i].resize(this->gridSize());//depends on saturation
- }
- capillaryPressure_.resize(this->gridSize());//depends on saturation
- volumecorrection_.resize(this->gridSize());//dS/dt for correction of pressure equation
-
- //initialise variables
- saturation_ = initialSat;
- }
-
-public:
- //Write saturation and pressure into file
- //! adds variables to output
- template<class MultiWriter>
- void addOutputVtkFields(MultiWriter &writer)
- {
- if (codim_ == 0)
- {
- ScalarSolutionType *pressure = writer.allocateManagedBuffer (this->gridSize());
- ScalarSolutionType *saturation = writer.allocateManagedBuffer (this->gridSize());
-
- *pressure = this->pressure();
- *saturation = saturation_;
-
- if (pressureType_ == pw)
- {
- writer.attachCellData(*pressure, "wetting pressure");
- }
- if (pressureType_ == pn)
- {
- writer.attachCellData(*pressure, "nonwetting pressure");
- }
- if (pressureType_ == pglobal)
- {
- writer.attachCellData(*pressure, "global pressure");
- }
-
- if (saturationType_ == Sw)
- {
- writer.attachCellData(*saturation, "wetting saturation");
- }
- if (saturationType_ == Sn)
- {
- writer.attachCellData(*saturation, "nonwetting saturation");
- }
-
- // output phase-dependent stuff
- ScalarSolutionType *pC = writer.allocateManagedBuffer (this->gridSize());
- *pC = capillaryPressure_;
- writer.attachCellData(*pC, "capillary pressure");
-
- ScalarSolutionType *densityWetting = writer.allocateManagedBuffer (this->gridSize());
- *densityWetting = density_[wPhaseIdx];
- writer.attachCellData(*densityWetting, "wetting density");
-
- ScalarSolutionType *densityNonwetting = writer.allocateManagedBuffer (this->gridSize());
- *densityNonwetting = density_[nPhaseIdx];
- writer.attachCellData(*densityNonwetting, "nonwetting density");
-
- ScalarSolutionType *viscosityWetting = writer.allocateManagedBuffer (this->gridSize());
- *viscosityWetting = viscosity_[wPhaseIdx];
- writer.attachCellData(*viscosityWetting, "wetting viscosity");
-
- ScalarSolutionType *viscosityNonwetting = writer.allocateManagedBuffer (this->gridSize());
- *viscosityNonwetting = viscosity_[nPhaseIdx];
- writer.attachCellData(*viscosityNonwetting, "nonwetting viscosity");
- }
- if (codim_ == dim)
- {
- ScalarSolutionType *pressure = writer.allocateManagedBuffer (this->gridSize());
- ScalarSolutionType *saturation = writer.allocateManagedBuffer (this->gridSize());
-
- *pressure = this->pressure();
- *saturation = saturation_;
-
- writer.attachVertexData(*pressure, "pressure");
- writer.attachVertexData(*saturation, "saturation");
- }
-
- return;
- }
-
- ////////////////////////////////////////
- // Adapt variable size of all variables
- ////////////////////////////////////////
-
- void adaptVariableSize(int size)
- {
- //resize to grid size
- this->setGridSize(size);
- for (int i=0; i<2; i++) //for both phases
- {
- density_[i].resize(size);
- viscosity_[i].resize(size);
- mobility_[i].resize(size);
- fracFlowFunc_[i].resize(size);
- }
- velocitySecondPhase_.resize(size);
- velocitySecondPhase_ = Dune::FieldVector<Scalar, dim>(0);
- saturation_.resize(size);
- capillaryPressure_.resize(size);
- volumecorrection_.resize(size);
- ParentType::adaptVariableSize(size); //Also adapt pressure, velocity and potential
- }
-
-
- ////////////////////////////////////////////////////////////
- // functions returning the vectors of the primary variables
- ////////////////////////////////////////////////////////////
-
- void communicateTransportedQuantity()
- {
-#if HAVE_MPI
- ElementHandleAssign<typename ScalarSolutionType::block_type, ScalarSolutionType, ElementMapper> elementHandle(saturation_, this->elementMapper());
- this->gridView().communicate(elementHandle,
- Dune::InteriorBorder_All_Interface,
- Dune::ForwardCommunication);
-#endif
- }
-
- //! Return the vector of the transported quantity, which is the saturation for an IMPES scheme
- ScalarSolutionType& transportedQuantity()
- {
- return saturation_;
- }
-
- //! Return saturation vector
- ScalarSolutionType& saturation()
- {
- return saturation_;
- }
-
- const ScalarSolutionType& saturation() const
- {
- return saturation_;
- }
-
- //////////////////////////////////////////////////////////////
- // functions returning the vectors of secondary variables
- //////////////////////////////////////////////////////////////
-
- //! Return velocity vector
- DimVecElemFaceType& velocitySecondPhase()
- {
- return velocitySecondPhase_;
- }
-
- const DimVecElemFaceType& velocitySecondPhase() const
- {
- return velocitySecondPhase_;
- }
-
- //! Return vector of wetting phase mobilities
- ScalarSolutionType& mobilityWetting()
- {
- return mobility_[wPhaseIdx];
- }
-
- const ScalarSolutionType& mobilityWetting() const
- {
- return mobility_[wPhaseIdx];
- }
-
- //! Return vector of non-wetting phase mobilities
- ScalarSolutionType& mobilityNonwetting()
- {
- return mobility_[nPhaseIdx];
- }
-
- const ScalarSolutionType& mobilityNonwetting() const
- {
- return mobility_[nPhaseIdx];
- }
-
- //! Return vector of wetting phase fractional flow functions
- ScalarSolutionType& fracFlowFuncWetting()
- {
- return fracFlowFunc_[wPhaseIdx];
- }
-
- const ScalarSolutionType& fracFlowFuncWetting() const
- {
- return fracFlowFunc_[wPhaseIdx];
- }
-
- //! Return vector of non-wetting phase fractional flow functions
- ScalarSolutionType& fracFlowFuncNonwetting()
- {
- return fracFlowFunc_[nPhaseIdx];
- }
-
- const ScalarSolutionType& fracFlowFuncNonwetting() const
- {
- return fracFlowFunc_[nPhaseIdx];
- }
-
- //! Return capillary pressure vector
- ScalarSolutionType& capillaryPressure()
- {
- return capillaryPressure_;
- }
-
- const ScalarSolutionType& capillaryPressure() const
- {
- return capillaryPressure_;
- }
-
- //! Return density vector
- ScalarSolutionType& densityWetting()
- {
- return density_[wPhaseIdx];
- }
- ScalarSolutionType& densityWetting() const
- {
- return density_[wPhaseIdx];
- }
-
- //! Return density vector
- ScalarSolutionType& densityNonwetting()
- {
- return density_[nPhaseIdx];
- }
-
- const ScalarSolutionType& densityNonwetting() const
- {
- return density_[nPhaseIdx];
- }
-
- //! Return density vector
- ScalarSolutionType& viscosityWetting()
- {
- return viscosity_[wPhaseIdx];
- }
-
- const ScalarSolutionType& viscosityWetting() const
- {
- return viscosity_[wPhaseIdx];
- }
-
- //! Return density vector
- ScalarSolutionType& viscosityNonwetting()
- {
- return viscosity_[nPhaseIdx];
- }
-
- const ScalarSolutionType& viscosityNonwetting() const
- {
- return viscosity_[nPhaseIdx];
- }
-
- ////////////////////////////////////////////////////////////////////////
- // functions returning entries of the vectors of secondary variables
- ////////////////////////////////////////////////////////////////////////
-
- //! Return vector of wetting phase potential gradients
- Scalar& potentialWetting(int Idx1, int Idx2)
- {
- return this->potential(Idx1, Idx2)[wPhaseIdx];
- }
-
- const Scalar& potentialWetting(int Idx1, int Idx2) const
- {
- return this->potential(Idx1, Idx2)[wPhaseIdx];
- }
-
-
- //! Return vector of non-wetting phase potential gradients
- Scalar& potentialNonwetting(int Idx1, int Idx2)
- {
- return this->potential(Idx1, Idx2)[nPhaseIdx];
- }
-
- const Scalar& potentialNonwetting(int Idx1, int Idx2) const
- {
- return this->potential(Idx1, Idx2)[nPhaseIdx];
- }
-
-
- //! Return vector of wetting phase mobilities
- Scalar& mobilityWetting(int Idx)
- {
- return mobility_[wPhaseIdx][Idx][0];
- }
-
- const Scalar& mobilityWetting(int Idx) const
- {
- return mobility_[wPhaseIdx][Idx][0];
- }
-
-
- //! Return vector of non-wetting phase mobilities
- Scalar& mobilityNonwetting(int Idx)
- {
- return mobility_[nPhaseIdx][Idx][0];
- }
-
- const Scalar& mobilityNonwetting(int Idx) const
- {
- return mobility_[nPhaseIdx][Idx][0];
- }
-
-
- //! Return vector of wetting phase fractional flow functions
- Scalar& fracFlowFuncWetting(int Idx)
- {
- return fracFlowFunc_[wPhaseIdx][Idx][0];
- }
-
- const Scalar& fracFlowFuncWetting(int Idx) const
- {
- return fracFlowFunc_[wPhaseIdx][Idx][0];
- }
-
-
- //! Return vector of non-wetting phase fractional flow functions
- Scalar& fracFlowFuncNonwetting(int Idx)
- {
- return fracFlowFunc_[nPhaseIdx][Idx][0];
- }
-
- const Scalar& fracFlowFuncNonwetting(int Idx) const
- {
- return fracFlowFunc_[nPhaseIdx][Idx][0];
- }
-
-
- //! Return capillary pressure vector
- Scalar& capillaryPressure(int Idx)
- {
- return capillaryPressure_[Idx][0];
- }
-
- const Scalar& capillaryPressure(int Idx) const
- {
- return capillaryPressure_[Idx][0];
- }
-
-
- //! Return density vector
- Scalar& densityWetting(int Idx)
- {
- return density_[wPhaseIdx][Idx][0];
- }
- const Scalar& densityWetting(int Idx) const
- {
- return density_[wPhaseIdx][Idx][0];
- }
-
-
- //! Return density vector
- Scalar& densityNonwetting(int Idx)
- {
- return density_[nPhaseIdx][Idx][0];
- }
-
- const Scalar& densityNonwetting(int Idx) const
- {
- return density_[nPhaseIdx][Idx][0];
- }
-
-
- //! Return density vector
- Scalar& viscosityWetting(int Idx)
- {
- return viscosity_[wPhaseIdx][Idx][0];
- }
-
- const Scalar& viscosityWetting(int Idx) const
- {
- return viscosity_[wPhaseIdx][Idx][0];
- }
-
-
- //! Return density vector
- Scalar& viscosityNonwetting(int Idx)
- {
- return viscosity_[nPhaseIdx][Idx][0];
- }
-
- const Scalar& viscosityNonwetting(int Idx) const
- {
- return viscosity_[nPhaseIdx][Idx][0];
- }
-
-
- ScalarSolutionType& volumecorrection()
- {
- return volumecorrection_;
- }
-
- const ScalarSolutionType& volumecorrection() const
- {
- return volumecorrection_;
- }
-
- Scalar& volumecorrection(int Idx)
- {
- return volumecorrection_[Idx][0];
- }
-
- const Scalar& volumecorrection(int Idx) const
- {
- return volumecorrection_[Idx][0];
- }
-
- //! Get saturation
- /*! evaluate saturation at given element
- @param element entity of codim 0
- \return value of saturation
- */
- Dune::FieldVector<Scalar, 1>& satElement(const Element& element)
- {
- return saturation_[this->elementMapper().map(element)];
- }
-
- const Dune::FieldVector<Scalar, 1>& satElement(const Element& element) const
- {
- return saturation_[this->elementMapper().map(element)];
- }
-
- //! Get velocity at given element face
- /*! evaluate velocity at given location
- @param element entity of codim 0
- @param indexInInside index in reference element
- \return vector of velocity
- */
- Dune::FieldVector<Scalar, dim>& velocitySecondPhaseElementFace(const Element& element, const int indexInInside)
- {
- int elemId = this->index(element);
-
- return (velocitySecondPhase_[elemId][indexInInside]);
- }
-
- const Dune::FieldVector<Scalar, dim>& velocitySecondPhaseElementFace(const Element& element, const int indexInInside) const
- {
- int elemId = this->index(element);
-
- return (velocitySecondPhase_[elemId][indexInInside]);
- }
-
-};
-}
-#endif
diff --git a/dumux/decoupled/common/decoupledproperties_old.hh b/dumux/decoupled/common/decoupledproperties_old.hh
deleted file mode 100644
index 7c90f18d630633830a01ac25dfbee3409b510374..0000000000000000000000000000000000000000
--- a/dumux/decoupled/common/decoupledproperties_old.hh
+++ /dev/null
@@ -1,234 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2009 by Markus Wolff *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-#ifndef DUMUX_DECOUPLED_PROPERTIES_HH
-#define DUMUX_DECOUPLED_PROPERTIES_HH
-
-#include <dumux/common/propertysystem.hh>
-#include <dumux/common/basicproperties.hh>
-#include <dumux/linear/linearsolverproperties.hh>
-
-/*!
- * \ingroup Sequential
- * \ingroup Properties
- */
-/*!
- * \file
- * \brief Base file for properties related to sequential (decoupled) models
- */
-namespace Dumux
-{
-namespace Properties
-{
-
-//////////////////////////////////////////////////////////////////
-// Type tags
-//////////////////////////////////////////////////////////////////
-
-//! Create a type tag for all decoupled models
-NEW_TYPE_TAG(DecoupledModel, INHERITS_FROM(ExplicitModel, LinearSolverTypeTag));
-
-//////////////////////////////////////////////////////////////////
-// Property tags
-//////////////////////////////////////////////////////////////////
-
-//! Property tag for types associated with the solution of the PDE.
-//! This means vectors of primary variables, solution functions on the
-//! grid, and elements, and shape functions.
-NEW_PROP_TAG( SolutionTypes);
-NEW_PROP_TAG( TransportSolutionType);
-NEW_PROP_TAG( PrimaryVariables);
-
-NEW_PROP_TAG( Grid); //!< The type of the DUNE grid
-NEW_PROP_TAG( GridView); //!< The type of the grid view
-NEW_PROP_TAG( AdaptiveGrid); //!< Defines if the grid is h-adaptive
-
-NEW_PROP_TAG( Problem); //!< The type of the problem
-NEW_PROP_TAG( Model); //!< The type of the discretizations
-NEW_PROP_TAG( NumEq ); //!< Number of equations in the system of PDEs
-NEW_PROP_TAG( NumPhases); //!< Number of phases in the system
-NEW_PROP_TAG( NumComponents); //!< Number of components in the system
-NEW_PROP_TAG( Variables); //!< The type of the container of global variables
-NEW_PROP_TAG(TimeManager); //!< Manages the simulation time
-NEW_PROP_TAG(BoundaryTypes); //!< Stores the boundary types of a single degree of freedom
-}
-}
-
-#include <dune/grid/common/mcmgmapper.hh>
-#include <dune/istl/bvector.hh>
-
-#include <dumux/common/timemanager.hh>
-#include <dumux/common/boundarytypes.hh>
-#include<dumux/common/boundaryconditions.hh>
-
-template<class TypeTag>
-class VariableClass;
-
-namespace Dumux
-{
-namespace Properties
-{
-//////////////////////////////////////////////////////////////////
-// Properties
-//////////////////////////////////////////////////////////////////
-
-//! Use the leaf grid view if not defined otherwise
-SET_PROP(DecoupledModel, GridView)
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
-
-public:
- typedef typename Grid::LeafGridView type;
-};
-
-//! Disables Grid Adaptivity as standard
-SET_BOOL_PROP(DecoupledModel, AdaptiveGrid, false);
-
-//! Use the parent VariableClass
-SET_TYPE_PROP(DecoupledModel, Variables, VariableClass<TypeTag>);
-
-NEW_PROP_TAG(MaxIntersections); //!< maximum number of intersections per element
-SET_PROP(DecoupledModel, MaxIntersections)
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
- static const int dim = Grid::dimension;
-public:
- typedef int type;
- static const int value = 2*dim;
-};
-
-/*!
- * \brief Specifies the types which are assoicated with a solution.
- *
- * This means shape functions, solution vectors, etc.
- */
-SET_PROP(DecoupledModel, SolutionTypes)
-{
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
- typedef typename GridView::Grid Grid;
- typedef typename Grid::ctype CoordScalar;
- typedef typename GET_PROP_TYPE(TypeTag, Variables) Variables;
-
- enum
- {
- dim = GridView::dimension,
- numEq = GET_PROP_VALUE(TypeTag, NumEq),
- numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
- numComponents = GET_PROP_VALUE(TypeTag, NumComponents),
- maxIntersections = GET_PROP_VALUE(TypeTag, MaxIntersections)
- };
-
- template<int dim>
- struct VertexLayout
- {
- bool contains (Dune::GeometryType gt) const
- { return gt.dim() == 0;}
- };
-
- template<int dim>
- struct ElementLayout
- {
- bool contains (Dune::GeometryType gt) const
- { return gt.dim() == dim;}
- };
-
-public:
- /*!
- * \brief Mapper for the grid view's vertices.
- */
- typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, VertexLayout> VertexMapper;
-
- /*!
- * \brief Mapper for the grid view's elements.
- */
- typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, ElementLayout> ElementMapper;
-
- /*!
- * \brief The type of a solution at a fixed time.
- *
- * This defines the primary and secondary variable vectors at each degree of freedom.
- */
- typedef Dune::FieldVector<Scalar, numEq> PrimaryVariables;
- typedef Dune::BlockVector<Dune::FieldVector<Scalar, 1> > ScalarSolution;//!<type for vector of scalars
- typedef Dune::FieldVector<Dune::BlockVector<Dune::FieldVector<Scalar,1> >, numComponents> ComponentProperty;//!<type for vector of phase properties
- typedef Dune::FieldVector<Dune::BlockVector<Dune::FieldVector<Scalar,1> >, numPhases> PhaseProperty;//!<type for vector of phase properties
- typedef Dune::FieldVector<Dune::BlockVector<Dune::FieldVector<Scalar,1> >, numPhases> FluidProperty;//!<type for vector of fluid properties: Vector[element][phase]
- typedef Dune::BlockVector<Dune::FieldVector<Dune::FieldVector<Scalar, numPhases>, maxIntersections > > PhasePropertyElemFace;//!<type for vector of vectors (of size 2 x dimension) of scalars
- typedef Dune::BlockVector<Dune::FieldVector<Dune::FieldVector<Scalar, dim>, maxIntersections > > DimVecElemFace;//!<type for vector of vectors (of size 2 x dimension) of vector (of size dimension) of scalars
-};
-
-SET_TYPE_PROP(DecoupledModel, PrimaryVariables, typename GET_PROP(TypeTag, SolutionTypes)::PrimaryVariables);
-
-/*!
- * \brief Default implementation for the Vector of the transportet quantity
- *
- * This type defines the data type of the transportet quantity. In case of a
- * immiscible 2p system, this would represent a vector holding the saturation
- * of one phase.
- */
-SET_PROP(DecoupledModel, TransportSolutionType)
-{
- private:
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionType;
-
- public:
- typedef typename SolutionType::ScalarSolution type;//!<type for vector of scalar properties
-};
-
-//! Set the default type for the time manager
-SET_TYPE_PROP(DecoupledModel, TimeManager, Dumux::TimeManager<TypeTag>);
-
-/*!
- * \brief Boundary types at a single degree of freedom.
- */
-SET_PROP(DecoupledModel, BoundaryTypes)
-{ private:
- enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
-public:
- typedef Dumux::BoundaryTypes<numEq> type;
-};
-
-//! set the default for the reduction of the initial residual
-SET_PROP(DecoupledModel, LinearSolverResidualReduction)
-{public:
- static constexpr double value = 1e-13;
-};
-
-//! set the default number of maximum iterations for the linear solver
-SET_PROP(DecoupledModel, LinearSolverMaxIterations)
-{public:
- static constexpr int value = 500;
-};
-
-//! set the default number of maximum iterations for the linear solver
-SET_PROP(DecoupledModel, LinearSolverBlockSize)
-{public:
- static constexpr int value = 1;
-};
-// \}
-
-}
-}
-
-#endif
diff --git a/dumux/decoupled/common/gridadapt_old.hh b/dumux/decoupled/common/gridadapt_old.hh
deleted file mode 100644
index 54539f362a5b9e7444a7798543a37614a28a43a1..0000000000000000000000000000000000000000
--- a/dumux/decoupled/common/gridadapt_old.hh
+++ /dev/null
@@ -1,430 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2010 by Benjamin Faigle *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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
- * \brief Base class for h-adaptive sequential models.
- */
-#ifndef DUMUX_GIRDADAPT_HH
-#define DUMUX_GIRDADAPT_HH
-
-#include "decoupledproperties.hh"
-
-namespace Dumux
-{
-
-/*!
- * @brief Standard Module for h-adaptive simulations
- *
- * This class is created by the problem class with the template
- * parameters <TypeTag, true> and provides basic functionality
- * for adaptive methods:
- *
- * A standard implementation adaptGrid() will prepare everything
- * to calculate the next pressure field on the new grid.
- */
-template<class TypeTag, bool adaptive>
-class GridAdapt
-{
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
-
-
- //*******************************
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
- typedef typename GridView::Grid Grid;
- typedef typename Grid::LeafGridView LeafGridView;
- typedef typename LeafGridView::template Codim<0>::Iterator LeafIterator;
- typedef typename GridView::IntersectionIterator LeafIntersectionIterator;
- typedef typename Grid::template Codim<0>::Entity Entity;
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
- typedef typename SolutionTypes::ScalarSolution ScalarSolutionType;
-
-public:
- /*!
- * Constructor for h-adaptive simulations (adaptive grids)
- * @param problem The problem
- * @param levelMin minimum refinement level
- * @param levelMax maximum refinement level
- */
- GridAdapt (Problem& problem, const int levelMin = 0,const int levelMax=1)
- : levelMin_(levelMin), levelMax_(levelMax), problem_(problem)
- {
- if (levelMin_ < 0)
- Dune::dgrave << __FILE__<< ":" <<__LINE__
- << " : Dune cannot coarsen to gridlevels smaller 0! "<< std::endl;
-
- standardIndicator_ = true;
- refinetol_ = 0.05;
- coarsentol_ = 0.001;
- }
-
- /*!
- * @brief Standard method to adapt the grid
- *
- * This method is called in preTimeStep() of the problem if
- * adaptive grids are used in the simulation.
- *
- * It uses a standard procedure for adaptivity:
- * 1) Determine the refinement indicator
- * 2) Mark the elements
- * 3) Store primary variables in a map
- * 4) Adapt the grid, adapt variables sizes, update mappers
- * 5) Reconstruct primary variables, regain secondary variables
- */
- void adaptGrid()
- {
- // reset internal counter for marked elements
- marked_ = coarsened_ = 0;
-
- // prepare an indicator for refinement
- if(indicatorVector_.size()!=problem_.variables().gridSize())
- {
- indicatorVector_.resize(problem_.variables().gridSize());
- indicatorVector_ = -1e100;
- }
-
- /**** 1) determine refining parameter if standard is used ***/
- // if not, the indicatorVector and refinement Bounds have to
- // specified by the problem through setIndicator()
- if(standardIndicator_)
- indicator();
-
- /**** 2) mark elements according to indicator *********/
- markElements(indicatorVector_, refineBound_, coarsenBound_);
-
- // abort if nothing in grid is marked
- if (marked_==0 && coarsened_ == 0)
- return;
- else
- Dune::dinfo << marked_ << " cells have been marked_ to be refined, "
- << coarsened_ << " to be coarsened." << std::endl;
-
- /**** 2b) Do pre-adaption step *****/
- problem_.grid().preAdapt();
- problem_.preAdapt();
-
- /**** 3) Put primary variables in a map *********/
- problem_.variables().storePrimVars(problem_);
-
- /**** 4) Adapt Grid and size of variable vectors *****/
-// problem_.grid().preAdapt();
- problem_.grid().adapt();
-
-// forceRefineRatio(1);
-
- // update mapper to new cell idice
- problem_.variables().elementMapper().update();
-
- // adapt size of vectors (
- problem_.variables().adaptVariableSize(problem_.variables().elementMapper().size());
-
- /**** 5) (Re-)construct primary variables to new grid **/
- problem_.variables().reconstructPrimVars(problem_);
- // delete markers in grid
- problem_.grid().postAdapt();
- // adapt secondary variables
- problem_.pressureModel().updateMaterialLaws();
-
- // write out new grid
-// Dune::VTKWriter<LeafGridView> vtkwriter(leafView);
-// vtkwriter.write("latestgrid",Dune::VTKOptions::binaryappended);
- return;
- };
-
- /*!
- * Mark Elements for grid refinement according to applied Indicator
- * @param indicator Vector where the refinement indicator is stored
- * @param refineThreshold lower threshold where to refine
- * @param coarsenThreshold upper threshold where to coarsen
- * @return Total ammount of marked cells
- */
- int markElements(const ScalarSolutionType &indicator, const double refineThreshold, const double coarsenThreshold)
- {
- for (LeafIterator eIt = problem_.gridView().template begin<0>();
- eIt!=problem_.gridView().template end<0>(); ++eIt)
- {
- // Verfeinern?
- if (indicator[problem_.variables().elementMapper().map(*eIt)] > refineThreshold
- && eIt.level()<levelMax_)
- {
- const Entity &entity =*eIt;
- problem_.grid().mark( 1, entity );
- ++marked_;
-
- // this also refines the neighbor elements
- checkNeighborsRefine_(entity);
- }
- }
- // coarsen
- for (LeafIterator eIt = problem_.gridView().template begin<0>();
- eIt!=problem_.gridView().template end<0>(); ++eIt)
- {
- if (indicator[problem_.variables().elementMapper().map(*eIt)] < coarsenThreshold
- && eIt.level()>levelMin_ && problem_.grid().getMark(*eIt) == 0)
- {
- // check if coarsening is possible
- bool coarsenPossible = true;
- LeafIntersectionIterator isend = problem_.gridView().iend(*eIt);
- for(LeafIntersectionIterator is = problem_.gridView().ibegin(*eIt); is != isend; ++is)
- {
- if(!is->neighbor())
- continue;
-
- const Entity &outside = *is->outside();
- if ((problem_.grid().getMark(outside) > 0)
- || (outside.level()>eIt.level()))
- coarsenPossible = false;
- }
-
- if(coarsenPossible)
- {
- problem_.grid().mark( -1, *eIt );
- ++coarsened_;
- }
- }
- }
-
- return marked_;
- }
-
- /*!
- * Sets minimum and maximum refinement levels
- *
- * @param levMin minimum level for coarsening
- * @param levMax maximum level for refinement
- */
- void setLevels(int levMin, int levMax)
- {
- if (levMin < 0)
- Dune::dgrave << __FILE__<< ":" <<__LINE__
- << " : Dune cannot coarsen to gridlevels smaller 0! "<< std::endl;
- levelMin_ = levMin;
- levelMax_ = levMax;
- }
- /*!
- * Sets minimum and maximum refinement tolerances
- *
- * @param coarsentol minimum tolerance when to coarsen
- * @param refinetol maximum tolerance when to refine
- */
- void setTolerance(Scalar coarsentol, Scalar refinetol)
- {
- if (coarsentol < 0. or refinetol < 0. )
- Dune::dgrave << __FILE__<< ":" <<__LINE__
- << " : Tolerance levels out of meaningful bounds! "<< std::endl;
- if (coarsentol > refinetol )
- Dune::dgrave << __FILE__<< ":" <<__LINE__
- << " : Check tolerance levels: coarsentol > refinetol! "<< std::endl;
-
- coarsentol_ = coarsentol;
- refinetol_ = refinetol;
- }
- /*!
- * Gets maximum refinement level
- *
- * @return levelMax_ maximum level for refinement
- */
- const int getMaxLevel() const
- {
- return levelMax_;
- }
- /*!
- * Gets minimum refinement level
- *
- * @return levelMin_ minimum level for coarsening
- */
- const int getMinLevel() const
- {
- return levelMin_;
- }
- /*!
- * @brief Adapter for external Refinement indicators.
- *
- * External indicators to refine/coarsen the grid can be set
- * from outside this container. Both indicator itsself as well as
- * the coarsening and refinement bounds have to be specified.
- * @param indicatorVector Vector holding indicator values
- * @param coarsenLowerBound bounding value where to be coarsened
- * @param refineUpperBound bounding value where to be refined
- */
- const void setIndicator(const ScalarSolutionType& indicatorVector,
- const Scalar& coarsenLowerBound, const Scalar& refineUpperBound)
- {
- // switch off usage of standard (saturation) indicator: by settting this,
- // the standard function indicator() called by adaptGrid() is disabled.
- standardIndicator_=false;
-
- indicatorVector_ = indicatorVector;
- refineBound_ = refineUpperBound;
- coarsenBound_ = coarsenLowerBound;
- }
-private:
- /*!
- * @brief Simple standard indicator.
- *
- * Mehod computes the refinement and coarsening bounds through a
- * standard refinement criteria.
- */
- void indicator()
- {
- Scalar globalMax = -1e100;
- Scalar globalMin = 1e100;
-
- /**** determine refining parameter *************/
- problem_.transportModel().indicatorSaturation(indicatorVector_, globalMin, globalMax);
- Scalar globaldelta = globalMax- globalMin;
-// globaldelta = std::max(globaldelta,0.1);
-
- refineBound_ = refinetol_*globaldelta;
- coarsenBound_ = coarsentol_*globaldelta;
-
- return;
- }
- /*!
- * @brief Method ensuring the refinement ratio of 2:1
- *
- * For any given entity, a loop over the neighbors checks weather the
- * entities refinement would require that any of the neighbors has
- * to be refined, too.
- * This is done recursively over all levels of the grid.
- *
- * @param entity Entity of interest that is to be refined
- * @param level level of the refined entity: it is at least 1
- * @return true if everything was successful
- */
- bool checkNeighborsRefine_(const Entity &entity, int level = 1)
- {
- // this also refines the neighbor elements
- LeafIntersectionIterator isend = problem_.gridView().iend(entity);
- for(LeafIntersectionIterator is = problem_.gridView().ibegin(entity); is != isend; ++is)
- {
- if(!is->neighbor())
- continue;
-
- const Entity &outside =*is->outside();
- if ((outside.level()<levelMax_)
- && (outside.level()<entity.level()))
- {
- problem_.grid().mark(1, outside);
- ++marked_;
-
- if(level != levelMax_)
- checkNeighborsRefine_(outside, ++level);
- }
- }
- return true;
- };
-
-
- /*!
- * \brief Enforces a given refine ratio after grid was adapted
- *
- * If the refine ratio is not taken into consideration during
- * marking, then this method ensures a certain ratio.
- *
- * @param maxLevelDelta The maximum level difference (refine ratio)
- * between neighbors.
- */
- void forceRefineRatio(int maxLevelDelta = 1)
- {
- LeafGridView leafView = problem_.gridView();
- // delete all existing marks
- problem_.grid().postAdapt();
- bool done;
- do
- {
- // run through all cells
- done=true;
- for (LeafIterator it = leafView.template begin<0>();
- it!=leafView.template end<0>(); ++it)
- {
- // run through all neighbor-cells (intersections)
- LeafIntersectionIterator isend = leafView.iend(*it);
- for (LeafIntersectionIterator is = leafView.ibegin(*it); is!= isend; ++is)
- {
- const typename LeafIntersectionIterator::Intersection intersection = *is;
- if(!intersection.neighbor())
- continue;
-
- const Entity &outside =intersection.outside();
- if (it.level()+maxLevelDelta<outside.level())
- {
- const Entity &entity =*it;
- problem_.grid().mark( 1, entity );
- done=false;
- }
- }
- }
- if (done==false)
- {
- // adapt the grid
- problem_.grid().adapt();
- // delete marks
- problem_.grid().postAdapt();
- }
- }
- while (done!=true);
- }
-
- // private Variables
- ScalarSolutionType indicatorVector_;
- bool standardIndicator_;
- Scalar refineBound_, coarsenBound_;
- int marked_, coarsened_;
- int levelMin_, levelMax_;
- Scalar refinetol_;
- Scalar coarsentol_;
-
- Problem& problem_;
-};
-
-/*!
- * @brief Class for NON-adaptive simulations
- *
- * This class provides empty methods for non-adaptive simulations
- * for compilation reasons. If adaptivity is desired, create the
- * class with template arguments <TypeTag, true> instead.
- */
-template<class TypeTag>
-class GridAdapt<TypeTag, false>
-{
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
- typedef typename SolutionTypes::ScalarSolution ScalarSolutionType;
-
-public:
- void adaptGrid()
- {};
- void setLevels(int, int)
- {};
- void setTolerance(int, int)
- {};
- const void setIndicator(const ScalarSolutionType&,
- const Scalar&, const Scalar&)
- {};
- GridAdapt (Problem& problem)
- {}
-};
-
-}
-#endif /* DUMUX_GIRDADAPT_HH */
diff --git a/dumux/decoupled/common/onemodelproblem_old.hh b/dumux/decoupled/common/onemodelproblem_old.hh
deleted file mode 100644
index 868e980265aeb0c5fec9dc9038692bde5fe6f988..0000000000000000000000000000000000000000
--- a/dumux/decoupled/common/onemodelproblem_old.hh
+++ /dev/null
@@ -1,663 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2010 by Markus Wolff *
- * Copyright (C) 2009 by Andreas Lauser *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-
-#ifndef DUMUX_ONE_MODEL_PROBLEM_HH
-#define DUMUX_ONE_MODEL_PROBLEM_HH
-
-#include <dumux/decoupled/common/decoupledproperties_old.hh>
-#include <dumux/io/vtkmultiwriter.hh>
-#include <dumux/io/restart.hh>
-
-
-/**
- * @file
- * @brief Base class for definition of an decoupled diffusion (pressure) or transport problem
- * @author Markus Wolff
- */
-
-namespace Dumux
-{
-
-/*! \ingroup IMPET
- *
- * @brief Base class for definition of an decoupled diffusion (pressure) or transport problem
- *
- * @tparam TypeTag The Type Tag
- */
-template<class TypeTag>
-class OneModelProblem
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
-
- typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
-
- typedef Dumux::VtkMultiWriter<GridView> VtkMultiWriter;
-
- typedef typename GET_PROP_TYPE(TypeTag, Variables) Variables;
-
- typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
-
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
-
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
- typedef typename SolutionTypes::VertexMapper VertexMapper;
- typedef typename SolutionTypes::ElementMapper ElementMapper;
- typedef typename SolutionTypes::ScalarSolution Solution;
-
- enum
- {
- dim = GridView::dimension,
- dimWorld = GridView::dimensionworld
- };
- enum
- {
- wetting = 0, nonwetting = 1
- };
-
- typedef Dune::FieldVector<Scalar,dimWorld> GlobalPosition;
- typedef typename GridView::template Codim<dim>::Iterator VertexIterator;
- typedef typename GridView::Traits::template Codim<0>::Entity Element;
- typedef typename GridView::Intersection Intersection;
-
- typedef typename SolutionTypes::PrimaryVariables PrimaryVariables;
- typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
-
- // private!! copy constructor
- OneModelProblem(const OneModelProblem&)
- {}
-
-public:
-
- //! Constructs an object of type OneModelProblemProblem
- /*!
- * \tparam TypeTag The TypeTag
- * \tparam verbose Output level for Dumux::TimeManager
- */
- OneModelProblem(const GridView &gridView, bool verbose = true)
- : gridView_(gridView),
- bboxMin_(std::numeric_limits<double>::max()),
- bboxMax_(-std::numeric_limits<double>::max()),
- deleteTimeManager_(true),
- variables_(gridView),
- outputInterval_(1)
- {
- // calculate the bounding box of the grid view
- VertexIterator vIt = gridView.template begin<dim>();
- const VertexIterator vEndIt = gridView.template end<dim>();
- for (; vIt!=vEndIt; ++vIt) {
- for (int i=0; i<dim; i++) {
- bboxMin_[i] = std::min(bboxMin_[i], vIt->geometry().center()[i]);
- bboxMax_[i] = std::max(bboxMax_[i], vIt->geometry().center()[i]);
- }
- }
-
- timeManager_ = new TimeManager(verbose);
-
- model_ = new Model(asImp_()) ;
-
- resultWriter_ = NULL;
- }
-
- //! Constructs an object of type OneModelProblemProblem
- /*!
- * \tparam TypeTag The TypeTag
- * \tparam verbose Output level for Dumux::TimeManager
- */
- OneModelProblem(TimeManager &timeManager, const GridView &gridView)
- : gridView_(gridView),
- bboxMin_(std::numeric_limits<double>::max()),
- bboxMax_(-std::numeric_limits<double>::max()),
- timeManager_(&timeManager),
- deleteTimeManager_(false),
- variables_(gridView),
- outputInterval_(1)
- {
- // calculate the bounding box of the grid view
- VertexIterator vIt = gridView.template begin<dim>();
- const VertexIterator vEndIt = gridView.template end<dim>();
- for (; vIt!=vEndIt; ++vIt) {
- for (int i=0; i<dim; i++) {
- bboxMin_[i] = std::min(bboxMin_[i], vIt->geometry().center()[i]);
- bboxMax_[i] = std::max(bboxMax_[i], vIt->geometry().center()[i]);
- }
- }
-
- model_ = new Model(asImp_()) ;
-
- resultWriter_ = NULL;
- }
-
- //! destructor
- virtual ~OneModelProblem ()
- {
- delete model_;
- delete resultWriter_;
- if (deleteTimeManager_)
- delete timeManager_;
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param bcTypes The boundary types for the conservation equations
- * \param intersection The intersection for which the boundary type is set
- */
- void boundaryTypes(BoundaryTypes &bcTypes,
- const Intersection &intersection) const
- {
- // forward it to the method which only takes the global coordinate
- asImp_().boundaryTypesAtPos(bcTypes, intersection.geometry().center());
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param bcTypes The boundary types for the conservation equations
- * \param globalPos The position of the center of the boundary intersection
- */
- void boundaryTypesAtPos(BoundaryTypes &bcTypes,
- const GlobalPosition &globalPos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a boundaryTypesAtPos() method.");
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a dirichlet
- * control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param intersection The boundary intersection
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void dirichlet(PrimaryVariables &values,
- const Intersection &intersection) const
- {
- // forward it to the method which only takes the global coordinate
- asImp_().dirichletAtPos(values, intersection.geometry().center());
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a dirichlet
- * control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param globalPos The position of the center of the boundary intersection
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void dirichletAtPos(PrimaryVariables &values,
- const GlobalPosition &globalPos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem specifies that some boundary "
- "segments are dirichlet, but does not provide "
- "a dirichletAtPos() method.");
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param intersection The boundary intersection
- *
- * For this method, the \a values parameter stores the mass flux
- * in normal direction of each phase. Negative values mean influx.
- */
- void neumann(PrimaryVariables &values,
- const Intersection &intersection) const
- {
- // forward it to the interface with only the global position
- asImp_().neumannAtPos(values, intersection.geometry().center());
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param globalPos The position of the center of the boundary intersection
- *
- * For this method, the \a values parameter stores the mass flux
- * in normal direction of each phase. Negative values mean influx.
- */
- void neumannAtPos(PrimaryVariables &values,
- const GlobalPosition &globalPos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Neumann conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem specifies that some boundary "
- "segments are neumann, but does not provide "
- "a neumannAtPos() method.");
- }
-
- /*!
- * \brief Evaluate the source term
- *
- * \param values The source and sink values for the conservation equations
- * \param element The element
- *
- * For this method, the \a values parameter stores the rate mass
- * generated or annihilate per volume unit. Positive values mean
- * that mass is created, negative ones mean that it vanishes.
- */
- void source(PrimaryVariables &values,
- const Element &element) const
- {
- // forward to generic interface
- asImp_().sourceAtPos(values, element.geometry().center());
- }
-
- /*!
- * \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
- * \param globalPos The position of the center of the finite volume
- * for which the source term ought to be
- * specified in global coordinates
- *
- * For this method, the \a values parameter stores the rate mass
- * generated or annihilate per volume unit. Positive values mean
- * that mass is created, negative ones mean that it vanishes.
- */
- void sourceAtPos(PrimaryVariables &values,
- const GlobalPosition &globalPos) const
- { // Throw an exception (there is no initial condition)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a sourceAtPos() method.");
- }
-
- /*!
- * \brief Evaluate the initial value for a control volume.
- *
- * \param values The initial values for the primary variables
- * \param element The element
- *
- * For this method, the \a values parameter stores primary
- * variables.
- */
- void initial(PrimaryVariables &values,
- const Element &element) const
- {
- // forward to generic interface
- asImp_().initialAtPos(values, element.geometry().center());
- }
-
- /*!
- * \brief Evaluate the initial value for a control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param globalPos The position of the center of the finite volume
- * for which the initial values ought to be
- * set (in global coordinates)
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void initialAtPos(PrimaryVariables &values,
- const GlobalPosition &globalPos) const
- {
- // Throw an exception (there is no initial condition)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a initialAtPos() method.");
- }
-
- /*!
- * \brief Called by the Dumux::TimeManager in order to
- * initialize the problem.
- */
- void init()
- {
- // set the initial condition of the model
- model().initialize();
- }
-
- /*!
- * \brief Called by Dumux::TimeManager just before the time
- * integration.
- */
- void preTimeStep()
- { };
-
- /*!
- * \brief Called by Dumux::TimeManager in order to do a time
- * integration on the model.
- */
- void timeIntegration()
- { };
-
- /*!
- * \brief Called by Dumux::TimeManager whenever a solution for a
- * timestep has been computed and the simulation time has
- * been updated.
- *
- * This is used to do some janitorial tasks like writing the
- * current solution to disk.
- */
- void postTimeStep()
- { };
-
- /*!
- * \brief Called by the time manager after everything which can be
- * done about the current time step is finished and the
- * model should be prepared to do the next time integration.
- */
- void advanceTimeLevel()
- {}
-
- /*!
- * \brief Returns the current time step size [seconds].
- */
- Scalar timeStepSize() const
- { return timeManager().timeStepSize(); }
-
- /*!
- * \brief Sets the current time step size [seconds].
- */
- void setTimeStepSize(Scalar dt)
- { timeManager().setTimeStepSize(dt); }
-
- /*!
- * \brief Called by Dumux::TimeManager whenever a solution for a
- * timestep has been computed and the simulation time has
- * been updated.
- */
- Scalar nextTimeStepSize(Scalar dt)
- { return timeManager().timeStepSize();}
-
- /*!
- * \brief Returns true if a restart file should be written to
- * disk.
- *
- * The default behaviour is to write one restart file every 5 time
- * steps. This file is intented to be overwritten by the
- * implementation.
- */
- bool shouldWriteRestartFile() const
- {
- return
- timeManager().timeStepIndex() > 0 &&
- (timeManager().timeStepIndex() % 5 == 0);
- }
-
- /*!
- * \brief Sets the interval for Output
- *
- * The default is 1 -> Output every time step
- */
- void setOutputInterval(int interval)
- { outputInterval_ = interval; }
-
- /*!
- * \brief Returns true if the current solution should be written to
- * disk (i.e. as a VTK file)
- *
- * The default behaviour is to write out every the solution for
- * very time step. This file is intented to be overwritten by the
- * implementation.
- */
- bool shouldWriteOutput() const
- {
- if (this->timeManager().timeStepIndex() % outputInterval_ == 0 || this->timeManager().willBeFinished())
- {
- return true;
- }
- return false;
- }
-
- void addOutputVtkFields()
- {}
-
- //! Write the fields current solution into an VTK output file.
- void writeOutput(bool verbose = true)
- {
- if (verbose && gridView().comm().rank() == 0)
- std::cout << "Writing result file for current time step\n";
- if (!resultWriter_)
- resultWriter_ = new VtkMultiWriter(gridView(), asImp_().name());
- resultWriter_->beginWrite(timeManager().time() + timeManager().timeStepSize());
- model().addOutputVtkFields(*resultWriter_);
- asImp_().addOutputVtkFields();
- resultWriter_->endWrite();
- }
-
- /*!
- * \brief Called when the end of an simulation episode is reached.
- */
- void episodeEnd()
- {
- std::cerr << "The end of an episode is reached, but the problem "
- << "does not override the episodeEnd() method. "
- << "Doing nothing!\n";
- };
-
- // \}
-
- /*!
- * \brief The problem name.
- *
- * This is used as a prefix for files generated by the simulation.
- * It could be either overwritten by the problem files, or simply
- * declared over the setName() function in the application file.
- */
- const char *name() const
- {
- return simname_.c_str();
- }
-
- /*!
- * \brief Set the problem name.
- *
- * This function sets the simulation name, which should be called before
- * the application problem is declared! If not, the default name "sim"
- * will be used.
- */
- void setName(const char *newName)
- {
- simname_ = newName;
- }
-
- /*!
- * \brief The GridView which used by the problem.
- */
- const GridView &gridView() const
- { return gridView_; }
-
- /*!
- * \brief Returns the mapper for vertices to indices.
- */
- const VertexMapper &vertexMapper() const
- { return variables_.vertexMapper(); }
-
- /*!
- * \brief Returns the mapper for elements to indices.
- */
- const ElementMapper &elementMapper() const
- { return variables_.elementMapper(); }
-
- /*!
- * \brief The coordinate of the corner of the GridView's bounding
- * box with the smallest values.
- */
- const GlobalPosition &bboxMin() const
- { return bboxMin_; }
-
- /*!
- * \brief The coordinate of the corner of the GridView's bounding
- * box with the largest values.
- */
- const GlobalPosition &bboxMax() const
- { return bboxMax_; }
-
- /*!
- * \brief Returns TimeManager object used by the simulation
- */
- TimeManager &timeManager()
- { return *timeManager_; }
-
- /*!
- * \brief \copybrief Dumux::OneModelProblem::timeManager()
- */
- const TimeManager &timeManager() const
- { return *timeManager_; }
-
- /*!
- * \brief Returns variables object.
- */
- Variables& variables ()
- { return variables_; }
-
- /*!
- * \brief \copybrief Dumux::OneModelProblem::variables()
- */
- const Variables& variables () const
- { return variables_; }
-
- /*!
- * \brief Returns numerical model used for the problem.
- */
- Model &model()
- { return *model_; }
-
- /*!
- * \brief \copybrief Dumux::OneModelProblem::model()
- */
- const Model &model() const
- { return *model_; }
- // \}
-
-
- /*!
- * \name Restart mechanism
- */
- // \{
-
- /*!
- * \brief This method writes the complete state of the problem
- * to the harddisk.
- *
- * The file will start with the prefix returned by the name()
- * method, has the current time of the simulation clock in it's
- * name and uses the extension <tt>.drs</tt>. (Dumux ReStart
- * file.) See Dumux::Restart for details.
- */
- void serialize()
- {
- typedef Dumux::Restart Restarter;
-
- Restarter res;
- res.serializeBegin(asImp_());
- std::cerr << "Serialize to file " << res.fileName() << "\n";
-
- timeManager().serialize(res);
- resultWriter().serialize(res);
- model().serialize(res);
-
- res.serializeEnd();
- }
-
- /*!
- * \brief This method restores the complete state of the problem
- * from disk.
- *
- * It is the inverse of the serialize() method.
- */
- void restart(double tRestart)
- {
- typedef Dumux::Restart Restarter;
-
- Restarter res;
- res.deserializeBegin(asImp_(), tRestart);
- std::cerr << "Deserialize from file " << res.fileName() << "\n";
-
- timeManager().deserialize(res);
- resultWriter().deserialize(res);
- model().deserialize(res);
-
- res.deserializeEnd();
- };
- //! Use restart() function instead!
- void deserialize(double tRestart) DUNE_DEPRECATED
- { restart(tRestart);}
-
- // \}
-
-protected:
- VtkMultiWriter& resultWriter()
- {
- if (!resultWriter_)
- resultWriter_ = new VtkMultiWriter(gridView_, asImp_().name());
- return *resultWriter_;
- }
-
- VtkMultiWriter& resultWriter() const
- {
- if (!resultWriter_)
- resultWriter_ = new VtkMultiWriter(gridView_, asImp_().name());
- return *resultWriter_;
- }
-
-private:
- //! Returns the implementation of the problem (i.e. static polymorphism)
- Implementation &asImp_()
- { return *static_cast<Implementation *>(this); }
-
- //! \brief \copybrief Dumux::OneModelProblem::asImp_()
- const Implementation &asImp_() const
- { return *static_cast<const Implementation *>(this); }
-
- std::string simname_; // a string for the name of the current simulation,
- // which could be set by means of an program argument,
- // for example.
- const GridView gridView_;
-
- GlobalPosition bboxMin_;
- GlobalPosition bboxMax_;
-
- TimeManager *timeManager_;
- bool deleteTimeManager_;
-
- Variables variables_;
-
- Model* model_;
-
- VtkMultiWriter *resultWriter_;
- int outputInterval_;
-};
-
-}
-#endif
diff --git a/dumux/decoupled/common/variableclass_old.hh b/dumux/decoupled/common/variableclass_old.hh
deleted file mode 100644
index 7d4b150df275680b4bc3a606e1f424f180a6b6ed..0000000000000000000000000000000000000000
--- a/dumux/decoupled/common/variableclass_old.hh
+++ /dev/null
@@ -1,375 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * Copyright (C) 2009 by Markus Wolff *
- * Institute for Modelling Hydraulic and Environmental Systems *
- * University of Stuttgart, Germany *
- * email: <givenname>.<name>@iws.uni-stuttgart.de *
- * *
- * 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/>. *
- *****************************************************************************/
-#ifndef DUMUX_VARIABLECLASS_HH
-#define DUMUX_VARIABLECLASS_HH
-
-//#define HACK_SINTEF_RESPROP
-
-#include <dune/istl/bvector.hh>
-#include <dumux/io/vtkmultiwriter.hh>
-#include "decoupledproperties_old.hh"
-
-// for parallelization
-#include <dumux/parallel/elementhandles.hh>
-
-/**
- * @file
- * @brief Base class holding the variables for sequential models.
- * @author Markus Wolff
- */
-
-namespace Dumux
-{
-/*!
- * \ingroup Sequential
- */
-//! Base class holding the variables and discretized data for sequential models.
-/*!
- * Stores global information and variables that are common for all sequential models and also functions needed to access these variables.
- * Can be directly used for a single phase model.
- *
- * @tparam TypeTag The Type Tag
- *
- */
-template<class TypeTag>
-class VariableClass
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
- typedef typename GET_PROP(TypeTag, SolutionTypes) SolutionTypes;
-
- enum
- {
- dim = GridView::dimension, dimWorld = GridView::dimensionworld, numPhase = GET_PROP_VALUE(TypeTag, PTAG(
- NumPhases))
- };
-
- typedef typename GridView::Grid Grid;
- typedef typename GridView::Traits::template Codim<0>::Entity Element;
- typedef typename GridView::Traits::template Codim<dim>::Entity Vertex;
- typedef typename GridView::template Codim<0>::Iterator ElementIterator;
- typedef typename GridView::IntersectionIterator IntersectionIterator;
-
- typedef typename SolutionTypes::VertexMapper VertexMapper;
- typedef typename SolutionTypes::ElementMapper ElementMapper;
-
- typedef Dune::FieldVector<Scalar, dim> LocalPosition;
- typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
-public:
- typedef typename SolutionTypes::ScalarSolution ScalarSolutionType;//!<type for vector of scalars
- typedef typename SolutionTypes::PhaseProperty PhasePropertyType;//!<type for vector of phase properties
- typedef typename SolutionTypes::FluidProperty FluidPropertyType;//!<type for vector of fluid properties
- typedef typename SolutionTypes::PhasePropertyElemFace PhasePropertyElemFaceType;//!<type for vector of vectors (of size 2 x dimension) of scalars
- typedef typename SolutionTypes::DimVecElemFace DimVecElemFaceType;//!<type for vector of vectors (of size 2 x dimension) of vector (of size dimension) of scalars
-
-private:
- const GridView& gridView_;
- ElementMapper elementMapper_;
- VertexMapper vertexMapper_;
- int gridSize_;
-
- const int codim_;
-
- DimVecElemFaceType velocity_;
- PhasePropertyElemFaceType potential_;
-protected:
- ScalarSolutionType pressure_;
-
-
-public:
- //! Constructs a VariableClass object
- /**
- * @param gridView a DUNE gridview object corresponding to diffusion and transport equation
- * @param initialVel initial value for the velocity (only necessary if only transport part is solved)
- */
-
- VariableClass(const GridView& gridView, Dune::FieldVector<Scalar, dim>& initialVel = *(new Dune::FieldVector<Scalar, dim>(0))) :
- gridView_(gridView), elementMapper_(gridView), vertexMapper_(gridView), gridSize_(gridView_.size(0)), codim_(0)
- {
- initializeGlobalVariables(initialVel);
- }
-
- //! Constructs a VariableClass object
- /**
- * @param gridView a DUNE gridview object corresponding to diffusion and transport equation
- * @param codim codimension of the entity of which data has to be strored
- * @param initialVel initial value for the velocity (only necessary if only transport part is solved)
- */
- VariableClass(const GridView& gridView, int codim, Dune::FieldVector<Scalar,
- dim>& initialVel = *(new Dune::FieldVector<Scalar, dim>(0))) :
- gridView_(gridView), elementMapper_(gridView), vertexMapper_(gridView), gridSize_(gridView_.size(codim)), codim_(codim)
- {
- initializeGlobalVariables(initialVel);
- }
-
- // serialization methods
- //! Function needed for restart option.
- template<class Restarter>
- void serialize(Restarter &res)
- {
- res.template serializeEntities<0> (*this, gridView_);
- }
-
- //! Function needed for restart option.
- template<class Restarter>
- void deserialize(Restarter &res)
- {
- res.template deserializeEntities<0> (*this, gridView_);
- }
-
- //! Function needed for restart option.
- void serializeEntity(std::ostream &outstream, const Element &element)
- {
- int globalIdx = elementMapper_.map(element);
- outstream << pressure_[globalIdx];
- }
-
- //! Function needed for restart option.
- void deserializeEntity(std::istream &instream, const Element &element)
- {
- int globalIdx = elementMapper_.map(element);
- instream >> pressure_[globalIdx];
- }
-
- //! initializes the potential differences stored at the element faces.
- void initializePotentials(Dune::FieldVector<Scalar, dim>& initialPot)
- {
- if (initialPot.two_norm())
- {
- // compute update vector
- ElementIterator eItEnd = gridView_.template end<0> ();
- for (ElementIterator eIt = gridView_.template begin<0> (); eIt != eItEnd; ++eIt)
- {
- // cell index
- int globalIdxI = elementMapper_.map(*eIt);
-
- // run through all intersections with neighbors and boundary
- IntersectionIterator isItEnd = gridView_.iend(*eIt);
- for (IntersectionIterator isIt = gridView_.ibegin(*eIt); isIt != isItEnd; ++isIt)
- {
- // local number of facet
- int indexInInside = isIt->indexInInside();
-
- const GlobalPosition& unitOuterNormal = isIt->centerUnitOuterNormal();
-
- for (int i = 0; i < numPhase; i++) {potential_[globalIdxI][indexInInside][i] = initialPot * unitOuterNormal;}
- }
- }
- }
- else
- {
- potential_ = Dune::FieldVector<Scalar, numPhase> (0);
- }
- return;
- }
-
- //! Resizes decoupled variable vectors
- /*! Method that change the size of the vectors for h-adaptive simulations.
- *
- *\param size Size of the current (refined and coarsened) grid
- */
- void adaptVariableSize(int size)
- {
- pressure_.resize(size);
- velocity_.resize(size);
- velocity_ = Dune::FieldVector<Scalar, dim>(0);
- potential_.resize(size);
- potential_ = Dune::FieldVector<Scalar, dim>(0);
- }
-private:
- void initializeGlobalVariables(Dune::FieldVector<Scalar, dim>& initialVel)
- {
- //resize to grid size
- pressure_.resize(gridSize_);
- velocity_.resize(gridSize_);//depends on pressure
- potential_.resize(gridSize_);//depends on pressure
-
- //initialise variables
- pressure_ = 0;
- velocity_ = initialVel;
- initializePotentials(initialVel);
- }
-
- //Write saturation and pressure into file
- template<class MultiWriter>
- void addOutputVtkFields(MultiWriter &writer)
- {
- if (codim_ == 0)
- {
- ScalarSolutionType *pressure = writer.allocateManagedBuffer (this->gridSize());
-
- *pressure = this->pressure();
-
- writer.attachCellData(*pressure, "pressure");
- }
- if (codim_ == dim)
- {
- ScalarSolutionType *pressure = writer.allocateManagedBuffer (this->gridSize());
-
- *pressure = this->pressure();
-
- writer.attachVertexData(*pressure, "pressure");
- }
-
- return;
- }
-public:
-
- void communicatePressure()
- {
-#if HAVE_MPI
- ElementHandleAssign<typename ScalarSolutionType::block_type, ScalarSolutionType, ElementMapper> elementHandle(pressure_, elementMapper_);
- gridView_.communicate(elementHandle,
- Dune::InteriorBorder_All_Interface,
- Dune::ForwardCommunication);
-#endif
- }
-
- //! Return pressure vector
- const ScalarSolutionType& pressure() const
- {
- return pressure_;
- }
-
- ScalarSolutionType& pressure()
- {
- return pressure_;
- }
-
- //! Return velocity vector
- const DimVecElemFaceType& velocity() const
- {
- return velocity_;
- }
-
- DimVecElemFaceType& velocity()
- {
- return velocity_;
- }
-
- const PhasePropertyElemFaceType& potential() const
- {
- return potential_;
- }
-
- PhasePropertyElemFaceType& potential()
- {
- return potential_;
- }
-
- //! Return vector of wetting phase potential gradients
- Dune::FieldVector<Scalar, numPhase>& potential(int Idx1, int Idx2)
- {
- return potential_[Idx1][Idx2];
- }
-
- //! Get index of element (codim 0 entity)
- /*! Get index of element (codim 0 entity).
- * @param element codim 0 entity
- * \return element index
- */
- int index(const Element& element) const
- {
- return elementMapper_.map(element);
- }
-
- //! Get index of vertex (codim dim entity)
- /*! Get index of vertex (codim dim entity).
- * @param vertex codim dim entity
- * \return vertex index
- */
- int index(const Vertex& vertex) const
- {
- return vertexMapper_.map(vertex);
- }
-
- //!Return the number of data elements
- int gridSize() const
- {
- return gridSize_;
- }
-
- void setGridSize(int size)
- {
- gridSize_=size;
- }
-
- //!Return gridView
- const GridView& gridView() const
- {
- return gridView_;
- }
-
- //! Return mapper for elements (for adaptive grids)
- ElementMapper& elementMapper()
- {
- return elementMapper_;
- }
- //! Return mapper for elements (for static grids)
- const ElementMapper& elementMapper() const
- {
- return elementMapper_;
- }
- //! Return mapper for vertices (for adaptive grids)
- VertexMapper& vertexMapper()
- {
- return vertexMapper_;
- }
- //! Return mapper for vertices (for static grids)
- const VertexMapper& vertexMapper() const
- {
- return vertexMapper_;
- }
-
- //! Get pressure
- /*! evaluate pressure at given element
- @param element entity of codim 0
- \return value of pressure
- */
- const Dune::FieldVector<Scalar, 1>& pressElement(const Element& element) const
- {
- return pressure_[elementMapper_.map(element)];
- }
-
- //! Get velocity at given element face
- /*! evaluate velocity at given location
- @param element entity of codim 0
- @param indexInInside index in reference element
- \return vector of velocity
- */
- Dune::FieldVector<Scalar, dim>& velocityElementFace(const Element& element, const int indexInInside)
- {
- int elemId = elementMapper_.map(element);
-
- return (velocity_[elemId][indexInInside]);
- }
-
- const Dune::FieldVector<Scalar, dim>& velocityElementFace(const Element& element, const int indexInInside) const
- {
- int elemId = elementMapper_.map(element);
-
- return (velocity_[elemId][indexInInside]);
- }
-};
-}
-#endif