From 4d3b52f0766a869b9bbfe002c19d2d8139adc24a Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt <kilian.weishaupt@iws.uni-stuttgart.de> Date: Fri, 24 Aug 2018 11:42:18 +0200 Subject: [PATCH] [2p1cni][volVars] Use function to get temperature * allows to use thermal non-equilibrium model (T_soil != T_fluid) --- .../porousmediumflow/2p1c/volumevariables.hh | 57 +++++++++++++------ 1 file changed, 41 insertions(+), 16 deletions(-) diff --git a/dumux/porousmediumflow/2p1c/volumevariables.hh b/dumux/porousmediumflow/2p1c/volumevariables.hh index 4f8cff6c9e..660029e45f 100644 --- a/dumux/porousmediumflow/2p1c/volumevariables.hh +++ b/dumux/porousmediumflow/2p1c/volumevariables.hh @@ -212,22 +212,8 @@ public: : priVars[pressureIdx] + pc_); } - // get temperature - Scalar temperature; - if (phasePresence == liquidPhaseOnly || phasePresence == gasPhaseOnly) - temperature = priVars[switchIdx]; - else if (phasePresence == twoPhases) - temperature = FluidSystem::vaporTemperature(fluidState, wPhaseIdx); - else - DUNE_THROW(Dune::InvalidStateException, "phasePresence: " << phasePresence << " is invalid."); - - Valgrind::CheckDefined(temperature); - - for(int phaseIdx=0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) - { - fluidState.setTemperature(phaseIdx, temperature); - } - solidState.setTemperature(temperature); + // set the temperature + updateTemperature(elemSol, problem, element, scv, fluidState, solidState); // set the densities for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) @@ -257,6 +243,45 @@ public: } } + //! Depending on the phase state, the fluid temperature is either obtained as a primary variable from the solution vector + //! or calculated from the liquid's vapor pressure. + template<class ElemSol, class Problem, class Element, class Scv> + void updateTemperature(const ElemSol& elemSol, + const Problem& problem, + const Element& element, + const Scv& scv, + FluidState& fluidState, + SolidState& solidState) + { + const auto& priVars = elemSol[scv.localDofIndex()]; + const auto phasePresence = priVars.state(); + const int wPhaseIdx = problem.spatialParams().template wettingPhase<FluidSystem>(element, scv, elemSol); + + // get temperature + Scalar fluidTemperature; + if (phasePresence == liquidPhaseOnly || phasePresence == gasPhaseOnly) + fluidTemperature = priVars[switchIdx]; + else if (phasePresence == twoPhases) + fluidTemperature = FluidSystem::vaporTemperature(fluidState, wPhaseIdx); + else + DUNE_THROW(Dune::InvalidStateException, "phasePresence: " << phasePresence << " is invalid."); + + Valgrind::CheckDefined(fluidTemperature); + + // the model assumes that all fluid phases have the same temperature + for (int phaseIdx=0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) + fluidState.setTemperature(phaseIdx, fluidTemperature); + + // the solid phase could have a different temperature + if (Traits::ModelTraits::numEnergyEq() == 1) + solidState.setTemperature(fluidTemperature); + else + { + const Scalar solidTemperature = elemSol[scv.localDofIndex()][Traits::ModelTraits::numEq()-1]; + solidState.setTemperature(solidTemperature); + } + } + /*! * \brief Returns the fluid state for the control-volume. */ -- GitLab