Dear Dumux users,
i want to set up a simulation similar to
"ex_interface_coupling_ff-pm.cc" from the dumux-courses. However i
struggle with implementing periodic boundary conditions in x direction
for the staggered freeflow grid, it seems to work for the pm-grid.
Following this hint
<https://www.mail-archive.com/dumux@listserv.uni-stuttgart.de/msg01926.html>
from the mailing list i changed the yasp grid to a spgrid and tried to
make the grid periodic:
.input file: Periodic true false or
.dgf file: PERIODICFACETRANSFORMATION: 1 0, 0 1 + 1 0
I had no success and made a MWE for the ff-region which you can find in
the attachment.
I am thankful for every hint.
Best regards,
Lars
# executables for Navierstokeseriodic
dune_add_test(NAME stokesperiodic
SOURCES stokesperiodic.cc)
# add tutorial to the common target
add_dependencies(test_exercises stokesperiodic)
# add a symlink for each input file
add_input_file_links()
// -*- 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 The free flow sub problem
*/
#ifndef DUMUX_STOKES_SUBPROBLEM_HH
#define DUMUX_STOKES_SUBPROBLEM_HH
/*#include <dune/grid/yaspgrid.hh> #doesn't support periodic!*/
#include <dune/grid/spgrid.hh>
//****** uncomment for the third task *****//
// #include <dumux/io/grid/subgridmanager.hh>
#include <dumux/material/fluidsystems/1pliquid.hh>
#include <dumux/material/components/simpleh2o.hh>
#include <dumux/freeflow/navierstokes/problem.hh>
#include <dumux/discretization/staggered/freeflow/properties.hh>
#include <dumux/freeflow/navierstokes/model.hh>
namespace Dumux
{
template <class TypeTag>
class StokesSubProblem;
namespace Properties
{
// Create new type tags
namespace TTag {
struct StokesOneP { using InheritsFrom = std::tuple<NavierStokes, StaggeredFreeFlowModel>; };
} // end namespace TTag
// the fluid system
template<class TypeTag>
struct FluidSystem<TypeTag, TTag::StokesOneP>
{
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::SimpleH2O<Scalar> > ;
};
// Set the grid type
template<class TypeTag>
struct Grid<TypeTag, TTag::StokesOneP>
{
static constexpr auto dim = 2;
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using SPGrid = Dune::SPGrid<Scalar, 2>;
//using TensorGrid = Dune::YaspGrid<2, Dune::TensorProductCoordinates<Scalar, dim> >;
//****** comment out for the third task *****//
using type = SPGrid;
//****** uncomment for the third task *****//
// using HostGrid = TensorGrid;
// using type = Dune::SubGrid<dim, HostGrid>;
};
// Set the problem property
template<class TypeTag>
struct Problem<TypeTag, TTag::StokesOneP> { using type = Dumux::StokesSubProblem<TypeTag> ; };
template<class TypeTag>
struct EnableFVGridGeometryCache<TypeTag, TTag::StokesOneP> { static constexpr bool value = true; };
template<class TypeTag>
struct EnableGridFluxVariablesCache<TypeTag, TTag::StokesOneP> { static constexpr bool value = true; };
template<class TypeTag>
struct EnableGridVolumeVariablesCache<TypeTag, TTag::StokesOneP> { static constexpr bool value = true; };
}
/*!
* \brief The free flow sub problem
*/
template <class TypeTag>
class StokesSubProblem : public NavierStokesProblem<TypeTag>
{
using ParentType = NavierStokesProblem<TypeTag>;
using GridView = GetPropType<TypeTag, Properties::GridView>;
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
using BoundaryTypes = GetPropType<TypeTag, Properties::BoundaryTypes>;
using FVGridGeometry = GetPropType<TypeTag, Properties::FVGridGeometry>;
using FVElementGeometry = typename FVGridGeometry::LocalView;
using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
using Element = typename GridView::template Codim<0>::Entity;
using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
using NumEqVector = GetPropType<TypeTag, Properties::NumEqVector>;
using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
//using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
public:
StokesSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry)//, std::shared_ptr<CouplingManager> couplingManager)
: ParentType(fvGridGeometry, "Stokes"), eps_(1e-6)//, couplingManager_(couplingManager)
{
deltaP_ = getParamFromGroup<Scalar>(this->paramGroup(), "Problem.PressureDifference");
}
/*!
* \name Problem parameters
*/
// \{
/*!
* \brief Return the temperature within the domain in [K].
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
Scalar temperature() const
{ return 273.15 + 10; } // 10°C
/*!
* \brief Return the sources within the domain.
*
* \param globalPos The global position
*/
NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
{ return NumEqVector(0.0); }
// \}
/*!
* \name Boundary conditions
*/
// \{
/*!
* \brief Specifies which kind of boundary condition should be
* used for which equation on a given boundary segment.
*
* \param element The finite element
* \param scvf The sub control volume face
*/
BoundaryTypes boundaryTypes(const Element& element,
const SubControlVolumeFace& scvf) const
{
BoundaryTypes values;
const auto& globalPos = scvf.dofPosition();
if(onUpperBoundary_(globalPos))
{
values.setDirichlet(Indices::velocityXIdx);
values.setDirichlet(Indices::velocityYIdx);
}
// left/right wall
if (onRightBoundary_(globalPos) || (onLeftBoundary_(globalPos)))
{
values.setDirichlet(Indices::velocityXIdx);
values.setDirichlet(Indices::velocityYIdx);
}
// coupling interface
// if(couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf))
// {
// values.setCouplingNeumann(Indices::conti0EqIdx);
// values.setCouplingNeumann(Indices::momentumYBalanceIdx);
// values.setDirichlet(Indices::velocityXIdx); // assume no slip on interface
// }
return values;
}
/*!
* \brief Evaluate the boundary conditions for a Dirichlet control volume.
*
* \param globalPos The global position
*/
PrimaryVariables dirichletAtPos(const GlobalPosition& globalPos) const
{
PrimaryVariables values(0.0);
if(onUpperBoundary_(globalPos))
{
values=initialAtPos(globalPos);
}
// left/right wall
if (onRightBoundary_(globalPos))// || (onLeftBoundary_(globalPos)))
{ //periodic?!
values[Indices::velocityYIdx]=-1e-6 * (globalPos[1]-1);
}
// coupling interface
if(onLowerBoundary_(globalPos))
{
values[Indices::velocityXIdx]=0; // assume no slip on interface
}
return values;
}
/*!
* \brief Evaluate the boundary conditions for a Neumann control volume.
*
* \param element The element for which the Neumann boundary condition is set
* \param fvGeomentry The fvGeometry
* \param elemVolVars The element volume variables
* \param elemFaceVars The element face variables
* \param scvf The boundary sub control volume face
*/
template<class ElementVolumeVariables, class ElementFaceVariables>
NumEqVector neumann(const Element& element,
const FVElementGeometry& fvGeometry,
const ElementVolumeVariables& elemVolVars,
const ElementFaceVariables& elemFaceVars,
const SubControlVolumeFace& scvf) const
{
NumEqVector values(0.0);
// if(couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf))
// {
// values[Indices::conti0EqIdx] = couplingManager().couplingData().massCouplingCondition(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
// values[Indices::momentumYBalanceIdx] = couplingManager().couplingData().momentumCouplingCondition(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
// }
return values;
}
// \}
// //! Set the coupling manager
// void setCouplingManager(std::shared_ptr<CouplingManager> cm)
// { couplingManager_ = cm; }
//
// //! Get the coupling manager
// const CouplingManager& couplingManager() const
// { return *couplingManager_; }
/*!
* \name Volume terms
*/
// \{
/*!
* \brief Evaluate the initial value for a control volume.
*
* \param globalPos The global position
*/
PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
{
PrimaryVariables values(0.0);
values[Indices::velocityYIdx] = -1e-6 * globalPos[0];// * (this->fvGridGeometry().bBoxMax()[0] - globalPos[0]);
return values;
}
/*!
* \brief Returns the intrinsic permeability of required as input parameter for the Beavers-Joseph-Saffman boundary condition
*/
// Scalar permeability(const Element& element, const SubControlVolumeFace& scvf) const
// {
// return couplingManager().couplingData().darcyPermeability(element, scvf);
// }
/*!
* \brief Returns the alpha value required as input parameter for the Beavers-Joseph-Saffman boundary condition
*/
// Scalar alphaBJ(const SubControlVolumeFace& scvf) const
// {
// return couplingManager().problem(CouplingManager::darcyIdx).spatialParams().beaversJosephCoeffAtPos(scvf.center());
// }
/*!
* \brief calculate the analytical velocity in x direction based on Beavers & Joseph (1967)
*/
void calculateAnalyticalVelocityX() const
{
analyticalVelocityX_.resize(this->fvGridGeometry().gridView().size(0));
using std::sqrt;
const Scalar dPdX = -deltaP_ / (this->fvGridGeometry().bBoxMax()[0] - this->fvGridGeometry().bBoxMin()[0]);
static const Scalar mu = FluidSystem::viscosity(temperature(), 1e5);
static const Scalar alpha = getParam<Scalar>("Darcy.SpatialParams.AlphaBeaversJoseph");
static const Scalar K = getParam<Scalar>("Darcy.SpatialParams.Permeability");
static const Scalar sqrtK = sqrt(K);
const Scalar sigma = (this->fvGridGeometry().bBoxMax()[1] - this->fvGridGeometry().bBoxMin()[1])/sqrtK;
const Scalar uB = -K/(2.0*mu) * ((sigma*sigma + 2.0*alpha*sigma) / (1.0 + alpha*sigma)) * dPdX;
for (const auto& element : elements(this->fvGridGeometry().gridView()))
{
const auto eIdx = this->fvGridGeometry().gridView().indexSet().index(element);
const Scalar y = element.geometry().center()[1] - this->fvGridGeometry().bBoxMin()[1];
const Scalar u = uB*(1.0 + alpha/sqrtK*y) + 1.0/(2.0*mu) * (y*y + 2*alpha*y*sqrtK) * dPdX;
analyticalVelocityX_[eIdx] = u;
}
}
/*!
* \brief Get the analytical velocity in x direction
*/
const std::vector<Scalar>& getAnalyticalVelocityX() const
{
if(analyticalVelocityX_.empty())
calculateAnalyticalVelocityX();
return analyticalVelocityX_;
}
// \}
private:
bool onLeftBoundary_(const GlobalPosition &globalPos) const
{ return globalPos[0] < this->fvGridGeometry().bBoxMin()[0] + eps_; }
bool onRightBoundary_(const GlobalPosition &globalPos) const
{ return globalPos[0] > this->fvGridGeometry().bBoxMax()[0] - eps_; }
bool onLowerBoundary_(const GlobalPosition &globalPos) const
{ return globalPos[1] < this->fvGridGeometry().bBoxMin()[1] + eps_; }
bool onUpperBoundary_(const GlobalPosition &globalPos) const
{ return globalPos[1] > this->fvGridGeometry().bBoxMax()[1] - eps_; }
Scalar eps_;
Scalar deltaP_;
// std::shared_ptr<CouplingManager> couplingManager_;
mutable std::vector<Scalar> analyticalVelocityX_;
};
} //end namespace
#endif // DUMUX_STOKES_SUBPROBLEM_HH
DGF
INTERVAL
0 0 % lower left corner
1 1 % upper right corner
100 100 % number of cells in each direction
#
PERIODICFACETRANSFORMATION
1 0, 0 1 + 1 0 % make periodic in x-direction
#
BOUNDARYDOMAIN
DEFAULT 1
#
GRIDPARAMETER
OVERLAP 1
#
#
// -*- 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 A test problem for the coupled Stokes/Darcy problem (1p)
*/
#include <config.h>
#include <iostream>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/istl/io.hh>
#include <dumux/common/properties.hh>
#include <dumux/common/parameters.hh>
#include <dumux/common/dumuxmessage.hh>
#include <dumux/common/partial.hh>
#include <dumux/linear/seqsolverbackend.hh>
//#include <dumux/assembly/fvassembler.hh>
#include <dumux/assembly/staggeredfvassembler.hh>
#include <dumux/assembly/diffmethod.hh>
#include <dumux/discretization/method.hh>
#include <dumux/io/vtkoutputmodule.hh>
#include <dumux/io/staggeredvtkoutputmodule.hh>
#include <dumux/io/grid/gridmanager_sp.hh>
#include <dumux/multidomain/staggeredtraits.hh>
#include <dumux/multidomain/fvassembler.hh>
#include <dumux/multidomain/newtonsolver.hh>
// #include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
// #include "pmproblem.hh"
#include "ffproblem.hh"
// namespace Dumux {
// namespace Properties {
//
// template<class TypeTag>
// struct CouplingManager<TypeTag, TTag::StokesOneP>
// {
// using Traits = StaggeredMultiDomainTraits<TypeTag, TypeTag, Properties::TTag::DarcyOneP>;
// using type = Dumux::StokesDarcyCouplingManager<Traits>;
// };
//
// template<class TypeTag>
// struct CouplingManager<TypeTag, TTag::DarcyOneP>
// {
// using Traits = StaggeredMultiDomainTraits<Properties::TTag::StokesOneP, Properties::TTag::StokesOneP, TypeTag>;
// using type = Dumux::StokesDarcyCouplingManager<Traits>;
// };
//
// } // end namespace Properties
// } // end namespace Dumux
int main(int argc, char** argv) try
{
using namespace Dumux;
// initialize MPI, finalize is done automatically on exit
const auto& mpiHelper = Dune::MPIHelper::instance(argc, argv);
// print dumux start message
if (mpiHelper.rank() == 0)
DumuxMessage::print(/*firstCall=*/true);
// parse command line arguments and input file
Parameters::init(argc, argv);
// Define the sub problem type tags
using StokesTypeTag = Properties::TTag::StokesOneP;
// using DarcyTypeTag = Properties::TTag::DarcyOneP;
// create two individual grids (from the given grid file or the input file)
// for both sub-domains
// using DarcyGridManager = Dumux::GridManager<GetPropType<DarcyTypeTag, Properties::Grid>>;
// DarcyGridManager darcyGridManager;
// darcyGridManager.init("Darcy"); // pass parameter group
using StokesGridManager = Dumux::GridManager<GetPropType<StokesTypeTag, Properties::Grid>>;
StokesGridManager stokesGridManager;
stokesGridManager.init("Stokes"); // pass parameter group
// we compute on the leaf grid view
// const auto& darcyGridView = darcyGridManager.grid().leafGridView();
const auto& stokesGridView = stokesGridManager.grid().leafGridView();
// create the finite volume grid geometry
using StokesFVGridGeometry = GetPropType<StokesTypeTag, Properties::FVGridGeometry>;
auto stokesFvGridGeometry = std::make_shared<StokesFVGridGeometry>(stokesGridView);
stokesFvGridGeometry->update();
// using DarcyFVGridGeometry = GetPropType<DarcyTypeTag, Properties::FVGridGeometry>;
// auto darcyFvGridGeometry = std::make_shared<DarcyFVGridGeometry>(darcyGridView);
// darcyFvGridGeometry->update();
// Check if a periodic mesh is used
//stokesFvGridGeometry->setPeriodic(true);
bool periodic = stokesFvGridGeometry->isPeriodic();// && darcyFvGridGeometry->isPeriodic();
//periodic = darcyFvGridGeometry->gridView().comm().max(periodic);
if (!periodic)
DUNE_THROW(Dune::GridError, "Your grid is not periodic. Maybe the grid manager doesn't support periodic boundaries.");
// using Traits = StaggeredMultiDomainTraits<StokesTypeTag, StokesTypeTag, DarcyTypeTag>;
// the coupling manager
// using CouplingManager = StokesDarcyCouplingManager<Traits>;
// auto couplingManager = std::make_shared<CouplingManager>(stokesFvGridGeometry, darcyFvGridGeometry);
// the indices
// constexpr auto stokesCellCenterIdx = CouplingManager::stokesCellCenterIdx;
// constexpr auto stokesFaceIdx = CouplingManager::stokesFaceIdx;
constexpr auto stokesCellCenterIdx = StokesFVGridGeometry::cellCenterIdx();
constexpr auto stokesFaceIdx = StokesFVGridGeometry::faceIdx();
// constexpr auto darcyIdx = CouplingManager::darcyIdx;
// the problem (initial and boundary conditions)
using StokesProblem = GetPropType<StokesTypeTag, Properties::Problem>;
auto stokesProblem = std::make_shared<StokesProblem>(stokesFvGridGeometry);//, couplingManager);
// using DarcyProblem = GetPropType<DarcyTypeTag, Properties::Problem>;
// auto darcyProblem = std::make_shared<DarcyProblem>(darcyFvGridGeometry, couplingManager);
// the solution vector
//Traits::SolutionVector sol;
using SolutionVector = GetPropType<StokesTypeTag, Properties::SolutionVector>;
SolutionVector sol;
sol[stokesCellCenterIdx].resize(stokesFvGridGeometry->numCellCenterDofs());
sol[stokesFaceIdx].resize(stokesFvGridGeometry->numFaceDofs());
// sol[darcyIdx].resize(darcyFvGridGeometry->numDofs());
// auto stokesSol = partial(sol, stokesCellCenterIdx, stokesFaceIdx);
stokesProblem->applyInitialSolution(sol);
// darcyProblem->applyInitialSolution(sol[darcyIdx]);
// couplingManager->init(stokesProblem, darcyProblem, sol);
// the grid variables
using StokesGridVariables = GetPropType<StokesTypeTag, Properties::GridVariables>;
auto stokesGridVariables = std::make_shared<StokesGridVariables>(stokesProblem, stokesFvGridGeometry);
stokesGridVariables->init(sol);
// using DarcyGridVariables = GetPropType<DarcyTypeTag, Properties::GridVariables>;
// auto darcyGridVariables = std::make_shared<DarcyGridVariables>(darcyProblem, darcyFvGridGeometry);
// darcyGridVariables->init(sol[darcyIdx]);
// intialize the vtk output module
const auto stokesName = getParam<std::string>("Problem.Name") + "_" + stokesProblem->name();
// const auto darcyName = getParam<std::string>("Problem.Name") + "_" + darcyProblem->name();
StaggeredVtkOutputModule<StokesGridVariables, decltype(sol)> stokesVtkWriter(*stokesGridVariables, sol, stokesName);
GetPropType<StokesTypeTag, Properties::IOFields>::initOutputModule(stokesVtkWriter);
//****** uncomment the add analytical solution of v_x *****//
// stokesVtkWriter.addField(stokesProblem->getAnalyticalVelocityX(), "analyticalV_x");
// stokesVtkWriter.write(0.0);
// VtkOutputModule<DarcyGridVariables, GetPropType<DarcyTypeTag, Properties::SolutionVector>> darcyVtkWriter(*darcyGridVariables, sol[darcyIdx], darcyName);
// using DarcyVelocityOutput = GetPropType<DarcyTypeTag, Properties::VelocityOutput>;
// darcyVtkWriter.addVelocityOutput(std::make_shared<DarcyVelocityOutput>(*darcyGridVariables));
// GetPropType<DarcyTypeTag, Properties::IOFields>::initOutputModule(darcyVtkWriter);
// darcyVtkWriter.write(0.0);
// the assembler for a stationary problem
// using Assembler = MultiDomainFVAssembler<Traits, CouplingManager, DiffMethod::numeric>;
// auto assembler = std::make_shared<Assembler>(std::make_tuple(stokesProblem, stokesProblem, darcyProblem),
// std::make_tuple(stokesFvGridGeometry->cellCenterFVGridGeometryPtr(),
// stokesFvGridGeometry->faceFVGridGeometryPtr(),
// darcyFvGridGeometry),
// std::make_tuple(stokesGridVariables->cellCenterGridVariablesPtr(),
// stokesGridVariables->faceGridVariablesPtr(),
// darcyGridVariables),
// couplingManager);
using Assembler = StaggeredFVAssembler<StokesTypeTag, DiffMethod::numeric>;
auto assembler = std::make_shared<Assembler>(stokesProblem, stokesFvGridGeometry, stokesGridVariables);
// the linear solver
using LinearSolver = UMFPackBackend;
auto linearSolver = std::make_shared<LinearSolver>();
// the non-linear solver
// using NewtonSolver = MultiDomainNewtonSolver<Assembler, LinearSolver>;
using NewtonSolver = Dumux::NewtonSolver<Assembler, LinearSolver>;
NewtonSolver nonLinearSolver(assembler, linearSolver);
// solve the non-linear system
nonLinearSolver.solve(sol);
// write vtk output
stokesVtkWriter.write(1.0);
// darcyVtkWriter.write(1.0);
////////////////////////////////////////////////////////////
// finalize, print dumux message to say goodbye
////////////////////////////////////////////////////////////
// print dumux end message
if (mpiHelper.rank() == 0)
{
Parameters::print();
DumuxMessage::print(/*firstCall=*/false);
}
return 0;
} // end main
catch (Dumux::ParameterException &e)
{
std::cerr << std::endl << e << " ---> Abort!" << std::endl;
return 1;
}
catch (Dune::DGFException & e)
{
std::cerr << "DGF exception thrown (" << e <<
"). Most likely, the DGF file name is wrong "
"or the DGF file is corrupted, "
"e.g. missing hash at end of file or wrong number (dimensions) of entries."
<< " ---> Abort!" << std::endl;
return 2;
}
catch (Dune::Exception &e)
{
std::cerr << "Dune reported error: " << e << " ---> Abort!" << std::endl;
return 3;
}
catch (...)
{
std::cerr << "Unknown exception thrown! ---> Abort!" << std::endl;
return 4;
}
[Stokes.Grid]
File = Stokes.dgf
# SpGrid
#LowerLeft = 0.0 1.0
#UpperRight = 1.0 2.0
#Cells = 100 100
#Verbosity = true
#Periodic = true false
[Stokes.Problem]
Name = stokes
PressureDifference = 1e-9
EnableInertiaTerms = false
[Vtk]
AddVelocity = 1
_______________________________________________
Dumux mailing list
Dumux@listserv.uni-stuttgart.de
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux