From 6429ab93cf1b557a2800e60aeea5614a71dc0ae2 Mon Sep 17 00:00:00 2001
From: Bernd Flemisch <bernd@iws.uni-stuttgart.de>
Date: Fri, 30 Nov 2012 14:54:39 +0000
Subject: [PATCH] implicit branch: unify Box- and CC-Problem to
ImplicitProblem. Delete ccproblem.hh. Deprecate BoxProblem in favor of
ImplicitProblem. Adapt includes and names.
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/branches/implicit@9744 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
dumux/boxmodels/common/boxproblem.hh | 66 +-
dumux/freeflow/stokes/stokesproblem.hh | 2 +-
dumux/implicit/box/boxproblem.hh | 835 -----------------
dumux/implicit/box/porousmediaboxproblem.hh | 6 +-
dumux/implicit/cellcentered/ccproblem.hh | 839 ------------------
.../cellcentered/porousmediaccproblem.hh | 6 +-
dumux/implicit/common/implicitproblem.hh | 67 +-
.../common/porousmediaimplicitproblem.hh | 6 +-
dumux/implicit/richards/richardsproblem.hh | 2 +-
test/freeflow/stokes/stokestestproblem.hh | 8 +-
test/freeflow/stokes2c/stokes2ctestproblem.hh | 10 +-
.../stokes2cni/stokes2cnitestproblem.hh | 10 +-
12 files changed, 146 insertions(+), 1711 deletions(-)
delete mode 100644 dumux/implicit/box/boxproblem.hh
delete mode 100644 dumux/implicit/cellcentered/ccproblem.hh
diff --git a/dumux/boxmodels/common/boxproblem.hh b/dumux/boxmodels/common/boxproblem.hh
index 2a70a0ccd8..7599784201 100644
--- a/dumux/boxmodels/common/boxproblem.hh
+++ b/dumux/boxmodels/common/boxproblem.hh
@@ -1,3 +1,65 @@
-#warning This file is deprecated. Include dumux/implicit/box/boxproblem.hh instead.
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ * See the file COPYING for full copying permissions. *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 2 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ *****************************************************************************/
+/*!
+ * \file
+ * \brief Base class for all problems which use the box scheme
+ */
+#ifndef DUMUX_BOX_PROBLEM_HH
+#define DUMUX_BOX_PROBLEM_HH
-#include <dumux/implicit/box/boxproblem.hh>
+#include <dumux/implicit/common/implicitproblem.hh>
+
+namespace Dumux
+{
+/*!
+ * \ingroup BoxModel
+ * \ingroup BoxBaseProblems
+ * \brief Base class for all problems which use the box scheme.
+ *
+ * \note All quantities are specified assuming a threedimensional
+ * world. Problems discretized using 2D grids are assumed to be
+ * extruded by \f$1 m\f$ and 1D grids are assumed to have a
+ * cross section of \f$1m \times 1m\f$.
+ */
+template<class TypeTag>
+class BoxProblem : public ImplicitProblem<TypeTag>
+{
+ typedef ImplicitProblem<TypeTag> ParentType;
+ typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+ typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
+
+ // copying a problem is not a good idea
+ BoxProblem(const BoxProblem &);
+
+public:
+ /*!
+ * \brief Constructor
+ *
+ * \param timeManager The TimeManager which is used by the simulation
+ * \param gridView The simulation's idea about physical space
+ */
+ DUNE_DEPRECATED_MSG("Use ImplicitProblem from dumux/implicit/common/implicitproblem.hh.")
+ BoxProblem(TimeManager &timeManager, const GridView &gridView)
+ : ParentType(timeManager, gridView)
+ {}
+};
+
+}
+
+#endif
diff --git a/dumux/freeflow/stokes/stokesproblem.hh b/dumux/freeflow/stokes/stokesproblem.hh
index afceb9474b..20636c89d4 100644
--- a/dumux/freeflow/stokes/stokesproblem.hh
+++ b/dumux/freeflow/stokes/stokesproblem.hh
@@ -23,7 +23,7 @@
#ifndef DUMUX_STOKES_PROBLEM_HH
#define DUMUX_STOKES_PROBLEM_HH
-#include <dumux/implicit/box/boxproblem.hh>
+#include <dumux/implicit/common/implicitproblem.hh>
#include "stokesproperties.hh"
diff --git a/dumux/implicit/box/boxproblem.hh b/dumux/implicit/box/boxproblem.hh
deleted file mode 100644
index c455443807..0000000000
--- a/dumux/implicit/box/boxproblem.hh
+++ /dev/null
@@ -1,835 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- * See the file COPYING for full copying permissions. *
- * *
- * This program is free software: you can redistribute it and/or modify *
- * it under the terms of the GNU General Public License as published by *
- * the Free Software Foundation, either version 2 of the License, or *
- * (at your option) any later version. *
- * *
- * This program is distributed in the hope that it will be useful, *
- * but WITHOUT ANY WARRANTY; without even the implied warranty of *
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
- * GNU General Public License for more details. *
- * *
- * You should have received a copy of the GNU General Public License *
- * along with this program. If not, see <http://www.gnu.org/licenses/>. *
- *****************************************************************************/
-/*!
- * \file
- * \brief Base class for all problems which use the box scheme
- */
-#ifndef DUMUX_BOX_PROBLEM_HH
-#define DUMUX_BOX_PROBLEM_HH
-
-#include "boxproperties.hh"
-
-#include <dumux/io/vtkmultiwriter.hh>
-#include <dumux/io/restart.hh>
-
-namespace Dumux
-{
-/*!
- * \ingroup BoxModel
- * \ingroup BoxBaseProblems
- * \brief Base class for all problems which use the box scheme.
- *
- * \note All quantities are specified assuming a threedimensional
- * world. Problems discretized using 2D grids are assumed to be
- * extruded by \f$1 m\f$ and 1D grids are assumed to have a
- * cross section of \f$1m \times 1m\f$.
- */
-template<class TypeTag>
-class BoxProblem
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
-
- typedef Dumux::VtkMultiWriter<GridView> VtkMultiWriter;
-
- typedef typename GET_PROP_TYPE(TypeTag, NewtonMethod) NewtonMethod;
- typedef typename GET_PROP_TYPE(TypeTag, NewtonController) NewtonController;
-
- typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
-
- typedef typename GET_PROP_TYPE(TypeTag, VertexMapper) VertexMapper;
- typedef typename GET_PROP_TYPE(TypeTag, ElementMapper) ElementMapper;
-
- typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
- typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
- typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
- typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
-
- enum {
- dim = GridView::dimension,
- dimWorld = GridView::dimensionworld
- };
-
- typedef typename GridView::template Codim<0>::Entity Element;
- typedef typename GridView::template Codim<dim>::Entity Vertex;
- typedef typename GridView::template Codim<dim>::Iterator VertexIterator;
- typedef typename GridView::Intersection Intersection;
-
- typedef typename GridView::Grid::ctype CoordScalar;
- typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
-
- // copying a problem is not a good idea
- BoxProblem(const BoxProblem &);
-
-public:
- /*!
- * \brief Constructor
- *
- * \param timeManager The TimeManager which is used by the simulation
- * \param gridView The simulation's idea about physical space
- */
- BoxProblem(TimeManager &timeManager, const GridView &gridView)
- : gridView_(gridView)
- , bboxMin_(std::numeric_limits<double>::max())
- , bboxMax_(-std::numeric_limits<double>::max())
- , elementMapper_(gridView)
- , vertexMapper_(gridView)
- , timeManager_(&timeManager)
- , newtonMethod_(asImp_())
- , newtonCtl_(asImp_())
- {
- // calculate the bounding box of the local partition 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().corner(0)[i]);
- bboxMax_[i] = std::max(bboxMax_[i], vIt->geometry().corner(0)[i]);
- }
- }
-
- // communicate to get the bounding box of the whole domain
- if (gridView.comm().size() > 1)
- for (int i = 0; i < dim; ++i) {
- bboxMin_[i] = gridView.comm().min(bboxMin_[i]);
- bboxMax_[i] = gridView.comm().max(bboxMax_[i]);
- }
-
- // set a default name for the problem
- simName_ = "sim";
-
- resultWriter_ = NULL;
- }
-
- ~BoxProblem()
- {
- delete resultWriter_;
- };
-
-
- /*!
- * \brief Called by the Dumux::TimeManager in order to
- * initialize the problem.
- *
- * If you overload this method don't forget to call
- * ParentType::init()
- */
- void init()
- {
- // set the initial condition of the model
- model().init(asImp_());
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param values The boundary types for the conservation equations
- * \param vertex The vertex for which the boundary type is set
- */
- void boundaryTypes(BoundaryTypes &values,
- const Vertex &vertex) const
- {
- // forward it to the method which only takes the global coordinate
- asImp_().boundaryTypesAtPos(values, vertex.geometry().center());
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param values The boundary types for the conservation equations
- * \param pos The position of the finite volume in global coordinates
- */
- void boundaryTypesAtPos(BoundaryTypes &values,
- const GlobalPosition &pos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a boundaryTypes() method.");
- }
-
-
- /*!
- * \brief Evaluate the boundary conditions for a dirichlet
- * control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param vertex The vertex representing the "half volume on the boundary"
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void dirichlet(PrimaryVariables &values,
- const Vertex &vertex) const
- {
- // forward it to the method which only takes the global coordinate
- asImp_().dirichletAtPos(values, vertex.geometry().center());
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a dirichlet
- * control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param pos The position of the center of the finite volume
- * for which the dirichlet condition ought to be
- * set in global coordinates
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void dirichletAtPos(PrimaryVariables &values,
- const GlobalPosition &pos) 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 dirichlet() method.");
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * This is the method for the case where the Neumann condition is
- * potentially solution dependent and requires some box method
- * specific things.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param element The finite element
- * \param fvGeometry The finite-volume geometry in the box scheme
- * \param is The intersection between element and boundary
- * \param scvIdx The local vertex index
- * \param boundaryFaceIdx The index of the boundary face
- * \param elemVolVars All volume variables for the element
- *
- * For this method, the \a values parameter stores the mass flux
- * in normal direction of each phase. Negative values mean influx.
- */
- void boxSDNeumann(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvGeometry,
- const Intersection &is,
- const int scvIdx,
- const int boundaryFaceIdx,
- const ElementVolumeVariables &elemVolVars) const
- {
- // forward it to the interface without the volume variables
- asImp_().neumann(values,
- element,
- fvGeometry,
- is,
- scvIdx,
- boundaryFaceIdx);
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param element The finite element
- * \param fvGeometry The finite-volume geometry in the box scheme
- * \param is The intersection between element and boundary
- * \param scvIdx The local vertex index
- * \param boundaryFaceIdx The index of the boundary face
- *
- * 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 Element &element,
- const FVElementGeometry &fvGeometry,
- const Intersection &is,
- const int scvIdx,
- const int boundaryFaceIdx) const
- {
- // forward it to the interface with only the global position
- asImp_().neumannAtPos(values, fvGeometry.boundaryFace[boundaryFaceIdx].ipGlobal);
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param pos The position of the boundary face's integration point in global coordinates
- *
- * 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 &pos) 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 for all phases within a given
- * sub-control-volume.
- *
- * This is the method for the case where the source term is
- * potentially solution dependent and requires some box method
- * specific things.
- *
- * \param values The source and sink values for the conservation equations
- * \param element The finite element
- * \param fvGeometry The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- * \param elemVolVars All volume variables for 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 boxSDSource(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvGeometry,
- const int scvIdx,
- const ElementVolumeVariables &elemVolVars) const
- {
- // forward to solution independent, box specific interface
- asImp_().source(values, element, fvGeometry, scvIdx);
- }
-
- /*!
- * \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 element The finite element
- * \param fvGeometry The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- *
- * 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 FVElementGeometry &fvGeometry,
- const int scvIdx) const
- {
- // forward to generic interface
- asImp_().sourceAtPos(values, fvGeometry.subContVol[scvIdx].global);
- }
-
- /*!
- * \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 pos 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 &pos) const
- {
- 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 finite element
- * \param fvGeometry The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- *
- * For this method, the \a values parameter stores primary
- * variables.
- */
- void initial(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvGeometry,
- const int scvIdx) const
- {
- // forward to generic interface
- asImp_().initialAtPos(values, fvGeometry.subContVol[scvIdx].global);
- }
-
- /*!
- * \brief Evaluate the initial value for a control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param pos 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 &pos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a initialAtPos() method.");
- }
-
- /*!
- * \brief Return how much the domain is extruded at a given sub-control volume.
- *
- * This means the factor by which a lower-dimensional (1D or 2D)
- * entity needs to be expanded to get a full dimensional cell. The
- * default is 1.0 which means that 1D problems are actually
- * thought as pipes with a cross section of 1 m^2 and 2D problems
- * are assumed to extend 1 m to the back.
- */
- Scalar boxExtrusionFactor(const Element &element,
- const FVElementGeometry &fvGeometry,
- const int scvIdx) const
- {
- // forward to generic interface
- return asImp_().extrusionFactorAtPos(fvGeometry.subContVol[scvIdx].global);
- }
-
- /*!
- * \brief Return how much the domain is extruded at a given position.
- *
- * This means the factor by which a lower-dimensional (1D or 2D)
- * entity needs to be expanded to get a full dimensional cell. The
- * default is 1.0 which means that 1D problems are actually
- * thought as pipes with a cross section of 1 m^2 and 2D problems
- * are assumed to extend 1 m to the back.
- */
- Scalar extrusionFactorAtPos(const GlobalPosition &pos) const
- { return 1.0; }
-
- /*!
- * \brief If model coupling is used, this updates the parameters
- * required to calculate the coupling fluxes between the
- * sub-models.
- *
- * By default it does nothing
- *
- * \param element The DUNE Codim<0> entity for which the coupling
- * parameters should be computed.
- */
- void updateCouplingParams(const Element &element) const
- {}
-
- /*!
- * \name Simulation steering
- */
- // \{
-
- /*!
- * \brief Called by the time manager before the time integration.
- */
- void preTimeStep()
- {}
-
- /*!
- * \brief Called by Dumux::TimeManager in order to do a time
- * integration on the model.
- */
- void timeIntegration()
- {
- const int maxFails = 10;
- for (int i = 0; i < maxFails; ++i) {
- if (model_.update(newtonMethod_, newtonCtl_))
- return;
-
- Scalar dt = timeManager().timeStepSize();
- Scalar nextDt = dt / 2;
- timeManager().setTimeStepSize(nextDt);
-
- // update failed
- std::cout << "Newton solver did not converge with dt="<<dt<<" seconds. Retrying with time step of "
- << nextDt << " seconds\n";
- }
-
- DUNE_THROW(Dune::MathError,
- "Newton solver didn't converge after "
- << maxFails
- << " time-step divisions. dt="
- << timeManager().timeStepSize());
- }
-
- /*!
- * \brief Returns the newton method object
- */
- NewtonMethod &newtonMethod()
- { return newtonMethod_; }
-
- /*!
- * \copydoc newtonMethod()
- */
- const NewtonMethod &newtonMethod() const
- { return newtonMethod_; }
-
- /*!
- * \brief Returns the newton contoller object
- */
- NewtonController &newtonController()
- { return newtonCtl_; }
-
- /*!
- * \copydoc newtonController()
- */
- const NewtonController &newtonController() const
- { return newtonCtl_; }
-
- /*!
- * \brief Called by Dumux::TimeManager whenever a solution for a
- * time step has been computed and the simulation time has
- * been updated.
- *
- * \param dt The current time-step size
- */
- Scalar nextTimeStepSize(const Scalar dt)
- {
- return std::min(GET_PARAM_FROM_GROUP(TypeTag, Scalar, TimeManager, MaxTimeStepSize),
- newtonCtl_.suggestTimeStepSize(dt));
- };
-
- /*!
- * \brief Returns true if a restart file should be written to
- * disk.
- *
- * The default behavior is to write one restart file every 5 time
- * steps. This file is intended to be overwritten by the
- * implementation.
- */
- bool shouldWriteRestartFile() const
- {
- return timeManager().timeStepIndex() > 0 &&
- (timeManager().timeStepIndex() % 10 == 0);
- }
-
- /*!
- * \brief Returns true if the current solution should be written to
- * disk (i.e. as a VTK file)
- *
- * The default behavior is to write out every the solution for
- * very time step. This file is intended to be overwritten by the
- * implementation.
- */
- bool shouldWriteOutput() const
- { return true; }
-
- /*!
- * \brief Called by the time manager after the time integration to
- * do some post processing on the solution.
- */
- 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()
- {
- model_.advanceTimeLevel();
- }
-
- /*!
- * \brief Called when the end of an simulation episode is reached.
- *
- * Typically a new episode should be started in this method.
- */
- 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 static method sets the simulation name, which should be
- * called before the application problem is declared! If not, the
- * default name "sim" will be used.
- *
- * \param newName The problem's name
- */
- void setName(const char *newName)
- {
- simName_ = newName;
- }
-
-
- /*!
- * \brief Returns the number of the current VTK file.
- */
- int currentVTKFileNumber()
- {
- createResultWriter_();
- return resultWriter_->curWriterNum();
- }
-
- /*!
- * \brief The GridView which used by the problem.
- */
- const GridView &gridView() const
- { return gridView_; }
-
- /*!
- * \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 the mapper for vertices to indices.
- */
- const VertexMapper &vertexMapper() const
- { return vertexMapper_; }
-
- /*!
- * \brief Returns the mapper for elements to indices.
- */
- const ElementMapper &elementMapper() const
- { return elementMapper_; }
-
- /*!
- * \brief Returns TimeManager object used by the simulation
- */
- TimeManager &timeManager()
- { return *timeManager_; }
-
- /*!
- * \copydoc timeManager()
- */
- const TimeManager &timeManager() const
- { return *timeManager_; }
-
- /*!
- * \brief Returns numerical model used for the problem.
- */
- Model &model()
- { return model_; }
-
- /*!
- * \copydoc model()
- */
- const Model &model() const
- { return model_; }
- // \}
-
- /*!
- * \name Restart mechanism
- */
- // \{
-
- /*!
- * \brief This method writes the complete state of the simulation
- * 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_());
- if (gridView().comm().rank() == 0)
- std::cout << "Serialize to file '" << res.fileName() << "'\n";
-
- timeManager().serialize(res);
- asImp_().serialize(res);
- res.serializeEnd();
- }
-
- /*!
- * \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.
- *
- * \tparam Restarter The serializer type
- *
- * \param res The serializer object
- */
- template <class Restarter>
- void serialize(Restarter &res)
- {
- createResultWriter_();
- resultWriter_->serialize(res);
- model().serialize(res);
- }
-
- /*!
- * \brief Load a previously saved state of the whole simulation
- * from disk.
- *
- * \param tRestart The simulation time on which the program was
- * written to disk.
- */
- void restart(const Scalar tRestart)
- {
- typedef Dumux::Restart Restarter;
-
- Restarter res;
-
- res.deserializeBegin(asImp_(), tRestart);
- if (gridView().comm().rank() == 0)
- std::cout << "Deserialize from file '" << res.fileName() << "'\n";
- timeManager().deserialize(res);
- asImp_().deserialize(res);
- res.deserializeEnd();
- }
-
- /*!
- * \brief This method restores the complete state of the problem
- * from disk.
- *
- * It is the inverse of the serialize() method.
- *
- * \tparam Restarter The deserializer type
- *
- * \param res The deserializer object
- */
- template <class Restarter>
- void deserialize(Restarter &res)
- {
- createResultWriter_();
- resultWriter_->deserialize(res);
- model().deserialize(res);
- }
-
- // \}
-
- /*!
- * \brief Adds additional VTK output data to the VTKWriter. Function is called by writeOutput().
- */
- void addOutputVtkFields()
- {}
-
- /*!
- * \brief Write the relevant secondary variables of the current
- * solution into an VTK output file.
- */
- void writeOutput(const bool verbose = true)
- {
- // write the current result to disk
- if (asImp_().shouldWriteOutput()) {
- if (verbose && gridView().comm().rank() == 0)
- std::cout << "Writing result file for \"" << asImp_().name() << "\"\n";
-
- // calculate the time _after_ the time was updated
- Scalar t = timeManager().time() + timeManager().timeStepSize();
- createResultWriter_();
- resultWriter_->beginWrite(t);
- model().addOutputVtkFields(model().curSol(), *resultWriter_);
- asImp_().addOutputVtkFields();
- resultWriter_->endWrite();
- }
- }
-
-protected:
- //! Returns the implementation of the problem (i.e. static polymorphism)
- Implementation &asImp_()
- { return *static_cast<Implementation *>(this); }
-
- //! \copydoc asImp_()
- const Implementation &asImp_() const
- { return *static_cast<const Implementation *>(this); }
-
- //! Returns the applied VTK-writer for the output
- VtkMultiWriter& resultWriter()
- {
- createResultWriter_();
- return *resultWriter_;
- }
- //! \copydoc Dumux::IMPETProblem::resultWriter()
- VtkMultiWriter& resultWriter() const
- {
- createResultWriter_();
- return *resultWriter_;
- }
-
-
-private:
- // makes sure that the result writer exists
- void createResultWriter_()
- { if (!resultWriter_) resultWriter_ = new VtkMultiWriter(gridView_, asImp_().name()); };
-
- std::string simName_;
- const GridView gridView_;
-
- GlobalPosition bboxMin_;
- GlobalPosition bboxMax_;
-
- ElementMapper elementMapper_;
- VertexMapper vertexMapper_;
-
- TimeManager *timeManager_;
-
- Model model_;
-
- NewtonMethod newtonMethod_;
- NewtonController newtonCtl_;
-
- VtkMultiWriter *resultWriter_;
-};
-
-}
-
-#endif
diff --git a/dumux/implicit/box/porousmediaboxproblem.hh b/dumux/implicit/box/porousmediaboxproblem.hh
index 3fc4d54d9c..b4107b35c4 100644
--- a/dumux/implicit/box/porousmediaboxproblem.hh
+++ b/dumux/implicit/box/porousmediaboxproblem.hh
@@ -26,7 +26,7 @@
#include "boxproperties.hh"
-#include <dumux/implicit/box/boxproblem.hh>
+#include <dumux/implicit/common/implicitproblem.hh>
namespace Dumux
{
@@ -41,9 +41,9 @@ NEW_PROP_TAG(ProblemEnableGravity); //!< Returns whether gravity is considered i
* \brief Base class for all porous media box problems
*/
template<class TypeTag>
-class PorousMediaBoxProblem : public BoxProblem<TypeTag>
+class PorousMediaBoxProblem : public ImplicitProblem<TypeTag>
{
- typedef BoxProblem<TypeTag> ParentType;
+ typedef ImplicitProblem<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
diff --git a/dumux/implicit/cellcentered/ccproblem.hh b/dumux/implicit/cellcentered/ccproblem.hh
deleted file mode 100644
index 716bc59f56..0000000000
--- a/dumux/implicit/cellcentered/ccproblem.hh
+++ /dev/null
@@ -1,839 +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 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/>. *
- *****************************************************************************/
-/*!
- * \file
- * \brief Base class for all problems which use the box scheme
- */
-#ifndef DUMUX_CC_PROBLEM_HH
-#define DUMUX_CC_PROBLEM_HH
-
-#include "ccproperties.hh"
-
-#include <dumux/io/vtkmultiwriter.hh>
-#include <dumux/io/restart.hh>
-
-namespace Dumux
-{
-/*!
- * \ingroup CCModel
- * \ingroup CCBaseProblems
- * \brief Base class for all problems which use the box scheme.
- *
- * \note All quantities are specified assuming a threedimensional
- * world. Problems discretized using 2D grids are assumed to be
- * extruded by \f$1 m\f$ and 1D grids are assumed to have a
- * cross section of \f$1m \times 1m\f$.
- */
-template<class TypeTag>
-class CCProblem
-{
-private:
- typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
- typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
-
- typedef Dumux::VtkMultiWriter<GridView> VtkMultiWriter;
-
- typedef typename GET_PROP_TYPE(TypeTag, NewtonMethod) NewtonMethod;
- typedef typename GET_PROP_TYPE(TypeTag, NewtonController) NewtonController;
-
- typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
- typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
- typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
-
- typedef typename GET_PROP_TYPE(TypeTag, VertexMapper) VertexMapper;
- typedef typename GET_PROP_TYPE(TypeTag, ElementMapper) ElementMapper;
-
- typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
- typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
- typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
- typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
-
- enum {
- dim = GridView::dimension,
- dimWorld = GridView::dimensionworld
- };
-
- typedef typename GridView::template Codim<0>::Entity Element;
- typedef typename GridView::template Codim<dim>::Entity Vertex;
- typedef typename GridView::template Codim<dim>::Iterator VertexIterator;
- typedef typename GridView::template Codim<0>::Iterator ElementIterator;
- typedef typename GridView::Intersection Intersection;
-
- typedef typename GridView::Grid::ctype CoordScalar;
- typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
-
- // copying a problem is not a good idea
- CCProblem(const CCProblem &);
-
-public:
- /*!
- * \brief Constructor
- *
- * \param timeManager The TimeManager which is used by the simulation
- * \param gridView The simulation's idea about physical space
- */
- CCProblem(TimeManager &timeManager, const GridView &gridView)
- : gridView_(gridView)
- , bboxMin_(std::numeric_limits<double>::max())
- , bboxMax_(-std::numeric_limits<double>::max())
- , elementMapper_(gridView)
- , vertexMapper_(gridView)
- , timeManager_(&timeManager)
- , newtonMethod_(asImp_())
- , newtonCtl_(asImp_())
- {
- // calculate the bounding box of the local partition 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().corner(0)[i]);
- bboxMax_[i] = std::max(bboxMax_[i], vIt->geometry().corner(0)[i]);
- }
- }
-
- // communicate to get the bounding box of the whole domain
- if (gridView.comm().size() > 1)
- for (int i = 0; i < dim; ++i) {
- bboxMin_[i] = gridView.comm().min(bboxMin_[i]);
- bboxMax_[i] = gridView.comm().max(bboxMax_[i]);
- }
-
- // set a default name for the problem
- simName_ = "sim";
-
- resultWriter_ = NULL;
- }
-
- ~CCProblem()
- {
- delete resultWriter_;
- };
-
-
- /*!
- * \brief Called by the Dumux::TimeManager in order to
- * initialize the problem.
- *
- * If you overload this method don't forget to call
- * ParentType::init()
- */
- void init()
- {
- // set the initial condition of the model
- model().init(asImp_());
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param values The boundary types for the conservation equations
- * \param intersection The intersection for which the boundary type is set
- */
- void boundaryTypes(BoundaryTypes &values,
- const Intersection &intersection) const
- {
- // forward it to the method which only takes the global coordinate
- asImp_().boundaryTypesAtPos(values, intersection.geometry().center());
- }
-
- /*!
- * \brief Specifies which kind of boundary condition should be
- * used for which equation on a given boundary segment.
- *
- * \param values The boundary types for the conservation equations
- * \param pos The position of the finite volume in global coordinates
- */
- void boundaryTypesAtPos(BoundaryTypes &values,
- const GlobalPosition &pos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a boundaryTypes() method.");
- }
-
-
- /*!
- * \brief Evaluate the boundary conditions for a dirichlet
- * control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param intersection The intersection representing the "half volume on the boundary"
- *
- * 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 pos The position of the center of the finite volume
- * for which the dirichlet condition ought to be
- * set in global coordinates
- *
- * For this method, the \a values parameter stores primary variables.
- */
- void dirichletAtPos(PrimaryVariables &values,
- const GlobalPosition &pos) 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 dirichlet() method.");
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * This is the method for the case where the Neumann condition is
- * potentially solution dependent and requires some box method
- * specific things.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param element The finite element
- * \param fvElemGeom The finite-volume geometry in the box scheme
- * \param is The intersection between element and boundary
- * \param scvIdx The local vertex index
- * \param boundaryFaceIdx The index of the boundary face
- * \param elemVolVars All volume variables for the element
- *
- * For this method, the \a values parameter stores the mass flux
- * in normal direction of each phase. Negative values mean influx.
- */
- void boxSDNeumann(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvElemGeom,
- const Intersection &is,
- int scvIdx,
- int boundaryFaceIdx,
- const ElementVolumeVariables &elemVolVars) const
- {
- // forward it to the interface without the volume variables
- asImp_().neumann(values,
- element,
- fvElemGeom,
- is,
- scvIdx,
- boundaryFaceIdx);
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param element The finite element
- * \param fvElemGeom The finite-volume geometry in the box scheme
- * \param is The intersection between element and boundary
- * \param scvIdx The local vertex index
- * \param boundaryFaceIdx The index of the boundary face
- *
- * 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 Element &element,
- const FVElementGeometry &fvElemGeom,
- const Intersection &is,
- int scvIdx,
- int boundaryFaceIdx) const
- {
- // forward it to the interface with only the global position
- asImp_().neumannAtPos(values, fvElemGeom.boundaryFace[boundaryFaceIdx].ipGlobal);
- }
-
- /*!
- * \brief Evaluate the boundary conditions for a neumann
- * boundary segment.
- *
- * \param values The neumann values for the conservation equations [kg / (m^2 *s )]
- * \param pos The position of the boundary face's integration point in global coordinates
- *
- * 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 &pos) 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 for all phases within a given
- * sub-control-volume.
- *
- * This is the method for the case where the source term is
- * potentially solution dependent and requires some box method
- * specific things.
- *
- * \param values The source and sink values for the conservation equations
- * \param element The finite element
- * \param fvElemGeom The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- * \param elemVolVars All volume variables for 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 boxSDSource(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvElemGeom,
- int scvIdx,
- const ElementVolumeVariables &elemVolVars) const
- {
- // forward to solution independent, box specific interface
- asImp_().source(values, element, fvElemGeom, scvIdx);
- }
-
- /*!
- * \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 element The finite element
- * \param fvElemGeom The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- *
- * 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 FVElementGeometry &fvElemGeom,
- int scvIdx) const
- {
- // forward to generic interface
- asImp_().sourceAtPos(values, fvElemGeom.subContVol[scvIdx].global);
- }
-
- /*!
- * \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 pos 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 &pos) const
- {
- 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 finite element
- * \param fvElemGeom The finite-volume geometry in the box scheme
- * \param scvIdx The local vertex index
- *
- * For this method, the \a values parameter stores primary
- * variables.
- */
- void initial(PrimaryVariables &values,
- const Element &element,
- const FVElementGeometry &fvElemGeom,
- int scvIdx) const
- {
- // forward to generic interface
- asImp_().initialAtPos(values, fvElemGeom.subContVol[scvIdx].global);
- }
-
- /*!
- * \brief Evaluate the initial value for a control volume.
- *
- * \param values The dirichlet values for the primary variables
- * \param pos 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 &pos) const
- {
- // Throw an exception (there is no reasonable default value
- // for Dirichlet conditions)
- DUNE_THROW(Dune::InvalidStateException,
- "The problem does not provide "
- "a initialAtPos() method.");
- }
-
- /*!
- * \brief Return how much the domain is extruded at a given sub-control volume.
- *
- * This means the factor by which a lower-dimensional (1D or 2D)
- * entity needs to be expanded to get a full dimensional cell. The
- * default is 1.0 which means that 1D problems are actually
- * thought as pipes with a cross section of 1 m^2 and 2D problems
- * are assumed to extend 1 m to the back.
- */
- Scalar boxExtrusionFactor(const Element &element,
- const FVElementGeometry &fvElemGeom,
- int scvIdx) const
- {
- // forward to generic interface
- return asImp_().extrusionFactorAtPos(fvElemGeom.subContVol[scvIdx].global);
- }
-
- /*!
- * \brief Return how much the domain is extruded at a given position.
- *
- * This means the factor by which a lower-dimensional (1D or 2D)
- * entity needs to be expanded to get a full dimensional cell. The
- * default is 1.0 which means that 1D problems are actually
- * thought as pipes with a cross section of 1 m^2 and 2D problems
- * are assumed to extend 1 m to the back.
- */
- Scalar extrusionFactorAtPos(const GlobalPosition &pos) const
- { return 1.0; }
-
- /*!
- * \brief If model coupling is used, this updates the parameters
- * required to calculate the coupling fluxes between the
- * sub-models.
- *
- * By default it does nothing
- *
- * \param element The DUNE Codim<0> entity for which the coupling
- * parameters should be computed.
- */
- void updateCouplingParams(const Element &element) const
- {}
-
- /*!
- * \name Simulation steering
- */
- // \{
-
- /*!
- * \brief Called by the time manager before the time integration.
- */
- void preTimeStep()
- {}
-
- /*!
- * \brief Called by Dumux::TimeManager in order to do a time
- * integration on the model.
- */
- void timeIntegration()
- {
- const int maxFails = 10;
- for (int i = 0; i < maxFails; ++i) {
- if (model_.update(newtonMethod_, newtonCtl_))
- return;
-
- Scalar dt = timeManager().timeStepSize();
- Scalar nextDt = dt / 2;
- timeManager().setTimeStepSize(nextDt);
-
- // update failed
- std::cout << "Newton solver did not converge with dt="<<dt<<" seconds. Retrying with time step of "
- << nextDt << " seconds\n";
- }
-
- DUNE_THROW(Dune::MathError,
- "Newton solver didn't converge after "
- << maxFails
- << " time-step divisions. dt="
- << timeManager().timeStepSize());
- }
-
- /*!
- * \brief Returns the newton method object
- */
- NewtonMethod &newtonMethod()
- { return newtonMethod_; }
-
- /*!
- * \copydoc newtonMethod()
- */
- const NewtonMethod &newtonMethod() const
- { return newtonMethod_; }
-
- /*!
- * \brief Returns the newton contoller object
- */
- NewtonController &newtonController()
- { return newtonCtl_; }
-
- /*!
- * \copydoc newtonController()
- */
- const NewtonController &newtonController() const
- { return newtonCtl_; }
-
- /*!
- * \brief Called by Dumux::TimeManager whenever a solution for a
- * time step has been computed and the simulation time has
- * been updated.
- *
- * \param dt The current time-step size
- */
- Scalar nextTimeStepSize(Scalar dt)
- {
- return std::min(GET_PARAM_FROM_GROUP(TypeTag, Scalar, TimeManager, MaxTimeStepSize),
- newtonCtl_.suggestTimeStepSize(dt));
- };
-
- /*!
- * \brief Returns true if a restart file should be written to
- * disk.
- *
- * The default behavior is to write one restart file every 5 time
- * steps. This file is intended to be overwritten by the
- * implementation.
- */
- bool shouldWriteRestartFile() const
- {
- return timeManager().timeStepIndex() > 0 &&
- (timeManager().timeStepIndex() % 10 == 0);
- }
-
- /*!
- * \brief Returns true if the current solution should be written to
- * disk (i.e. as a VTK file)
- *
- * The default behavior is to write out every the solution for
- * very time step. This file is intended to be overwritten by the
- * implementation.
- */
- bool shouldWriteOutput() const
- { return true; }
-
- /*!
- * \brief Called by the time manager after the time integration to
- * do some post processing on the solution.
- */
- 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()
- {
- model_.advanceTimeLevel();
- }
-
- /*!
- * \brief Called when the end of an simulation episode is reached.
- *
- * Typically a new episode should be started in this method.
- */
- 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 static method sets the simulation name, which should be
- * called before the application problem is declared! If not, the
- * default name "sim" will be used.
- *
- * \param newName The problem's name
- */
- void setName(const char *newName)
- {
- simName_ = newName;
- }
-
-
- /*!
- * \brief Returns the number of the current VTK file.
- */
- int currentVTKFileNumber()
- {
- createResultWriter_();
- return resultWriter_->curWriterNum();
- }
-
- /*!
- * \brief The GridView which used by the problem.
- */
- const GridView &gridView() const
- { return gridView_; }
-
- /*!
- * \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 the mapper for vertices to indices.
- */
- const VertexMapper &vertexMapper() const
- { return vertexMapper_; }
-
- /*!
- * \brief Returns the mapper for elements to indices.
- */
- const ElementMapper &elementMapper() const
- { return elementMapper_; }
-
- /*!
- * \brief Returns TimeManager object used by the simulation
- */
- TimeManager &timeManager()
- { return *timeManager_; }
-
- /*!
- * \copydoc timeManager()
- */
- const TimeManager &timeManager() const
- { return *timeManager_; }
-
- /*!
- * \brief Returns numerical model used for the problem.
- */
- Model &model()
- { return model_; }
-
- /*!
- * \copydoc model()
- */
- const Model &model() const
- { return model_; }
- // \}
-
- /*!
- * \name Restart mechanism
- */
- // \{
-
- /*!
- * \brief This method writes the complete state of the simulation
- * 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_());
- if (gridView().comm().rank() == 0)
- std::cout << "Serialize to file '" << res.fileName() << "'\n";
-
- timeManager().serialize(res);
- asImp_().serialize(res);
- res.serializeEnd();
- }
-
- /*!
- * \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.
- *
- * \tparam Restarter The serializer type
- *
- * \param res The serializer object
- */
- template <class Restarter>
- void serialize(Restarter &res)
- {
- createResultWriter_();
- resultWriter_->serialize(res);
- model().serialize(res);
- }
-
- /*!
- * \brief Load a previously saved state of the whole simulation
- * from disk.
- *
- * \param tRestart The simulation time on which the program was
- * written to disk.
- */
- void restart(Scalar tRestart)
- {
- typedef Dumux::Restart Restarter;
-
- Restarter res;
-
- res.deserializeBegin(asImp_(), tRestart);
- if (gridView().comm().rank() == 0)
- std::cout << "Deserialize from file '" << res.fileName() << "'\n";
- timeManager().deserialize(res);
- asImp_().deserialize(res);
- res.deserializeEnd();
- }
-
- /*!
- * \brief This method restores the complete state of the problem
- * from disk.
- *
- * It is the inverse of the serialize() method.
- *
- * \tparam Restarter The deserializer type
- *
- * \param res The deserializer object
- */
- template <class Restarter>
- void deserialize(Restarter &res)
- {
- createResultWriter_();
- resultWriter_->deserialize(res);
- model().deserialize(res);
- }
-
- // \}
-
- /*!
- * \brief Adds additional VTK output data to the VTKWriter. Function is called by writeOutput().
- */
- void addOutputVtkFields()
- {}
-
- /*!
- * \brief Write the relevant secondary variables of the current
- * solution into an VTK output file.
- */
- void writeOutput(bool verbose = true)
- {
- // write the current result to disk
- if (asImp_().shouldWriteOutput()) {
- if (verbose && gridView().comm().rank() == 0)
- std::cout << "Writing result file for \"" << asImp_().name() << "\"\n";
-
- // calculate the time _after_ the time was updated
- Scalar t = timeManager().time() + timeManager().timeStepSize();
- createResultWriter_();
- resultWriter_->beginWrite(t);
- model().addOutputVtkFields(model().curSol(), *resultWriter_);
- asImp_().addOutputVtkFields();
- resultWriter_->endWrite();
- }
- }
-
-protected:
- //! Returns the implementation of the problem (i.e. static polymorphism)
- Implementation &asImp_()
- { return *static_cast<Implementation *>(this); }
-
- //! \copydoc asImp_()
- const Implementation &asImp_() const
- { return *static_cast<const Implementation *>(this); }
-
- //! Returns the applied VTK-writer for the output
- VtkMultiWriter& resultWriter()
- {
- createResultWriter_();
- return *resultWriter_;
- }
- //! \copydoc Dumux::IMPETProblem::resultWriter()
- VtkMultiWriter& resultWriter() const
- {
- createResultWriter_();
- return *resultWriter_;
- }
-
-
-private:
- // makes sure that the result writer exists
- void createResultWriter_()
- { if (!resultWriter_) resultWriter_ = new VtkMultiWriter(gridView_, asImp_().name()); };
-
- std::string simName_;
- const GridView gridView_;
-
- GlobalPosition bboxMin_;
- GlobalPosition bboxMax_;
-
- ElementMapper elementMapper_;
- VertexMapper vertexMapper_;
-
- TimeManager *timeManager_;
-
- Model model_;
-
- NewtonMethod newtonMethod_;
- NewtonController newtonCtl_;
-
- VtkMultiWriter *resultWriter_;
-};
-
-}
-
-#endif
diff --git a/dumux/implicit/cellcentered/porousmediaccproblem.hh b/dumux/implicit/cellcentered/porousmediaccproblem.hh
index 5556f53468..3e32676486 100644
--- a/dumux/implicit/cellcentered/porousmediaccproblem.hh
+++ b/dumux/implicit/cellcentered/porousmediaccproblem.hh
@@ -26,7 +26,7 @@
#ifndef DUMUX_POROUS_MEDIA_CC_PROBLEM_HH
#define DUMUX_POROUS_MEDIA_CC_PROBLEM_HH
-#include <dumux/implicit/cellcentered/ccproblem.hh>
+#include <dumux/implicit/common/implicitproblem.hh>
#include "ccproperties.hh"
namespace Dumux
@@ -43,9 +43,9 @@ namespace Properties
* \brief Base class for all fully implicit cell centered porous media problems
*/
template<class TypeTag>
-class PorousMediaCCProblem : public CCProblem<TypeTag>
+class PorousMediaCCProblem : public ImplicitProblem<TypeTag>
{
- typedef CCProblem<TypeTag> ParentType;
+ typedef ImplicitProblem<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
diff --git a/dumux/implicit/common/implicitproblem.hh b/dumux/implicit/common/implicitproblem.hh
index 2c96e0b458..492ffd3d34 100644
--- a/dumux/implicit/common/implicitproblem.hh
+++ b/dumux/implicit/common/implicitproblem.hh
@@ -18,10 +18,10 @@
*****************************************************************************/
/*!
* \file
- * \brief Base class for all problems which use the box scheme
+ * \brief Base class for all fully implicit problems
*/
-#ifndef DUMUX_BOX_PROBLEM_HH
-#define DUMUX_BOX_PROBLEM_HH
+#ifndef DUMUX_IMPLICIT_PROBLEM_HH
+#define DUMUX_IMPLICIT_PROBLEM_HH
#include "implicitproperties.hh"
@@ -31,8 +31,8 @@
namespace Dumux
{
/*!
- * \ingroup BoxModel
- * \ingroup BoxBaseProblems
+ * \ingroup ImplicitModel
+ * \ingroup ImplicitBaseProblems
* \brief Base class for all problems which use the box scheme.
*
* \note All quantities are specified assuming a threedimensional
@@ -41,7 +41,7 @@ namespace Dumux
* cross section of \f$1m \times 1m\f$.
*/
template<class TypeTag>
-class BoxProblem
+class ImplicitProblem
{
private:
typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
@@ -77,8 +77,10 @@ private:
typedef typename GridView::Grid::ctype CoordScalar;
typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
+ enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
+
// copying a problem is not a good idea
- BoxProblem(const BoxProblem &);
+ ImplicitProblem(const ImplicitProblem &);
public:
/*!
@@ -87,7 +89,7 @@ public:
* \param timeManager The TimeManager which is used by the simulation
* \param gridView The simulation's idea about physical space
*/
- BoxProblem(TimeManager &timeManager, const GridView &gridView)
+ ImplicitProblem(TimeManager &timeManager, const GridView &gridView)
: gridView_(gridView)
, bboxMin_(std::numeric_limits<double>::max())
, bboxMax_(-std::numeric_limits<double>::max())
@@ -120,7 +122,7 @@ public:
resultWriter_ = NULL;
}
- ~BoxProblem()
+ ~ImplicitProblem()
{
delete resultWriter_;
};
@@ -149,10 +151,32 @@ public:
void boundaryTypes(BoundaryTypes &values,
const Vertex &vertex) const
{
+ if (!isBox)
+ DUNE_THROW(Dune::InvalidStateException,
+ "boundaryTypes(..., vertex) called for cell-centered method.");
+
// forward it to the method which only takes the global coordinate
asImp_().boundaryTypesAtPos(values, vertex.geometry().center());
}
+ /*!
+ * \brief Specifies which kind of boundary condition should be
+ * used for which equation on a given boundary segment.
+ *
+ * \param values The boundary types for the conservation equations
+ * \param intersection The intersection for which the boundary type is set
+ */
+ void boundaryTypes(BoundaryTypes &values,
+ const Intersection &intersection) const
+ {
+ if (isBox)
+ DUNE_THROW(Dune::InvalidStateException,
+ "boundaryTypes(..., intersection) called for box method.");
+
+ // forward it to the method which only takes the global coordinate
+ asImp_().boundaryTypesAtPos(values, intersection.geometry().center());
+ }
+
/*!
* \brief Specifies which kind of boundary condition should be
* used for which equation on a given boundary segment.
@@ -170,7 +194,6 @@ public:
"a boundaryTypes() method.");
}
-
/*!
* \brief Evaluate the boundary conditions for a dirichlet
* control volume.
@@ -183,10 +206,34 @@ public:
void dirichlet(PrimaryVariables &values,
const Vertex &vertex) const
{
+ if (!isBox)
+ DUNE_THROW(Dune::InvalidStateException,
+ "dirichlet(..., vertex) called for cell-centered method.");
+
// forward it to the method which only takes the global coordinate
asImp_().dirichletAtPos(values, vertex.geometry().center());
}
+ /*!
+ * \brief Evaluate the boundary conditions for a dirichlet
+ * control volume.
+ *
+ * \param values The dirichlet values for the primary variables
+ * \param intersection The intersection for which the condition is evaluated
+ *
+ * For this method, the \a values parameter stores primary variables.
+ */
+ void dirichlet(PrimaryVariables &values,
+ const Intersection &intersection) const
+ {
+ if (isBox)
+ DUNE_THROW(Dune::InvalidStateException,
+ "dirichlet(..., intersection) called for box method.");
+
+ // 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.
diff --git a/dumux/implicit/common/porousmediaimplicitproblem.hh b/dumux/implicit/common/porousmediaimplicitproblem.hh
index 80e0465358..277edf2de9 100644
--- a/dumux/implicit/common/porousmediaimplicitproblem.hh
+++ b/dumux/implicit/common/porousmediaimplicitproblem.hh
@@ -26,7 +26,7 @@
#include "implicitproperties.hh"
-#include <dumux/implicit/box/boxproblem.hh>
+#include <dumux/implicit/common/implicitproblem.hh>
namespace Dumux
{
@@ -41,9 +41,9 @@ NEW_PROP_TAG(ProblemEnableGravity); //!< Returns whether gravity is considered i
* \brief Base class for all porous media box problems
*/
template<class TypeTag>
-class PorousMediaBoxProblem : public BoxProblem<TypeTag>
+class PorousMediaBoxProblem : public ImplicitProblem<TypeTag>
{
- typedef BoxProblem<TypeTag> ParentType;
+ typedef ImplicitProblem<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
diff --git a/dumux/implicit/richards/richardsproblem.hh b/dumux/implicit/richards/richardsproblem.hh
index 59c080ce90..0a991ec208 100644
--- a/dumux/implicit/richards/richardsproblem.hh
+++ b/dumux/implicit/richards/richardsproblem.hh
@@ -56,7 +56,7 @@ public:
* jacobian assembler, etc inside the constructor.
*
* If the problem requires information from these, the
- * BoxProblem::init() method be overloaded.
+ * ImplicitProblem::init() method be overloaded.
*
* \param timeManager The TimeManager which keeps track of time
* \param gridView The GridView used by the problem.
diff --git a/test/freeflow/stokes/stokestestproblem.hh b/test/freeflow/stokes/stokestestproblem.hh
index 5b5bfe4f0a..a104e300a1 100644
--- a/test/freeflow/stokes/stokestestproblem.hh
+++ b/test/freeflow/stokes/stokestestproblem.hh
@@ -158,7 +158,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::boundaryTypes()
+ //! \copydoc ImplicitProblem::boundaryTypes()
void boundaryTypes(BoundaryTypes &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -179,7 +179,7 @@ public:
values.setDirichlet(massBalanceIdx);
}
- //! \copydoc BoxProblem::dirichlet()
+ //! \copydoc ImplicitProblem::dirichlet()
void dirichlet(PrimaryVariables &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -228,7 +228,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::source()
+ //! \copydoc ImplicitProblem::source()
void source(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
@@ -239,7 +239,7 @@ public:
values = Scalar(0.0);
}
- //! \copydoc BoxProblem::initial()
+ //! \copydoc ImplicitProblem::initial()
void initial(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
diff --git a/test/freeflow/stokes2c/stokes2ctestproblem.hh b/test/freeflow/stokes2c/stokes2ctestproblem.hh
index 47d039202b..6cee5fdc6f 100644
--- a/test/freeflow/stokes2c/stokes2ctestproblem.hh
+++ b/test/freeflow/stokes2c/stokes2ctestproblem.hh
@@ -166,7 +166,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::boundaryTypes()
+ //! \copydoc ImplicitProblem::boundaryTypes()
void boundaryTypes(BoundaryTypes &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -187,7 +187,7 @@ public:
values.setDirichlet(massBalanceIdx);
}
- //! \copydoc BoxProblem::dirichlet()
+ //! \copydoc ImplicitProblem::dirichlet()
void dirichlet(PrimaryVariables &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -199,7 +199,7 @@ public:
}
}
- //! \copydoc BoxProblem::neumann()
+ //! \copydoc ImplicitProblem::neumann()
void neumann(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
@@ -216,7 +216,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::source()
+ //! \copydoc ImplicitProblem::source()
void source(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
@@ -227,7 +227,7 @@ public:
values = Scalar(0.0);
}
- //! \copydoc BoxProblem::initial()
+ //! \copydoc ImplicitProblem::initial()
void initial(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
diff --git a/test/freeflow/stokes2cni/stokes2cnitestproblem.hh b/test/freeflow/stokes2cni/stokes2cnitestproblem.hh
index 92a67cee16..4f9a60ad31 100644
--- a/test/freeflow/stokes2cni/stokes2cnitestproblem.hh
+++ b/test/freeflow/stokes2cni/stokes2cnitestproblem.hh
@@ -157,7 +157,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::boundaryTypes()
+ //! \copydoc ImplicitProblem::boundaryTypes()
void boundaryTypes(BoundaryTypes &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -177,7 +177,7 @@ public:
values.setDirichlet(massBalanceIdx);
}
- //! \copydoc BoxProblem::dirichlet()
+ //! \copydoc ImplicitProblem::dirichlet()
void dirichlet(PrimaryVariables &values, const Vertex &vertex) const
{
const GlobalPosition globalPos = vertex.geometry().center();
@@ -185,7 +185,7 @@ public:
initial_(values, globalPos);
}
- //! \copydoc BoxProblem::neumann()
+ //! \copydoc ImplicitProblem::neumann()
void neumann(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
@@ -203,7 +203,7 @@ public:
*/
// \{
- //! \copydoc BoxProblem::source()
+ //! \copydoc ImplicitProblem::source()
void source(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
@@ -214,7 +214,7 @@ public:
values = Scalar(0.0);
}
- //! \copydoc BoxProblem::initial()
+ //! \copydoc ImplicitProblem::initial()
void initial(PrimaryVariables &values,
const Element &element,
const FVElementGeometry &fvGeometry,
--
GitLab