FV discretizations: add checkConservativeness() debugging method

this method checks that the difference in the storage terms before and
after a time step is the same as the accumulated fluxes over the
domain boundary plus the source terms.
This commit is contained in:
Andreas Lauser 2014-07-17 17:18:28 +02:00
parent f7982e15ce
commit eefecea33f

View File

@ -200,6 +200,7 @@ class LensProblem : public GET_PROP_TYPE(TypeTag, BaseProblem)
dimWorld = GridView::dimensionworld
};
typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
typedef typename GET_PROP_TYPE(TypeTag, BoundaryRateVector) BoundaryRateVector;
@ -363,14 +364,18 @@ public:
*/
void endTimeStep()
{
#ifndef NDEBUG
this->model().checkConservativeness();
// Calculate storage terms
PrimaryVariables storage;
EqVector storage;
this->model().globalStorage(storage);
// Write mass balance information for rank 0
if (this->gridView().comm().rank() == 0) {
std::cout << "Storage: " << storage << std::endl << std::flush;
}
#endif // NDEBUG
}
//! \}
@ -384,8 +389,8 @@ public:
* \copydoc FvBaseProblem::boundary
*/
template <class Context>
void boundary(BoundaryRateVector &values, const Context &context,
int spaceIdx, int timeIdx) const
void boundary(BoundaryRateVector &values,
const Context &context, int spaceIdx, int timeIdx) const
{
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx);
@ -456,8 +461,7 @@ public:
* \copydoc FvBaseProblem::initial
*/
template <class Context>
void initial(PrimaryVariables &values, const Context &context, int spaceIdx,
int timeIdx) const
void initial(PrimaryVariables &values, const Context &context, int spaceIdx, int timeIdx) const
{
const GlobalPosition &pos = context.pos(spaceIdx, timeIdx);
Scalar depth = this->boundingBoxMax()[1] - pos[1];