Skip to content

Instantly share code, notes, and snippets.

@jwpeterson
Created August 6, 2013 20:52
Show Gist options
  • Save jwpeterson/6168524 to your computer and use it in GitHub Desktop.
Save jwpeterson/6168524 to your computer and use it in GitHub Desktop.
Sample code that demonstrates write_element_data not working with ParallelMesh.
// This utility lists all nodes on a given boundary ID.
#include "libmesh/libmesh.h"
#include "libmesh/serial_mesh.h"
#include "libmesh/parallel_mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/equation_systems.h"
#include "libmesh/explicit_system.h"
#include "libmesh/dof_map.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/exodusII_io.h"
int main (int argc, char** argv)
{
LibMeshInit init (argc, argv);
{
// Switch back and forth to see the difference in output... should work
// fine with SerialMesh and fail with ParallelMesh.
ParallelMesh mesh;
// SerialMesh mesh;
MeshTools::Generation::build_square (mesh,
/*nx*/10,
/*ny*/100,
/*xmin*/0., /*xmax*/1.,
/*ymin*/0., /*ymax*/10.,
QUAD4);
EquationSystems equation_systems (mesh);
// Create an explicit system for writing out mesh "data"
ExplicitSystem& system =
equation_systems.add_system<ExplicitSystem> ("Dummy Data");
// Default is linear Lagrange
system.add_variable("proc_id_nodal");
// Cell-based data
system.add_variable("proc_id_elemental", FEType(CONSTANT, MONOMIAL));
// Initialize the data structures for the equation system.
equation_systems.init();
// Print information about the system to the screen.
mesh.print_info();
equation_systems.print_info();
// Now, ready to set "solution" values
{
NumericVector<Number>& solution = *system.solution;
// Vectors to store various dof indices
std::vector<unsigned int> dof_indices_proc_id_nodal, dof_indices_proc_id_elemental;
// Get variable numbers for the variables on the system
const unsigned int proc_id_nodal_var_num = system.variable_number ("proc_id_nodal");
const unsigned int proc_id_elemental_var_num = system.variable_number ("proc_id_elemental");
const DofMap & dof_map = system.get_dof_map();
MeshBase::const_element_iterator el = mesh.active_local_elements_begin();
const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
for ( ; el != end_el; ++el)
{
const Elem* elem = *el;
// Get dof indices for the current element
dof_map.dof_indices (elem, dof_indices_proc_id_nodal, proc_id_nodal_var_num);
dof_map.dof_indices (elem, dof_indices_proc_id_elemental, proc_id_elemental_var_num);
// Loop over each dof index. If it's local to this processor, set the value.
// For linear LAGRANGE, there should be 1 DOF per node.
for (unsigned int i=0; i<dof_indices_proc_id_nodal.size(); i++)
if ((dof_indices_proc_id_nodal[i] >= solution.first_local_index()) &&
(dof_indices_proc_id_nodal[i] < solution.last_local_index()))
{
solution.set (dof_indices_proc_id_nodal[i], elem->get_node(i)->processor_id());
} // end if dof index is local
// Maybe assert this instead...
if (dof_indices_proc_id_elemental.size() > 1)
libmesh_error();
// Set the one MONOMIAL degree of freedom.
solution.set(dof_indices_proc_id_elemental[0], elem->processor_id());
} // end loop over elements
// system.update() doesn't automatically close the solution vector (did it do that before?)
// so we do it ourselves.
solution.close();
system.update();
} // done setting solution values
// Write the equation systems to file
ExodusII_IO writer(mesh);
// The Exodus writer writes constant monomials as nodal data
// unless you provide an explicit list of variables to include in
// the output by calling set_output_variables().
std::vector<std::string> output_vars(1, "proc_id_nodal");
writer.set_output_variables(output_vars);
writer.write_equation_systems("write_element_data_test.e", equation_systems);
// Remarks:
// * Can only call _after_ calling write_timestep() or write_equation_systems()
// * Doesn't take a filename.
// * Seems to be buggy with ParallelMesh enabled.
writer.write_element_data(equation_systems);
}
return 0;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment