Skip to content
Snippets Groups Projects
Commit 9d5451e5 authored by Santiago Ospina De Los Ríos's avatar Santiago Ospina De Los Ríos
Browse files

Move linear solver out & rename get grid operator method from model obj 

* Is easier to reason about the solver out of the model
* For the testing of the jacobian we only need the stationary grid operator so renaming is necessary
parent 845d7195
No related branches found
No related tags found
1 merge request!39Resolve "Improve testing on jacobians"
Hide whitespace changes
Inline Side-by-side
dune/copasi/common/enum.hh
+ 1
− 6
View file @ 9d5451e5
......@@ -21,8 +21,7 @@ enum class Stages
Constraints = 1 << 4,
LocalOperator = 1 << 5,
GridOperator = 1 << 6,
JacobianOperator = 1 << 7,
Writer = 1 << 8
Writer = 1 << 7
};
//! Policy that setup the grid function spaces in a model
......@@ -48,10 +47,6 @@ static constexpr Stages setup_local_operator =
static constexpr Stages setup_grid_operator =
setup_local_operator | Stages::GridOperator;
//! Policy that setup the jacobian operator in a model
static constexpr Stages setup_jacobian_operator =
setup_grid_operator | Stages::JacobianOperator;
//! Policy that setup the writer in a model
static constexpr Stages setup_writer =
setup_coefficient_vector | Stages::Writer;
......
dune/copasi/common/stepper.hh
+ 13
− 13
View file @ 9d5451e5
......@@ -299,45 +299,45 @@ private:
auto& get_solver(const System& system) const
{
using Coefficients = typename System::State::Coefficients;
using GridOperator = typename System::GridOperator;
using JacobianOperator = typename System::JacobianOperator;
using InstationaryGridOperator = typename System::InstationaryGridOperator;
using LinearSolver = Dune::PDELab::ISTLBackend_NOVLP_BCGS_SSORk<InstationaryGridOperator>;
using NonLinearOperator =
PDELab::Newton<GridOperator, JacobianOperator, Coefficients>;
PDELab::Newton<InstationaryGridOperator, LinearSolver, Coefficients>;
using OneStepOperator = PDELab::
OneStepMethod<Time, GridOperator, NonLinearOperator, Coefficients>;
OneStepMethod<Time, InstationaryGridOperator, NonLinearOperator, Coefficients>;
using InternalState = std::tuple<System const*,
GridOperator const*,
JacobianOperator const*,
InstationaryGridOperator const*,
std::shared_ptr<LinearSolver>,
std::shared_ptr<NonLinearOperator>,
std::shared_ptr<OneStepOperator>>;
GridOperator& grid_operator = *system.get_grid_operator();
JacobianOperator& jacobian_operator = *system.get_jacobian_operator();
InstationaryGridOperator& grid_operator = *system.get_instationary_grid_operator();
auto linear_solver = std::make_unique<LinearSolver>(grid_operator);
// If internal data correspond to input, return cached one-step-operator.
if (_internal_state.has_value() and
_internal_state.type() == typeid(InternalState)) {
auto& internal_state = *std::any_cast<InternalState>(&_internal_state);
if (std::get<0>(internal_state) == &system and
std::get<1>(internal_state) == &grid_operator and
std::get<2>(internal_state) == &jacobian_operator)
std::get<1>(internal_state) == &grid_operator)
return *std::get<4>(internal_state);
}
_logger.trace("Get non-linear operator"_fmt);
auto non_linear_operator =
std::make_shared<NonLinearOperator>(grid_operator, jacobian_operator);
std::make_unique<NonLinearOperator>(grid_operator, *linear_solver);
_logger.trace("Get one step operator"_fmt);
auto one_step_operator = std::make_shared<OneStepOperator>(
auto one_step_operator = std::make_unique<OneStepOperator>(
*_rk_method, grid_operator, *non_linear_operator);
one_step_operator->setVerbosityLevel(5);
_internal_state =
std::make_any<InternalState>(&system,
&grid_operator,
&jacobian_operator,
std::move(linear_solver),
std::move(non_linear_operator),
std::move(one_step_operator));
return *std::get<4>(*std::any_cast<InternalState>(&_internal_state));
......
dune/copasi/model/diffusion_reaction.cc
+ 0
− 10
View file @ 9d5451e5
......@@ -267,14 +267,6 @@ ModelDiffusionReaction<Traits>::setup_grid_operator()
std::make_shared<GOI>(*_spatial_grid_operator, *_temporal_grid_operator);
}
template<class Traits>
void
ModelDiffusionReaction<Traits>::setup_jacobian_operator()
{
_logger.debug("Create jacobian operator"_fmt);
_jacobian_operator = std::make_shared<JacobianOperator>(*_grid_operator);
}
template<class Traits>
void
ModelDiffusionReaction<Traits>::setup_vtk_writer()
......@@ -366,8 +358,6 @@ ModelDiffusionReaction<Traits>::setup(BitFlags<ModelSetup::Stages> setup_policy)
setup_local_operator();
if (setup_policy.test(ModelSetup::Stages::GridOperator))
setup_grid_operator();
if (setup_policy.test(ModelSetup::Stages::JacobianOperator))
setup_jacobian_operator();
if (setup_policy.test(ModelSetup::Stages::Writer))
setup_vtk_writer();
}
......
dune/copasi/model/diffusion_reaction.hh
+ 7
− 15
View file @ 9d5451e5
......@@ -269,9 +269,6 @@ public:
//! Instationary grid operator
using GOI = Dune::PDELab::OneStepGridOperator<GOS, GOT>;
//! Jacobian operator
using JO = Dune::PDELab::ISTLBackend_NOVLP_BCGS_SSORk<GOI>;
using DataHandler = PDELab::vtk::
DGFTreeCommonData<const GFS, const X, PDELab::vtk::DefaultPredicate, GV>;
......@@ -289,10 +286,10 @@ public:
using ConstState = Dune::Copasi::ConstModelState<Grid, GFS, X>;
//! Grid operator type
using GridOperator = GOI;
using StationaryGridOperator = GOS;
//! Jacobian operator type
using JacobianOperator = JO;
//! Grid operator type
using InstationaryGridOperator = GOI;
/**
* @brief Constructs the model
......@@ -398,8 +395,6 @@ public:
setup_local_operator();
if (_grid_operator)
setup_grid_operator();
if (_jacobian_operator)
setup_jacobian_operator();
}
}
......@@ -488,16 +483,15 @@ public:
std::vector<std::shared_ptr<ComponentGridFunction>> get_grid_functions()
const;
//! warning this is not completely const correct. grid operators may modify local operators and thus the model on certain calls
std::shared_ptr<GridOperator> get_grid_operator() const
std::shared_ptr<StationaryGridOperator> get_stationary_grid_operator() const
{
return _grid_operator;
return _spatial_grid_operator;
}
std::shared_ptr<JacobianOperator> get_jacobian_operator() const
std::shared_ptr<InstationaryGridOperator> get_instationary_grid_operator() const
{
return _jacobian_operator;
return _grid_operator;
}
protected:
......@@ -520,7 +514,6 @@ private:
void setup_constraints();
void setup_local_operator();
void setup_grid_operator();
void setup_jacobian_operator();
void setup_vtk_writer();
static auto get_data_handler(const ConstState& state);
......@@ -540,7 +533,6 @@ private:
std::shared_ptr<GOS> _spatial_grid_operator;
std::shared_ptr<GOT> _temporal_grid_operator;
std::shared_ptr<GOI> _grid_operator;
std::shared_ptr<JO> _jacobian_operator;
};
} // namespace Dune::Copasi
......
dune/copasi/model/multidomain_diffusion_reaction.cc
+ 0
− 10
View file @ 9d5451e5
......@@ -247,14 +247,6 @@ ModelMultiDomainDiffusionReaction<Traits>::setup_grid_operator()
*_temporal_grid_operator);
}
template<class Traits>
void
ModelMultiDomainDiffusionReaction<Traits>::setup_jacobian_operator()
{
_logger.debug("Create jacobian operator"_fmt);
_jacobian_operator = std::make_shared<JacobianOperator>(*_grid_operator);
}
template<class Traits>
void
ModelMultiDomainDiffusionReaction<Traits>::setup_vtk_writer()
......@@ -351,8 +343,6 @@ ModelMultiDomainDiffusionReaction<Traits>::setup(
setup_local_operator();
if (setup_policy.test(ModelSetup::Stages::GridOperator))
setup_grid_operator();
if (setup_policy.test(ModelSetup::Stages::JacobianOperator))
setup_jacobian_operator();
if (setup_policy.test(ModelSetup::Stages::Writer))
setup_vtk_writer();
}
......
dune/copasi/model/multidomain_diffusion_reaction.hh
+ 8
− 14
View file @ 9d5451e5
......@@ -160,9 +160,6 @@ private:
//! Instationary grid operator
using GOI = Dune::PDELab::OneStepGridOperator<GOS, GOT>;
//! Linear solver backend
using JO = Dune::PDELab::ISTLBackend_NOVLP_BCGS_SSORk<GOI>;
//! Entity transformation between grids
using EntityTransformation =
Dune::Copasi::MultiDomainEntityTransformation<Grid>;
......@@ -189,10 +186,10 @@ public:
using ConstState = Dune::Copasi::ConstModelState<Grid, GFS, X>;
//! Grid operator type
using GridOperator = GOI;
using StationaryGridOperator = GOS;
//! Jacobian operator type
using JacobianOperator = JO;
//! Grid operator type
using InstationaryGridOperator = GOI;
/**
* @brief Constructs a new instance.
......@@ -252,8 +249,6 @@ public:
setup_local_operator();
if (_grid_operator)
setup_grid_operator();
if (_jacobian_operator)
setup_jacobian_operator();
}
}
......@@ -352,14 +347,15 @@ public:
get_grid_functions() const;
//! warning this is not completely const correct. grid operators may modify local operators and thus the model on certain calls
std::shared_ptr<GridOperator> get_grid_operator() const
std::shared_ptr<StationaryGridOperator> get_stationary_grid_operator() const
{
return _grid_operator;
return _spatial_grid_operator;
}
std::shared_ptr<JacobianOperator> get_jacobian_operator() const
//! warning this is not completely const correct. grid operators may modify local operators and thus the model on certain calls
std::shared_ptr<InstationaryGridOperator> get_instationary_grid_operator() const
{
return _jacobian_operator;
return _grid_operator;
}
protected:
......@@ -379,7 +375,6 @@ private:
auto setup_local_operator(std::size_t) const;
void setup_local_operator();
void setup_grid_operator();
void setup_jacobian_operator();
void setup_solvers();
void setup_vtk_writer();
......@@ -400,7 +395,6 @@ private:
std::shared_ptr<GOS> _spatial_grid_operator;
std::shared_ptr<GOT> _temporal_grid_operator;
std::shared_ptr<GOI> _grid_operator;
std::shared_ptr<JO> _jacobian_operator;
std::size_t _domains;
};
......
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