From e4e349eb7e525be8d0eb12f06428296d8525005f Mon Sep 17 00:00:00 2001 From: hnil Date: Thu, 4 Jun 2020 15:37:44 +0200 Subject: [PATCH] Add LinearizationType, which will be used for sequential methods. --- opm/models/blackoil/blackoilenergymodules.hh | 2 +- opm/models/blackoil/blackoilpolymermodules.hh | 6 ++- opm/models/blackoil/blackoilsolventmodules.hh | 7 +-- .../common/fvbasediscretization.hh | 9 ++++ .../common/fvbaseelementcontext.hh | 9 ++++ .../discretization/common/fvbaselinearizer.hh | 22 +++++++++- .../common/fvbaseprimaryvariables.hh | 6 +-- .../common/linearizationtype.hh | 44 +++++++++++++++++++ 8 files changed, 95 insertions(+), 10 deletions(-) create mode 100644 opm/models/discretization/common/linearizationtype.hh diff --git a/opm/models/blackoil/blackoilenergymodules.hh b/opm/models/blackoil/blackoilenergymodules.hh index aa6e74724..9158cb05d 100644 --- a/opm/models/blackoil/blackoilenergymodules.hh +++ b/opm/models/blackoil/blackoilenergymodules.hh @@ -358,7 +358,7 @@ public: const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx); // set temperature - fs.setTemperature(priVars.makeEvaluation(temperatureIdx, timeIdx)); + fs.setTemperature(priVars.makeEvaluation(temperatureIdx, timeIdx, elemCtx.linearizationType())); } /*! diff --git a/opm/models/blackoil/blackoilpolymermodules.hh b/opm/models/blackoil/blackoilpolymermodules.hh index 6cc874941..ffb37e8bf 100644 --- a/opm/models/blackoil/blackoilpolymermodules.hh +++ b/opm/models/blackoil/blackoilpolymermodules.hh @@ -884,6 +884,7 @@ public: // Use log(v0) as initial value for u auto u = v0AbsLog; bool converged = false; + // TODO make this into parameters for (int i = 0; i < 20; ++i) { auto f = F(u); auto df = dF(u); @@ -1044,10 +1045,11 @@ public: unsigned dofIdx, unsigned timeIdx) { + const auto linearizationType = elemCtx.linearizationType(); const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx); - polymerConcentration_ = priVars.makeEvaluation(polymerConcentrationIdx, timeIdx); + polymerConcentration_ = priVars.makeEvaluation(polymerConcentrationIdx, timeIdx, linearizationType); if (enablePolymerMolarWeight) { - polymerMoleWeight_ = priVars.makeEvaluation(polymerMoleWeightIdx, timeIdx); + polymerMoleWeight_ = priVars.makeEvaluation(polymerMoleWeightIdx, timeIdx, linearizationType); } const Scalar cmax = PolymerModule::plymaxMaxConcentration(elemCtx, dofIdx, timeIdx); diff --git a/opm/models/blackoil/blackoilsolventmodules.hh b/opm/models/blackoil/blackoilsolventmodules.hh index cc5c9bc01..5e349847d 100644 --- a/opm/models/blackoil/blackoilsolventmodules.hh +++ b/opm/models/blackoil/blackoilsolventmodules.hh @@ -902,7 +902,7 @@ public: { const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx); auto& fs = asImp_().fluidState_; - solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx); + solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, elemCtx.linearizationType()); hydrocarbonSaturation_ = fs.saturation(gasPhaseIdx); // apply a cut-off. Don't waste calculations if no solvent @@ -951,10 +951,11 @@ public: MaterialLaw::capillaryPressures(pC, materialParams, fs); //oil is the reference phase for pressure + const auto linearizationType = elemCtx.linearizationType(); if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_pg_Rv) - pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx); + pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx, linearizationType); else { - const Evaluation& po = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx); + const Evaluation& po = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx, linearizationType); pgMisc = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]); } diff --git a/opm/models/discretization/common/fvbasediscretization.hh b/opm/models/discretization/common/fvbasediscretization.hh index dfb5b92f2..0f7e001ee 100644 --- a/opm/models/discretization/common/fvbasediscretization.hh +++ b/opm/models/discretization/common/fvbasediscretization.hh @@ -806,6 +806,15 @@ public: bool enableStorageCache() const { return enableStorageCache_; } + /*! + * \brief Set the value of enable storage cache + * + * Be aware that calling the *CachedStorage() methods if the storage cache is + * disabled will crash the program. + */ + void setEnableStorageCache(bool enableStorageCache) + { enableStorageCache_= enableStorageCache; } + /*! * \brief Retrieve an entry of the cache for the storage term. * diff --git a/opm/models/discretization/common/fvbaseelementcontext.hh b/opm/models/discretization/common/fvbaseelementcontext.hh index 77fc8b2f7..e899a5d31 100644 --- a/opm/models/discretization/common/fvbaseelementcontext.hh +++ b/opm/models/discretization/common/fvbaseelementcontext.hh @@ -30,6 +30,7 @@ #include "fvbaseproperties.hh" +#include #include #include @@ -268,6 +269,14 @@ public: unsigned focusDofIndex() const { return focusDofIdx_; } + /*! + * \brief Returns the linearization type. + * + * \copydetails setLinearizationType() + */ + LinearizationType linearizationType() const + { return this->model().linearizer().getLinearizationType(); } + /*! * \brief Return a reference to the simulator. */ diff --git a/opm/models/discretization/common/fvbaselinearizer.hh b/opm/models/discretization/common/fvbaselinearizer.hh index 1e21f7b5f..bb0c9ccfb 100644 --- a/opm/models/discretization/common/fvbaselinearizer.hh +++ b/opm/models/discretization/common/fvbaselinearizer.hh @@ -29,6 +29,7 @@ #define EWOMS_FV_BASE_LINEARIZER_HH #include "fvbaseproperties.hh" +#include "linearizationtype.hh" #include #include @@ -143,6 +144,12 @@ public: { simulatorPtr_ = &simulator; eraseMatrix(); + auto it = elementCtx_.begin(); + const auto& endIt = elementCtx_.end(); + for (; it != endIt; ++it){ + delete *it; + } + elementCtx_.resize(0); } /*! @@ -160,7 +167,11 @@ public: /*! * \brief Linearize the full system of non-linear equations. * - * This means the spatial domain plus all auxiliary equations. + * The linearizationType() controls the scheme used and the focus + * time index. The default is fully implicit scheme, and focus index + * equal to 0, i.e. current time (end of step). + * + * This linearizes the spatial domain and all auxiliary equations. */ void linearize() { @@ -263,6 +274,14 @@ public: GlobalEqVector& residual() { return residual_; } + void setLinearizationType(LinearizationType linearizationType){ + linearizationType_ = linearizationType; + }; + + const LinearizationType& getLinearizationType() const{ + return linearizationType_; + }; + /*! * \brief Returns the map of constraint degrees of freedom. * @@ -573,6 +592,7 @@ private: // the right-hand side GlobalEqVector residual_; + LinearizationType linearizationType_; std::mutex globalMatrixMutex_; }; diff --git a/opm/models/discretization/common/fvbaseprimaryvariables.hh b/opm/models/discretization/common/fvbaseprimaryvariables.hh index 5b97f80bd..5b4c13cbc 100644 --- a/opm/models/discretization/common/fvbaseprimaryvariables.hh +++ b/opm/models/discretization/common/fvbaseprimaryvariables.hh @@ -31,7 +31,7 @@ #include #include "fvbaseproperties.hh" - +#include "linearizationtype.hh" #include #include #include @@ -87,13 +87,13 @@ public: * it represents the a constant f = x_i. (the difference is that in the first case, * the derivative w.r.t. x_i is 1, while it is 0 in the second case. */ - Evaluation makeEvaluation(unsigned varIdx, unsigned timeIdx) const + Evaluation makeEvaluation(unsigned varIdx, unsigned timeIdx, LinearizationType linearizationType = LinearizationType()) const { if (std::is_same::value) return (*this)[varIdx]; // finite differences else { // automatic differentiation - if (timeIdx == 0) + if (timeIdx == linearizationType.time) return Toolbox::createVariable((*this)[varIdx], varIdx); else return Toolbox::createConstant((*this)[varIdx]); diff --git a/opm/models/discretization/common/linearizationtype.hh b/opm/models/discretization/common/linearizationtype.hh new file mode 100644 index 000000000..b4f0cbea4 --- /dev/null +++ b/opm/models/discretization/common/linearizationtype.hh @@ -0,0 +1,44 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/* + This file is part of the Open Porous Media project (OPM). + + OPM 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. + + OPM 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 OPM. If not, see . + + Consult the COPYING file in the top-level source directory of this + module for the precise wording of the license and the list of + copyright holders. +*/ +/*! + * \file + * + * \copydoc Opm::FvBaseLinearizer + */ +#ifndef EWOMS_FV_BASE_LINEARIZATIONTYPE_HH +#define EWOMS_FV_BASE_LINEARIZATIONTYPE_HH + +namespace Opm +{ + +struct LinearizationType +{ + enum VarType { implicit, pressure, seqtransport }; + VarType type = implicit; + unsigned int time = 0; +}; + +} // namespace Opm + + +#endif