20 #include "libmesh/distributed_mesh.h" 21 #include "libmesh/equation_systems.h" 22 #include "libmesh/libmesh_logging.h" 23 #include "libmesh/mesh_output.h" 24 #include "libmesh/unstructured_mesh.h" 25 #include "libmesh/numeric_vector.h" 33 const std::set<std::string> * system_names)
35 LOG_SCOPE(
"write_equation_systems()",
"MeshOutput");
39 MT & my_mesh =
const_cast<MT &
>(*_obj);
43 libmesh_assert_equal_to(&es.
get_mesh(), _obj);
47 if (!_is_parallel_format &&
48 (my_mesh.max_elem_id() != my_mesh.n_elem() ||
49 my_mesh.max_node_id() != my_mesh.n_nodes()))
56 "Warning: This MeshOutput subclass only supports meshes which are contiguously renumbered!" 59 my_mesh.allow_renumbering(
true);
61 my_mesh.renumber_nodes_and_elements();
65 my_mesh.allow_renumbering(
false);
68 if (!_is_parallel_format)
70 MeshSerializer serialize(const_cast<MT &>(*_obj), !_is_parallel_format, _serial_only_needed_on_proc_0);
73 std::vector<std::string> names;
78 std::vector<Number> soln;
80 this->get_add_sides());
82 this->write_nodal_data (fname, soln, names);
85 this->write_nodal_data (fname, es, system_names);
91 const std::set<std::string> * system_names)
93 LOG_SCOPE(
"write_discontinuous_equation_systems()",
"MeshOutput");
97 MT & my_mesh =
const_cast<MT &
>(*_obj);
101 libmesh_assert_equal_to(&es.
get_mesh(), _obj);
105 if (my_mesh.max_elem_id() != my_mesh.n_elem() ||
106 my_mesh.max_node_id() != my_mesh.n_nodes())
113 "Warning: This MeshOutput subclass only supports meshes which are contiguously renumbered!" 116 my_mesh.allow_renumbering(
true);
118 my_mesh.renumber_nodes_and_elements();
122 my_mesh.allow_renumbering(
false);
125 MeshSerializer serialize(const_cast<MT &>(*_obj), !_is_parallel_format, _serial_only_needed_on_proc_0);
128 std::vector<std::string> names;
131 if (!_is_parallel_format)
135 std::vector<Number> soln;
138 this->get_add_sides());
140 this->write_nodal_data_discontinuous (fname, soln, names);
144 libmesh_not_implemented();
151 const std::vector<std::string> & names)
156 std::vector<Number> soln;
158 this->write_nodal_data(fname, soln, names);
164 const std::set<std::string> * system_names)
166 std::vector<std::string> names;
169 std::unique_ptr<NumericVector<Number>> parallel_soln =
172 this->write_nodal_data (fname, *parallel_soln, names);
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
This is the EquationSystems class.
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.
virtual void write_equation_systems(const std::string &, const EquationSystems &, const std::set< std::string > *system_names=nullptr)
This method implements writing a mesh with data to a specified file where the data is taken from the ...
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.
The libMesh namespace provides an interface to certain functionality in the library.
virtual void write_discontinuous_equation_systems(const std::string &, const EquationSystems &, const std::set< std::string > *system_names=nullptr)
This method implements writing a mesh with discontinuous data to a specified file where the data is t...
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.
Temporarily serialize a DistributedMesh for non-distributed-mesh capable code paths.
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...
const MeshBase & get_mesh() const
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.