remove unused ebos/eclnewtonmethod.hh

This commit is contained in:
Arne Morten Kvarving 2024-02-05 19:22:15 +01:00
parent fb7c4f6fd2
commit 2a91fb4c2f
6 changed files with 0 additions and 384 deletions

View File

@ -410,7 +410,6 @@ list (APPEND TEST_DATA_FILES
# originally generated with the command:
# find opm -name '*.h*' -a ! -name '*-pch.hpp' -printf '\t%p\n' | sort
list (APPEND PUBLIC_HEADER_FILES
ebos/eclnewtonmethod.hh
opm/simulators/flow/ActionHandler.hpp
opm/simulators/flow/AluGridCartesianIndexMapper.hpp
opm/simulators/flow/AluGridVanguard.hpp

View File

@ -1,269 +0,0 @@
// -*- 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 <http://www.gnu.org/licenses/>.
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::EclNewtonMethod
*/
#ifndef EWOMS_ECL_NEWTON_METHOD_HH
#define EWOMS_ECL_NEWTON_METHOD_HH
#include <opm/common/Exceptions.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/models/blackoil/blackoilnewtonmethod.hh>
#include <opm/models/utils/signum.hh>
namespace Opm::Properties {
template<class TypeTag, class MyTypeTag>
struct EclNewtonSumTolerance {
using type = UndefinedProperty;
};
template<class TypeTag, class MyTypeTag>
struct EclNewtonStrictIterations {
using type = UndefinedProperty;
};
template<class TypeTag, class MyTypeTag>
struct EclNewtonRelaxedVolumeFraction {
using type = UndefinedProperty;
};
template<class TypeTag, class MyTypeTag>
struct EclNewtonSumToleranceExponent {
using type = UndefinedProperty;
};
template<class TypeTag, class MyTypeTag>
struct EclNewtonRelaxedTolerance {
using type = UndefinedProperty;
};
} // namespace Opm::Properties
namespace Opm {
/*!
* \brief A newton solver which is ebos specific.
*/
template <class TypeTag>
class EclNewtonMethod : public BlackOilNewtonMethod<TypeTag>
{
using ParentType = BlackOilNewtonMethod<TypeTag>;
using DiscNewtonMethod = GetPropType<TypeTag, Properties::DiscNewtonMethod>;
using Simulator = GetPropType<TypeTag, Properties::Simulator>;
using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
using EqVector = GetPropType<TypeTag, Properties::EqVector>;
using Indices = GetPropType<TypeTag, Properties::Indices>;
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using Linearizer = GetPropType<TypeTag, Properties::Linearizer>;
using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
static constexpr int contiSolventEqIdx = Indices::contiSolventEqIdx;
static constexpr int contiPolymerEqIdx = Indices::contiPolymerEqIdx;
static constexpr int contiEnergyEqIdx = Indices::contiEnergyEqIdx;
friend NewtonMethod<TypeTag>;
friend DiscNewtonMethod;
friend ParentType;
public:
EclNewtonMethod(Simulator& simulator) : ParentType(simulator)
{
errorPvFraction_ = 1.0;
relaxedMaxPvFraction_ = EWOMS_GET_PARAM(TypeTag, Scalar, EclNewtonRelaxedVolumeFraction);
sumTolerance_ = 0.0; // this gets determined in the error calculation proceedure
relaxedTolerance_ = EWOMS_GET_PARAM(TypeTag, Scalar, EclNewtonRelaxedTolerance);
numStrictIterations_ = EWOMS_GET_PARAM(TypeTag, int, EclNewtonStrictIterations);
}
/*!
* \brief Register all run-time parameters for the Newton method.
*/
static void registerParameters()
{
ParentType::registerParameters();
EWOMS_REGISTER_PARAM(TypeTag, Scalar, EclNewtonSumTolerance,
"The maximum error tolerated by the Newton"
"method for considering a solution to be "
"converged");
EWOMS_REGISTER_PARAM(TypeTag, int, EclNewtonStrictIterations,
"The number of Newton iterations where the"
" volumetric error is considered.");
EWOMS_REGISTER_PARAM(TypeTag, Scalar, EclNewtonRelaxedVolumeFraction,
"The fraction of the pore volume of the reservoir "
"where the volumetric error may be voilated during "
"strict Newton iterations.");
EWOMS_REGISTER_PARAM(TypeTag, Scalar, EclNewtonSumToleranceExponent,
"The the exponent used to scale the sum tolerance by "
"the total pore volume of the reservoir.");
EWOMS_REGISTER_PARAM(TypeTag, Scalar, EclNewtonRelaxedTolerance,
"The maximum error which the volumetric residual "
"may exhibit if it is in a 'relaxed' "
"region during a strict iteration.");
}
/*!
* \brief Returns true if the error of the solution is below the
* tolerance.
*/
bool converged() const
{
if (errorPvFraction_ < relaxedMaxPvFraction_)
return (this->error_ < relaxedTolerance_ && errorSum_ < sumTolerance_) ;
else if (this->numIterations() > numStrictIterations_)
return (this->error_ < relaxedTolerance_ && errorSum_ < sumTolerance_) ;
return this->error_ <= this->tolerance() && errorSum_ <= sumTolerance_;
}
void preSolve_(const SolutionVector&,
const GlobalEqVector& currentResidual)
{
const auto& constraintsMap = this->model().linearizer().constraintsMap();
this->lastError_ = this->error_;
Scalar newtonMaxError = EWOMS_GET_PARAM(TypeTag, Scalar, NewtonMaxError);
// calculate the error as the maximum weighted tolerance of
// the solution's residual
this->error_ = 0.0;
Dune::FieldVector<Scalar, numEq> componentSumError;
std::fill(componentSumError.begin(), componentSumError.end(), 0.0);
Scalar sumPv = 0.0;
errorPvFraction_ = 0.0;
const Scalar dt = this->simulator_.timeStepSize();
for (unsigned dofIdx = 0; dofIdx < currentResidual.size(); ++dofIdx) {
// do not consider auxiliary DOFs for the error
if (dofIdx >= this->model().numGridDof()
|| this->model().dofTotalVolume(dofIdx) <= 0.0)
continue;
if (!this->model().isLocalDof(dofIdx))
continue;
// also do not consider DOFs which are constraint
if (this->enableConstraints_()) {
if (constraintsMap.count(dofIdx) > 0)
continue;
}
const auto& r = currentResidual[dofIdx];
Scalar pvValue =
this->simulator_.problem().referencePorosity(dofIdx, /*timeIdx=*/0)
* this->model().dofTotalVolume(dofIdx);
sumPv += pvValue;
bool cnvViolated = false;
Scalar dofVolume = this->model().dofTotalVolume(dofIdx);
for (unsigned eqIdx = 0; eqIdx < r.size(); ++eqIdx) {
Scalar tmpError = r[eqIdx] * dt * this->model().eqWeight(dofIdx, eqIdx) / pvValue;
Scalar tmpError2 = r[eqIdx] * this->model().eqWeight(dofIdx, eqIdx);
// in the case of a volumetric formulation, the residual in the above is
// per cubic meter
if (getPropValue<TypeTag, Properties::UseVolumetricResidual>()) {
tmpError *= dofVolume;
tmpError2 *= dofVolume;
}
this->error_ = max(std::abs(tmpError), this->error_);
if (std::abs(tmpError) > this->tolerance_)
cnvViolated = true;
componentSumError[eqIdx] += std::abs(tmpError2);
}
if (cnvViolated)
errorPvFraction_ += pvValue;
}
// take the other processes into account
this->error_ = this->comm_.max(this->error_);
componentSumError = this->comm_.sum(componentSumError);
sumPv = this->comm_.sum(sumPv);
errorPvFraction_ = this->comm_.sum(errorPvFraction_);
componentSumError /= sumPv;
componentSumError *= dt;
errorPvFraction_ /= sumPv;
errorSum_ = 0;
for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx)
errorSum_ = std::max(std::abs(componentSumError[eqIdx]), errorSum_);
// scale the tolerance for the total error with the pore volume. by default, the
// exponent is 1/3, i.e., cubic root.
Scalar x = EWOMS_GET_PARAM(TypeTag, Scalar, EclNewtonSumTolerance);
Scalar y = EWOMS_GET_PARAM(TypeTag, Scalar, EclNewtonSumToleranceExponent);
sumTolerance_ = x*std::pow(sumPv, y);
this->endIterMsg() << " (max: " << this->tolerance_ << ", violated for " << errorPvFraction_*100 << "% of the pore volume), aggegate error: " << errorSum_ << " (max: " << sumTolerance_ << ")";
// make sure that the error never grows beyond the maximum
// allowed one
if (this->error_ > newtonMaxError)
throw NumericalProblem("Newton: Error "+std::to_string(double(this->error_))
+ " is larger than maximum allowed error of "
+ std::to_string(double(newtonMaxError)));
// make sure that the error never grows beyond the maximum
// allowed one
if (errorSum_ > newtonMaxError)
throw NumericalProblem("Newton: Sum of the error "+std::to_string(double(errorSum_))
+ " is larger than maximum allowed error of "
+ std::to_string(double(newtonMaxError)));
}
void endIteration_(SolutionVector& nextSolution,
const SolutionVector& currentSolution)
{
ParentType::endIteration_(nextSolution, currentSolution);
OpmLog::debug( "Newton iteration " + std::to_string(this->numIterations_) + ""
+ " error: " + std::to_string(double(this->error_))
+ this->endIterMsg().str());
this->endIterMsg().str("");
}
private:
Scalar errorPvFraction_;
Scalar errorSum_;
Scalar relaxedTolerance_;
Scalar relaxedMaxPvFraction_;
Scalar sumTolerance_;
int numStrictIterations_;
};
} // namespace Opm
#endif

View File

@ -155,11 +155,6 @@ namespace Opm {
EWOMS_HIDE_PARAM(TypeTag, NewtonTargetIterations);
EWOMS_HIDE_PARAM(TypeTag, NewtonVerbose);
EWOMS_HIDE_PARAM(TypeTag, NewtonWriteConvergence);
EWOMS_HIDE_PARAM(TypeTag, EclNewtonSumTolerance);
EWOMS_HIDE_PARAM(TypeTag, EclNewtonSumToleranceExponent);
EWOMS_HIDE_PARAM(TypeTag, EclNewtonStrictIterations);
EWOMS_HIDE_PARAM(TypeTag, EclNewtonRelaxedVolumeFraction);
EWOMS_HIDE_PARAM(TypeTag, EclNewtonRelaxedTolerance);
// the default eWoms checkpoint/restart mechanism does not work with flow
EWOMS_HIDE_PARAM(TypeTag, RestartTime);

View File

@ -34,8 +34,6 @@
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <ebos/eclnewtonmethod.hh>
#include <opm/common/utility/TimeService.hpp>
#include <opm/core/props/satfunc/RelpermDiagnostics.hpp>

View File

@ -28,8 +28,6 @@
#ifndef OPM_FLOW_PROBLEM_PROPERTIES_HPP
#define OPM_FLOW_PROBLEM_PROPERTIES_HPP
#include <ebos/eclnewtonmethod.hh>
#include <opm/input/eclipse/Parser/ParserKeywords/E.hpp>
#include <opm/material/fluidmatrixinteractions/EclMaterialLawManager.hpp>
@ -306,71 +304,6 @@ struct NewtonTolerance<TypeTag, TTag::FlowBaseProblem> {
static constexpr type value = 1e-2;
};
// the tolerated amount of "incorrect" amount of oil per time step for the complete
// reservoir. this is scaled by the pore volume of the reservoir, i.e., larger reservoirs
// will tolerate larger residuals.
template<class TypeTag>
struct EclNewtonSumTolerance<TypeTag, TTag::FlowBaseProblem> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 1e-4;
};
// set the exponent for the volume scaling of the sum tolerance: larger reservoirs can
// tolerate a higher amount of mass lost per time step than smaller ones! since this is
// not linear, we use the cube root of the overall pore volume by default, i.e., the
// value specified by the NewtonSumTolerance parameter is the "incorrect" mass per
// timestep for an reservoir that exhibits 1 m^3 of pore volume. A reservoir with a total
// pore volume of 10^3 m^3 will tolerate 10 times as much.
template<class TypeTag>
struct EclNewtonSumToleranceExponent<TypeTag, TTag::FlowBaseProblem> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 1.0/3.0;
};
// set number of Newton iterations where the volumetric residual is considered for
// convergence
template<class TypeTag>
struct EclNewtonStrictIterations<TypeTag, TTag::FlowBaseProblem> {
static constexpr int value = 8;
};
// set fraction of the pore volume where the volumetric residual may be violated during
// strict Newton iterations
template<class TypeTag>
struct EclNewtonRelaxedVolumeFraction<TypeTag, TTag::FlowBaseProblem> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 0.03;
};
// the maximum volumetric error of a cell in the relaxed region
template<class TypeTag>
struct EclNewtonRelaxedTolerance<TypeTag, TTag::FlowBaseProblem> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 1e9;
};
// Ignore the maximum error mass for early termination of the newton method.
template<class TypeTag>
struct NewtonMaxError<TypeTag, TTag::FlowBaseProblem> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 10e9;
};
// set the maximum number of Newton iterations to 14 because the likelyhood that a time
// step succeeds at more than 14 Newton iteration is rather small
template<class TypeTag>
struct NewtonMaxIterations<TypeTag, TTag::FlowBaseProblem> {
static constexpr int value = 14;
};
// also, reduce the target for the "optimum" number of Newton iterations to 6. Note that
// this is only relevant if the time step is reduced from the report step size for some
// reason. (because ebos first tries to do a report step using a single time step.)
template<class TypeTag>
struct NewtonTargetIterations<TypeTag, TTag::FlowBaseProblem> {
static constexpr int value = 6;
};
// Disable the VTK output by default for this problem ...
template<class TypeTag>
struct EnableVtkOutput<TypeTag, TTag::FlowBaseProblem> {
@ -504,13 +437,6 @@ struct GradientCalculator<TypeTag, TTag::FlowBaseProblem> {
using type = DummyGradientCalculator<TypeTag>;
};
// Use a custom Newton-Raphson method class for ebos in order to attain more
// sophisticated update and error computation mechanisms
template<class TypeTag>
struct NewtonMethod<TypeTag, TTag::FlowBaseProblem> {
using type = EclNewtonMethod<TypeTag>;
};
// The frequency of writing restart (*.ers) files. This is the number of time steps
// between writing restart files
template<class TypeTag>
@ -524,7 +450,6 @@ struct RestartWritingInterval<TypeTag, TTag::FlowBaseProblem> {
template<class TypeTag>
struct EnableDriftCompensation<TypeTag, TTag::FlowBaseProblem> {
static constexpr bool value = true;
};
// By default, we enable the debugging checks if we're compiled in debug mode

View File

@ -119,38 +119,6 @@ struct NewtonTolerance<TypeTag, TTag::TestTypeTag> {
static constexpr type value = 1e-1;
};
// set fraction of the pore volume where the volumetric residual may be violated during
// strict Newton iterations
template<class TypeTag>
struct EclNewtonRelaxedVolumeFraction<TypeTag, TTag::TestTypeTag> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 0.05;
};
// the maximum volumetric error of a cell in the relaxed region
template<class TypeTag>
struct EclNewtonRelaxedTolerance<TypeTag, TTag::TestTypeTag> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 1e6*getPropValue<TypeTag, Properties::NewtonTolerance>();
};
// the tolerated amount of "incorrect" amount of oil per time step for the complete
// reservoir. this is scaled by the pore volume of the reservoir, i.e., larger reservoirs
// will tolerate larger residuals.
template<class TypeTag>
struct EclNewtonSumTolerance<TypeTag, TTag::TestTypeTag> {
using type = GetPropType<TypeTag, Scalar>;
static constexpr type value = 1e-5;
};
// make all Newton iterations strict, i.e., the volumetric Newton tolerance must be
// always be upheld in the majority of the spatial domain. In this context, "majority"
// means 1 - EclNewtonRelaxedVolumeFraction.
template<class TypeTag>
struct EclNewtonStrictIterations<TypeTag, TTag::TestTypeTag> {
static constexpr int value = 100;
};
// set the maximum number of Newton iterations to 8 so that we fail quickly (albeit
// relatively often)
template<class TypeTag>