Commit c8d099f1 authored by Thomas Fetzer's avatar Thomas Fetzer
Browse files

[boundarytypes,non-isothermal] Revert boundarytypes to almost original state...

[boundarytypes,non-isothermal] Revert boundarytypes to almost original state (see below), update the non-isothermal models and reference solutions.

Problems occured if the boundary types are herited or more precise the problem was to reset out of the parent class the boundary types from the childs classes.
parent b4fd8324
......@@ -40,29 +40,41 @@ public:
{ reset(); }
/*!
* \brief Reset the boundary types.
* \brief Reset the boundary types for all equations.
*
* After this method no equations will be disabled and neither
* Neumann nor Dirichlet conditions will be evaluated. This
* corrosponds to a Neumann zero boundary.
* corresponds to a Neumann zero boundary.
*/
void reset()
{
for (int i=0; i < numEq; ++i) {
boundaryInfo_[i].visited = 0;
boundaryInfo_[i].isDirichlet = 0;
boundaryInfo_[i].isNeumann = 0;
boundaryInfo_[i].isOutflow = 0;
boundaryInfo_[i].isCouplingInflow = 0;
boundaryInfo_[i].isCouplingOutflow = 0;
boundaryInfo_[i].isMortarCoupling = 0;
eq2pvIdx_[i] = i;
pv2eqIdx_[i] = i;
for (int eqIdx=0; eqIdx < numEq; ++eqIdx)
{
resetEq(eqIdx);
}
}
/*!
* \brief Reset the boundary types for one equation.
*/
void resetEq(int eqIdx)
{
boundaryInfo_[eqIdx].visited = 0;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
boundaryInfo_[eqIdx].isCouplingDirichlet = 0;
boundaryInfo_[eqIdx].isCouplingNeumann = 0;
boundaryInfo_[eqIdx].isCouplingMortar = 0;
eq2pvIdx_[eqIdx] = eqIdx;
pv2eqIdx_[eqIdx] = eqIdx;
}
/*!
* \brief Returns true if the boundary types for a given equation
* has been specified.
......@@ -75,7 +87,7 @@ public:
/*!
* \brief Make sure the boundary conditions are well-posed.
*
* If they are not, an assertation fails and the program aborts!
* If they are not, an assertion fails and the program aborts!
* (if the NDEBUG macro is not defined)
*/
void checkWellPosed() const
......@@ -120,22 +132,47 @@ public:
}
}
/*!
* \brief Set all boundary conditions to Dirichlet-like coupling
*/
void setAllCouplingDirichlet()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingDirichlet(eqIdx);
}
}
/*!
* \brief Set all boundary conditions to Neumann-like coupling.
*/
void setAllCouplingNeumann()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingNeumann(eqIdx);
}
}
/*!
* \brief Set all boundary conditions to mortar coupling.
*/
void setAllCouplingMortar()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingDirichlet(eqIdx);
}
}
/*!
* \brief Set all boundary conditions to coupling inflow.
*/
DUNE_DEPRECATED_MSG("setAllCouplingInflow() is deprecated")
void setAllCouplingInflow()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 1;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingInflow(eqIdx);
}
}
......@@ -145,16 +182,9 @@ public:
DUNE_DEPRECATED_MSG("setAllCouplingOutflow() is deprecated")
void setAllCouplingOutflow()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 1;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingOutflow(eqIdx);
}
}
......@@ -164,16 +194,9 @@ public:
DUNE_DEPRECATED_MSG("setAllMortarCoupling() is deprecated")
void setAllMortarCoupling()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 1;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingMortar(eqIdx);
}
}
......@@ -185,13 +208,9 @@ public:
*/
void setNeumann(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 1;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
......@@ -207,13 +226,9 @@ public:
*/
void setDirichlet(int pvIdx, int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 1;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
// update the equation <-> primary variable mapping
eq2pvIdx_[eqIdx] = pvIdx;
......@@ -231,13 +246,48 @@ public:
*/
void setOutflow(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 1;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
/*!
* \brief Set a boundary condition for a single equation to
* a Dirichlet-like coupling condition.
*/
void setCouplingDirichlet(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isCouplingDirichlet = 1;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
/*!
* \brief Set a boundary condition for a single equation to
* a Neumann-like coupling condition.
*/
void setCouplingNeumann(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isCouplingNeumann = 1;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
/*!
* \brief Set a boundary condition for a single equation to
* a mortar coupling condition.
*/
void setCouplingMortar(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isCouplingMortar = 1;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
......@@ -248,13 +298,9 @@ public:
DUNE_DEPRECATED_MSG("setCouplingInflow() is deprecated")
void setCouplingInflow(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 1;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
......@@ -265,13 +311,9 @@ public:
DUNE_DEPRECATED_MSG("setCouplingOutflow() is deprecated")
void setCouplingOutflow(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 1;
boundaryInfo_[eqIdx].isMortarCoupling = 0;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
}
......@@ -282,12 +324,8 @@ public:
DUNE_DEPRECATED_MSG("setMortarCoupling() is deprecated")
void setMortarCoupling(int eqIdx)
{
resetEq(eqIdx);
boundaryInfo_[eqIdx].visited = 1;
boundaryInfo_[eqIdx].isDirichlet = 0;
boundaryInfo_[eqIdx].isNeumann = 0;
boundaryInfo_[eqIdx].isOutflow = 0;
boundaryInfo_[eqIdx].isCouplingInflow = 0;
boundaryInfo_[eqIdx].isCouplingOutflow = 0;
boundaryInfo_[eqIdx].isMortarCoupling = 1;
Valgrind::SetDefined(boundaryInfo_[eqIdx]);
......@@ -442,17 +480,82 @@ public:
return false;
}
/*!
* \brief Returns true if an equation is used to specify an
* Dirichlet coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingDirichlet(unsigned eqIdx) const
{ return boundaryInfo_[eqIdx].isCouplingDirichlet; }
/*!
* \brief Returns true if some equation is used to specify an
* Dirichlet coupling condition.
*/
bool hasCouplingDirichlet() const
{
for (int i = 0; i < numEq; ++i)
if (boundaryInfo_[i].isCouplingDirichlet)
return true;
return false;
}
/*!
* \brief Returns true if an equation is used to specify an
* Neumann coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingNeumann(unsigned eqIdx) const
{ return boundaryInfo_[eqIdx].isCouplingNeumann; }
/*!
* \brief Returns true if some equation is used to specify an
* Neumann coupling condition.
*/
bool hasCouplingNeumann() const
{
for (int i = 0; i < numEq; ++i)
if (boundaryInfo_[i].isCouplingNeumann)
return true;
return false;
}
/*!
* \brief Returns true if an equation is used to specify a
* coupling condition.
* Mortar coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingMortar(unsigned eqIdx) const
{
return boundaryInfo_[eqIdx].isCouplingMortar;
}
/*!
* \brief Returns true if some equation is used to specify an
* Mortar coupling condition.
*/
bool hasCouplingMortar() const
{
for (int i = 0; i < numEq; ++i)
if (boundaryInfo_[i].isCouplingMortar)
return true;
return false;
}
/*!
* \brief Returns true if an equation is used to specify a
* coupling condition.
*
* This only used for correcting the pressure in the Stokes equation.
* \param eqIdx The index of the equation
*/
bool isCoupling(unsigned eqIdx) const
{
return false;
return boundaryInfo_[eqIdx].isCouplingDirichlet
|| boundaryInfo_[eqIdx].isCouplingNeumann
|| boundaryInfo_[eqIdx].isCouplingMortar;
}
/*!
......@@ -461,6 +564,9 @@ public:
*/
bool hasCoupling() const
{
for (int i = 0; i < numEq; ++i)
if (isCoupling(i))
return true;
return false;
}
......@@ -486,7 +592,7 @@ public:
unsigned eqToDirichletIndex(unsigned eqIdx) const
{ return eq2pvIdx_[eqIdx]; }
private:
protected:
// this is a bitfield structure!
struct __attribute__((__packed__)) {
unsigned char visited : 1;
......@@ -496,6 +602,9 @@ private:
unsigned char isCouplingInflow : 1;
unsigned char isCouplingOutflow : 1;
unsigned char isMortarCoupling : 1;
unsigned char isCouplingDirichlet : 1;
unsigned char isCouplingNeumann : 1;
unsigned char isCouplingMortar : 1;
} boundaryInfo_[numEq];
unsigned char eq2pvIdx_[numEq];
......
......@@ -220,7 +220,7 @@ public:
cParams,
couplingRes1, couplingRes2);
if (cParams.boundaryTypes2.isCouplingInflow(energyEqIdx2))
if (cParams.boundaryTypes2.isCouplingNeumann(energyEqIdx2))
{
unsigned int blModel = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, BoundaryLayer, Model);
// only enter here, if a boundary layer model is used for the computation of the diffusive fluxes
......@@ -342,7 +342,7 @@ public:
globalProblem.localResidual1().residual(vertInElem1)[energyEqIdx1]);
}
}
if (cParams.boundaryTypes2.isCouplingOutflow(energyEqIdx2))
if (cParams.boundaryTypes2.isCouplingDirichlet(energyEqIdx2))
{
// set residualDarcy[energyEqIdx2] = T in 2p2cnilocalresidual.hh
couplingRes2.accumulate(lfsu_n.child(energyEqIdx2), vertInElem2,
......@@ -378,13 +378,13 @@ public:
normalMassFlux2[phaseIdx] = -boundaryVars2.volumeFlux(phaseIdx)*
cParams.elemVolVarsCur2[vertInElem2].density(phaseIdx);
if (cParams.boundaryTypes1.isCouplingOutflow(energyEqIdx1))
if (cParams.boundaryTypes1.isCouplingDirichlet(energyEqIdx1))
{
// set residualStokes[energyIdx1] = T in stokes2cnilocalresidual.hh
couplingRes1.accumulate(lfsu_s.child(energyEqIdx1), vertInElem1,
-cParams.elemVolVarsCur2[vertInElem2].temperature());
}
if (cParams.boundaryTypes1.isCouplingInflow(energyEqIdx1))
if (cParams.boundaryTypes1.isCouplingNeumann(energyEqIdx1))
{
if (globalProblem.sdProblem2().isCornerPoint(globalPos2))
{
......
// -*- 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 Class to specify the type of a boundary for multidomain problems.
*/
#ifndef MULTIDOMAIN_BOUNDARY_TYPES_HH
#define MULTIDOMAIN_BOUNDARY_TYPES_HH
#include <dumux/common/boundarytypes.hh>
#include <dumux/common/valgrind.hh>
namespace Dumux
{
/*!
* \ingroup BC
* \brief Class to specify the type of a boundary for multidomain problems.
*/
template <int numEq>
class MultidomainBoundaryTypes : public BoundaryTypes<numEq>
{
typedef BoundaryTypes<numEq> ParentType;
public:
MultidomainBoundaryTypes()
{ reset(); }
/*!
* \brief Reset the boundary types.
*
* After this method no equations will be disabled and neither
* neumann nor dirichlet conditions will be evaluated. This
* corrosponds to a neumann zero boundary.
*/
void reset()
{
ParentType::reset();
for (int i=0; i < numEq; ++i) {
boundaryCouplingInfo_[i].visited = 0;
boundaryCouplingInfo_[i].isCouplingDirichlet = 0;
boundaryCouplingInfo_[i].isCouplingNeumann = 0;
}
}
/*!
* \brief Returns true if the boundary types for a given equation
* has been specified.
*
* \param eqIdx The index of the equation
*/
bool isSet(int eqIdx) const
{
return ParentType::boundaryInfo_[eqIdx].visited
|| boundaryCouplingInfo_[eqIdx].visited;
}
/*!
* \brief Set all boundary conditions to coupling Dirichlet.
*/
void setAllCouplingDirichlet()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingDirichlet(eqIdx);
}
}
/*!
* \brief Set all boundary conditions to coupling Neumann.
*/
void setAllCouplingNeumann()
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
setCouplingNeumann(eqIdx);
}
}
/*!
* \brief Set a boundary condition for a single equation to coupling inflow.
*/
void setCouplingDirichlet(int eqIdx)
{
ParentType::reset();
boundaryCouplingInfo_[eqIdx].visited = 1;
boundaryCouplingInfo_[eqIdx].isCouplingDirichlet = 1;
boundaryCouplingInfo_[eqIdx].isCouplingNeumann = 0;
Valgrind::SetDefined(boundaryCouplingInfo_[eqIdx]);
}
/*!
* \brief Set a boundary condition for a single equation to coupling outflow.
*/
void setCouplingNeumann(int eqIdx)
{
ParentType::reset();
boundaryCouplingInfo_[eqIdx].visited = 1;
boundaryCouplingInfo_[eqIdx].isCouplingDirichlet = 0;
boundaryCouplingInfo_[eqIdx].isCouplingNeumann = 1;
Valgrind::SetDefined(boundaryCouplingInfo_[eqIdx]);
}
/*!
* \brief Returns true if an equation is used to specify an
* Dirichlet coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingDirichlet(unsigned eqIdx) const
{ return boundaryCouplingInfo_[eqIdx].isCouplingDirichlet; }
/*!
* \brief Returns true if some equation is used to specify an
* Dirichlet coupling condition.
*/
bool hasCouplingDirichlet() const
{
for (int i = 0; i < numEq; ++i)
if (boundaryCouplingInfo_[i].isCouplingDirichlet)
return true;
return false;
}
/*!
* \brief Returns true if an equation is used to specify an
* Neumann coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingNeumann(unsigned eqIdx) const
{ return boundaryCouplingInfo_[eqIdx].isCouplingNeumann; }
/*!
* \brief Returns true if some equation is used to specify an
* Neumann coupling condition.
*/
bool hasCouplingNeumann() const
{
for (int i = 0; i < numEq; ++i)
if (boundaryCouplingInfo_[i].isCouplingNeumann)
return true;
return false;
}
/*!
* \brief Returns true if an equation is used to specify a
* Mortar coupling condition.
*
* \param eqIdx The index of the equation
*/
bool isCouplingMortar(unsigned eqIdx) const
{
return false;
}
/*!
* \brief Returns true if some equation is used to specify an
* Mortar coupling condition.
*/ </