20 #include "libmesh/default_coupling.h" 21 #include "libmesh/dof_map.h" 22 #include "libmesh/eigen_system.h" 23 #include "libmesh/elem.h" 24 #include "libmesh/explicit_system.h" 25 #include "libmesh/fe_interface.h" 26 #include "libmesh/frequency_system.h" 27 #include "libmesh/int_range.h" 28 #include "libmesh/libmesh_logging.h" 29 #include "libmesh/linear_implicit_system.h" 30 #include "libmesh/mesh_base.h" 31 #include "libmesh/mesh_refinement.h" 32 #include "libmesh/newmark_system.h" 33 #include "libmesh/nonlinear_implicit_system.h" 34 #include "libmesh/parallel.h" 35 #include "libmesh/rb_construction.h" 36 #include "libmesh/remote_elem.h" 37 #include "libmesh/transient_rb_construction.h" 38 #include "libmesh/transient_system.h" 47 #include "libmesh/equation_systems.h" 55 _refine_in_reinit(true),
56 _enable_default_ghosting(true)
60 this->
parameters.
set<
unsigned int>(
"linear solver maximum iterations") = 5000;
104 const unsigned int n_sys = this->
n_systems();
106 libmesh_assert_not_equal_to (n_sys, 0);
110 for (
auto & node :
_mesh.node_ptr_range())
111 node->set_n_systems(n_sys);
113 for (
auto & elem :
_mesh.element_ptr_range())
114 elem->set_n_systems(n_sys);
119 #ifdef LIBMESH_ENABLE_AMR 133 parallel_object_only();
135 const unsigned int n_sys = this->
n_systems();
136 libmesh_assert_not_equal_to (n_sys, 0);
139 for (
unsigned int i=0; i != n_sys; ++i)
152 for (
auto & node :
_mesh.node_ptr_range())
153 node->set_n_systems(n_sys);
156 for (
auto & elem :
_mesh.element_ptr_range())
157 elem->set_n_systems(n_sys);
161 for (
unsigned int i=0; i != n_sys; ++i)
164 #ifdef LIBMESH_ENABLE_AMR 166 bool mesh_changed =
false;
172 for (
unsigned int i=0; i != n_sys; ++i)
250 #endif // #ifdef LIBMESH_ENABLE_AMR 271 const unsigned int n_sys = this->
n_systems();
273 libmesh_assert_not_equal_to (n_sys, 0);
280 for (
auto & node :
_mesh.node_ptr_range())
281 node->set_n_systems(n_sys);
283 for (
auto & elem :
_mesh.element_ptr_range())
284 elem->set_n_systems(n_sys);
311 mesh.add_ghosting_functor(
mesh.default_ghosting());
313 mesh.remove_ghosting_functor(
mesh.default_ghosting());
329 LOG_SCOPE(
"update()",
"EquationSystems");
339 std::string_view
name)
350 else if (sys_type ==
"Basic")
351 this->add_system<System> (
name);
354 else if (sys_type ==
"Newmark")
355 this->add_system<NewmarkSystem> (
name);
358 else if ((sys_type ==
"Explicit"))
359 this->add_system<ExplicitSystem> (
name);
362 else if ((sys_type ==
"Implicit") ||
363 (sys_type ==
"Steady" ))
364 this->add_system<ImplicitSystem> (
name);
367 else if ((sys_type ==
"Transient") ||
368 (sys_type ==
"TransientImplicit") ||
369 (sys_type ==
"TransientLinearImplicit"))
370 this->add_system<TransientLinearImplicitSystem> (
name);
373 else if (sys_type ==
"TransientNonlinearImplicit")
374 this->add_system<TransientNonlinearImplicitSystem> (
name);
377 else if (sys_type ==
"TransientExplicit")
378 this->add_system<TransientExplicitSystem> (
name);
381 else if (sys_type ==
"LinearImplicit")
382 this->add_system<LinearImplicitSystem> (
name);
385 else if (sys_type ==
"NonlinearImplicit")
386 this->add_system<NonlinearImplicitSystem> (
name);
389 else if (sys_type ==
"RBConstruction")
390 this->add_system<RBConstruction> (
name);
393 else if (sys_type ==
"TransientRBConstruction")
394 this->add_system<TransientRBConstruction> (
name);
396 #ifdef LIBMESH_HAVE_SLEPC 398 else if (sys_type ==
"Eigen")
399 this->add_system<EigenSystem> (
name);
400 else if (sys_type ==
"TransientEigenSystem")
401 this->add_system<TransientEigenSystem> (
name);
404 #if defined(LIBMESH_USE_COMPLEX_NUMBERS) 406 else if (sys_type ==
"Frequency")
407 this->add_system<FrequencySystem> (
name);
411 libmesh_error_msg(
"ERROR: Unknown system type: " << sys_type);
435 this->
get_system(i).sensitivity_solve(parameters_in);
444 for (
unsigned int i=this->
n_systems(); i != 0; --i)
445 this->
get_system(i-1).adjoint_solve(qoi_indices);
452 const std::set<std::string> * system_names)
const 455 unsigned int var_num = var_names.size();
459 std::unordered_multiset<std::string> seen_names;
465 unsigned int n_scalar_vars = 0;
466 unsigned int n_vector_vars = 0;
468 for (
const auto & [sys_name, sys_ptr] :
_systems)
471 bool use_current_system = (system_names ==
nullptr);
472 if (!use_current_system)
473 use_current_system = system_names->count(sys_name);
474 if (!use_current_system || sys_ptr->hide_output())
477 seen_names.insert(sys_ptr->variable_name(vn));
483 seen_names.insert(sys_ptr->variable_name(vn));
494 unsigned int nv = n_scalar_vars +
dim*n_vector_vars;
498 libmesh_assert_less_equal ( nv, (
dim > 0 ?
dim : 1)*this->
n_vars() );
503 var_names.resize( nv );
506 for (
const auto & [sys_name, sys_ptr] :
_systems)
509 bool use_current_system = (system_names ==
nullptr);
510 if (!use_current_system)
511 use_current_system = system_names->count(sys_name);
512 if (!use_current_system || sys_ptr->hide_output())
517 const std::string & var_name = sys_ptr->variable_name(vn);
518 const FEType & fe_type = sys_ptr->variable_type(vn);
523 if (type ==
nullptr || (type && *type == fe_type))
531 var_names[var_num++] = var_name;
532 libmesh_error_msg_if(seen_names.count(var_name) > 1,
533 "Duplicate variable name "+var_name);
536 var_names[var_num++] = var_name+
"_x";
537 var_names[var_num++] = var_name+
"_y";
538 libmesh_error_msg_if(seen_names.count(var_name+
"_x"),
539 "Duplicate variable name "+var_name+
"_x");
540 libmesh_error_msg_if(seen_names.count(var_name+
"_y"),
541 "Duplicate variable name "+var_name+
"_y");
544 var_names[var_num++] = var_name+
"_x";
545 var_names[var_num++] = var_name+
"_y";
546 var_names[var_num++] = var_name+
"_z";
547 libmesh_error_msg_if(seen_names.count(var_name+
"_x"),
548 "Duplicate variable name "+var_name+
"_x");
549 libmesh_error_msg_if(seen_names.count(var_name+
"_y"),
550 "Duplicate variable name "+var_name+
"_y");
551 libmesh_error_msg_if(seen_names.count(var_name+
"_z"),
552 "Duplicate variable name "+var_name+
"_z");
555 libmesh_error_msg(
"Invalid dim in build_variable_names");
559 var_names[var_num++] = var_name;
564 var_names.resize(var_num);
571 std::string_view)
const 574 libmesh_not_implemented();
580 std::unique_ptr<NumericVector<Number>>
582 bool add_sides)
const 584 LOG_SCOPE(
"build_parallel_solution_vector()",
"EquationSystems");
587 parallel_object_only();
607 unsigned int n_scalar_vars = 0;
608 unsigned int n_vector_vars = 0;
609 for (
const auto & [sys_name, sys_ptr] :
_systems)
612 bool use_current_system = (system_names ==
nullptr);
613 if (!use_current_system)
614 use_current_system = system_names->count(sys_name);
615 if (!use_current_system || sys_ptr->hide_output())
628 nv = n_scalar_vars +
dim*n_vector_vars;
634 _mesh.local_nodes_end()));
645 const dof_id_type gaps_per_processor = n_gaps / n_proc;
646 const dof_id_type remainder_gaps = n_gaps % n_proc;
648 n_local_nodes = n_local_nodes +
650 (my_pid < remainder_gaps);
655 added_side_nodes = 0;
658 std::vector<dof_id_type> others_added_side_nodes;
662 std::map<std::tuple<dof_id_type, unsigned short, unsigned short>,
670 std::vector<dof_id_type> true_node_offsets;
672 std::vector<dof_id_type> added_node_offsets;
674 auto node_id_to_vec_id =
675 [&true_node_offsets, &added_node_offsets]
678 if (true_node_offsets.empty())
682 const auto lb = std::upper_bound(true_node_offsets.begin(),
683 true_node_offsets.end(), node_id);
687 return node_id + added_node_offsets[p];
692 true_node_offsets.resize(n_proc);
693 added_node_offsets.resize(n_proc);
696 for (
const auto & elem :
_mesh.active_element_ptr_range())
698 for (
auto s : elem->side_index_range())
703 const std::vector<unsigned int> side_nodes =
704 elem->nodes_on_side(s);
707 local_added_side_nodes += side_nodes.size();
711 others_added_side_nodes.resize(n_proc);
713 others_added_side_nodes);
715 added_side_nodes = std::accumulate(others_added_side_nodes.begin(),
716 others_added_side_nodes.end(), 0,
721 true_node_offsets[p+1] += true_node_offsets[p];
726 added_node_offsets[0] = 0;
728 added_node_offsets[p+1] =
729 added_node_offsets[p] + others_added_side_nodes[p];
734 dof_id_type node_counter = true_node_offsets[my_pid];
736 node_counter += others_added_side_nodes[p];
739 for (
const auto & elem :
_mesh.active_local_element_ptr_range())
741 for (
auto s : elem->side_index_range())
746 const std::vector<unsigned int> side_nodes =
747 elem->nodes_on_side(s);
750 discontinuous_node_indices
751 [std::make_tuple(elem->id(),s,n)] = node_counter++;
757 n_global_vals = (max_nn + added_side_nodes) * nv,
758 n_local_vals = (n_local_nodes + local_added_side_nodes) * nv;
763 parallel_soln.
init(n_global_vals, n_local_vals,
false,
PARALLEL);
769 repeat_count.
init(n_global_vals, n_local_vals,
false,
PARALLEL);
771 repeat_count.
close();
773 unsigned int var_num=0;
780 for (
const auto & [sys_name, sys_ptr] :
_systems)
783 bool use_current_system = (system_names ==
nullptr);
784 if (!use_current_system)
785 use_current_system = system_names->count(sys_name);
786 if (!use_current_system || sys_ptr->hide_output())
789 const System & system = *sys_ptr;
790 const unsigned int nv_sys = system.
n_vars();
791 const unsigned int sys_num = system.
number();
794 unsigned int n_scalar_vars = 0;
795 unsigned int n_vector_vars = 0;
806 unsigned int nv_sys_split = n_scalar_vars +
dim*n_vector_vars;
817 if (!non_const_sys.
solution->closed())
826 std::vector<Number> elem_soln;
827 std::vector<Number> nodal_soln;
828 std::vector<dof_id_type> dof_indices;
830 unsigned var_inc = 0;
831 for (
unsigned int var=0; var<nv_sys; var++)
839 for (
const auto & elem :
_mesh.active_local_element_ptr_range())
844 sys_soln.
get(dof_indices, elem_soln);
855 if (!elem->infinite())
857 libmesh_assert_equal_to (nodal_soln.size(), n_vec_dim*elem->n_nodes());
859 for (
auto n : elem->node_index_range())
861 const Node & node = elem->node_ref(n);
864 nv * node_id_to_vec_id(node.
id());
866 for (
unsigned int d=0; d < n_vec_dim; d++)
870 parallel_soln.
add(node_idx + (var_inc+d + var_num), nodal_soln[n_vec_dim*n+d]);
873 repeat_count.
add(node_idx + (var_inc+d + var_num), 1);
879 for (
auto s : elem->side_index_range())
887 (fe_type, elem, s, elem_soln,
888 nodal_soln, add_p_level, n_vec_dim);
891 const std::vector<unsigned int> side_nodes =
892 elem->nodes_on_side(s);
894 libmesh_assert_equal_to
902 std::size_t node_index =
903 nv * libmesh_map_find(discontinuous_node_indices,
904 std::make_tuple(elem->id(), s, n));
906 for (
unsigned int d=0; d < n_vec_dim; d++)
908 parallel_soln.
add(node_index + (var_inc+d + var_num), nodal_soln[n_vec_dim*n+d]);
909 repeat_count.
add(node_index + (var_inc+d + var_num), 1);
917 for (
auto n : elem->node_index_range())
919 const Node & node = elem->node_ref(n);
923 if (!node.
n_dofs(sys_num, var))
926 nv * node_id_to_vec_id(node.
id());
928 for (
unsigned int d=0; d < n_vec_dim; d++)
929 repeat_count.
add(node_idx + (var_inc+d + var_num), 1);
934 var_inc += n_vec_dim;
937 var_num += nv_sys_split;
941 parallel_soln.
close();
942 repeat_count.
close();
956 if (std::abs(repeat_count(i)) <
TOLERANCE)
957 repeat_count.
set(i, 1.);
962 repeat_count.
close();
966 parallel_soln /= repeat_count;
968 return parallel_soln_ptr;
974 const std::set<std::string> * system_names,
975 bool add_sides)
const 977 LOG_SCOPE(
"build_solution_vector()",
"EquationSystems");
980 std::unique_ptr<NumericVector<Number>> parallel_soln =
984 parallel_soln->localize_to_one(soln);
990 std::vector<std::set<subdomain_id_type>> & vars_active_subdomains)
const 992 vars_active_subdomains.clear();
993 vars_active_subdomains.resize(names.size());
997 const auto & sys_ptr = pr.second;
1000 const std::string & var_name = sys_ptr->variable_name(vn);
1002 auto names_it = std::find(names.begin(), names.end(), var_name);
1003 if(names_it != names.end())
1005 const Variable & variable = sys_ptr->variable(vn);
1006 const std::set<subdomain_id_type> & active_subdomains = variable.
active_subdomains();
1007 vars_active_subdomains[
std::distance(names.begin(), names_it)] = active_subdomains;
1015 #ifdef LIBMESH_ENABLE_DEPRECATED 1017 std::vector<std::string> & names)
const 1019 libmesh_deprecated();
1022 #endif // LIBMESH_ENABLE_DEPRECATED 1028 std::vector<std::string> & names)
const 1031 std::unique_ptr<NumericVector<Number>> parallel_soln =
1041 parallel_soln->localize_to_one(soln);
1044 std::vector<std::pair<unsigned int, unsigned int>>
1046 (std::vector<std::string> & names,
const FEType * type,
const std::vector<FEType> * types)
const 1049 parallel_object_only();
1054 libmesh_assert_msg(!type || !types,
1055 "Input 'type', 'types', or neither in find_variable_numbers, but not both.");
1057 std::vector<FEType> type_filter;
1059 type_filter.push_back(*type);
1061 type_filter = *types;
1067 std::vector<std::string> name_filter = names;
1068 bool is_names_empty = name_filter.empty();
1075 const std::vector<std::string> component_suffix = {
"_x",
"_y",
"_z"};
1076 unsigned int dim = _mesh.spatial_dimension();
1077 libmesh_error_msg_if(
dim > 3,
"Invalid dim in find_variable_numbers");
1081 std::vector<std::pair<unsigned int, unsigned int>> var_nums;
1083 for (
const auto & pr : _systems)
1085 const System & system = *(pr.second);
1091 if (type_filter.size() &&
1092 std::find(type_filter.begin(), type_filter.end(), var_type) == type_filter.end())
1098 std::vector<std::string> component_names;
1099 for (
unsigned int comp = 0; comp <
dim; ++comp)
1102 if (is_names_empty ||
1103 (std::find(name_filter.begin(), name_filter.end(),
name) != name_filter.end()))
1104 component_names.push_back(
name);
1107 if (! component_names.empty())
1108 names.insert(names.end(), component_names.begin(), component_names.end());
1115 if (is_names_empty ||
1116 (std::find(name_filter.begin(), name_filter.end(),
name) != name_filter.end()))
1117 names.push_back(
name);
1123 var_nums.emplace_back(system.
number(), var);
1128 std::vector<unsigned int> sort_index(var_nums.size());
1129 std::iota(sort_index.begin(), sort_index.end(), 0);
1130 std::sort(sort_index.begin(), sort_index.end(),
1131 [&](
const unsigned int & lhs,
const unsigned int & rhs)
1132 {
return this->get_system(var_nums[lhs].first).variable_name(var_nums[lhs].second) <
1133 this->get_system(var_nums[rhs].first).variable_name(var_nums[rhs].second);});
1135 std::vector<std::pair<unsigned int, unsigned int>> var_nums_sorted(var_nums.size());
1138 var_nums_sorted[i].first = var_nums[sort_index[i]].first;
1139 var_nums_sorted[i].second = var_nums[sort_index[i]].second;
1143 std::sort(names.begin(), names.end());
1146 return var_nums_sorted;
1150 std::unique_ptr<NumericVector<Number>>
1157 std::vector<std::pair<unsigned int, unsigned int>> var_nums =
1160 const std::size_t nv = names.size();
1166 return std::unique_ptr<NumericVector<Number>>(
nullptr);
1180 parallel_soln_local_size = div+1;
1185 parallel_soln.
init(parallel_soln_global_size,
1186 parallel_soln_local_size,
1190 unsigned int sys_ctr = 0;
1191 unsigned int var_ctr = 0;
1194 std::pair<unsigned int, unsigned int> var_num = var_nums[i];
1198 if (sys_ctr != var_num.first)
1207 if (!non_const_sys.
solution->closed())
1210 sys_ctr = var_num.first;
1216 std::vector<dof_id_type> dof_indices;
1218 const unsigned int var = var_num.second;
1227 const unsigned int n_comps =
1231 for (
const auto & elem :
_mesh.active_local_element_ptr_range())
1238 libmesh_assert_equal_to(dof_indices.size(), n_comps);
1240 for (
unsigned int comp = 0; comp < n_comps; comp++)
1241 parallel_soln.
set(ne * (var_ctr + comp) + elem->id(), sys_soln(dof_indices[comp]));
1252 libmesh_assert_equal_to(names.size(), var_ctr);
1254 parallel_soln.
close();
1255 return parallel_soln_ptr;
1262 (std::vector<Number> & soln,
1263 const std::set<std::string> * system_names,
1264 const std::vector<std::string> * var_names,
1266 bool add_sides)
const 1268 LOG_SCOPE(
"build_discontinuous_solution_vector()",
"EquationSystems");
1274 unsigned int nv = 0;
1276 for (
const auto & [sys_name, sys_ptr] : _systems)
1279 bool use_current_system = (system_names ==
nullptr);
1280 if (!use_current_system)
1281 use_current_system = system_names->count(sys_name);
1282 if (!use_current_system || sys_ptr->hide_output())
1287 for (
auto var_id :
make_range(sys_ptr->n_vars()))
1289 bool use_current_var = (var_names ==
nullptr);
1290 if (!use_current_var)
1291 use_current_var = std::count(var_names->begin(),
1293 sys_ptr->variable_name(var_id));
1297 if (use_current_var)
1305 for (
const auto & elem : _mesh.active_element_ptr_range())
1307 tw += vertices_only ? elem->n_vertices() : elem->n_nodes();
1311 for (
auto s : elem->side_index_range())
1313 if (redundant_added_side(*elem,s))
1316 const std::vector<unsigned int> side_nodes =
1317 elem->nodes_on_side(s);
1320 tw += side_nodes.size();
1323 if (elem->is_vertex(side_nodes[n]))
1331 if (_mesh.processor_id() == 0)
1334 std::vector<Number> sys_soln;
1338 unsigned int var_offset = 0;
1345 for (
const auto & [sys_name, system] : _systems)
1348 bool use_current_system = (system_names ==
nullptr);
1349 if (!use_current_system)
1350 use_current_system = system_names->count(sys_name);
1351 if (!use_current_system || system->hide_output())
1354 const unsigned int nv_sys = system->n_vars();
1355 const auto & dof_map = system->get_dof_map();
1357 system->update_global_solution (sys_soln, 0);
1360 unsigned int n_vars_written_current_system = 0;
1362 if (_mesh.processor_id() == 0)
1364 std::vector<Number> soln_coeffs;
1365 std::vector<Number> nodal_soln;
1366 std::vector<dof_id_type> dof_indices;
1373 for (
unsigned int var=0; var<nv_sys; var++)
1375 bool use_current_var = (var_names ==
nullptr);
1376 if (!use_current_var)
1377 use_current_var = std::count(var_names->begin(),
1379 system->variable_name(var));
1383 if (!use_current_var)
1386 const FEType & fe_type = system->variable_type(var);
1387 const Variable & var_description = system->variable(var);
1388 const auto vg = dof_map.var_group_from_var_number(var);
1389 const bool add_p_level = dof_map.should_p_refine(vg);
1393 for (
auto & elem : _mesh.active_element_ptr_range())
1397 system->get_dof_map().dof_indices (elem, dof_indices, var);
1399 soln_coeffs.resize(dof_indices.size());
1402 soln_coeffs[i] = sys_soln[dof_indices[i]];
1416 if (!elem->infinite())
1418 libmesh_assert_equal_to (nodal_soln.size(), elem->n_nodes());
1420 const unsigned int n_vals =
1421 vertices_only ? elem->n_vertices() : elem->n_nodes();
1423 for (
unsigned int n=0; n<n_vals; n++)
1427 nv * (nn++) + (n_vars_written_current_system + var_offset);
1429 soln[index] += nodal_soln[n];
1434 nn += vertices_only ? elem->n_vertices() : elem->n_nodes();
1445 std::vector<std::vector<const Elem *>>
1446 elems_by_pid(_mesh.n_processors());
1448 for (
const auto & elem : _mesh.active_element_ptr_range())
1449 elems_by_pid[elem->processor_id()].push_back(elem);
1452 for (
const Elem * elem : elems_by_pid[p])
1456 system->get_dof_map().dof_indices (elem, dof_indices, var);
1458 soln_coeffs.resize(dof_indices.size());
1461 soln_coeffs[i] = sys_soln[dof_indices[i]];
1463 for (
auto s : elem->side_index_range())
1465 if (redundant_added_side(*elem,s))
1468 const std::vector<unsigned int> side_nodes =
1469 elem->nodes_on_side(s);
1476 (fe_type, elem, s, soln_coeffs,
1477 nodal_soln, add_p_level,
1480 libmesh_assert_equal_to
1496 neigh->
level() == elem->level() &&
1499 std::vector<dof_id_type> neigh_indices;
1500 system->get_dof_map().dof_indices (neigh, neigh_indices, var);
1501 std::vector<Number> neigh_coeffs(neigh_indices.size());
1504 neigh_coeffs[i] = sys_soln[neigh_indices[i]];
1506 const unsigned int s_neigh =
1508 std::vector<Number> neigh_soln;
1510 (fe_type, neigh, s_neigh,
1511 neigh_coeffs, neigh_soln, add_p_level,
1514 const std::vector<unsigned int> neigh_nodes =
1518 if (neigh->
node_ptr(neigh_nodes[neigh_n])
1519 == elem->node_ptr(side_nodes[n]))
1521 nodal_soln[n] += neigh_soln[neigh_n];
1528 if (vertices_only &&
1529 !elem->is_vertex(n))
1534 nv * (nn++) + (n_vars_written_current_system + var_offset);
1536 soln[index] += nodal_soln[n];
1542 nn += vertices_only ? elem->n_vertices() : elem->n_nodes();
1544 for (
auto s : elem->side_index_range())
1546 if (redundant_added_side(*elem,s))
1549 const std::vector<unsigned int> side_nodes =
1550 elem->nodes_on_side(s);
1554 if (vertices_only &&
1555 !elem->is_vertex(n))
1565 n_vars_written_current_system++;
1571 var_offset += n_vars_written_current_system;
1592 if (!neigh->active())
1597 return (neigh->id() < elem.
id());
1603 const Real threshold,
1604 const bool verbose)
const 1608 std::vector<bool> os_result;
1614 libMesh::out <<
" Fatal difference. This system handles " 1615 << this->
n_systems() <<
" systems," << std::endl
1616 <<
" while the other system handles " 1618 <<
" systems." << std::endl
1619 <<
" Aborting comparison." << std::endl;
1626 for (
const auto & [sys_name, sys_ptr] :
_systems)
1631 os_result.push_back (sys_ptr->compare (other_system, threshold, verbose));
1639 if (os_result.size()==0)
1647 os_identical = os_result[n];
1650 while (os_identical && n<os_result.size());
1651 return os_identical;
1659 std::ostringstream oss;
1661 oss <<
" EquationSystems\n" 1662 <<
" n_systems()=" << this->
n_systems() <<
'\n';
1666 oss << pr.second->get_info();
1712 tot += pr.second->n_vars();
1724 tot += pr.second->n_dofs();
1737 tot += pr.second->n_active_dofs();
1746 for (
auto & node :
_mesh.node_ptr_range())
1750 for (
auto & elem :
_mesh.element_ptr_range())
1756 this->
get_system(sys_num).get_dof_map().remove_default_ghosting();
EquationSystems(MeshBase &mesh)
Constructor.
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
This is the EquationSystems class.
unsigned int n_vars() const
void build_variable_names(std::vector< std::string > &var_names, const FEType *type=nullptr, const std::set< std::string > *system_names=nullptr) const
Fill the input vector var_names with the names of the variables for each system.
A Node is like a Point, but with more information.
unsigned int n_systems() const
const Variable & variable(unsigned int var) const
Return a constant reference to Variable var.
bool _enable_default_ghosting
Flag for whether to enable default ghosting on newly added Systems.
virtual void get(const std::vector< numeric_index_type > &index, T *values) const
Access multiple components at once.
Data structure for specifying which Parameters should be independent variables in a parameter sensiti...
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
std::unique_ptr< NumericVector< Number > > build_parallel_elemental_solution_vector(std::vector< std::string > &names) const
Builds a parallel vector of CONSTANT MONOMIAL and/or components of CONSTANT MONOMIAL_VEC solution val...
std::size_t distribute_dofs(MeshBase &)
Distribute dofs on the current mesh.
static constexpr Real TOLERANCE
std::size_t n_dofs() const
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
void build_discontinuous_solution_vector(std::vector< Number > &soln, const std::set< std::string > *system_names=nullptr, const std::vector< std::string > *var_names=nullptr, bool vertices_only=false, bool add_sides=false) const
Fill the input vector soln with solution values.
void add_default_ghosting()
Add the default functor(s) for coupling and algebraic ghosting.
bool refine_elements()
Only refines the user-requested elements.
bool reinit_solutions()
Handle any mesh changes and project any solutions onto the updated mesh.
virtual void clear()
Restores the data structure to a pristine state.
virtual void allgather()
Gathers all elements and nodes of the mesh onto every processor.
virtual ~EquationSystems()
Destructor.
This is the base class from which all geometric element types are derived.
void print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
static FEFieldType field_type(const FEType &fe_type)
std::ostream & operator<<(std::ostream &os, const OrderWrapper &order)
Overload stream operators.
processor_id_type rank() const
const Parallel::Communicator & comm() const
virtual void enable_default_ghosting(bool enable)
Enable or disable default ghosting functors on the Mesh and on all Systems.
void clean_refinement_flags()
Sets the refinement flag to Elem::DO_NOTHING for each element in the mesh.
void should_p_refine(unsigned int g, bool p_refine)
Describe whether the given variable group should be p-refined.
virtual void adjoint_solve(const QoISet &qoi_indices=QoISet())
Call adjoint_solve on all the individual equation systems.
const Parallel::Communicator & _communicator
The libMesh namespace provides an interface to certain functionality in the library.
const T_sys & get_system(std::string_view name) const
Real distance(const Point &p)
virtual void init(const numeric_index_type n, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC)=0
Change the dimension of the vector to n.
void iota(ForwardIter first, ForwardIter last, T value)
Utility::iota was created back when std::iota was just an SGI STL extension.
static void side_nodal_soln(const FEType &fe_t, const Elem *elem, const unsigned int side, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln on one side from the (full) element soln.
uint8_t processor_id_type
This is the MeshBase class.
virtual std::string get_info() const
virtual void prolong_vectors()
Prolong vectors after the mesh has refined.
MeshBase & _mesh
The mesh data structure.
bool coarsen_elements()
Only coarsens the user-requested elements.
unsigned int var_group_from_var_number(unsigned int var_num) const
signed char & underrefined_boundary_limit()
If underrefined_boundary_limit is set to a nonnegative value, then refinement and coarsening will pro...
unsigned int n_dofs(const unsigned int s, const unsigned int var=libMesh::invalid_uint) const
This class handles the numbering of degrees of freedom on a mesh.
void allgather()
Serializes a distributed mesh and its associated degree of freedom numbering for all systems...
void get_vars_active_subdomains(const std::vector< std::string > &names, std::vector< std::set< subdomain_id_type >> &vars_active_subdomains) const
Retrieve vars_active_subdomains, which indicates the active subdomains for each variable in names...
processor_id_type size() const
virtual void sensitivity_solve(const ParameterVector ¶meters)
Call sensitivity_solve on all the individual equation systems.
virtual void reinit_systems()
Reinitialize all systems on the current mesh.
processor_id_type n_processors() const
virtual bool is_serial() const
unsigned int number() const
This class defines the notion of a variable in the system.
const std::set< subdomain_id_type > & active_subdomains() const
Implements (adaptive) mesh refinement algorithms for a MeshBase.
dof_id_type numeric_index_type
unsigned int which_neighbor_am_i(const Elem *e) const
This function tells you which neighbor e is.
Manages consistently variables, degrees of freedom, and coefficient vectors.
virtual void reinit()
Handle any mesh changes and reinitialize all the systems on the updated mesh.
virtual dof_id_type max_elem_id() const =0
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
const std::string & variable_name(const unsigned int i) const
void build_solution_vector(std::vector< Number > &soln, std::string_view system_name, std::string_view variable_name="all_vars") const
Fill the input vector soln with the solution values for the system named name.
virtual void restrict_vectors()
Restrict vectors after the mesh has coarsened.
bool active_on_subdomain(subdomain_id_type sid) const
bool _refine_in_reinit
Flag for whether to call coarsen/refine in reinit().
virtual void reinit_mesh()
Handle the association of a completely new mesh with the EquationSystem and all the Systems assigned ...
unsigned char & face_level_mismatch_limit()
If face_level_mismatch_limit is set to a nonzero value, then refinement and coarsening will produce m...
An object whose state is distributed along a set of processors.
static unsigned int n_vec_dim(const MeshBase &mesh, const FEType &fe_type)
virtual void reinit_constraints()
Reinitializes the constraints for this system.
virtual void close()=0
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
static FEContinuity get_continuity(const FEType &fe_type)
Returns the input FEType's FEContinuity based on the underlying FEFamily and potentially the Order...
virtual void solve()
Call solve on all the individual equation systems.
virtual numeric_index_type first_local_index() const =0
const Elem * neighbor_ptr(unsigned int i) const
void get_solution(std::vector< Number > &soln, std::vector< std::string > &names) const
Retrieve the solution data for CONSTANT MONOMIALs and/or components of CONSTANT MONOMIAL_VECs.
virtual void update()
Update the local values to reflect the solution on neighboring processors.
unsigned int level() const
const FEType & variable_type(const unsigned int i) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool contract()=0
Delete subactive (i.e.
subdomain_id_type subdomain_id() const
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
const Node * node_ptr(const unsigned int i) const
unsigned int spatial_dimension() const
T & set(const std::string &)
std::unique_ptr< NumericVector< Number > > build_parallel_solution_vector(const std::set< std::string > *system_names=nullptr, bool add_sides=false) const
A version of build_solution_vector which is appropriate for "parallel" output formats like Nemesis...
virtual void clear()
Clears internal data structures & frees any allocated memory.
signed char & overrefined_boundary_limit()
If overrefined_boundary_limit is set to a nonnegative value, then refinement and coarsening will prod...
void remove_default_ghosting()
Remove any default ghosting functor(s).
void _remove_default_ghosting(unsigned int sys_num)
This just calls DofMap::remove_default_ghosting() but using a shim lets us forward-declare DofMap...
const MeshBase & get_mesh() const
std::size_t n_active_dofs() const
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
static std::unique_ptr< NumericVector< T > > build(const Parallel::Communicator &comm, SolverPackage solver_package=libMesh::default_solver_package(), ParallelType parallel_type=AUTOMATIC)
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...
std::vector< std::pair< unsigned int, unsigned int > > find_variable_numbers(std::vector< std::string > &names, const FEType *type=nullptr, const std::vector< FEType > *types=nullptr) const
Finds system and variable numbers for any variables of 'type' or of 'types' corresponding to the entr...
Parameters parameters
Data structure holding arbitrary parameters.
std::unique_ptr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
virtual void init()
Initialize all the systems.
virtual bool compare(const EquationSystems &other_es, const Real threshold, const bool verbose) const
unsigned int n_vars() const
virtual dof_id_type max_node_id() const =0
virtual dof_id_type n_elem() const =0
virtual void add(const numeric_index_type i, const T value)=0
Adds value to the vector entry specified by i.
virtual System & add_system(std::string_view system_type, std::string_view name)
Add the system of type system_type named name to the systems array.
processor_id_type processor_id() const
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln from the element soln.
void _add_system_to_nodes_and_elems()
This function is used in the implementation of add_system, it loops over the nodes and elements of th...
const DofMap & get_dof_map() const
virtual numeric_index_type last_local_index() const =0
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
void update()
Updates local values for all the systems.
virtual std::vector< unsigned int > nodes_on_side(const unsigned int) const =0
void build_elemental_solution_vector(std::vector< Number > &soln, std::vector< std::string > &names) const
Retrieve the solution data for CONSTANT MONOMIALs and/or components of CONSTANT MONOMIAL_VECs.
std::map< std::string, std::unique_ptr< System >, std::less<> > _systems
Data structure holding the systems.
static bool redundant_added_side(const Elem &elem, unsigned int side)
const RemoteElem * remote_elem