diff --git a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc index 0aa3d1647a7b7db4674ff8072453fdf14b67b2e4..9dfd0d76d2c3be725e8a12a171a8711481f5f013 100644 --- a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc +++ b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc @@ -55,7 +55,7 @@ namespace Dumux { namespace Properties { -SET_PROP(StokesOnePTwoCTypeTag, CouplingManager) +SET_PROP(ZeroEqTypeTag, CouplingManager) { using Traits = StaggeredMultiDomainTraits; using type = Dumux::StokesDarcyCouplingManager; @@ -63,7 +63,7 @@ SET_PROP(StokesOnePTwoCTypeTag, CouplingManager) SET_PROP(DarcyTwoPTwoCTypeTag, CouplingManager) { - using Traits = StaggeredMultiDomainTraits; + using Traits = StaggeredMultiDomainTraits; using type = Dumux::StokesDarcyCouplingManager; }; @@ -85,7 +85,7 @@ int main(int argc, char** argv) try Parameters::init(argc, argv); // Define the sub problem type tags - using StokesTypeTag = TTAG(StokesOnePTwoCTypeTag); + using StokesTypeTag = TTAG(ZeroEqTypeTag); using DarcyTypeTag = TTAG(DarcyTwoPTwoCTypeTag); // try to create a grid (from the given grid file or the input file) @@ -160,6 +160,8 @@ int main(int argc, char** argv) try std::get<0>(stokesSol) = sol[stokesCellCenterIdx]; std::get<1>(stokesSol) = sol[stokesFaceIdx]; stokesProblem->applyInitialSolution(stokesSol); + // TODO: update static wall properties + // TODO: update dynamic wall properties auto solStokesOld = stokesSol; sol[stokesCellCenterIdx] = stokesSol[stokesCellCenterIdx]; sol[stokesFaceIdx] = stokesSol[stokesFaceIdx]; @@ -226,7 +228,7 @@ int main(int argc, char** argv) try // make the new solution the old solution solOld = sol; - // update dynamic wall properties + // TODO: update dynamic wall properties // post time step treatment of Darcy problem darcyProblem->postTimeStep(sol[darcyIdx], *darcyGridVariables, timeLoop->timeStepSize()); diff --git a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input index 515dfd2c988e234ff4ef277ca2fef5ba3be0e69d..514cf156d55211777c2ec0980f37d1a95124867e 100644 --- a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input +++ b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input @@ -23,7 +23,7 @@ Verbosity = true [Stokes.Problem] Name = stokes -RefVelocity = 0.1 #3.5 # [m/s] +RefVelocity = 3.5 # [m/s] RefPressure = 1e5 # [Pa] refMoleFrac = 0 # [-] RefTemperature = 298.15 # [K] diff --git a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh index a5ec83fb49e6a3dfca0e7ba0df091aa5f60f8e65..9260ad6f463b526b7540d576c59fd8aca0d1484b 100644 --- a/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh +++ b/exercises/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh @@ -29,8 +29,8 @@ #include #include -#include #include +#include namespace Dumux { @@ -39,30 +39,30 @@ class StokesSubProblem; namespace Properties { -NEW_TYPE_TAG(StokesOnePTwoCTypeTag, INHERITS_FROM(StaggeredFreeFlowModel, NavierStokesNCNI)); +NEW_TYPE_TAG(ZeroEqTypeTag, INHERITS_FROM(StaggeredFreeFlowModel, NavierStokesNCNI)); // Set the grid type -SET_TYPE_PROP(StokesOnePTwoCTypeTag, Grid, Dune::YaspGrid<2, Dune::TensorProductCoordinates >); +SET_TYPE_PROP(ZeroEqTypeTag, Grid, Dune::YaspGrid<2, Dune::TensorProductCoordinates >); // set the fluid system -SET_TYPE_PROP(StokesOnePTwoCTypeTag, FluidSystem, FluidSystems::H2OAir); +SET_TYPE_PROP(ZeroEqTypeTag, FluidSystem, FluidSystems::H2OAir); // set phase index (air) -SET_INT_PROP(StokesOnePTwoCTypeTag, PhaseIdx, GET_PROP_TYPE(TypeTag, FluidSystem)::gasPhaseIdx); +SET_INT_PROP(ZeroEqTypeTag, PhaseIdx, GET_PROP_TYPE(TypeTag, FluidSystem)::gasPhaseIdx); -SET_INT_PROP(StokesOnePTwoCTypeTag, ReplaceCompEqIdx, 3); +SET_INT_PROP(ZeroEqTypeTag, ReplaceCompEqIdx, 3); // Use formulation based on mass fractions -SET_BOOL_PROP(StokesOnePTwoCTypeTag, UseMoles, true); +SET_BOOL_PROP(ZeroEqTypeTag, UseMoles, true); // Set the problem property -SET_TYPE_PROP(StokesOnePTwoCTypeTag, Problem, Dumux::StokesSubProblem ); +SET_TYPE_PROP(ZeroEqTypeTag, Problem, Dumux::StokesSubProblem ); -SET_BOOL_PROP(StokesOnePTwoCTypeTag, EnableFVGridGeometryCache, true); -SET_BOOL_PROP(StokesOnePTwoCTypeTag, EnableGridFluxVariablesCache, true); -SET_BOOL_PROP(StokesOnePTwoCTypeTag, EnableGridVolumeVariablesCache, true); +SET_BOOL_PROP(ZeroEqTypeTag, EnableFVGridGeometryCache, true); +SET_BOOL_PROP(ZeroEqTypeTag, EnableGridFluxVariablesCache, true); +SET_BOOL_PROP(ZeroEqTypeTag, EnableGridVolumeVariablesCache, true); -SET_BOOL_PROP(StokesOnePTwoCTypeTag, EnableInertiaTerms, true); +SET_BOOL_PROP(ZeroEqTypeTag, EnableInertiaTerms, true); } /*! @@ -161,22 +161,30 @@ public: const auto& globalPos = scvf.center(); - values.setNeumann(Indices::energyBalanceIdx); + if (onLeftBoundary_(globalPos)) + { + values.setDirichlet(Indices::velocityXIdx); + values.setDirichlet(Indices::velocityYIdx); + values.setDirichlet(Indices::conti0EqIdx); + values.setDirichlet(Indices::energyBalanceIdx); + } - if (onLowerBoundary_(globalPos) || onLeftBoundary_(globalPos)) + if (onLowerBoundary_(globalPos)) { values.setDirichlet(Indices::velocityXIdx); values.setDirichlet(Indices::velocityYIdx); values.setNeumann(Indices::conti0EqIdx); values.setNeumann(Indices::conti0EqIdx + 1); + values.setNeumann(Indices::energyBalanceIdx); } if (onUpperBoundary_(globalPos)) { - values.setDirichlet(Indices::velocityXIdx); - values.setDirichlet(Indices::velocityYIdx); - values.setNeumann(Indices::conti0EqIdx); - values.setNeumann(Indices::conti0EqIdx + 1); + values.setDirichlet(Indices::velocityXIdx); + values.setDirichlet(Indices::velocityYIdx); + values.setNeumann(Indices::conti0EqIdx); + values.setNeumann(Indices::conti0EqIdx + 1); + values.setNeumann(Indices::energyBalanceIdx); } if (onRightBoundary_(globalPos)) @@ -227,25 +235,6 @@ public: const SubControlVolumeFace& scvf) const { PrimaryVariables values(0.0); - const auto& globalPos = scvf.dofPosition(); - const auto& scv = fvGeometry.scv(scvf.insideScvIdx()); - - FluidState fluidState; - updateFluidStateForBC_(fluidState, elemVolVars[scv].temperature(), - elemVolVars[scv].pressure(), - elemVolVars[scv].moleFraction(transportCompIdx)); - - const Scalar density = useMoles ? fluidState.molarDensity(phaseIdx) : fluidState.density(phaseIdx); - const Scalar xVelocity = refVelocity(); - - if (onLeftBoundary_(globalPos)) - { - // rho*v*X at inflow - values[Indices::conti0EqIdx] = -xVelocity * density * refMoleFrac(); - values[Indices::conti0EqIdx + 1] = -xVelocity * density * (1.0 - refMoleFrac()); - values[Indices::energyBalanceIdx] = -xVelocity * fluidState.density(phaseIdx) * fluidState.enthalpy(phaseIdx); - } - if(couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf)) { values[Indices::momentumYBalanceIdx] = couplingManager().couplingData().momentumCouplingCondition(fvGeometry, elemVolVars, elemFaceVars, scvf); diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/CMakeLists.txt b/exercises/solution/exercise-coupling-ff-pm/turbulence/CMakeLists.txt index 30a98ad93c4e9a82d47e832c5821c0cb18ea251d..a4a6689dbca1c74f886e3915ce68875a09f7ceee 100644 --- a/exercises/solution/exercise-coupling-ff-pm/turbulence/CMakeLists.txt +++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/CMakeLists.txt @@ -1,5 +1,10 @@ add_input_file_links() +dune_add_test(NAME orig_turbulence_coupling_ff-pm + SOURCES ex_turbulence_coupling_ff-pm.cc + CMD_ARGS ex_turbulence_coupling_ff-pm.input) +target_compile_definitions(orig_turbulence_coupling_ff-pm PUBLIC "EXNUMBER=0") + dune_add_test(NAME sol1_turbulence_coupling_ff-pm SOURCES ex_turbulence_coupling_ff-pm.cc CMD_ARGS ex_turbulence_coupling_ff-pm.input) diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc index 51a4b8c2a8ae8af29d28f31f1aed89995511879e..5a1d6eb07fde7e819f881e81c1d67cda325039f0 100644 --- a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc +++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.cc @@ -160,6 +160,8 @@ int main(int argc, char** argv) try std::get<0>(stokesSol) = sol[stokesCellCenterIdx]; std::get<1>(stokesSol) = sol[stokesFaceIdx]; stokesProblem->applyInitialSolution(stokesSol); + // TODO: update static wall properties + // TODO: update dynamic wall properties #if EXNUMBER >= 1 stokesProblem->updateStaticWallProperties(); stokesProblem->updateDynamicWallProperties(stokesSol); @@ -230,7 +232,7 @@ int main(int argc, char** argv) try // make the new solution the old solution solOld = sol; - // update dynamic wall properties + // TODO: update dynamic wall properties #if EXNUMBER >= 1 std::get<0>(stokesSol) = sol[stokesCellCenterIdx]; std::get<1>(stokesSol) = sol[stokesFaceIdx]; diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input index 515dfd2c988e234ff4ef277ca2fef5ba3be0e69d..514cf156d55211777c2ec0980f37d1a95124867e 100644 --- a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input +++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_coupling_ff-pm.input @@ -23,7 +23,7 @@ Verbosity = true [Stokes.Problem] Name = stokes -RefVelocity = 0.1 #3.5 # [m/s] +RefVelocity = 3.5 # [m/s] RefPressure = 1e5 # [Pa] refMoleFrac = 0 # [-] RefTemperature = 298.15 # [K] diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh index 3a912917a614fd0ba2109d296468e016978974b9..435cd8ae357ac15a163cc4ef6872f7f57e1a014a 100644 --- a/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh +++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/ex_turbulence_ffproblem.hh @@ -29,9 +29,13 @@ #include #include -#include +#if EXNUMBER >= 1 #include #include +#else +#include +#include +#endif namespace Dumux { @@ -40,7 +44,11 @@ class StokesSubProblem; namespace Properties { +#if EXNUMBER >= 1 NEW_TYPE_TAG(ZeroEqTypeTag, INHERITS_FROM(StaggeredFreeFlowModel, ZeroEqNCNI)); +#else +NEW_TYPE_TAG(ZeroEqTypeTag, INHERITS_FROM(StaggeredFreeFlowModel, NavierStokesNCNI)); +#endif // Set the grid type SET_TYPE_PROP(ZeroEqTypeTag, Grid, Dune::YaspGrid<2, Dune::TensorProductCoordinates >); @@ -73,9 +81,15 @@ SET_BOOL_PROP(ZeroEqTypeTag, EnableInertiaTerms, true); * Horizontal flow from left to right with a parabolic velocity profile. */ template +#if EXNUMBER >= 1 class StokesSubProblem : public ZeroEqProblem { using ParentType = ZeroEqProblem; +#else +class StokesSubProblem : public NavierStokesProblem +{ + using ParentType = NavierStokesProblem; +#endif using GridView = typename GET_PROP_TYPE(TypeTag, GridView); using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar); @@ -162,24 +176,33 @@ public: const auto& globalPos = scvf.center(); - values.setNeumann(Indices::energyBalanceIdx); + if (onLeftBoundary_(globalPos)) + { + values.setDirichlet(Indices::velocityXIdx); + values.setDirichlet(Indices::velocityYIdx); + values.setDirichlet(Indices::conti0EqIdx); + values.setDirichlet(Indices::energyBalanceIdx); + } - if (onLowerBoundary_(globalPos) || onLeftBoundary_(globalPos)) + if (onLowerBoundary_(globalPos)) { values.setDirichlet(Indices::velocityXIdx); values.setDirichlet(Indices::velocityYIdx); values.setNeumann(Indices::conti0EqIdx); values.setNeumann(Indices::conti0EqIdx + 1); + values.setNeumann(Indices::energyBalanceIdx); } if (onUpperBoundary_(globalPos)) { - values.setDirichlet(Indices::velocityXIdx); - values.setDirichlet(Indices::velocityYIdx); - values.setNeumann(Indices::conti0EqIdx); - values.setNeumann(Indices::conti0EqIdx + 1); #if EXNUMBER >=2 - values.setAllSymmetry(); + values.setAllSymmetry(); +#else + values.setDirichlet(Indices::velocityXIdx); + values.setDirichlet(Indices::velocityYIdx); + values.setNeumann(Indices::conti0EqIdx); + values.setNeumann(Indices::conti0EqIdx + 1); + values.setNeumann(Indices::energyBalanceIdx); #endif } @@ -231,25 +254,6 @@ public: const SubControlVolumeFace& scvf) const { PrimaryVariables values(0.0); - const auto& globalPos = scvf.dofPosition(); - const auto& scv = fvGeometry.scv(scvf.insideScvIdx()); - - FluidState fluidState; - updateFluidStateForBC_(fluidState, elemVolVars[scv].temperature(), - elemVolVars[scv].pressure(), - elemVolVars[scv].moleFraction(transportCompIdx)); - - const Scalar density = useMoles ? fluidState.molarDensity(phaseIdx) : fluidState.density(phaseIdx); - const Scalar xVelocity = refVelocity(); - - if (onLeftBoundary_(globalPos)) - { - // rho*v*X at inflow - values[Indices::conti0EqIdx] = -xVelocity * density * refMoleFrac(); - values[Indices::conti0EqIdx + 1] = -xVelocity * density * (1.0 - refMoleFrac()); - values[Indices::energyBalanceIdx] = -xVelocity * fluidState.density(phaseIdx) * fluidState.enthalpy(phaseIdx); - } - if(couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf)) { values[Indices::momentumYBalanceIdx] = couplingManager().couplingData().momentumCouplingCondition(fvGeometry, elemVolVars, elemFaceVars, scvf); @@ -276,11 +280,16 @@ public: const CouplingManager& couplingManager() const { return *couplingManager_; } -#if EXNUMBER >= 1 +#if EXNUMBER >= 2 bool isOnWall(const GlobalPosition& globalPos) const { return (onLowerBoundary_(globalPos)); } +#elif EXNUMBER >= 1 + bool isOnWall(const GlobalPosition& globalPos) const + { + return (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos)); + } #endif /*!