Skip to content

GitLab

Commit 225f3f71 authored by Timo Koch's avatar Timo Koch
Browse files

Merge branch 'cleanup/poromechanics' into 'master' 

Cleanup/poromechanics

See merge request !1675
parents 97ea5b46 13f68da4
Hide whitespace changes
Inline Side-by-side
dumux/geomechanics/poroelastic/couplingmanager.hh
View file @ 225f3f71
......@@ -66,6 +66,7 @@ class PoroMechanicsCouplingManager : public virtual CouplingManager< MDTraits >
template<std::size_t id> using PrimaryVariables = typename GridVariables<id>::PrimaryVariables;
template<std::size_t id> using GridVolumeVariables = typename GridVariables<id>::GridVolumeVariables;
template<std::size_t id> using ElementVolumeVariables = typename GridVolumeVariables<id>::LocalView;
template<std::size_t id> using VolumeVariables = typename GridVolumeVariables<id>::VolumeVariables;
template<std::size_t id> using FVGridGeometry = typename GridVariables<id>::GridGeometry;
template<std::size_t id> using FVElementGeometry = typename FVGridGeometry<id>::LocalView;
template<std::size_t id> using GridView = typename FVGridGeometry<id>::GridView;
......@@ -138,8 +139,10 @@ public:
std::shared_ptr< Problem<PoroMechId> > poroMechanicalProblem,
const SolutionVector& curSol)
{
// set up tuple containing the sub-problems
problemTuple_ = std::make_tuple(pmFlowProblem, poroMechanicalProblem);
// set the sub problems
this->setSubProblem(pmFlowProblem, pmFlowId);
this->setSubProblem(poroMechanicalProblem, poroMechId);
// copy the solution vector
ParentType::updateSolution(curSol);
// set up the coupling map pmfow -> poromechanics
......@@ -153,7 +156,7 @@ public:
const Element<PMFlowId>& element,
Dune::index_constant<PoroMechId> poroMechDomainId) const
{
return pmFlowCouplingMap_[ problem<PMFlowId>().fvGridGeometry().elementMapper().index(element) ];
return pmFlowCouplingMap_[ this->problem(pmFlowId).fvGridGeometry().elementMapper().index(element) ];
}
/*!
......@@ -163,7 +166,7 @@ public:
const Element<PoroMechId>& element,
Dune::index_constant<PMFlowId> pmFlowDomainId) const
{
const auto eIdx = problem<PMFlowId>().fvGridGeometry().elementMapper().index(element);
const auto eIdx = this->problem(pmFlowId).fvGridGeometry().elementMapper().index(element);
return CouplingStencilType<PoroMechId>{ {eIdx} };
}
......@@ -185,11 +188,12 @@ public:
// prepare the fvGeometry and the element volume variables
// these quantities will be used later to obtain the effective pressure
auto fvGeometry = localView( problem<PMFlowId>().fvGridGeometry() );
auto fvGeometry = localView( this->problem(pmFlowId).fvGridGeometry() );
auto elemVolVars = localView( assembler.gridVariables(Dune::index_constant<PMFlowId>()).curGridVolVars() );
fvGeometry.bindElement(element);
elemVolVars.bindElement(element, fvGeometry, this->curSol()[Dune::index_constant<PMFlowId>()]);
poroMechCouplingContext_.pmFlowElement = element;
poroMechCouplingContext_.pmFlowFvGeometry = std::make_unique< FVElementGeometry<PMFlowId> >(fvGeometry);
poroMechCouplingContext_.pmFlowElemVolVars = std::make_unique< ElementVolumeVariables<PMFlowId> >(elemVolVars);
......@@ -346,20 +350,20 @@ public:
poroMechLocalAssembler.elemBcTypes());
}
//! Return a const reference to one of the problems
template<std::size_t id, std::enable_if_t<(id == PMFlowId || id == PoroMechId), int> = 0>
const Problem<id>& problem() const
{ return *std::get<(id == PMFlowId ? 0 : 1)>(problemTuple_); }
//! Return reference to one of the problems
template<std::size_t id, std::enable_if_t<(id == PMFlowId || id == PoroMechId), int> = 0>
Problem<id>& problem()
{ return *std::get<(id == PMFlowId ? 0 : 1)>(problemTuple_); }
//! Return the coupling context (used in mechanical sub-problem to compute effective pressure)
[[deprecated("Obtain the volume variables directly calling getPMFlowVolVars(element). Will be removed after 3.1!")]]
const PoroMechanicsCouplingContext& poroMechanicsCouplingContext() const
{ return poroMechCouplingContext_; }
//! Return the porous medium flow variables an element/scv of the poromech domain couples to
const VolumeVariables<PMFlowId>& getPMFlowVolVars(const Element<PoroMechId>& element) const
{
//! If we do not yet have the queried object, build it first
const auto eIdx = this->problem(poroMechId).fvGridGeometry().elementMapper().index(element);
return (*poroMechCouplingContext_.pmFlowElemVolVars)[eIdx];
}
/*!
* \brief the solution vector of the coupled problem
* \note in case of numeric differentiation the solution vector always carries the deflected solution
......@@ -376,8 +380,8 @@ private:
void initializeCouplingMap_()
{
// some references for convenience
const auto& pmFlowGridGeom = problem<PMFlowId>().fvGridGeometry();
const auto& poroMechGridGeom = problem<PoroMechId>().fvGridGeometry();
const auto& pmFlowGridGeom = this->problem(pmFlowId).fvGridGeometry();
const auto& poroMechGridGeom = this->problem(poroMechId).fvGridGeometry();
// make sure the two grids are really the same. Note that if the two grids
// happen to have equal number of elements by chance, we don't detect this source of error.
......@@ -411,11 +415,6 @@ private:
}
}
// tuple for storing pointers to the sub-problems
using PMFlowProblemPtr = std::shared_ptr< Problem<PMFlowId> >;
using PoroMechanicalProblemPtr = std::shared_ptr< Problem<PoroMechId> >;
std::tuple<PMFlowProblemPtr, PoroMechanicalProblemPtr> problemTuple_;
// Container for storing the coupling element stencils for the pm flow domain
std::vector< CouplingStencilType<PMFlowId> > pmFlowCouplingMap_;
......
test/multidomain/poromechanics/el1p/problem_poroelastic.hh
View file @ 225f3f71
......@@ -138,11 +138,9 @@ public:
Scalar effectiveFluidDensity(const Element& element,
const SubControlVolume& scv) const
{
// get context from coupling manager
// here, we know that the flow problem uses cell-centered finite volumes,
// thus, we simply take the volume variables of the element and return the density
const auto& context = couplingManager().poroMechanicsCouplingContext();
return (*context.pmFlowElemVolVars)[scv.elementIndex()].density();
// get porous medium flow volume variables from coupling manager
const auto pmFlowVolVars = couplingManager().getPMFlowVolVars(element);
return pmFlowVolVars.density();
}
/*!
......@@ -154,12 +152,9 @@ public:
const ElementVolumeVariables& elemVolVars,
const FluxVarsCache& fluxVarsCache) const
{
// get context from coupling manager
// here, we know that the flow problem uses cell-centered finite volumes,
// thus, we simply take the volume variables of the element and return the pressure
const auto& context = couplingManager().poroMechanicsCouplingContext();
const auto eIdx = this->fvGridGeometry().elementMapper().index(element);
return (*context.pmFlowElemVolVars)[eIdx].pressure();
// get porous medium flow volume variables from coupling manager
const auto pmFlowVolVars = couplingManager().getPMFlowVolVars(element);
return pmFlowVolVars.pressure();
}
/*!
......
test/multidomain/poromechanics/el1p/spatialparams_1p.hh
View file @ 225f3f71
......@@ -75,7 +75,7 @@ public:
{
static constexpr auto poroMechId = CouplingManager::poroMechId;
const auto& poroMechGridGeom = couplingManagerPtr_->template problem<poroMechId>().fvGridGeometry();
const auto& poroMechGridGeom = couplingManagerPtr_->problem(poroMechId).fvGridGeometry();
const auto poroMechElemSol = elementSolution(element, couplingManagerPtr_->curSol()[poroMechId], poroMechGridGeom);
// evaluate the deformation-dependent porosity at the scv center
......
test/multidomain/poromechanics/el2p/problem_poroelastic.hh
View file @ 225f3f71
......@@ -139,16 +139,13 @@ public:
*/
Scalar effectiveFluidDensity(const Element& element, const SubControlVolume& scv) const
{
// get context from coupling manager
const auto& context = couplingManager().poroMechanicsCouplingContext();
// here, we know that the flow problem uses cell-centered finite volumes, thus,
// we simply take the volume variables of the scv (i.e. element) to obtain fluid properties
const auto& facetVolVars = (*context.pmFlowElemVolVars)[scv.elementIndex()];
Scalar wPhaseDensity = facetVolVars.density(FluidSystem::phase0Idx);
Scalar nPhaseDensity = facetVolVars.density(FluidSystem::phase1Idx);
Scalar Sw = facetVolVars.saturation(FluidSystem::phase0Idx);
Scalar Sn = facetVolVars.saturation(FluidSystem::phase1Idx);
// get porous medium flow volume variables from coupling manager
const auto pmFlowVolVars = couplingManager().getPMFlowVolVars(element);
Scalar wPhaseDensity = pmFlowVolVars.density(FluidSystem::phase0Idx);
Scalar nPhaseDensity = pmFlowVolVars.density(FluidSystem::phase1Idx);
Scalar Sw = pmFlowVolVars.saturation(FluidSystem::phase0Idx);
Scalar Sn = pmFlowVolVars.saturation(FluidSystem::phase1Idx);
return (wPhaseDensity * Sw + nPhaseDensity * Sn);
}
......@@ -161,17 +158,13 @@ public:
const ElementVolumeVariables& elemVolVars,
const FluxVarsCache& fluxVarsCache) const
{
// get context from coupling manager
const auto& context = couplingManager().poroMechanicsCouplingContext();
// here, we know that the flow problem uses cell-centered finite volumes, thus,
// we simply take the volume variables of the element to obtain fluid properties
const auto eIdx = this->fvGridGeometry().elementMapper().index(element);
const auto& facetVolVars = (*context.pmFlowElemVolVars)[eIdx];
Scalar pw = facetVolVars.pressure(FluidSystem::phase0Idx);
Scalar pn = facetVolVars.pressure(FluidSystem::phase1Idx);
Scalar Sw = facetVolVars.saturation(FluidSystem::phase0Idx);
Scalar Sn = facetVolVars.saturation(FluidSystem::phase1Idx);
// get porous medium flow volume variables from coupling manager
const auto pmFlowVolVars = couplingManager().getPMFlowVolVars(element);
Scalar pw = pmFlowVolVars.pressure(FluidSystem::phase0Idx);
Scalar pn = pmFlowVolVars.pressure(FluidSystem::phase1Idx);
Scalar Sw = pmFlowVolVars.saturation(FluidSystem::phase0Idx);
Scalar Sn = pmFlowVolVars.saturation(FluidSystem::phase1Idx);
return (pw * Sw + pn * Sn);
}
......
test/multidomain/poromechanics/el2p/spatialparams_2p.hh
View file @ 225f3f71
......@@ -89,7 +89,7 @@ public:
{
static constexpr auto poroMechId = CouplingManager::poroMechId;
const auto& poroMechGridGeom = couplingManagerPtr_->template problem<poroMechId>().fvGridGeometry();
const auto& poroMechGridGeom = couplingManagerPtr_->problem(poroMechId).fvGridGeometry();
const auto poroMechElemSol = elementSolution(element, couplingManagerPtr_->curSol()[poroMechId], poroMechGridGeom);
// evaluate the deformation-dependent porosity at the scv center
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment