13#include <range/v3/algorithm/any_of.hpp>
14#include <range/v3/algorithm/set_algorithm.hpp>
15#include <range/v3/view/drop.hpp>
16#include <range/v3/view/transform.hpp>
27 std::vector<GlobalVector*>
const& x_prev)
29 for (
auto*
const solution : x)
33 for (
auto*
const solution : x_prev)
47 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
48 std::vector<std::unique_ptr<ParameterLib::ParameterBase>>
const& parameters,
49 unsigned const integration_order,
50 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
53 const bool use_monolithic_scheme)
54 : name(std::move(name_)),
56 _secondary_variables(std::move(secondary_variables)),
57 cell_average_data_(mesh),
58 _jacobian_assembler(std::move(jacobian_assembler)),
59 _global_assembler(*_jacobian_assembler),
60 _use_monolithic_scheme(use_monolithic_scheme),
61 _integration_order(integration_order),
62 _process_variables(std::move(process_variables)),
64 [&](const std::size_t number_of_process_variables)
67 std::vector<BoundaryConditionCollection> pcs_BCs;
68 pcs_BCs.reserve(number_of_process_variables);
69 for (std::size_t i = 0; i < number_of_process_variables; i++)
74 }(_process_variables.size())),
75 _source_term_collections(
76 [&](
const std::size_t number_of_processes)
77 -> std::vector<SourceTermCollection>
79 std::vector<SourceTermCollection> pcs_sts;
80 pcs_sts.reserve(number_of_processes);
81 for (std::size_t i = 0; i < number_of_processes; i++)
86 }(_process_variables.size()))
92 std::map<
int, std::shared_ptr<MaterialPropertyLib::Medium>>
const& media)
97 per_process_BCs.addBCsForProcessVariables(per_process_variables, dof_table,
101 per_process_sts.addSourceTermsForProcessVariables(
106 std::map<
int, std::shared_ptr<MaterialPropertyLib::Medium>>
const& media)
112 for (std::size_t pcs_id = 0; pcs_id < number_of_processes; pcs_id++)
120 std::map<
int, std::shared_ptr<MaterialPropertyLib::Medium>>
const& media)
122 DBUG(
"Initialize process.");
124 DBUG(
"Construct dof mappings.");
127 DBUG(
"Compute sparsity pattern");
130 DBUG(
"Initialize the extrapolator");
136 DBUG(
"Initialize boundary conditions.");
141 std::vector<GlobalVector*>& process_solutions,
142 std::vector<GlobalVector*>
const& process_solutions_prev,
144 int const process_id)
146 auto& x = *process_solutions[process_id];
147 auto& x_prev = *process_solutions_prev[process_id];
150 auto const& dof_table_of_process =
getDOFTable(process_id);
153 for (std::size_t variable_id = 0;
154 variable_id < per_process_variables.size();
159 auto const& pv = per_process_variables[variable_id];
160 DBUG(
"Set the initial condition of variable {:s} of process {:d}.",
161 pv.get().getName().data(), process_id);
163 auto const& ic = pv.get().getInitialCondition();
165 auto const num_comp = pv.get().getNumberOfGlobalComponents();
167 for (
int component_id = 0; component_id < num_comp; ++component_id)
169 auto const& mesh_subset =
170 dof_table_of_process.getMeshSubset(variable_id, component_id);
171 auto const mesh_id = mesh_subset.getMeshID();
172 for (
auto const* node : mesh_subset.getNodes())
178 node->getID(), {}, *node};
179 auto const& ic_value = ic(t, pos);
182 std::abs(dof_table_of_process.getGlobalIndex(l, variable_id,
193 if (global_index == x.size())
196 x.set(global_index, ic_value[component_id]);
214 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
219 const int process_id)
222 for (
auto const& variable : variables_per_process)
224 variable.get().updateDeactivatedSubdomains(time);
229 auto active_elements_ids = ranges::views::transform(
230 [](
auto const& variable)
231 {
return variable.get().getActiveElementIDs(); });
234 if (ranges::any_of(variables_per_process | active_elements_ids,
235 [](
auto const& vector) {
return vector.empty(); }))
244 variables_per_process[0].get().getActiveElementIDs();
246 for (
auto const& pv_active_element_ids :
247 variables_per_process | ranges::views::drop(1) | active_elements_ids)
249 std::vector<std::size_t> new_active_elements;
251 pv_active_element_ids.size());
253 std::back_inserter(new_active_elements));
266 std::vector<GlobalVector*>
const& x,
267 std::vector<GlobalVector*>
const& x_prev,
271 assert(x.size() == x_prev.size());
273 setLocalAccessibleVectors(x, x_prev);
287 std::vector<GlobalVector*>
const& x,
288 std::vector<GlobalVector*>
const& x_prev,
292 assert(x.size() == x_prev.size());
294 setLocalAccessibleVectors(x, x_prev);
316 const int specified_process_id = 0;
327 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
333 std::generate_n(std::back_inserter(all_mesh_subsets),
334 pv.getNumberOfGlobalComponents(),
335 [&]() { return *_mesh_subset_all_nodes; });
339 std::vector<int> vec_var_n_components;
341 back_inserter(vec_var_n_components),
346 std::make_unique<NumLib::LocalToGlobalIndexMap>(
347 std::move(all_mesh_subsets), vec_var_n_components,
354 const int specified_process_id)
361 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
364 std::vector<int> vec_var_n_components;
366 std::generate_n(std::back_inserter(all_mesh_subsets),
369 .getNumberOfGlobalComponents(),
375 .getNumberOfGlobalComponents());
377 std::make_unique<NumLib::LocalToGlobalIndexMap>(
378 std::move(all_mesh_subsets), vec_var_n_components,
385 int const number_of_processes)
const
387 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
388 dof_tables.reserve(number_of_processes);
389 std::generate_n(std::back_inserter(dof_tables), number_of_processes,
390 [&]() {
return &
getDOFTable(dof_tables.size()); });
394std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
401 const bool manage_storage =
false;
407 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component;
410 const bool manage_storage =
true;
413 std::move(all_mesh_subsets_single_component),
424 std::tie(dof_table_single_component, manage_storage) =
427 std::unique_ptr<NumLib::Extrapolator> extrapolator(
429 *dof_table_single_component));
433 std::move(extrapolator), dof_table_single_component, manage_storage);
443 const double delta_t,
const int process_id)
445 for (
auto*
const solution : x)
453 std::vector<GlobalVector*>
const& x_prev,
454 const double t,
const double delta_t,
455 int const process_id)
457 setLocalAccessibleVectors(x, x_prev);
465 std::vector<GlobalVector*>
const& x_prev,
466 const double t,
double const dt,
467 int const process_id)
469 setLocalAccessibleVectors(x, x_prev);
476 std::vector<GlobalVector*>
const& x,
478 int const process_id)
480 for (
auto const* solution : x)
500 std::vector<GlobalVector*>
const& x,
501 std::vector<GlobalVector*>
const& x_prev,
502 int const process_id)
504 for (
auto const* solution : x)
510std::vector<GlobalIndexType>
513 std::vector<GlobalIndexType> indices;
518 auto const& dof_table_of_process =
getDOFTable(process_id);
521 for (std::size_t variable_id = 0;
522 variable_id < per_process_variables.size();
525 auto const& pv = per_process_variables[variable_id];
526 DBUG(
"Set the initial condition of variable {:s} of process {:d}.",
527 pv.get().getName().data(), process_id);
529 if ((pv.get().compensateNonEquilibriumInitialResiduum()))
534 auto const num_comp = pv.get().getNumberOfGlobalComponents();
536 for (
int component_id = 0; component_id < num_comp; ++component_id)
538 auto const& mesh_subset = dof_table_of_process.getMeshSubset(
539 variable_id, component_id);
544 std::abs(dof_table_of_process.getGlobalIndex(
545 l, variable_id, component_id));
547 indices.push_back(global_index);
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Global vector based on Eigen vector.
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
int getNumberOfGlobalComponents() const
Returns the number of components of the process variable.
virtual NumLib::IterationResult postIterationConcreteProcess(GlobalVector const &)
virtual void preTimestepConcreteProcess(std::vector< GlobalVector * > const &, const double, const double, const int)
std::vector< BoundaryConditionCollection > _boundary_conditions
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables(int const number_of_processes) const
virtual void preIterationConcreteProcess(const unsigned, GlobalVector const &)
void constructDofTableOfSpecifiedProcessStaggeredScheme(const int specified_process_id)
static PROCESSLIB_EXPORT const std::string constant_one_parameter_name
void preAssemble(const double t, double const dt, GlobalVector const &x) final
virtual void initializeBoundaryConditions(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
void initialize(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
virtual void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order)=0
Process specific initialization called by initialize().
std::vector< std::size_t > _ids_of_active_elements
Union of active element ids per process variable.
void assemble(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) final
virtual void postTimestepConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, const double, int const)
Process(std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const ¶meters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true)
virtual void computeSecondaryVariableConcrete(double const, double const, std::vector< GlobalVector * > const &, GlobalVector const &, int const)
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
void preOutput(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)
virtual void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)=0
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
virtual void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &, double const, int const)
NumLib::IterationResult postIteration(GlobalVector const &x) final
void postTimestep(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id)
Postprocessing after a complete timestep.
virtual void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalVector &b, GlobalMatrix &Jac)=0
void preTimestep(std::vector< GlobalVector * > const &x, const double t, const double delta_t, const int process_id)
Preprocessing before starting assembly for new timestep.
void computeSecondaryVariable(double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
compute secondary variables for the coupled equations or for output.
void initializeExtrapolator()
virtual void postNonLinearSolverConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, double const, int const)
virtual void constructDofTable()
void constructMonolithicProcessDofTable()
virtual void preAssembleConcreteProcess(const double, double const, GlobalVector const &)
void postNonLinearSolver(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id)
std::vector< GlobalIndexType > getIndicesOfResiduumWithoutInitialCompensation() const override
unsigned const _integration_order
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
std::vector< SourceTermCollection > _source_term_collections
void updateDeactivatedSubdomains(double const time, const int process_id)
void computeSparsityPattern()
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
void setInitialConditions(std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id)
void assembleWithJacobian(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalVector &b, GlobalMatrix &Jac) final
virtual void preOutputConcreteProcess(const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const)
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > const & getProcessVariables() const
ExtrapolatorData _extrapolator_data
GlobalSparsityPattern _sparsity_pattern
void preIteration(const unsigned iter, GlobalVector const &x) final
const bool _use_monolithic_scheme
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const
Handles configuration of several secondary variables from the project file.
IterationResult
Status flags telling the NonlinearSolver if an iteration succeeded.
void finalizeAssembly(PETScMatrix &A)
void copy(PETScVector const &x, PETScVector &y)
void setLocalAccessibleVector(PETScVector const &x)
auto meshLocations(Mesh const &mesh, MeshItemType const item_type)
@ BY_LOCATION
Ordering data by spatial location.
GlobalSparsityPattern computeSparsityPattern(LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
Computes a sparsity pattern for the given inputs.
void setLocalAccessibleVectors(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev)