From 92c1467a1cf14d055475762510a545184730d27c Mon Sep 17 00:00:00 2001
From: Edward 'Ned' Coltman <edward.coltman@iws.uni-stuttgart.de>
Date: Sun, 9 Jan 2022 05:21:14 +0100
Subject: [PATCH] [test][navierstokesnc] Port density driven flow test
---
.../densitydrivenflow/CMakeLists.txt | 11 +-
.../navierstokesnc/densitydrivenflow/main.cc | 106 +++++---
.../densitydrivenflow/problem.hh | 243 ++++++++++--------
.../densitydrivenflow/properties.hh | 52 +++-
4 files changed, 250 insertions(+), 162 deletions(-)
diff --git a/test/freeflow/navierstokesnc/densitydrivenflow/CMakeLists.txt b/test/freeflow/navierstokesnc/densitydrivenflow/CMakeLists.txt
index d439f49b32..acbf3779a1 100644
--- a/test/freeflow/navierstokesnc/densitydrivenflow/CMakeLists.txt
+++ b/test/freeflow/navierstokesnc/densitydrivenflow/CMakeLists.txt
@@ -1,9 +1,10 @@
# SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
# SPDX-License-Identifier: GPL-3.0-or-later
+dune_symlink_to_source_files(FILES "params.input")
dumux_add_test(NAME test_ff_stokes2c_densitydrivenflow
- LABELS freeflow navierstokes
+ LABELS freeflow navierstokesnc navierstokes
SOURCES main.cc
CMAKE_GUARD HAVE_UMFPACK
COMMAND ${CMAKE_SOURCE_DIR}/bin/testing/runtest.py
@@ -11,10 +12,4 @@ dumux_add_test(NAME test_ff_stokes2c_densitydrivenflow
--files ${CMAKE_SOURCE_DIR}/test/references/test_ff_stokes2c_densitydrivenflow-reference.vtu
${CMAKE_CURRENT_BINARY_DIR}/test_ff_stokes2c_densitydrivenflow-00021.vtu
--command "${CMAKE_CURRENT_BINARY_DIR}/test_ff_stokes2c_densitydrivenflow params.input
- -Problem.Name test_ff_stokes2c_densitydrivenflow"
- --zeroThreshold {"X^Air_liq":1e-12}
- --zeroThreshold {"x^Air_liq":1e-12}
- --zeroThreshold {"velocity_liq \(m/s\)":1e-10}
- --zeroThreshold {"deltaRho":1e-7})
-
-dune_symlink_to_source_files(FILES "params.input")
+ -Problem.Name test_ff_stokes2c_densitydrivenflow")
diff --git a/test/freeflow/navierstokesnc/densitydrivenflow/main.cc b/test/freeflow/navierstokesnc/densitydrivenflow/main.cc
index 14f759c687..1ae9410047 100644
--- a/test/freeflow/navierstokesnc/densitydrivenflow/main.cc
+++ b/test/freeflow/navierstokesnc/densitydrivenflow/main.cc
@@ -7,7 +7,7 @@
/*!
* \file
* \ingroup NavierStokesNCTests
- * \brief Test for the staggered grid multi-component (Navier-)Stokes model.
+ * \brief Density Driven Test for the staggered grid multi-component (Navier-)Stokes model
*/
#include <config.h>
@@ -17,21 +17,25 @@
#include <dune/common/parallel/mpihelper.hh>
#include <dune/common/timer.hh>
-#include <dune/grid/io/file/vtk.hh>
#include <dune/istl/io.hh>
-#include <dumux/assembly/staggeredfvassembler.hh>
#include <dumux/assembly/diffmethod.hh>
#include <dumux/common/initialize.hh>
#include <dumux/common/dumuxmessage.hh>
#include <dumux/common/parameters.hh>
#include <dumux/common/properties.hh>
+
#include <dumux/io/grid/gridmanager_yasp.hh>
-#include <dumux/io/staggeredvtkoutputmodule.hh>
#include <dumux/linear/istlsolvers.hh>
#include <dumux/linear/linearsolvertraits.hh>
#include <dumux/linear/linearalgebratraits.hh>
#include <dumux/nonlinear/newtonsolver.hh>
+#include <dumux/io/vtkoutputmodule.hh>
+#include <dumux/freeflow/navierstokes/velocityoutput.hh>
+
+#include <dumux/multidomain/fvassembler.hh>
+#include <dumux/multidomain/traits.hh>
+#include <dumux/multidomain/newtonsolver.hh>
#include "properties.hh"
@@ -40,7 +44,8 @@ int main(int argc, char** argv)
using namespace Dumux;
// define the type tag for this problem
- using TypeTag = Properties::TTag::DensityDrivenFlow;
+ using MomentumTypeTag = Properties::TTag::DensityDrivenFlowMomentum;
+ using MassTypeTag = Properties::TTag::DensityDrivenFlowMass;
// maybe initialize MPI and/or multithreading backend
initialize(argc, argv);
@@ -54,7 +59,7 @@ int main(int argc, char** argv)
Parameters::init(argc, argv);
// try to create a grid (from the given grid file or the input file)
- GridManager<GetPropType<TypeTag, Properties::Grid>> gridManager;
+ GridManager<GetPropType<MomentumTypeTag, Properties::Grid>> gridManager;
gridManager.init();
////////////////////////////////////////////////////////////
@@ -65,54 +70,84 @@ int main(int argc, char** argv)
const auto& leafGridView = gridManager.grid().leafGridView();
// create the finite volume grid geometry
- using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
- auto gridGeometry = std::make_shared<GridGeometry>(leafGridView);
+ using MomentumGridGeometry = GetPropType<MomentumTypeTag, Properties::GridGeometry>;
+ auto momentumGridGeometry = std::make_shared<MomentumGridGeometry>(leafGridView);
+ using MassGridGeometry = GetPropType<MassTypeTag, Properties::GridGeometry>;
+ auto massGridGeometry = std::make_shared<MassGridGeometry>(leafGridView);
+
+ // the coupling manager
+ using CouplingManager = GetPropType<MomentumTypeTag, Properties::CouplingManager>;
+ auto couplingManager = std::make_shared<CouplingManager>();
+
+ // the problem (boundary conditions)
+ using MomentumProblem = GetPropType<MomentumTypeTag, Properties::Problem>;
+ auto momentumProblem = std::make_shared<MomentumProblem>(momentumGridGeometry, couplingManager);
+ using MassProblem = GetPropType<MassTypeTag, Properties::Problem>;
+ auto massProblem = std::make_shared<MassProblem>(massGridGeometry, couplingManager);
- // the problem (initial and boundary conditions)
- using Problem = GetPropType<TypeTag, Properties::Problem>;
- auto problem = std::make_shared<Problem>(gridGeometry);
+ // Initialize the fluid system
+ GetPropType<MassTypeTag, Properties::FluidSystem>::init();
// get some time loop parameters
- using Scalar = GetPropType<TypeTag, Properties::Scalar>;
+ using Scalar = GetPropType<MomentumTypeTag, Properties::Scalar>;
const auto tEnd = getParam<Scalar>("TimeLoop.TEnd");
const auto maxDt = getParam<Scalar>("TimeLoop.MaxTimeStepSize");
auto dt = getParam<Scalar>("TimeLoop.DtInitial");
- // instantiate time loop
- auto timeLoop = std::make_shared<TimeLoop<Scalar>>(0, dt, tEnd);
- timeLoop->setMaxTimeStepSize(maxDt);
+ // check if we are about to restart a previously interrupted simulation
+ Scalar restartTime = getParam<Scalar>("Restart.Time", 0);
// the solution vector
- using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
+ constexpr auto momentumIdx = CouplingManager::freeFlowMomentumIndex;
+ constexpr auto massIdx = CouplingManager::freeFlowMassIndex;
+ using Traits = MultiDomainTraits<MomentumTypeTag, MassTypeTag>;
+ using SolutionVector = typename Traits::SolutionVector;
SolutionVector x;
- x[GridGeometry::cellCenterIdx()].resize(gridGeometry->numCellCenterDofs());
- x[GridGeometry::faceIdx()].resize(gridGeometry->numFaceDofs());
- problem->applyInitialSolution(x);
+ x[momentumIdx].resize(momentumGridGeometry->numDofs());
+ x[massIdx].resize(massGridGeometry->numDofs());
+ momentumProblem->applyInitialSolution(x[momentumIdx]);
+ massProblem->applyInitialSolution(x[massIdx]);
auto xOld = x;
+ // instantiate time loop
+ auto timeLoop = std::make_shared<CheckPointTimeLoop<Scalar>>(restartTime, dt, tEnd);
+ timeLoop->setMaxTimeStepSize(maxDt);
+
+ massProblem->setTimeLoop(timeLoop);
+ momentumProblem->setTimeLoop(timeLoop);
+
// the grid variables
- using GridVariables = GetPropType<TypeTag, Properties::GridVariables>;
- auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry);
- gridVariables->init(x);
+ using MomentumGridVariables = GetPropType<MomentumTypeTag, Properties::GridVariables>;
+ auto momentumGridVariables = std::make_shared<MomentumGridVariables>(momentumProblem, momentumGridGeometry);
- // initialize the vtk output module
- using IOFields = GetPropType<TypeTag, Properties::IOFields>;
- StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name());
- IOFields::initOutputModule(vtkWriter); // Add model specific output fields
- vtkWriter.addField(problem->getDeltaRho(), "deltaRho");
- vtkWriter.write(0.0);
+ using MassGridVariables = GetPropType<MassTypeTag, Properties::GridVariables>;
+ auto massGridVariables = std::make_shared<MassGridVariables>(massProblem, massGridGeometry);
- // the assembler with time loop for instationary problem
- using Assembler = StaggeredFVAssembler<TypeTag, DiffMethod::numeric>;
- auto assembler = std::make_shared<Assembler>(problem, gridGeometry, gridVariables, timeLoop, xOld);
+ couplingManager->init(momentumProblem, massProblem, std::make_tuple(momentumGridVariables, massGridVariables), x, xOld);
+ momentumGridVariables->init(x[momentumIdx]);
+ massGridVariables->init(x[massIdx]);
+ // initialize the vtk output module
+ using IOFields = GetPropType<MassTypeTag, Properties::IOFields>;
+ VtkOutputModule vtkWriter(*massGridVariables, x[massIdx], massProblem->name());
+ IOFields::initOutputModule(vtkWriter); // Add model specific output fields
+ vtkWriter.addVelocityOutput(std::make_shared<NavierStokesVelocityOutput<MassGridVariables>>());
+ vtkWriter.addVolumeVariable([](const auto& volVars){ return volVars.density() - 999.694; }, "deltaRho");
+ vtkWriter.write(0);
+
+ using Assembler = MultiDomainFVAssembler<Traits, CouplingManager, DiffMethod::numeric>;
+ auto assembler = std::make_shared<Assembler>(std::make_tuple(momentumProblem, massProblem),
+ std::make_tuple(momentumGridGeometry, massGridGeometry),
+ std::make_tuple(momentumGridVariables, massGridVariables),
+ couplingManager,
+ timeLoop, xOld);
// the linear solver
using LinearSolver = Dumux::UMFPackIstlSolver<SeqLinearSolverTraits, LinearAlgebraTraitsFromAssembler<Assembler>>;
auto linearSolver = std::make_shared<LinearSolver>();
// the non-linear solver
- using NewtonSolver = Dumux::NewtonSolver<Assembler, LinearSolver>;
- NewtonSolver nonLinearSolver(assembler, linearSolver);
+ using NewtonSolver = MultiDomainNewtonSolver<Assembler, LinearSolver, CouplingManager>;
+ NewtonSolver nonLinearSolver(assembler, linearSolver, couplingManager);
// time loop
timeLoop->start(); do
@@ -122,13 +157,12 @@ int main(int argc, char** argv)
// make the new solution the old solution
xOld = x;
- gridVariables->advanceTimeStep();
+ momentumGridVariables->advanceTimeStep();
+ massGridVariables->advanceTimeStep();
// advance to the time loop to the next step
timeLoop->advanceTimeStep();
- problem->calculateDeltaRho(*gridVariables, x);
-
// write vtk output
vtkWriter.write(timeLoop->time());
diff --git a/test/freeflow/navierstokesnc/densitydrivenflow/problem.hh b/test/freeflow/navierstokesnc/densitydrivenflow/problem.hh
index b900d90f4a..e11b64e639 100644
--- a/test/freeflow/navierstokesnc/densitydrivenflow/problem.hh
+++ b/test/freeflow/navierstokesnc/densitydrivenflow/problem.hh
@@ -7,195 +7,224 @@
/*!
* \file
* \ingroup NavierStokesNCTests
- * \brief Test for the compositional staggered grid (Navier-)Stokes model.
+ * \brief Density driven flow test for the multi-component staggered grid (Navier-)Stokes model.
*/
-#ifndef DUMUX_DENSITY_FLOW_NC_TEST_PROBLEM_HH
-#define DUMUX_DENSITY_FLOW_NC_TEST_PROBLEM_HH
+#ifndef DUMUX_DENSITY_DRIVEN_NC_TEST_PROBLEM_HH
+#define DUMUX_DENSITY_DRIVEN_NC_TEST_PROBLEM_HH
#include <dumux/common/parameters.hh>
#include <dumux/common/properties.hh>
-#include <dumux/common/numeqvector.hh>
#include <dumux/freeflow/navierstokes/boundarytypes.hh>
-#include <dumux/freeflow/navierstokes/staggered/problem.hh>
+
+#include <dumux/freeflow/navierstokes/momentum/fluxhelper.hh>
+#include <dumux/freeflow/navierstokes/scalarfluxhelper.hh>
namespace Dumux {
/*!
* \ingroup NavierStokesNCTests
- * \brief Test problem for the one-phase model.
+ * \brief Test problem for the one-phase (Navier-)Stokes model.
*
- * Here, a quadratic two-dimensional domain with closed and non-moving walls at
- * all sides is considered. Initially, the domain is filled with pure water.
- * At the top, a fixed concentration of the air component is set.
- * The air slowly dissolves in the water which leads to a local increase of density.
- * Due to the influence of gravity and
- * small numerical instabilities, fingers of denser water will form and sink downwards.
+ * Density driven flow test for the multi-component staggered grid (Navier-)Stokes model.
*/
-template <class TypeTag>
-class DensityDrivenFlowProblem : public NavierStokesStaggeredProblem<TypeTag>
+template <class TypeTag, class BaseProblem>
+class DensityDrivenFlowProblem : public BaseProblem
{
- using ParentType = NavierStokesStaggeredProblem<TypeTag>;
-
- using BoundaryTypes = Dumux::NavierStokesBoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
+ using ParentType = BaseProblem;
+ using BoundaryTypes = typename ParentType::BoundaryTypes;
using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
+ using FVElementGeometry = typename GridGeometry::LocalView;
+ using SubControlVolume = typename FVElementGeometry::SubControlVolume;
+ using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
- using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
- using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
+ using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>;
+ using InitialValues = typename ParentType::InitialValues;
+ using Sources = typename ParentType::Sources;
+ using DirichletValues = typename ParentType::DirichletValues;
+ using BoundaryFluxes = typename ParentType::BoundaryFluxes;
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
+ using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
+ static constexpr auto dimWorld = GridGeometry::GridView::dimensionworld;
using Element = typename GridGeometry::GridView::template Codim<0>::Entity;
using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+ using VelocityVector = Dune::FieldVector<Scalar, dimWorld>;
- static constexpr auto dimWorld = GridGeometry::GridView::dimensionworld;
+ using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
- static constexpr auto transportCompIdx = Indices::conti0EqIdx + 1;
- static constexpr auto transportEqIdx = Indices::conti0EqIdx + 1;
+ using TimeLoopPtr = std::shared_ptr<CheckPointTimeLoop<Scalar>>;
+
+ static constexpr auto compIdx = 1;
public:
- DensityDrivenFlowProblem(std::shared_ptr<const GridGeometry> gridGeometry)
- : ParentType(gridGeometry), eps_(1e-6)
+ DensityDrivenFlowProblem(std::shared_ptr<const GridGeometry> gridGeometry,
+ std::shared_ptr<CouplingManager> couplingManager)
+ : ParentType(gridGeometry, couplingManager)
+ , eps_(1e-6)
{
useWholeLength_ = getParam<bool>("Problem.UseWholeLength");
- FluidSystem::init();
- deltaRho_.resize(this->gridGeometry().numCellCenterDofs());
}
- /*!
- * \name Boundary conditions
+ /*!
+ * \brief Returns a reference pressure at a given sub control volume face.
+ * This pressure is subtracted from the actual pressure for the momentum balance
+ * which potentially helps to improve numerical accuracy by avoiding issues related to floating point arithmetic.
*/
- // \{
+ Scalar referencePressure(const Element& element,
+ const FVElementGeometry& fvGeometry,
+ const SubControlVolumeFace& scvf) const
+ { return 1.1e5; }
- /*!
+ /*!
* \brief Specifies which kind of boundary condition should be
* used for which equation on a given boundary control volume.
*
* \param globalPos The position of the center of the finite volume
*/
- BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
+ BoundaryTypes boundaryTypesAtPos(const GlobalPosition& globalPos) const
{
BoundaryTypes values;
- // set Dirichlet values for the velocity everywhere
- values.setDirichlet(Indices::velocityXIdx);
- values.setDirichlet(Indices::velocityYIdx);
- values.setNeumann(Indices::conti0EqIdx);
- values.setNeumann(transportEqIdx);
-
- if(globalPos[1] > this->gridGeometry().bBoxMax()[1] - eps_)
+ if constexpr (ParentType::isMomentumProblem())
+ values.setAllDirichlet();
+ else
{
- if(useWholeLength_)
- values.setDirichlet(transportCompIdx);
- else
- if(globalPos[0] > 0.4 && globalPos[0] < 0.6)
- values.setDirichlet(transportCompIdx);
+ values.setAllNeumann();
+ if (globalPos[1] > this->gridGeometry().bBoxMax()[1] - eps_)
+ {
+ if (useWholeLength_ || (globalPos[0] > 0.4 && globalPos[0] < 0.6))
+ values.setAllDirichlet();
+ }
+
}
return values;
}
/*!
- * \brief Returns whether a fixed Dirichlet value shall be used at a given cell.
+ * \brief Evaluates the boundary conditions for a Dirichlet control volume.
*
- * \param element The finite element
- * \param fvGeometry The finite-volume geometry
- * \param scv The sub control volume
- * \param pvIdx The primary variable index in the solution vector
+ * \param globalPos The center of the finite volume which ought to be set.
*/
- template<class FVElementGeometry, class SubControlVolume>
- bool isDirichletCell(const Element& element,
- const FVElementGeometry& fvGeometry,
- const SubControlVolume& scv,
- int pvIdx) const
+ DirichletValues dirichlet(const Element& element, const SubControlVolumeFace& scvf) const
{
- // set a fixed pressure in one cell
- return (isLowerLeftCell_(scv) && pvIdx == Indices::pressureIdx);
+ const auto& globalPos = scvf.ipGlobal();
+ DirichletValues values(0.0);
+
+ if constexpr (!ParentType::isMomentumProblem())
+ {
+ values[Indices::pressureIdx] = 1.1e+5;
+ if (useWholeLength_ || (globalPos[0] > 0.4 && globalPos[0] < 0.6))
+ values[Indices::conti0EqIdx + compIdx] = 1e-3;
+ }
+
+ return values;
}
- /*!
- * \brief Evaluates the boundary conditions for a Dirichlet control volume.
+ /*!
+ * \brief Evaluates the boundary conditions for a Neumann control volume.
*
- * \param globalPos The center of the finite volume which ought to be set.
+ * \param element The element for which the Neumann boundary condition is set
+ * \param fvGeometry The fvGeometry
+ * \param elemVolVars The element volume variables
+ * \param elemFaceVars The element face variables
+ * \param scvf The boundary sub control volume face
*/
- PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
+ template<class ElementVolumeVariables, class ElementFluxVariablesCache>
+ BoundaryFluxes neumann(const Element& element,
+ const FVElementGeometry& fvGeometry,
+ const ElementVolumeVariables& elemVolVars,
+ const ElementFluxVariablesCache& elemFluxVarsCache,
+ const SubControlVolumeFace& scvf) const
{
- PrimaryVariables values;
+ BoundaryFluxes values(0.0);
- values[Indices::pressureIdx] = 1.1e+5;
- values[transportCompIdx] = 1e-3;
- values[Indices::velocityXIdx] = 0.0;
- values[Indices::velocityYIdx] = 0.0;
+ if constexpr (!ParentType::isMomentumProblem())
+ {
+ const auto insideDensity = elemVolVars[scvf.insideScvIdx()].density();
+
+ // The resulting flux over the boundary is zero anyway (velocity is zero), but this will add some non-zero derivatives to the
+ // Jacobian and makes the BC more general.
+ values[Indices::conti0EqIdx] = this->faceVelocity(element, fvGeometry, scvf) * insideDensity * scvf.unitOuterNormal();
+ }
return values;
}
- // \}
-
- /*!
- * \name Volume terms
- */
- // \{
-
- /*!
+ /*!
* \brief Evaluates the initial value for a control volume.
*
* \param globalPos The global position
*/
- PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
+ InitialValues initialAtPos(const GlobalPosition& globalPos) const
{
- PrimaryVariables values;
- values[Indices::pressureIdx] = 1.1e+5;
- values[transportCompIdx] = 0.0;
- values[Indices::velocityXIdx] = 0.0;
- values[Indices::velocityYIdx] = 0.0;
+ InitialValues values(0.0);
+
+ if constexpr (!ParentType::isMomentumProblem())
+ values[Indices::pressureIdx] = 1.1e5;
return values;
}
- /*!
- * \brief Adds additional VTK output data to the VTKWriter.
- *
- * Function is called by the output module on every write.
+ void setTimeLoop(TimeLoopPtr timeLoop)
+ { timeLoop_ = timeLoop; }
+
+ Scalar time() const
+ { return timeLoop_->time(); }
+
+ //! Enable internal Dirichlet constraints
+ static constexpr bool enableInternalDirichletConstraints()
+ { return !ParentType::isMomentumProblem(); }
+
+ /*!
+ * \brief Tag a degree of freedom to carry internal Dirichlet constraints.
+ * If true is returned for a dof, the equation for this dof is replaced
+ * by the constraint that its primary variable values must match the
+ * user-defined values obtained from the function internalDirichlet(),
+ * which must be defined in the problem.
*
- * \param gridVariables The grid variables
- * \param sol The solution vector
+ * \param element The finite element
+ * \param scv The sub-control volume
*/
- template<class GridVariables, class SolutionVector>
- void calculateDeltaRho(const GridVariables& gridVariables, const SolutionVector& sol)
+ std::bitset<DirichletValues::dimension> hasInternalDirichletConstraint(const Element& element, const SubControlVolume& scv) const
{
- auto fvGeometry = localView(this->gridGeometry());
- auto elemVolVars = localView(gridVariables.curGridVolVars());
- for (const auto& element : elements(this->gridGeometry().gridView()))
- {
- fvGeometry.bindElement(element);
- elemVolVars.bind(element, fvGeometry, sol);
- for (auto&& scv : scvs(fvGeometry))
- {
- const auto ccDofIdx = scv.dofIndex();
- deltaRho_[ccDofIdx] = elemVolVars[scv].density() - 999.694;
- }
- }
- }
+ std::bitset<DirichletValues::dimension> values;
- const std::vector<Scalar>& getDeltaRho() const
- { return deltaRho_; }
+ // the pure Neumann problem is only defined up to a constant
+ // we create a well-posed problem by fixing the pressure at one dof
+ if constexpr (!ParentType::isMomentumProblem())
+ {
+ const bool isLowerLeftCell = (scv.dofIndex() == 0);
+ if (isLowerLeftCell)
+ values.set(0);
+ }
+ return values;
+ }
- // \}
+ /*!
+ * \brief Define the values of internal Dirichlet constraints for a degree of freedom.
+ * \param element The finite element
+ * \param scv The sub-control volume
+ */
+ DirichletValues internalDirichlet(const Element& element, const SubControlVolume& scv) const
+ { return DirichletValues(1.1e5); }
private:
+ bool isInlet_(const GlobalPosition& globalPos) const
+ { return globalPos[0] < eps_; }
- template<class SubControlVolume>
- bool isLowerLeftCell_(const SubControlVolume& scv) const
- { return scv.dofIndex() == 0; }
+ bool isOutlet_(const GlobalPosition& globalPos) const
+ { return globalPos[0] > this->gridGeometry().bBoxMax()[0] - eps_; }
const Scalar eps_;
bool useWholeLength_;
- std::vector<Scalar> deltaRho_;
+ TimeLoopPtr timeLoop_;
};
-} // end namespace Dumux
+} // end namespace Dumux
#endif
diff --git a/test/freeflow/navierstokesnc/densitydrivenflow/properties.hh b/test/freeflow/navierstokesnc/densitydrivenflow/properties.hh
index bb9cd4eacb..b012200401 100644
--- a/test/freeflow/navierstokesnc/densitydrivenflow/properties.hh
+++ b/test/freeflow/navierstokesnc/densitydrivenflow/properties.hh
@@ -14,25 +14,40 @@
#include <dune/grid/yaspgrid.hh>
-#include <dumux/discretization/staggered/freeflow/properties.hh>
+#include <dumux/discretization/fcstaggered.hh>
+#include <dumux/discretization/cctpfa.hh>
+
+#include <dumux/freeflow/navierstokes/mass/1pnc/model.hh>
+#include <dumux/freeflow/navierstokes/mass/problem.hh>
+
+#include <dumux/freeflow/navierstokes/momentum/problem.hh>
+#include <dumux/freeflow/navierstokes/momentum/model.hh>
-#include <dumux/freeflow/compositional/navierstokesncmodel.hh>
#include <dumux/material/components/simpleh2o.hh>
#include <dumux/material/fluidsystems/1padapter.hh>
#include <dumux/material/fluidsystems/h2oair.hh>
+#include <dumux/multidomain/traits.hh>
+#include <dumux/multidomain/freeflow/couplingmanager.hh>
+
#include "problem.hh"
namespace Dumux::Properties {
// Create new type tags
namespace TTag {
-struct DensityDrivenFlow { using InheritsFrom = std::tuple<NavierStokesNC, StaggeredFreeFlowModel>; };
+struct DensityDrivenFlowTest {};
+struct DensityDrivenFlowMomentum { using InheritsFrom = std::tuple<DensityDrivenFlowTest, NavierStokesMomentum, FaceCenteredStaggeredModel>; };
+#if !NONISOTHERMAL
+struct DensityDrivenFlowMass { using InheritsFrom = std::tuple<DensityDrivenFlowTest, NavierStokesMassOnePNC, CCTpfaModel>; };
+#else
+struct DensityDrivenFlowMass { using InheritsFrom = std::tuple<DensityDrivenFlowTest, NavierStokesMassOnePNCNI, CCTpfaModel>; };
+#endif
} // end namespace TTag
// Select the fluid system
template<class TypeTag>
-struct FluidSystem<TypeTag, TTag::DensityDrivenFlow>
+struct FluidSystem<TypeTag, TTag::DensityDrivenFlowTest>
{
using H2OAir = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>;
static constexpr int phaseIdx = H2OAir::liquidPhaseIdx;
@@ -40,25 +55,40 @@ struct FluidSystem<TypeTag, TTag::DensityDrivenFlow>
};
template<class TypeTag>
-struct ReplaceCompEqIdx<TypeTag, TTag::DensityDrivenFlow> { static constexpr int value = 0; };
+struct ReplaceCompEqIdx<TypeTag, TTag::DensityDrivenFlowTest> { static constexpr int value = 0; };
// Set the grid type
template<class TypeTag>
-struct Grid<TypeTag, TTag::DensityDrivenFlow> { using type = Dune::YaspGrid<2>; };
+struct Grid<TypeTag, TTag::DensityDrivenFlowTest> { using type = Dune::YaspGrid<2>; };
// Set the problem property
template<class TypeTag>
-struct Problem<TypeTag, TTag::DensityDrivenFlow> { using type = Dumux::DensityDrivenFlowProblem<TypeTag> ; };
+struct Problem<TypeTag, TTag::DensityDrivenFlowMomentum>
+{ using type = DensityDrivenFlowProblem<TypeTag, Dumux::NavierStokesMomentumProblem<TypeTag>>; };
+
+template<class TypeTag>
+struct Problem<TypeTag, TTag::DensityDrivenFlowMass>
+{ using type = DensityDrivenFlowProblem<TypeTag, Dumux::NavierStokesMassProblem<TypeTag>>; };
template<class TypeTag>
-struct EnableGridGeometryCache<TypeTag, TTag::DensityDrivenFlow> { static constexpr bool value = true; };
+struct EnableGridGeometryCache<TypeTag, TTag::DensityDrivenFlowTest> { static constexpr bool value = true; };
+template<class TypeTag>
+struct EnableGridFluxVariablesCache<TypeTag, TTag::DensityDrivenFlowTest> { static constexpr bool value = true; };
template<class TypeTag>
-struct EnableGridFluxVariablesCache<TypeTag, TTag::DensityDrivenFlow> { static constexpr bool value = true; };
+struct EnableGridVolumeVariablesCache<TypeTag, TTag::DensityDrivenFlowTest> { static constexpr bool value = true; };
+
template<class TypeTag>
-struct EnableGridVolumeVariablesCache<TypeTag, TTag::DensityDrivenFlow> { static constexpr bool value = true; };
+struct UseMoles<TypeTag, TTag::DensityDrivenFlowTest> { static constexpr bool value = true; };
+// Set the coupling manager
template<class TypeTag>
-struct UseMoles<TypeTag, TTag::DensityDrivenFlow> { static constexpr bool value = true; };
+struct CouplingManager<TypeTag, TTag::DensityDrivenFlowTest>
+{
+private:
+ using Traits = MultiDomainTraits<TTag::DensityDrivenFlowMomentum, TTag::DensityDrivenFlowMass>;
+public:
+ using type = FreeFlowCouplingManager<Traits>;
+};
} // end namespace Dumux::Properties
--
GitLab