Commit 9f61cb18 authored by Bernd Flemisch's avatar Bernd Flemisch
Browse files

[implicit] complete and test dim < dimWorld capability

All implicit porous-media models except 2pdfm are now able to run on
grids with dim < dimWorld. This still required some replacements of dim
by dimWorld, for example, in the velocity output of all models. In
implicit/1p, four new tests are added that run the 1p test problem on
1d-2d and 1d-3d Alberta grids with box and cell-centered, respectively.
Compilation has been tested also for all other models, but no runtime
testing has been performed.

This completes FS#183. Based upon preliminary work and a patch by
Natalie, thank you.

Reviewed by Natalie.


git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@14110 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent d51d1c55
......@@ -65,6 +65,7 @@ class OnePModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
enum { dofCodim = isBox ? dim : 0 };
......@@ -81,13 +82,13 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// create the required scalar fields
unsigned numDofs = this->numDofs();
ScalarField *p = writer.allocateManagedBuffer(numDofs);
ScalarField *K = writer.allocateManagedBuffer(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput())
......
......@@ -81,9 +81,8 @@ class OnePTwoCModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef typename GET_PROP_TYPE(TypeTag, SolutionVector) SolutionVector;
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
enum {
dim = GridView::dimension
};
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
......@@ -106,7 +105,7 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// create the required scalar fields
unsigned numDofs = this->numDofs();
......@@ -118,7 +117,7 @@ public:
ScalarField &massFraction1 = *writer.allocateManagedBuffer(numDofs);
ScalarField &rho = *writer.allocateManagedBuffer(numDofs);
ScalarField &mu = *writer.allocateManagedBuffer(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput())
......
......@@ -85,9 +85,8 @@ class TwoPModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
enum {
dim = GridView::dimension
};
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
......@@ -108,7 +107,7 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// get the number of degrees of freedom
unsigned numDofs = this->numDofs();
......@@ -125,8 +124,8 @@ public:
ScalarField *mobN = writer.allocateManagedBuffer(numDofs);
ScalarField *poro = writer.allocateManagedBuffer(numDofs);
ScalarField *Te = writer.allocateManagedBuffer(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput()) // check if velocity output is demanded
......
......@@ -296,7 +296,7 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// get the number of degrees of freedom
unsigned numDofs = this->numDofs();
......@@ -322,8 +322,8 @@ public:
moleFrac[phaseIdx][compIdx] = writer.allocateManagedBuffer(numDofs);
ScalarField *temperature = writer.allocateManagedBuffer(numDofs);
ScalarField *poro = writer.allocateManagedBuffer(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput()) // check if velocity output is demanded
......
......@@ -72,6 +72,7 @@ class ThreePModel: public GET_PROP_TYPE(TypeTag, BaseModel)
enum {
dim = GridView::dimension,
dimWorld = GridView::dimensionworld,
numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
......@@ -100,7 +101,7 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// get the number of degrees of freedom
unsigned numDofs = this->numDofs();
......@@ -119,9 +120,9 @@ public:
ScalarField *temperature = writer.allocateManagedBuffer (numDofs);
ScalarField *poro = writer.allocateManagedBuffer(numDofs);
ScalarField *perm = writer.allocateManagedBuffer(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityG = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityG = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput()) // check if velocity output is demanded
......
......@@ -297,7 +297,7 @@ public:
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// get the number of degrees of freedom
unsigned numDofs = this->numDofs();
......@@ -320,9 +320,9 @@ public:
moleFraction[i][j] = writer.allocateManagedBuffer (numDofs);
ScalarField *temperature = writer.allocateManagedBuffer (numDofs);
ScalarField *poro = writer.allocateManagedBuffer(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityG = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocityN = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityW = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
VectorField *velocityG = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput()) // check if velocity output is demanded
......
......@@ -92,10 +92,18 @@ public:
eg.subContVol[2].volume = eg.elementVolume/3;
break;
case 4: // 2D, quadrilinear
eg.subContVol[0].volume = eg.quadrilateralArea(eg.subContVol[0].global, edgeCoord[2], eg.elementGlobal, edgeCoord[0]);
eg.subContVol[1].volume = eg.quadrilateralArea(eg.subContVol[1].global, edgeCoord[1], eg.elementGlobal, edgeCoord[2]);
eg.subContVol[2].volume = eg.quadrilateralArea(eg.subContVol[2].global, edgeCoord[0], eg.elementGlobal, edgeCoord[3]);
eg.subContVol[3].volume = eg.quadrilateralArea(eg.subContVol[3].global, edgeCoord[3], eg.elementGlobal, edgeCoord[1]);
if (GlobalPosition::dimension == dim) {
eg.subContVol[0].volume = eg.quadrilateralArea(eg.subContVol[0].global, edgeCoord[2], eg.elementGlobal, edgeCoord[0]);
eg.subContVol[1].volume = eg.quadrilateralArea(eg.subContVol[1].global, edgeCoord[1], eg.elementGlobal, edgeCoord[2]);
eg.subContVol[2].volume = eg.quadrilateralArea(eg.subContVol[2].global, edgeCoord[0], eg.elementGlobal, edgeCoord[3]);
eg.subContVol[3].volume = eg.quadrilateralArea(eg.subContVol[3].global, edgeCoord[3], eg.elementGlobal, edgeCoord[1]);
}
else {
eg.subContVol[0].volume = eg.quadrilateralArea3D(eg.subContVol[0].global, edgeCoord[2], eg.elementGlobal, edgeCoord[0]);
eg.subContVol[1].volume = eg.quadrilateralArea3D(eg.subContVol[1].global, edgeCoord[1], eg.elementGlobal, edgeCoord[2]);
eg.subContVol[2].volume = eg.quadrilateralArea3D(eg.subContVol[2].global, edgeCoord[0], eg.elementGlobal, edgeCoord[3]);
eg.subContVol[3].volume = eg.quadrilateralArea3D(eg.subContVol[3].global, edgeCoord[3], eg.elementGlobal, edgeCoord[1]);
}
break;
default:
DUNE_THROW(Dune::NotImplemented, "_BoxFVElemGeomHelper::fillSubContVolData dim = " << dim << ", numVertices = " << numVertices);
......@@ -655,17 +663,17 @@ public:
// the compiler can optimize away all if
// cases which don't apply.
LocalPosition ipLocal;
GlobalPosition diffVec;
if (dim==1) {
subContVolFace[k].ipLocal = 0.5;
subContVolFace[k].normal = 1.0;
ipLocal = subContVolFace[k].ipLocal;
ipLocal = 0.5;
subContVolFace[k].ipLocal = ipLocal;
subContVolFace[k].normal = subContVol[j].global - subContVol[i].global;
subContVolFace[k].normal /= subContVolFace[k].normal.two_norm();
}
else if (dim==2) {
ipLocal = referenceElement.position(k, dim-1) + elementLocal;
ipLocal *= 0.5;
subContVolFace[k].ipLocal = ipLocal;
diffVec = elementGlobal - edgeCoordinates[k];
GlobalPosition diffVec = elementGlobal - edgeCoordinates[k];
subContVolFace[k].normal[0] = diffVec[1];
subContVolFace[k].normal[1] = -diffVec[0];
......
......@@ -103,7 +103,7 @@ public:
VertexIterator vIt = gridView.template begin<dim>();
const VertexIterator vEndIt = gridView.template end<dim>();
for (; vIt!=vEndIt; ++vIt) {
for (int i=0; i<dim; i++) {
for (int i=0; i<dimWorld; i++) {
bBoxMin_[i] = std::min(bBoxMin_[i], vIt->geometry().corner(0)[i]);
bBoxMax_[i] = std::max(bBoxMax_[i], vIt->geometry().corner(0)[i]);
}
......@@ -111,7 +111,7 @@ public:
// communicate to get the bounding box of the whole domain
if (gridView.comm().size() > 1)
for (int i = 0; i < dim; ++i) {
for (int i = 0; i < dimWorld; ++i) {
bBoxMin_[i] = gridView.comm().min(bBoxMin_[i]);
bBoxMax_[i] = gridView.comm().max(bBoxMax_[i]);
}
......
......@@ -210,7 +210,7 @@ public:
int vIdxGlobal = problem_.vertexMapper().map(element, scvIdx, dofCodim);
#endif
// calculate the subcontrolvolume velocity by the Piola transformation
Dune::FieldVector<CoordScalar, dim> scvVelocity(0);
Dune::FieldVector<CoordScalar, dimWorld> scvVelocity(0);
jacobianT2.mtv(scvVelocities[scvIdx], scvVelocity);
scvVelocity /= element.geometry().integrationElement(localPos)*cellNum_[vIdxGlobal];
......
......@@ -138,12 +138,13 @@ class MPNCFluxVariablesEnergy<TypeTag, /*enableEnergy=*/true, /*numEnergyEquatio
typedef typename GridView::template Codim<0>::Entity Element;
enum {dim = GridView::dimension};
enum {dimWorld = GridView::dimensionworld};
enum {numEnergyEqs = Indices::numPrimaryEnergyVars};
enum {wPhaseIdx = FluidSystem::wPhaseIdx};
enum {nPhaseIdx = FluidSystem::nPhaseIdx};
typedef Dune::FieldVector<CoordScalar, dim> GlobalPosition;
typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
typedef typename FVElementGeometry::SubControlVolume SCV;
typedef typename FVElementGeometry::SubControlVolumeFace SCVFace;
......
......@@ -58,6 +58,7 @@ class MPNCVtkWriterEnergy<TypeTag, /*enableEnergy = */ true, /* numEnergyEquatio
typedef typename GridView::template Codim<0>::Entity Element;
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numEnergyEqs = Indices::numPrimaryEnergyVars};
enum { velocityAveragingInModel = GET_PROP_VALUE(TypeTag, VelocityAveragingInModel) };
......@@ -76,9 +77,9 @@ class MPNCVtkWriterEnergy<TypeTag, /*enableEnergy = */ true, /* numEnergyEquatio
typedef typename ParentType::ComponentVector ComponentVector;
typedef Dune::array<ScalarVector, numEnergyEqs> EnergyEqVector;
typedef Dune::FieldVector<Scalar, dim> DimVector;
typedef Dune::BlockVector<DimVector> DimField;
typedef Dune::array<DimField, numPhases> PhaseDimField;
typedef Dune::FieldVector<Scalar, dimWorld> DimWorldVector;
typedef Dune::BlockVector<DimWorldVector> DimWorldField;
typedef Dune::array<DimWorldField, numPhases> PhaseDimWorldField;
public:
......@@ -266,7 +267,7 @@ private:
PhaseVector prandtlNumber_ ;
PhaseVector nusseltNumber_ ;
PhaseDimField velocity_;
PhaseDimWorldField velocity_;
ScalarVector awn_;
ScalarVector aws_;
......@@ -300,6 +301,7 @@ class MPNCVtkWriterEnergy<TypeTag, /*enableEnergy = */ true, /* numEnergyEquatio
typedef typename GridView::template Codim<0>::Entity Element;
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { wPhaseIdx = FluidSystem::wPhaseIdx};
......@@ -320,9 +322,9 @@ class MPNCVtkWriterEnergy<TypeTag, /*enableEnergy = */ true, /* numEnergyEquatio
typedef typename ParentType::ComponentVector ComponentVector;
typedef Dune::array<ScalarVector, numEnergyEqs> EnergyEqVector;
typedef Dune::FieldVector<Scalar, dim> DimVector;
typedef Dune::BlockVector<DimVector> DimField;
typedef Dune::array<DimField, numPhases> PhaseDimField;
typedef Dune::FieldVector<Scalar, dimWorld> DimWorldVector;
typedef Dune::BlockVector<DimWorldVector> DimWorldField;
typedef Dune::array<DimWorldField, numPhases> PhaseDimWorldField;
public:
......@@ -508,7 +510,7 @@ private:
PhaseVector nusseltNumber_ ;
ScalarVector qBoil_ ;
ScalarVector qsf_ ;
PhaseDimField velocity_;
PhaseDimWorldField velocity_;
};
......
......@@ -56,13 +56,14 @@ class MPNCVtkWriterCommon : public MPNCVtkWriterModule<TypeTag>
typedef typename ParentType::PhaseComponentMatrix PhaseComponentMatrix;
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
typedef Dune::FieldVector<Scalar, dim> DimVector;
typedef Dune::BlockVector<DimVector> DimField;
typedef Dune::array<DimField, numPhases> PhaseDimField;
typedef Dune::FieldVector<Scalar, dimWorld> DimWorldVector;
typedef Dune::BlockVector<DimWorldVector> DimWorldField;
typedef Dune::array<DimWorldField, numPhases> PhaseDimWorldField;
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
enum { dofCodim = isBox ? dim : 0 };
......@@ -239,7 +240,7 @@ private:
PhaseVector mobility_;
PhaseVector averageMolarMass_;
PhaseDimField velocity_;
PhaseDimWorldField velocity_;
ScalarVector porosity_;
ScalarVector permeability_;
......
......@@ -121,6 +121,7 @@ class RichardsModel : public GET_PROP_TYPE(TypeTag, BaseModel)
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
enum { dim = GridView::dimension };
enum { dimWorld = GridView::dimensionworld };
enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
enum { dofCodim = isBox ? dim : 0 };
......@@ -137,7 +138,7 @@ public:
void addOutputVtkFields(const SolutionVector &sol, MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<double, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<double, dim> > VectorField;
typedef Dune::BlockVector<Dune::FieldVector<double, dimWorld> > VectorField;
// get the number of degrees of freedom
unsigned numDofs = this->numDofs();
......@@ -154,7 +155,7 @@ public:
ScalarField *mobN = writer.allocateManagedBuffer(numDofs);
ScalarField *poro = writer.allocateManagedBuffer(numDofs);
ScalarField *Te = writer.allocateManagedBuffer(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dim>(numDofs);
VectorField *velocity = writer.template allocateManagedBuffer<double, dimWorld>(numDofs);
ImplicitVelocityOutput<TypeTag> velocityOutput(this->problem_());
if (velocityOutput.enableOutput())
......
......@@ -87,6 +87,19 @@ NEW_TYPE_TAG(OnePTestCCProblemWithAMG, INHERITS_FROM(OnePTestCCProblem));
SET_TYPE_PROP(OnePTestBoxProblemWithAMG, LinearSolver, Dumux::AMGBackend<TypeTag> );
SET_TYPE_PROP(OnePTestCCProblemWithAMG, LinearSolver, Dumux::AMGBackend<TypeTag> );
// if AlbertaGrid is available, test for dim < dimWorld
#if HAVE_ALBERTA
NEW_TYPE_TAG(OnePOneDThreeDTestProblem, INHERITS_FROM(OnePTestProblem));
NEW_TYPE_TAG(OnePOneDThreeDTestBoxProblem, INHERITS_FROM(BoxModel, OnePOneDThreeDTestProblem));
NEW_TYPE_TAG(OnePOneDThreeDTestCCProblem, INHERITS_FROM(CCModel, OnePOneDThreeDTestProblem));
SET_TYPE_PROP(OnePOneDThreeDTestProblem, Grid, Dune::AlbertaGrid<1, 3>);
NEW_TYPE_TAG(OnePTwoDThreeDTestProblem, INHERITS_FROM(OnePTestProblem));
NEW_TYPE_TAG(OnePTwoDThreeDTestBoxProblem, INHERITS_FROM(BoxModel, OnePTwoDThreeDTestProblem));
NEW_TYPE_TAG(OnePTwoDThreeDTestCCProblem, INHERITS_FROM(CCModel, OnePTwoDThreeDTestProblem));
SET_TYPE_PROP(OnePTwoDThreeDTestProblem, Grid, Dune::AlbertaGrid<2, 3>);
#endif
// Enable gravity
SET_BOOL_PROP(OnePTestProblem, ProblemEnableGravity, true);
}
......@@ -206,7 +219,7 @@ public:
const GlobalPosition &globalPos) const
{
double eps = 1.0e-3;
if (globalPos[dim-1] < eps || globalPos[dim-1] > this->bBoxMax()[dim-1] - eps)
if (globalPos[dimWorld-1] < eps || globalPos[dimWorld-1] > this->bBoxMax()[dimWorld-1] - eps)
values.setAllDirichlet();
else
values.setAllNeumann();
......@@ -224,7 +237,7 @@ public:
void dirichletAtPos(PrimaryVariables &values,
const GlobalPosition &globalPos) const
{
values[pressureIdx] = 1.0e+5*(2.0 - globalPos[dim-1]);
values[pressureIdx] = 1.0e+5*(2.0 - globalPos[dimWorld-1]);
}
/*!
......
......@@ -61,15 +61,15 @@ public:
try
{
lensLowerLeft_[0] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensLowerLeftX);
if (dim > 1)
if (dimWorld > 1)
lensLowerLeft_[1] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensLowerLeftY);
if (dim > 2)
if (dimWorld > 2)
lensLowerLeft_[2] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensLowerLeftZ);
lensUpperRight_[0] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensUpperRightX);
if (dim > 1)
if (dimWorld > 1)
lensUpperRight_[1] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensUpperRightY);
if (dim > 2)
if (dimWorld > 2)
lensUpperRight_[2] = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensUpperRightZ);
permeability_ = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.Permeability);
......@@ -119,7 +119,7 @@ public:
private:
bool isInLens_(const GlobalPosition &globalPos) const
{
for (int i = 0; i < dim; ++i) {
for (int i = 0; i < dimWorld; ++i) {
if (globalPos[i] < lensLowerLeft_[i] || globalPos[i] > lensUpperRight_[i])
return false;
}
......
......@@ -55,3 +55,37 @@ add_dumux_test(test_cc1pniconvection test_cc1pniconvection test_cc1pniconvection
${CMAKE_SOURCE_DIR}/test/references/1pniccconvection-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/cc1pniconvection-00011.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_cc1pniconvection)
# dim < dimWorld tests
add_dumux_test(test_box1p1d3d test_box1p1d3d test_box1p1d3d.cc
${CMAKE_SOURCE_DIR}/bin/runTest.sh
${CMAKE_SOURCE_DIR}/bin/fuzzycomparevtu.py
${CMAKE_SOURCE_DIR}/test/references/1ptestbox1d3d-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/1ptestbox1d3d-00002.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_box1p1d3d)
add_dune_alberta_flags(test_box1p1d3d GRIDDIM 1 WORLDDIM 3)
add_dumux_test(test_box1p2d3d test_box1p2d3d test_box1p2d3d.cc
${CMAKE_SOURCE_DIR}/bin/runTest.sh
${CMAKE_SOURCE_DIR}/bin/fuzzycomparevtu.py
${CMAKE_SOURCE_DIR}/test/references/1ptestbox2d3d-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/1ptestbox2d3d-00002.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_box1p2d3d)
add_dune_alberta_flags(test_box1p2d3d GRIDDIM 2 WORLDDIM 3)
add_dumux_test(test_cc1p1d3d test_cc1p1d3d test_cc1p1d3d.cc
${CMAKE_SOURCE_DIR}/bin/runTest.sh
${CMAKE_SOURCE_DIR}/bin/fuzzycomparevtu.py
${CMAKE_SOURCE_DIR}/test/references/1ptestcc1d3d-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/1ptestcc1d3d-00002.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_cc1p1d3d)
add_dune_alberta_flags(test_cc1p1d3d GRIDDIM 1 WORLDDIM 3)
add_dumux_test(test_cc1p2d3d test_cc1p2d3d test_cc1p2d3d.cc
${CMAKE_SOURCE_DIR}/bin/runTest.sh
${CMAKE_SOURCE_DIR}/bin/fuzzycomparevtu.py
${CMAKE_SOURCE_DIR}/test/references/1ptestcc2d3d-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/1ptestcc2d3d-00002.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_cc1p2d3d)
add_dune_alberta_flags(test_cc1p2d3d GRIDDIM 2 WORLDDIM 3)
check_PROGRAMS = test_box1p \
test_box1p1d3d \
test_box1p2d3d \
test_box1pniconduction \
test_box1pniconvection \
test_box1pwithamg \
test_cc1p \
test_cc1p1d3d \
test_cc1p2d3d \
test_cc1pniconduction \
test_cc1pniconvection \
test_cc1pwithamg
test_cc1pwithamg
noinst_HEADERS := $(wildcard *.hh)
EXTRA_DIST:=$(wildcard *.input) $(wildcard grids/*.dgf) CMakeLists.txt
test_box1p_SOURCES = test_box1p.cc
test_box1p1d3d_SOURCES = test_box1p1d3d.cc
test_box1p2d3d_SOURCES = test_box1p2d3d.cc
test_box1pniconduction_SOURCES = test_box1pniconduction.cc
test_box1pniconvection_SOURCES = test_box1pniconvection.cc
test_box1pwithamg_SOURCES = test_box1pwithamg.cc
test_cc1p_SOURCES = test_cc1p.cc
test_cc1p1d3d_SOURCES = test_cc1p1d3d.cc
test_cc1p2d3d_SOURCES = test_cc1p2d3d.cc
test_cc1pniconduction_SOURCES = test_cc1pniconduction.cc
test_cc1pniconvection_SOURCES = test_cc1pniconvection.cc
test_cc1pwithamg_SOURCES = test_cc1pwithamg.cc
test_cc1p1d3d_LDADD = $(AM_LDFLAGS) $(ALBERTA3D_LIBS)
test_cc1p1d3d_LDFLAGS= $(AM_LDFLAGS) $(ALBERTA3D_LDFLAGS)
test_cc1p1d3d_CPPFLAGS=$(AM_CPPFLAGS) $(ALBERTA3D_CPPFLAGS)
test_box1p1d3d_LDADD = $(AM_LDFLAGS) $(ALBERTA3D_LIBS)
test_box1p1d3d_LDFLAGS= $(AM_LDFLAGS) $(ALBERTA3D_LDFLAGS)
test_box1p1d3d_CPPFLAGS=$(AM_CPPFLAGS) $(ALBERTA3D_CPPFLAGS)
test_cc1p2d3d_LDADD = $(AM_LDFLAGS) $(ALBERTA3D_LIBS)
test_cc1p2d3d_LDFLAGS= $(AM_LDFLAGS) $(ALBERTA3D_LDFLAGS)
test_cc1p2d3d_CPPFLAGS=$(AM_CPPFLAGS) $(ALBERTA3D_CPPFLAGS)
test_box1p2d3d_LDADD = $(AM_LDFLAGS) $(ALBERTA3D_LIBS)
test_box1p2d3d_LDFLAGS= $(AM_LDFLAGS) $(ALBERTA3D_LDFLAGS)
test_box1p2d3d_CPPFLAGS=$(AM_CPPFLAGS) $(ALBERTA3D_CPPFLAGS)
include $(top_srcdir)/am/global-rules
DGF
Vertex
0 0 0
0.25 0.25 0.25
0.5 0.5 0.5
0.75 0.75 0.75
1 1 1
#
Simplex
0 1
1 2
2 3
3 4
#
Boundarydomain
default 1 % all boundaries have id 1
#
# unitcube.dgf
DGF
Vertex
0 0 0
0.5 0 0
1 0 0
0 0.5 0.5
0.5 0.5 0.5
1 0.5 0.5
0 1 1
0.5 1 1
1 1 1
#
Simplex
0 1 4
0 4 3
1 2 4
2 5 4
3 4 6
4 7 6
4 5 8
4 8 7
#
Boundarydomain
default 1 % all boundaries have id 1
#
# unitcube.dgf
// -*- 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/>. *
*****************************************************************************/