From 69c996ef13d6fdea3d2e64ec88130810c2610887 Mon Sep 17 00:00:00 2001 From: Alex Lindsay Date: Mon, 28 Oct 2024 11:23:52 -0700 Subject: [PATCH] Apply suggestions from code review - Remove commented method - `getPossiblyNullTimeIntegrator` -> `queryTimeIntegrator` - `const` addition to `from` vector in `copyVector` - Add method name to `mooseError` in `TimeIntegrator::solve` - Add comment regarding handling of null vectors - Create computeDuDotDu and duDotDuCoeff methods Co-authored-by: Logan Harbour --- framework/include/systems/SystemBase.h | 4 +--- framework/include/timeintegrators/BDF2.h | 2 ++ .../timeintegrators/CentralDifference.h | 2 ++ .../include/timeintegrators/CrankNicolson.h | 2 ++ .../timeintegrators/ExplicitSSPRungeKutta.h | 2 ++ .../include/timeintegrators/NewmarkBeta.h | 2 ++ .../include/timeintegrators/TimeIntegrator.h | 16 +++++++++++++- framework/src/restart/DataIO.C | 3 +++ framework/src/systems/SystemBase.C | 4 ++-- framework/src/timeintegrators/AStableDirk4.C | 4 +--- .../timeintegrators/ActuallyExplicitEuler.C | 5 +---- framework/src/timeintegrators/BDF2.C | 20 +++++++++--------- .../src/timeintegrators/CentralDifference.C | 10 ++++++--- framework/src/timeintegrators/CrankNicolson.C | 10 ++++++--- framework/src/timeintegrators/ExplicitEuler.C | 5 +---- framework/src/timeintegrators/ExplicitRK2.C | 6 +----- .../timeintegrators/ExplicitSSPRungeKutta.C | 10 ++++++--- .../src/timeintegrators/ExplicitTVDRK2.C | 4 +--- framework/src/timeintegrators/ImplicitEuler.C | 4 +--- .../src/timeintegrators/ImplicitMidpoint.C | 4 +--- framework/src/timeintegrators/LStableDirk2.C | 4 +--- framework/src/timeintegrators/LStableDirk3.C | 4 +--- framework/src/timeintegrators/LStableDirk4.C | 4 +--- framework/src/timeintegrators/NewmarkBeta.C | 10 ++++++--- .../src/timeintegrators/TimeIntegrator.C | 21 +++++++++++++------ framework/src/variables/MooseVariableData.C | 2 +- framework/src/variables/MooseVariableDataFV.C | 2 +- .../src/variables/MooseVariableDataLinearFV.C | 2 +- framework/src/variables/MooseVariableField.C | 2 +- 29 files changed, 98 insertions(+), 72 deletions(-) diff --git a/framework/include/systems/SystemBase.h b/framework/include/systems/SystemBase.h index de914eb7c288..b3be2f184dd6 100644 --- a/framework/include/systems/SystemBase.h +++ b/framework/include/systems/SystemBase.h @@ -882,8 +882,6 @@ class SystemBase : public libMesh::ParallelObject, public ConsoleStreamInterface const std::string & name, InputParameters & parameters); - // virtual void addTimeIntegrator(std::shared_ptr /*ti*/) {} - /// Whether or not there are variables to be restarted from an Exodus mesh file bool hasVarCopy() const { return _var_to_copy.size() > 0; } @@ -958,7 +956,7 @@ class SystemBase : public libMesh::ParallelObject, public ConsoleStreamInterface * integrator is found (this could happen for instance if we're solving a non-transient problem), * then a nullptr will be returned */ - const TimeIntegrator * getPossiblyNullTimeIntegrator(const unsigned int var_num) const; + const TimeIntegrator * queryTimeIntegrator(const unsigned int var_num) const; /** * @returns All the time integrators owned by this system diff --git a/framework/include/timeintegrators/BDF2.h b/framework/include/timeintegrators/BDF2.h index e5956adbb89f..5893756dbaae 100644 --- a/framework/include/timeintegrators/BDF2.h +++ b/framework/include/timeintegrators/BDF2.h @@ -39,6 +39,8 @@ class BDF2 : public TimeIntegrator void computeTimeDerivativeHelper(T & u_dot, const T2 & u, const T3 & u_old, const T4 & u_older) const; + virtual Real duDotDuCoeff() const override; + std::vector & _weight; /// The older solution diff --git a/framework/include/timeintegrators/CentralDifference.h b/framework/include/timeintegrators/CentralDifference.h index 175c600c804c..f5f5d40e4af0 100644 --- a/framework/include/timeintegrators/CentralDifference.h +++ b/framework/include/timeintegrators/CentralDifference.h @@ -31,6 +31,8 @@ class CentralDifference : public ActuallyExplicitEuler ADReal & ad_u_dotdot) const override; protected: + virtual Real duDotDuCoeff() const override; + /// solution vector for \f$ {du^dotdot}\over{du} \f$ Real & _du_dotdot_du; diff --git a/framework/include/timeintegrators/CrankNicolson.h b/framework/include/timeintegrators/CrankNicolson.h index a79b17f7e6ea..90e403707191 100644 --- a/framework/include/timeintegrators/CrankNicolson.h +++ b/framework/include/timeintegrators/CrankNicolson.h @@ -43,6 +43,8 @@ class CrankNicolson : public TimeIntegrator template void computeTimeDerivativeHelper(T & u_dot, const T2 & u_old) const; + virtual Real duDotDuCoeff() const override; + NumericVector & _residual_old; }; diff --git a/framework/include/timeintegrators/ExplicitSSPRungeKutta.h b/framework/include/timeintegrators/ExplicitSSPRungeKutta.h index 9423eb76acb3..5176f9130d90 100644 --- a/framework/include/timeintegrators/ExplicitSSPRungeKutta.h +++ b/framework/include/timeintegrators/ExplicitSSPRungeKutta.h @@ -38,6 +38,8 @@ class ExplicitSSPRungeKutta : public ExplicitTimeIntegrator */ bool solveStage(); + virtual Real duDotDuCoeff() const override; + /// Order of time integration const MooseEnum & _order; diff --git a/framework/include/timeintegrators/NewmarkBeta.h b/framework/include/timeintegrators/NewmarkBeta.h index e2bec22e5d55..aabfc15376b9 100644 --- a/framework/include/timeintegrators/NewmarkBeta.h +++ b/framework/include/timeintegrators/NewmarkBeta.h @@ -42,6 +42,8 @@ class NewmarkBeta : public TimeIntegrator T4 & u_dotdot, const T5 & u_dotdot_old) const; + virtual Real duDotDuCoeff() const override; + /// Newmark time integration parameter-beta Real _beta; diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index af14a4d2e60d..21c2896d322e 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -183,7 +183,18 @@ class TimeIntegrator : public MooseObject, public Restartable * Copy from one vector into another. If the time integrator has been restricted to a subset of * variables, then this will selectively copy their dofs */ - void copyVector(NumericVector & from, NumericVector & to); + void copyVector(const NumericVector & from, NumericVector & to); + + /** + * @returns The \p _du_dot_du multiplicative coefficient, e.g. if \p _du_dot_du is equivalent to + * 2/dt, then this method returns 2 + */ + virtual Real duDotDuCoeff() const { return 1; } + + /** + * Compute \p _du_dot_du + */ + void computeDuDotDu(); FEProblemBase & _fe_problem; SystemBase & _sys; @@ -235,4 +246,7 @@ class TimeIntegrator : public MooseObject, public Restartable /// The variables that this time integrator integrates std::unordered_set & _vars; + + /// A vector that is used for creating 'from' subvectors in the \p copyVector() method + std::unique_ptr> _from_subvector; }; diff --git a/framework/src/restart/DataIO.C b/framework/src/restart/DataIO.C index b47646c7d788..c5dd8abdbdc4 100644 --- a/framework/src/restart/DataIO.C +++ b/framework/src/restart/DataIO.C @@ -326,6 +326,9 @@ dataStore(std::ostream & stream, std::unique_ptr> & v, void * context) { + // Classes may declare unique pointers to vectors as restartable data and never actually create + // vector instances. This happens for example in the `TimeIntegrator` class where subvector + // instances are only created if multiple time integrators are present bool have_vector = v.get(); dataStore(stream, have_vector, context); if (!have_vector) diff --git a/framework/src/systems/SystemBase.C b/framework/src/systems/SystemBase.C index 845ee6525aed..b42544c57568 100644 --- a/framework/src/systems/SystemBase.C +++ b/framework/src/systems/SystemBase.C @@ -1596,7 +1596,7 @@ SystemBase::copyTimeIntegrators(const SystemBase & other_sys) } const TimeIntegrator * -SystemBase::getPossiblyNullTimeIntegrator(const unsigned int var_num) const +SystemBase::queryTimeIntegrator(const unsigned int var_num) const { for (auto & ti : _time_integrators) if (ti->integratesVar(var_num)) @@ -1608,7 +1608,7 @@ SystemBase::getPossiblyNullTimeIntegrator(const unsigned int var_num) const const TimeIntegrator & SystemBase::getTimeIntegrator(const unsigned int var_num) const { - const auto * const ti = getPossiblyNullTimeIntegrator(var_num); + const auto * const ti = queryTimeIntegrator(var_num); if (ti) return *ti; diff --git a/framework/src/timeintegrators/AStableDirk4.C b/framework/src/timeintegrators/AStableDirk4.C index 4a23f5580092..08bbbe5e9ac3 100644 --- a/framework/src/timeintegrators/AStableDirk4.C +++ b/framework/src/timeintegrators/AStableDirk4.C @@ -92,9 +92,7 @@ AStableDirk4::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/ActuallyExplicitEuler.C b/framework/src/timeintegrators/ActuallyExplicitEuler.C index 5b4a154da796..a02ef4f262f0 100644 --- a/framework/src/timeintegrators/ActuallyExplicitEuler.C +++ b/framework/src/timeintegrators/ActuallyExplicitEuler.C @@ -49,10 +49,7 @@ ActuallyExplicitEuler::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1.0 / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/BDF2.C b/framework/src/timeintegrators/BDF2.C index 43a3f6bb15ab..87f5123df0a1 100644 --- a/framework/src/timeintegrators/BDF2.C +++ b/framework/src/timeintegrators/BDF2.C @@ -50,21 +50,12 @@ BDF2::computeTimeDerivatives() NumericVector & u_dot = *_sys.solutionUDot(); if (_t_step == 1) - { u_dot = *_solution; - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; - } else - { u_dot.zero(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = _weight[0] / _dt; - } computeTimeDerivativeHelper(u_dot, *_solution, _solution_old, _solution_older); u_dot.close(); + computeDuDotDu(); } void @@ -85,3 +76,12 @@ BDF2::postResidual(NumericVector & residual) residual += _Re_non_time; residual.close(); } + +Real +BDF2::duDotDuCoeff() const +{ + if (_t_step == 1) + return 1; + else + return _weight[0]; +} diff --git a/framework/src/timeintegrators/CentralDifference.C b/framework/src/timeintegrators/CentralDifference.C index 869df4988321..1e77b4957930 100644 --- a/framework/src/timeintegrators/CentralDifference.C +++ b/framework/src/timeintegrators/CentralDifference.C @@ -87,9 +87,7 @@ CentralDifference::computeTimeDerivatives() u_dot.close(); // used for Jacobian calculations - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1.0 / (2 * _dt); + computeDuDotDu(); _du_dotdot_du = 1.0 / (_dt * _dt); // Computing udotdot "factor" @@ -116,3 +114,9 @@ CentralDifference::computeTimeDerivatives() u_dot_factor.close(); } } + +Real +CentralDifference::duDotDuCoeff() const +{ + return Real(1) / Real(2); +} diff --git a/framework/src/timeintegrators/CrankNicolson.C b/framework/src/timeintegrators/CrankNicolson.C index f5dbbc9c9b43..78cc9a5efa81 100644 --- a/framework/src/timeintegrators/CrankNicolson.C +++ b/framework/src/timeintegrators/CrankNicolson.C @@ -76,9 +76,7 @@ CrankNicolson::init() // time derivative is assumed to be zero on initial NumericVector & u_dot = *_sys.solutionUDot(); u_dot.zero(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 0; + computeDuDotDu(); // compute residual for the initial time step // Note: we can not directly pass _residual_old in computeResidualTag because @@ -141,3 +139,9 @@ CrankNicolson::postStep() { copyVector(_Re_non_time, _residual_old); } + +Real +CrankNicolson::duDotDuCoeff() const +{ + return 2; +} diff --git a/framework/src/timeintegrators/ExplicitEuler.C b/framework/src/timeintegrators/ExplicitEuler.C index bf56516eae59..26f0cd76d98e 100644 --- a/framework/src/timeintegrators/ExplicitEuler.C +++ b/framework/src/timeintegrators/ExplicitEuler.C @@ -43,10 +43,7 @@ ExplicitEuler::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1.0 / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/ExplicitRK2.C b/framework/src/timeintegrators/ExplicitRK2.C index b99af6045e1c..d1ae060345a8 100644 --- a/framework/src/timeintegrators/ExplicitRK2.C +++ b/framework/src/timeintegrators/ExplicitRK2.C @@ -53,11 +53,7 @@ ExplicitRK2::computeTimeDerivatives() NumericVector & u_dot = *_sys.solutionUDot(); u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old, _solution_older); - - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; - u_dot.close(); + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C index f238db646c4b..1dc48694012b 100644 --- a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C +++ b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C @@ -69,9 +69,7 @@ void ExplicitSSPRungeKutta::computeTimeDerivatives() { // Only the Jacobian needs to be computed, since the mass matrix needs it - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1.0 / (_b[_stage] * _dt); + computeDuDotDu(); } void @@ -209,3 +207,9 @@ ExplicitSSPRungeKutta::postResidual(NumericVector & residual) // Set time at which to evaluate nodal BCs _fe_problem.time() = _current_time; } + +Real +ExplicitSSPRungeKutta::duDotDuCoeff() const +{ + return Real(1) / _b[_stage]; +} diff --git a/framework/src/timeintegrators/ExplicitTVDRK2.C b/framework/src/timeintegrators/ExplicitTVDRK2.C index aed081626516..01d5a747b924 100644 --- a/framework/src/timeintegrators/ExplicitTVDRK2.C +++ b/framework/src/timeintegrators/ExplicitTVDRK2.C @@ -56,9 +56,7 @@ ExplicitTVDRK2::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old, _solution_older); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); u_dot.close(); } diff --git a/framework/src/timeintegrators/ImplicitEuler.C b/framework/src/timeintegrators/ImplicitEuler.C index 899865ccf8dc..5a2a6a447781 100644 --- a/framework/src/timeintegrators/ImplicitEuler.C +++ b/framework/src/timeintegrators/ImplicitEuler.C @@ -49,9 +49,7 @@ ImplicitEuler::computeTimeDerivatives() u_dot.close(); } - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1.0 / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/ImplicitMidpoint.C b/framework/src/timeintegrators/ImplicitMidpoint.C index 8bb834dec8de..f72077c24538 100644 --- a/framework/src/timeintegrators/ImplicitMidpoint.C +++ b/framework/src/timeintegrators/ImplicitMidpoint.C @@ -45,9 +45,7 @@ ImplicitMidpoint::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/LStableDirk2.C b/framework/src/timeintegrators/LStableDirk2.C index cfd22dc8d25b..06cc76457950 100644 --- a/framework/src/timeintegrators/LStableDirk2.C +++ b/framework/src/timeintegrators/LStableDirk2.C @@ -48,9 +48,7 @@ LStableDirk2::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/LStableDirk3.C b/framework/src/timeintegrators/LStableDirk3.C index 0412e0ae680b..006c5e1748d4 100644 --- a/framework/src/timeintegrators/LStableDirk3.C +++ b/framework/src/timeintegrators/LStableDirk3.C @@ -67,9 +67,7 @@ LStableDirk3::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/LStableDirk4.C b/framework/src/timeintegrators/LStableDirk4.C index ea9af65567da..df5a8f88f788 100644 --- a/framework/src/timeintegrators/LStableDirk4.C +++ b/framework/src/timeintegrators/LStableDirk4.C @@ -62,9 +62,7 @@ LStableDirk4::computeTimeDerivatives() u_dot = *_solution; computeTimeDerivativeHelper(u_dot, _solution_old); u_dot.close(); - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = 1. / _dt; + computeDuDotDu(); } void diff --git a/framework/src/timeintegrators/NewmarkBeta.C b/framework/src/timeintegrators/NewmarkBeta.C index 0fc03741a4f5..4ff7475fe254 100644 --- a/framework/src/timeintegrators/NewmarkBeta.C +++ b/framework/src/timeintegrators/NewmarkBeta.C @@ -92,9 +92,7 @@ NewmarkBeta::computeTimeDerivatives() // used for Jacobian calculations _du_dotdot_du = 1.0 / _beta / _dt / _dt; - for (const auto i : index_range(_du_dot_du)) - if (integratesVar(i)) - _du_dot_du[i] = _gamma / _beta / _dt; + computeDuDotDu(); } void @@ -120,3 +118,9 @@ NewmarkBeta::postResidual(NumericVector & residual) residual += _Re_non_time; residual.close(); } + +Real +NewmarkBeta::duDotDuCoeff() const +{ + return _gamma / _beta; +} diff --git a/framework/src/timeintegrators/TimeIntegrator.C b/framework/src/timeintegrators/TimeIntegrator.C index 3ea40728227e..ebf2b63f1a99 100644 --- a/framework/src/timeintegrators/TimeIntegrator.C +++ b/framework/src/timeintegrators/TimeIntegrator.C @@ -53,7 +53,8 @@ TimeIntegrator::TimeIntegrator(const InputParameters & parameters) _var_restriction(declareRestartableData( "var_restriction", !getParam>("variables").empty())), _local_indices(declareRestartableData>("local_indices")), - _vars(declareRestartableData>("vars")) + _vars(declareRestartableData>("vars")), + _from_subvector(NumericVector::build(this->comm())) { _fe_problem.setUDotRequested(true); } @@ -104,7 +105,7 @@ TimeIntegrator::init() void TimeIntegrator::solve() { - mooseError("Calling this method is no longer supported"); + mooseError("Calling TimeIntegrator::solve() is no longer supported"); } void @@ -130,17 +131,16 @@ TimeIntegrator::getNumLinearIterationsLastSolve() const } void -TimeIntegrator::copyVector(NumericVector & from, NumericVector & to) +TimeIntegrator::copyVector(const NumericVector & from, NumericVector & to) { if (!_var_restriction) to = from; else { auto to_sub = to.get_subvector(_local_indices); - auto from_sub = from.get_subvector(_local_indices); - *to_sub = *from_sub; + from.create_subvector(*_from_subvector, _local_indices, false); + *to_sub = *_from_subvector; to.restore_subvector(std::move(to_sub), _local_indices); - from.restore_subvector(std::move(from_sub), _local_indices); } } @@ -152,3 +152,12 @@ TimeIntegrator::integratesVar(const unsigned int var_num) const return _vars.count(var_num); } + +void +TimeIntegrator::computeDuDotDu() +{ + const auto coeff = duDotDuCoeff(); + for (const auto i : index_range(_du_dot_du)) + if (integratesVar(i)) + _du_dot_du[i] = coeff / _dt; +} diff --git a/framework/src/variables/MooseVariableData.C b/framework/src/variables/MooseVariableData.C index 14adab99773b..28bf58d72e7a 100644 --- a/framework/src/variables/MooseVariableData.C +++ b/framework/src/variables/MooseVariableData.C @@ -77,7 +77,7 @@ MooseVariableData::MooseVariableData(const MooseVariableField::MooseVariableDataFV(const MooseVariableFV(&_sys) ? true : false), _qrule(nullptr) diff --git a/framework/src/variables/MooseVariableDataLinearFV.C b/framework/src/variables/MooseVariableDataLinearFV.C index 94eac24b1361..0ab886a2e23e 100644 --- a/framework/src/variables/MooseVariableDataLinearFV.C +++ b/framework/src/variables/MooseVariableDataLinearFV.C @@ -41,7 +41,7 @@ MooseVariableDataLinearFV::MooseVariableDataLinearFV( _var_num(_var.number()), _assembly(_subproblem.assembly(_tid, var.kind() == Moose::VAR_SOLVER ? sys.number() : 0)), _element_type(element_type), - _time_integrator(_sys.getPossiblyNullTimeIntegrator(_var_num)), + _time_integrator(_sys.queryTimeIntegrator(_var_num)), _elem(elem), _displaced(dynamic_cast(&_sys) ? true : false), _qrule(nullptr) diff --git a/framework/src/variables/MooseVariableField.C b/framework/src/variables/MooseVariableField.C index 784c128521a6..979181084650 100644 --- a/framework/src/variables/MooseVariableField.C +++ b/framework/src/variables/MooseVariableField.C @@ -27,7 +27,7 @@ MooseVariableField::MooseVariableField(const InputParameters & param : MooseVariableFieldBase(parameters), Moose::FunctorBase::type>(name()), MeshChangedInterface(parameters), - _time_integrator(_sys.getPossiblyNullTimeIntegrator(_var_num)) + _time_integrator(_sys.queryTimeIntegrator(_var_num)) { }