libMesh::EquationSystems Class Reference

Manages multiples systems of equations. More...

#include <equation_systems.h>

Inheritance diagram for libMesh::EquationSystems:

Public Types

enum  ReadFlags {
  READ_HEADER = 1, READ_DATA = 2, READ_ADDITIONAL_DATA = 4, READ_LEGACY_FORMAT = 8,
  TRY_READ_IFEMS = 16, READ_BASIC_ONLY = 32
}
 
enum  WriteFlags { WRITE_DATA = 1, WRITE_ADDITIONAL_DATA = 2, WRITE_PARALLEL_FILES = 4, WRITE_SERIAL_FILES = 8 }
 

Public Member Functions

 EquationSystems (MeshBase &mesh)
 
virtual ~EquationSystems ()
 
virtual void clear ()
 
virtual void init ()
 
virtual void reinit ()
 
void update ()
 
unsigned int n_systems () const
 
bool has_system (const std::string &name) const
 
template<typename T_sys >
const T_sys & get_system (const std::string &name) const
 
template<typename T_sys >
T_sys & get_system (const std::string &name)
 
template<typename T_sys >
const T_sys & get_system (const unsigned int num) const
 
template<typename T_sys >
T_sys & get_system (const unsigned int num)
 
const Systemget_system (const std::string &name) const
 
Systemget_system (const std::string &name)
 
const Systemget_system (const unsigned int num) const
 
Systemget_system (const unsigned int num)
 
virtual Systemadd_system (const std::string &system_type, const std::string &name)
 
template<typename T_sys >
T_sys & add_system (const std::string &name)
 
void delete_system (const std::string &name)
 
unsigned int n_vars () const
 
std::size_t n_dofs () const
 
std::size_t n_active_dofs () const
 
virtual void solve ()
 
virtual void adjoint_solve (const QoISet &qoi_indices=QoISet())
 
virtual void sensitivity_solve (const ParameterVector &parameters)
 
void build_variable_names (std::vector< std::string > &var_names, const FEType *type=libmesh_nullptr, const std::set< std::string > *system_names=libmesh_nullptr) const
 
void build_solution_vector (std::vector< Number > &soln, const std::string &system_name, const std::string &variable_name="all_vars") const
 
void build_solution_vector (std::vector< Number > &soln, const std::set< std::string > *system_names=libmesh_nullptr) const
 
UniquePtr< NumericVector< Number > > build_parallel_solution_vector (const std::set< std::string > *system_names=libmesh_nullptr) const
 
void get_solution (std::vector< Number > &soln, std::vector< std::string > &names) const
 
void build_discontinuous_solution_vector (std::vector< Number > &soln, const std::set< std::string > *system_names=libmesh_nullptr) const
 
template<typename InValType >
void read (const std::string &name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
 
void read (const std::string &name, const XdrMODE mode, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
 
template<typename InValType >
void read (const std::string &name, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
 
void read (const std::string &name, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
 
void write (const std::string &name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
 
void write (const std::string &name, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
 
virtual bool compare (const EquationSystems &other_es, const Real threshold, const bool verbose) const
 
virtual std::string get_info () const
 
void print_info (std::ostream &os=libMesh::out) const
 
const MeshBaseget_mesh () const
 
MeshBaseget_mesh ()
 
void allgather ()
 
void enable_refine_in_reinit ()
 
void disable_refine_in_reinit ()
 
bool refine_in_reinit_flag ()
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Public Attributes

Parameters parameters
 

Protected Types

typedef std::map< std::string, System * >::iterator system_iterator
 
typedef std::map< std::string, System * >::const_iterator const_system_iterator
 
typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 

Protected Member Functions

void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

MeshBase_mesh
 
std::map< std::string, System * > _systems
 
bool _refine_in_reinit
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic< unsigned int > _n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

Private Member Functions

template<typename InValType >
void _read_impl (const std::string &name, const XdrMODE, const unsigned int read_flags, bool partition_agnostic=true)
 
void _add_system_to_nodes_and_elems ()
 

Friends

std::ostream & operator<< (std::ostream &os, const EquationSystems &es)
 

Detailed Description

Manages multiples systems of equations.

This is the EquationSystems class. It is in charge of handling all the various equation systems defined for a MeshBase. It may have multiple systems, which may be active or inactive, so that at different solution stages only a sub-set may be solved for. Also, through the templated access, different types of systems may be handled. Also other features, like flags, parameters, I/O etc are provided.

Author
Benjamin S. Kirk
Date
2002-2007

Definition at line 66 of file equation_systems.h.

Member Typedef Documentation

typedef std::map<std::string, System *>::const_iterator libMesh::EquationSystems::const_system_iterator
protected

Typedef for constatnt system iterators

Definition at line 486 of file equation_systems.h.

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information. The log is identified by the class name.

Definition at line 110 of file reference_counter.h.

typedef std::map<std::string, System *>::iterator libMesh::EquationSystems::system_iterator
protected

Typedef for system iterators

Definition at line 481 of file equation_systems.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

libMesh::EquationSystems::EquationSystems ( MeshBase mesh)

Constructor.

Definition at line 54 of file equation_systems.C.

References _refine_in_reinit, parameters, libMesh::Real, libMesh::Parameters::set(), and libMesh::TOLERANCE.

54  :
55  ParallelObject (m),
56  _mesh (m)
57 {
58  // Set default parameters
59  this->parameters.set<Real> ("linear solver tolerance") = TOLERANCE * TOLERANCE;
60  this->parameters.set<unsigned int>("linear solver maximum iterations") = 5000;
61  this->_refine_in_reinit = true; // default value
62 }
ParallelObject(const Parallel::Communicator &comm_in)
static const Real TOLERANCE
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
T & set(const std::string &)
Definition: parameters.h:469
libMesh::EquationSystems::~EquationSystems ( )
virtual

Destructor. Should be virtual, since the user may want to derive subclasses of EquationSystems.

Definition at line 66 of file equation_systems.C.

References clear().

67 {
68  this->clear ();
69 }

Member Function Documentation

void libMesh::EquationSystems::_add_system_to_nodes_and_elems ( )
private

This function is used in the implementation of add_system, it loops over the nodes and elements of the Mesh, adding the system to each one. The main reason to separate this part is to avoid coupling this header file to mesh.h, and elem.h.

Definition at line 1337 of file equation_systems.C.

References _mesh, libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), libMesh::MeshBase::nodes_begin(), and libMesh::MeshBase::nodes_end().

Referenced by add_system().

1338 {
1339  // All the nodes
1340  MeshBase::node_iterator node_it = _mesh.nodes_begin();
1341  const MeshBase::node_iterator node_end = _mesh.nodes_end();
1342 
1343  for ( ; node_it != node_end; ++node_it)
1344  (*node_it)->add_system();
1345 
1346  // All the elements
1347  MeshBase::element_iterator elem_it = _mesh.elements_begin();
1348  const MeshBase::element_iterator elem_end = _mesh.elements_end();
1349 
1350  for ( ; elem_it != elem_end; ++elem_it)
1351  (*elem_it)->add_system();
1352 }
virtual node_iterator nodes_begin()=0
virtual element_iterator elements_begin()=0
virtual element_iterator elements_end()=0
virtual node_iterator nodes_end()=0
template<typename InValType >
template void libMesh::EquationSystems::_read_impl< Real > ( const std::string &  name,
const XdrMODE  ,
const unsigned int  read_flags,
bool  partition_agnostic = true 
)
private

Actual read implementation. This can be called repeatedly inside a try-catch block in an attempt to read broken files.

Parameters
nameName of the file to be read.
read_flagsSingle flag created by bitwise-OR'ing several flags together.
partition_agnosticIf true then the mesh and degrees of freedom will be temporarily renumbered in a partition agnostic way so that files written using "n" mpi processes can be re-read on "m" mpi processes. Note that this renumbering is not compatible with meshes that have two nodes in exactly the same position!

This program implements the output of an EquationSystems object. This warrants some documentation. The output file essentially consists of 11 sections:

1.) A version header (for non-'legacy' formats, libMesh-0.7.0 and greater).
2.) The number of individual equation systems (unsigned int)

for each system

3.)  The name of the system (string)
4.)  The type of the system (string)

handled through System::read():

+-------------------------------------------------------------+
|  5.) The number of variables in the system (unsigned int)   |
|                                                             |
|   for each variable in the system                           |
|                                                             |
|    6.) The name of the variable (string)                    |
|                                                             |
|    7.) Combined in an FEType:                               |
|         - The approximation order(s) of the variable (Order |
|           Enum, cast to int/s)                              |
|         - The finite element family/ies of the variable     |
|           (FEFamily Enum, cast to int/s)                    |
|                                                             |
|   end variable loop                                         |
|                                                             |
| 8.) The number of additional vectors (unsigned int),        |
|                                                             |
|    for each additional vector in the equation system object |
|                                                             |
|    9.) the name of the additional vector  (string)          |
+-------------------------------------------------------------+

end system loop


for each system, handled through System::read_{serialized,parallel}_data():

+--------------------------------------------------------------+
| 10.) The global solution vector, re-ordered to be node-major |
|     (More on this later.)                                    |
|                                                              |
|    for each additional vector in the equation system object  |
|                                                              |
|    11.) The global additional vector, re-ordered to be       |
|         node-major (More on this later.)                     |
+--------------------------------------------------------------+

end system loop

Note that the actual IO is handled through the Xdr class (to be renamed later?) which provides a uniform interface to both the XDR (eXternal Data Representation) interface and standard ASCII output. Thus this one section of code will read XDR or ASCII files with no changes.

Definition at line 143 of file equation_systems_io.C.

References _mesh, add_system(), libMesh::Parallel::Communicator::broadcast(), libMesh::Xdr::close(), libMesh::ParallelObject::comm(), libMesh::Xdr::data(), libMesh::MeshBase::fix_broken_node_and_element_numbering(), get_mesh(), get_system(), libMesh::MeshTools::Private::globally_renumber_nodes_and_elements(), init(), libMesh::libmesh_assert(), mesh, libMesh::ParallelObject::processor_id(), read(), READ_ADDITIONAL_DATA, READ_BASIC_ONLY, READ_DATA, READ_HEADER, libMesh::System::read_header(), READ_LEGACY_FORMAT, libMesh::Xdr::reading(), libMesh::System::set_basic_system_only(), libMesh::Xdr::set_version(), libMesh::sys, TRY_READ_IFEMS, and update().

147 {
211  // Set booleans from the read_flags argument
212  const bool read_header = read_flags & EquationSystems::READ_HEADER;
213  const bool read_data = read_flags & EquationSystems::READ_DATA;
214  const bool read_additional_data = read_flags & EquationSystems::READ_ADDITIONAL_DATA;
215  const bool read_legacy_format = read_flags & EquationSystems::READ_LEGACY_FORMAT;
216  const bool try_read_ifems = read_flags & EquationSystems::TRY_READ_IFEMS;
217  const bool read_basic_only = read_flags & EquationSystems::READ_BASIC_ONLY;
218  bool read_parallel_files = false;
219 
220  std::vector<std::pair<std::string, System *> > xda_systems;
221 
222  // This will unzip a file with .bz2 as the extension, otherwise it
223  // simply returns the name if the file need not be unzipped.
224  Xdr io ((this->processor_id() == 0) ? name : "", mode);
225  libmesh_assert (io.reading());
226 
227  {
228  // 1.)
229  // Read the version header.
230  std::string version = "legacy";
231  if (!read_legacy_format)
232  {
233  if (this->processor_id() == 0) io.data(version);
234  this->comm().broadcast(version);
235 
236  // All processors have the version header, if it does not contain
237  // "libMesh" something then it is a legacy file.
238  std::string::size_type lm_pos = version.find("libMesh");
239  if (!(lm_pos < version.size()))
240  {
241  io.close();
242 
243  // Recursively call this read() function but with the
244  // EquationSystems::READ_LEGACY_FORMAT bit set.
245  this->read (name, mode, (read_flags | EquationSystems::READ_LEGACY_FORMAT), partition_agnostic);
246  return;
247  }
248 
249  // Figure out the libMesh version that created this file
250  std::istringstream iss(version.substr(lm_pos + 8));
251  int ver_major = 0, ver_minor = 0, ver_patch = 0;
252  char dot;
253  iss >> ver_major >> dot >> ver_minor >> dot >> ver_patch;
254  io.set_version(LIBMESH_VERSION_ID(ver_major, ver_minor, ver_patch));
255 
256 
257  read_parallel_files = (version.rfind(" parallel") < version.size());
258 
259  // If requested that we try to read infinite element information,
260  // and the string " with infinite elements" is not in the version,
261  // then tack it on. This is for compatibility reading ifem
262  // files written prior to 11/10/2008 - BSK
263  if (try_read_ifems)
264  if (!(version.rfind(" with infinite elements") < version.size()))
265  version += " with infinite elements";
266 
267  }
268  else
269  libmesh_deprecated();
270 
271  START_LOG("read()","EquationSystems");
272 
273  // 2.)
274  // Read the number of equation systems
275  unsigned int n_sys=0;
276  if (this->processor_id() == 0) io.data (n_sys);
277  this->comm().broadcast(n_sys);
278 
279  for (unsigned int sys=0; sys<n_sys; sys++)
280  {
281  // 3.)
282  // Read the name of the sys-th equation system
283  std::string sys_name;
284  if (this->processor_id() == 0) io.data (sys_name);
285  this->comm().broadcast(sys_name);
286 
287  // 4.)
288  // Read the type of the sys-th equation system
289  std::string sys_type;
290  if (this->processor_id() == 0) io.data (sys_type);
291  this->comm().broadcast(sys_type);
292 
293  if (read_header)
294  this->add_system (sys_type, sys_name);
295 
296  // 5.) - 9.)
297  // Let System::read_header() do the job
298  System & new_system = this->get_system(sys_name);
299  new_system.read_header (io,
300  version,
301  read_header,
302  read_additional_data,
303  read_legacy_format);
304 
305  xda_systems.push_back(std::make_pair(sys_name, &new_system));
306 
307  // If we're only creating "basic" systems, we need to tell
308  // each system that before we call init() later.
309  if (read_basic_only)
310  new_system.set_basic_system_only();
311  }
312  }
313 
314 
315 
316  // Now we are ready to initialize the underlying data
317  // structures. This will initialize the vectors for
318  // storage, the dof_map, etc...
319  if (read_header)
320  this->init();
321 
322  // 10.) & 11.)
323  // Read and set the numeric vector values
324  if (read_data)
325  {
326  // the EquationSystems::read() method should look constant from the mesh
327  // perspective, but we need to assign a temporary numbering to the nodes
328  // and elements in the mesh, which requires that we abuse const_cast
329  if (!read_legacy_format && partition_agnostic)
330  {
331  MeshBase & mesh = const_cast<MeshBase &>(this->get_mesh());
333  }
334 
335  Xdr local_io (read_parallel_files ? local_file_name(this->processor_id(),name) : "", mode);
336 
337  std::vector<std::pair<std::string, System *> >::iterator
338  pos = xda_systems.begin();
339 
340  for (; pos != xda_systems.end(); ++pos)
341  if (read_legacy_format)
342  {
343  libmesh_deprecated();
344  pos->second->read_legacy_data (io, read_additional_data);
345  }
346  else
347  if (read_parallel_files)
348  pos->second->read_parallel_data<InValType> (local_io, read_additional_data);
349  else
350  pos->second->read_serialized_data<InValType> (io, read_additional_data);
351 
352 
353  // Undo the temporary numbering.
354  if (!read_legacy_format && partition_agnostic)
356  }
357 
358  STOP_LOG("read()","EquationSystems");
359 
360  // Localize each system's data
361  this->update();
362 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
ImplicitSystem & sys
MeshBase & mesh
virtual void fix_broken_node_and_element_numbering()=0
void read(const std::string &name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
libmesh_assert(j)
virtual System & add_system(const std::string &system_type, const std::string &name)
void broadcast(T &data, const unsigned int root_id=0) const
void globally_renumber_nodes_and_elements(MeshBase &)
Definition: mesh_tools.C:1948
const Parallel::Communicator & comm() const
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
processor_id_type processor_id() const
System & libMesh::EquationSystems::add_system ( const std::string &  system_type,
const std::string &  name 
)
virtual

Add the system of type system_type named name to the systems array.

Definition at line 349 of file equation_systems.C.

References _systems, get_system(), and libMesh::Quality::name().

Referenced by _read_impl(), and libMesh::ErrorVector::plot_error().

351 {
352  // If the user already built a system with this name, we'll
353  // trust them and we'll use it. That way they can pre-add
354  // non-standard derived system classes, and if their restart file
355  // has some non-standard sys_type we won't throw an error.
356  if (_systems.count(name))
357  {
358  return this->get_system(name);
359  }
360  // Build a basic System
361  else if (sys_type == "Basic")
362  this->add_system<System> (name);
363 
364  // Build a Newmark system
365  else if (sys_type == "Newmark")
366  this->add_system<NewmarkSystem> (name);
367 
368  // Build an Explicit system
369  else if ((sys_type == "Explicit"))
370  this->add_system<ExplicitSystem> (name);
371 
372  // Build an Implicit system
373  else if ((sys_type == "Implicit") ||
374  (sys_type == "Steady" ))
375  this->add_system<ImplicitSystem> (name);
376 
377  // build a transient implicit linear system
378  else if ((sys_type == "Transient") ||
379  (sys_type == "TransientImplicit") ||
380  (sys_type == "TransientLinearImplicit"))
381  this->add_system<TransientLinearImplicitSystem> (name);
382 
383  // build a transient implicit nonlinear system
384  else if (sys_type == "TransientNonlinearImplicit")
385  this->add_system<TransientNonlinearImplicitSystem> (name);
386 
387  // build a transient explicit system
388  else if (sys_type == "TransientExplicit")
389  this->add_system<TransientExplicitSystem> (name);
390 
391  // build a linear implicit system
392  else if (sys_type == "LinearImplicit")
393  this->add_system<LinearImplicitSystem> (name);
394 
395  // build a nonlinear implicit system
396  else if (sys_type == "NonlinearImplicit")
397  this->add_system<NonlinearImplicitSystem> (name);
398 
399  // build a Reduced Basis Construction system
400  else if (sys_type == "RBConstruction")
401  this->add_system<RBConstruction> (name);
402 
403  // build a transient Reduced Basis Construction system
404  else if (sys_type == "TransientRBConstruction")
405  this->add_system<TransientRBConstruction> (name);
406 
407 #ifdef LIBMESH_HAVE_SLEPC
408  // build an eigen system
409  else if (sys_type == "Eigen")
410  this->add_system<EigenSystem> (name);
411  else if (sys_type == "TransientEigenSystem")
412  this->add_system<TransientEigenSystem> (name);
413 #endif
414 
415 #if defined(LIBMESH_USE_COMPLEX_NUMBERS)
416  // build a frequency system
417  else if (sys_type == "Frequency")
418  this->add_system<FrequencySystem> (name);
419 #endif
420 
421  else
422  libmesh_error_msg("ERROR: Unknown system type: " << sys_type);
423 
424  // Return a reference to the new system
425  //return (*this)(name);
426  return this->get_system(name);
427 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
const T_sys & get_system(const std::string &name) const
std::map< std::string, System * > _systems
template<typename T_sys >
T_sys & libMesh::EquationSystems::add_system ( const std::string &  name)
inline

Add the system named name to the systems array.

Definition at line 553 of file equation_systems.h.

References _add_system_to_nodes_and_elems(), _systems, libmesh_nullptr, n_systems(), and libMesh::Quality::name().

554 {
555  T_sys * ptr = libmesh_nullptr;
556 
557  if (!_systems.count(name))
558  {
559  ptr = new T_sys(*this, name, this->n_systems());
560 
561  _systems.insert (std::make_pair(name, ptr));
562 
563  // Tell all the \p DofObject entities to add a system.
565  }
566  else
567  {
568  // We now allow redundant add_system calls, to make it
569  // easier to load data from files for user-derived system
570  // subclasses
571  ptr = &(this->get_system<T_sys>(name));
572  }
573 
574  // Return a dynamically casted reference to the newly added System.
575  return *ptr;
576 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
const class libmesh_nullptr_t libmesh_nullptr
unsigned int n_systems() const
std::map< std::string, System * > _systems
void libMesh::EquationSystems::adjoint_solve ( const QoISet qoi_indices = QoISet())
virtual

Call adjoint_solve on all the individual equation systems.

By default this function solves each system's adjoint once, in the reverse order from that in which they were added. For more sophisticated decoupled problems the user may with to override this behavior in a derived class.

Definition at line 468 of file equation_systems.C.

References get_system(), libMesh::libmesh_assert(), and n_systems().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error().

469 {
470  libmesh_assert (this->n_systems());
471 
472  for (unsigned int i=this->n_systems(); i != 0; --i)
473  this->get_system(i-1).adjoint_solve(qoi_indices);
474 }
libmesh_assert(j)
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::allgather ( )

Serializes a distributed mesh and its associated degree of freedom numbering for all systems

Definition at line 291 of file equation_systems.C.

References _mesh, libMesh::MeshBase::allgather(), libMesh::DofMap::distribute_dofs(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), libMesh::System::get_dof_map(), get_system(), libMesh::MeshBase::is_serial(), n_systems(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), libMesh::DofMap::prepare_send_list(), libMesh::System::reinit_constraints(), and libMesh::sys.

Referenced by read().

292 {
293  // A serial mesh means nothing needs to be done
294  if (_mesh.is_serial())
295  return;
296 
297  const unsigned int n_sys = this->n_systems();
298 
299  libmesh_assert_not_equal_to (n_sys, 0);
300 
301  // Gather the mesh
302  _mesh.allgather();
303 
304  // Tell all the \p DofObject entities how many systems
305  // there are.
306  {
307  MeshBase::node_iterator node_it = _mesh.nodes_begin();
308  const MeshBase::node_iterator node_end = _mesh.nodes_end();
309 
310  for ( ; node_it != node_end; ++node_it)
311  (*node_it)->set_n_systems(n_sys);
312 
313  MeshBase::element_iterator elem_it = _mesh.elements_begin();
314  const MeshBase::element_iterator elem_end = _mesh.elements_end();
315 
316  for ( ; elem_it != elem_end; ++elem_it)
317  (*elem_it)->set_n_systems(n_sys);
318  }
319 
320  // And distribute each system's dofs
321  for (unsigned int i=0; i != this->n_systems(); ++i)
322  {
323  System & sys = this->get_system(i);
324  DofMap & dof_map = sys.get_dof_map();
325  dof_map.distribute_dofs(_mesh);
326 
327  // The user probably won't need constraint equations or the
328  // send_list after an allgather, but let's keep it in consistent
329  // shape just in case.
330  sys.reinit_constraints();
331  dof_map.prepare_send_list();
332  }
333 }
virtual bool is_serial() const
Definition: mesh_base.h:134
ImplicitSystem & sys
virtual void allgather()
Definition: mesh_base.h:148
virtual node_iterator nodes_begin()=0
virtual element_iterator elements_begin()=0
virtual element_iterator elements_end()=0
unsigned int n_systems() const
virtual node_iterator nodes_end()=0
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::build_discontinuous_solution_vector ( std::vector< Number > &  soln,
const std::set< std::string > *  system_names = libmesh_nullptr 
) const

Fill the input vector soln with solution values. The entries will be in variable-major format (corresponding to the names from build_variable_names()). If systems_names!=NULL, only include data from the specified systems.

Definition at line 1045 of file equation_systems.C.

References _mesh, _systems, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::DofMap::dof_indices(), end, libMesh::System::get_dof_map(), libMesh::Elem::infinite(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshBase::mesh_dimension(), libMesh::Elem::n_nodes(), n_systems(), libMesh::System::n_vars(), libMesh::FEInterface::nodal_soln(), libMesh::ParallelObject::processor_id(), libMesh::System::update_global_solution(), and libMesh::System::variable_type().

Referenced by libMesh::ExodusII_IO::write_discontinuous_exodusII(), and libMesh::GMVIO::write_discontinuous_gmv().

1047 {
1048  LOG_SCOPE("build_discontinuous_solution_vector()", "EquationSystems");
1049 
1050  libmesh_assert (this->n_systems());
1051 
1052  const unsigned int dim = _mesh.mesh_dimension();
1053 
1054  // Get the number of variables (nv) by counting the number of variables
1055  // in each system listed in system_names
1056  unsigned int nv = 0;
1057 
1058  {
1059  const_system_iterator pos = _systems.begin();
1060  const const_system_iterator end = _systems.end();
1061 
1062  for (; pos != end; ++pos)
1063  {
1064  // Check current system is listed in system_names, and skip pos if not
1065  bool use_current_system = (system_names == libmesh_nullptr);
1066  if (!use_current_system)
1067  use_current_system = system_names->count(pos->first);
1068  if (!use_current_system)
1069  continue;
1070 
1071  const System & system = *(pos->second);
1072  nv += system.n_vars();
1073  }
1074  }
1075 
1076  unsigned int tw=0;
1077 
1078  // get the total weight
1079  {
1080  MeshBase::element_iterator it = _mesh.active_elements_begin();
1081  const MeshBase::element_iterator end = _mesh.active_elements_end();
1082 
1083  for ( ; it != end; ++it)
1084  tw += (*it)->n_nodes();
1085  }
1086 
1087 
1088  // Only if we are on processor zero, allocate the storage
1089  // to hold (number_of_nodes)*(number_of_variables) entries.
1090  if (_mesh.processor_id() == 0)
1091  soln.resize(tw*nv);
1092 
1093  std::vector<Number> sys_soln;
1094 
1095 
1096  unsigned int var_num=0;
1097 
1098  // For each system in this EquationSystems object,
1099  // update the global solution and if we are on processor 0,
1100  // loop over the elements and build the nodal solution
1101  // from the element solution. Then insert this nodal solution
1102  // into the vector passed to build_solution_vector.
1103  {
1104  const_system_iterator pos = _systems.begin();
1105  const const_system_iterator end = _systems.end();
1106 
1107  for (; pos != end; ++pos)
1108  {
1109  // Check current system is listed in system_names, and skip pos if not
1110  bool use_current_system = (system_names == libmesh_nullptr);
1111  if (!use_current_system)
1112  use_current_system = system_names->count(pos->first);
1113  if (!use_current_system)
1114  continue;
1115 
1116  const System & system = *(pos->second);
1117  const unsigned int nv_sys = system.n_vars();
1118 
1119  system.update_global_solution (sys_soln, 0);
1120 
1121  if (_mesh.processor_id() == 0)
1122  {
1123  std::vector<Number> elem_soln; // The finite element solution
1124  std::vector<Number> nodal_soln; // The FE solution interpolated to the nodes
1125  std::vector<dof_id_type> dof_indices; // The DOF indices for the finite element
1126 
1127  for (unsigned int var=0; var<nv_sys; var++)
1128  {
1129  const FEType & fe_type = system.variable_type(var);
1130 
1131  MeshBase::element_iterator it = _mesh.active_elements_begin();
1132  const MeshBase::element_iterator end_elem = _mesh.active_elements_end();
1133 
1134  unsigned int nn=0;
1135 
1136  for ( ; it != end_elem; ++it)
1137  {
1138  const Elem * elem = *it;
1139  system.get_dof_map().dof_indices (elem, dof_indices, var);
1140 
1141  elem_soln.resize(dof_indices.size());
1142 
1143  for (std::size_t i=0; i<dof_indices.size(); i++)
1144  elem_soln[i] = sys_soln[dof_indices[i]];
1145 
1147  fe_type,
1148  elem,
1149  elem_soln,
1150  nodal_soln);
1151 
1152 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1153  // infinite elements should be skipped...
1154  if (!elem->infinite())
1155 #endif
1156  {
1157  libmesh_assert_equal_to (nodal_soln.size(), elem->n_nodes());
1158 
1159  for (unsigned int n=0; n<elem->n_nodes(); n++)
1160  {
1161  soln[nv*(nn++) + (var + var_num)] +=
1162  nodal_soln[n];
1163  }
1164  }
1165  }
1166  }
1167  }
1168 
1169  var_num += nv_sys;
1170  }
1171  }
1172 }
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
libmesh_assert(j)
std::map< std::string, System * >::const_iterator const_system_iterator
virtual element_iterator active_elements_begin()=0
unsigned int n_systems() const
virtual element_iterator active_elements_end()=0
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
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)
Definition: fe_interface.C:513
processor_id_type processor_id() const
std::map< std::string, System * > _systems
UniquePtr< NumericVector< Number > > libMesh::EquationSystems::build_parallel_solution_vector ( const std::set< std::string > *  system_names = libmesh_nullptr) const

A version of build_solution_vector which is appropriate for "parallel" output formats like Nemesis. Returns a UniquePtr to a node-major NumericVector of total length n_nodes*n_vars that various I/O classes can then use to get the local values they need to write on each processor.

Definition at line 669 of file equation_systems.C.

References libMesh::ParallelObject::_communicator, _mesh, _systems, libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::Variable::active_on_subdomain(), libMesh::NumericVector< T >::add(), libMesh::NumericVector< T >::build(), libMesh::NumericVector< T >::close(), libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::DofMap::dof_indices(), end, libMesh::FEInterface::field_type(), libMesh::System::get_dof_map(), libMesh::Elem::infinite(), libMesh::NumericVector< T >::init(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshBase::local_nodes_begin(), libMesh::MeshBase::local_nodes_end(), libMesh::MeshBase::max_node_id(), libMesh::MeshBase::mesh_dimension(), libMesh::DofObject::n_dofs(), libMesh::MeshBase::n_nodes(), libMesh::Elem::n_nodes(), libMesh::System::n_vars(), libMesh::FEInterface::n_vec_dim(), libMesh::FEInterface::nodal_soln(), libMesh::Elem::node_id(), libMesh::Elem::node_ptr(), libMesh::System::number(), libMesh::PARALLEL, libMesh::System::solution, libMesh::TYPE_VECTOR, libMesh::System::update(), libMesh::System::variable(), libMesh::System::variable_type(), and libMesh::Parallel::verify().

Referenced by build_solution_vector(), and libMesh::MeshOutput< MT >::write_equation_systems().

670 {
671  LOG_SCOPE("build_parallel_solution_vector()", "EquationSystems");
672 
673  // This function must be run on all processors at once
674  parallel_object_only();
675 
676  const unsigned int dim = _mesh.mesh_dimension();
677  const dof_id_type nn = _mesh.n_nodes();
678 
679  // We'd better have a contiguous node numbering
680  libmesh_assert_equal_to (nn, _mesh.max_node_id());
681 
682  // allocate storage to hold
683  // (number_of_nodes)*(number_of_variables) entries.
684  // We have to differentiate between between scalar and vector
685  // variables. We intercept vector variables and treat each
686  // component as a scalar variable (consistently with build_solution_names).
687 
688  unsigned int nv = 0;
689 
690  //Could this be replaced by a/some convenience methods?[PB]
691  {
692  unsigned int n_scalar_vars = 0;
693  unsigned int n_vector_vars = 0;
694  const_system_iterator pos = _systems.begin();
695  const const_system_iterator end = _systems.end();
696 
697  for (; pos != end; ++pos)
698  {
699  // Check current system is listed in system_names, and skip pos if not
700  bool use_current_system = (system_names == libmesh_nullptr);
701  if (!use_current_system)
702  use_current_system = system_names->count(pos->first);
703  if (!use_current_system)
704  continue;
705 
706  for (unsigned int vn=0; vn<pos->second->n_vars(); vn++)
707  {
708  if( FEInterface::field_type(pos->second->variable_type(vn)) ==
709  TYPE_VECTOR )
710  n_vector_vars++;
711  else
712  n_scalar_vars++;
713  }
714  }
715  // Here, we're assuming the number of vector components is the same
716  // as the mesh dimension. Will break for mixed dimension meshes.
717  nv = n_scalar_vars + dim*n_vector_vars;
718  }
719 
720  // Get the number of elements that share each node. We will
721  // compute the average value at each node. This is particularly
722  // useful for plotting discontinuous data.
723  MeshBase::element_iterator e_it = _mesh.active_local_elements_begin();
724  const MeshBase::element_iterator e_end = _mesh.active_local_elements_end();
725 
726  // Get the number of local nodes
727  dof_id_type n_local_nodes = cast_int<dof_id_type>
728  (std::distance(_mesh.local_nodes_begin(),
730 
731  // Create a NumericVector to hold the parallel solution
732  UniquePtr<NumericVector<Number> > parallel_soln_ptr = NumericVector<Number>::build(_communicator);
733  NumericVector<Number> & parallel_soln = *parallel_soln_ptr;
734  parallel_soln.init(nn*nv, n_local_nodes*nv, false, PARALLEL);
735 
736  // Create a NumericVector to hold the "repeat_count" for each node - this is essentially
737  // the number of elements contributing to that node's value
738  UniquePtr<NumericVector<Number> > repeat_count_ptr = NumericVector<Number>::build(_communicator);
739  NumericVector<Number> & repeat_count = *repeat_count_ptr;
740  repeat_count.init(nn*nv, n_local_nodes*nv, false, PARALLEL);
741 
742  repeat_count.close();
743 
744  unsigned int var_num=0;
745 
746  // For each system in this EquationSystems object,
747  // update the global solution and if we are on processor 0,
748  // loop over the elements and build the nodal solution
749  // from the element solution. Then insert this nodal solution
750  // into the vector passed to build_solution_vector.
751  const_system_iterator pos = _systems.begin();
752  const const_system_iterator end = _systems.end();
753 
754  for (; pos != end; ++pos)
755  {
756  // Check current system is listed in system_names, and skip pos if not
757  bool use_current_system = (system_names == libmesh_nullptr);
758  if (!use_current_system)
759  use_current_system = system_names->count(pos->first);
760  if (!use_current_system)
761  continue;
762 
763  const System & system = *(pos->second);
764  const unsigned int nv_sys = system.n_vars();
765  const unsigned int sys_num = system.number();
766 
767  //Could this be replaced by a/some convenience methods?[PB]
768  unsigned int n_scalar_vars = 0;
769  unsigned int n_vector_vars = 0;
770  for (unsigned int vn=0; vn<pos->second->n_vars(); vn++)
771  {
772  if( FEInterface::field_type(pos->second->variable_type(vn)) ==
773  TYPE_VECTOR )
774  n_vector_vars++;
775  else
776  n_scalar_vars++;
777  }
778 
779  // Here, we're assuming the number of vector components is the same
780  // as the mesh dimension. Will break for mixed dimension meshes.
781  unsigned int nv_sys_split = n_scalar_vars + dim*n_vector_vars;
782 
783  // Update the current_local_solution
784  {
785  System & non_const_sys = const_cast<System &>(system);
786  // We used to simply call non_const_sys.solution->close()
787  // here, but that is not allowed when the solution vector is
788  // locked read-only, for example when printing the solution
789  // during during the middle of a solve... So try to be a bit
790  // more careful about calling close() unnecessarily.
791  libmesh_assert(this->comm().verify(non_const_sys.solution->closed()));
792  if (!non_const_sys.solution->closed())
793  non_const_sys.solution->close();
794  non_const_sys.update();
795  }
796 
797  NumericVector<Number> & sys_soln(*system.current_local_solution);
798 
799  std::vector<Number> elem_soln; // The finite element solution
800  std::vector<Number> nodal_soln; // The FE solution interpolated to the nodes
801  std::vector<dof_id_type> dof_indices; // The DOF indices for the finite element
802 
803  for (unsigned int var=0; var<nv_sys; var++)
804  {
805  const FEType & fe_type = system.variable_type(var);
806  const Variable & var_description = system.variable(var);
807  const DofMap & dof_map = system.get_dof_map();
808 
809  unsigned int n_vec_dim = FEInterface::n_vec_dim( pos->second->get_mesh(), fe_type );
810 
811  MeshBase::element_iterator it = _mesh.active_local_elements_begin();
812  const MeshBase::element_iterator end_elem = _mesh.active_local_elements_end();
813 
814  for ( ; it != end_elem; ++it)
815  {
816  const Elem * elem = *it;
817 
818  if (var_description.active_on_subdomain((*it)->subdomain_id()))
819  {
820  dof_map.dof_indices (elem, dof_indices, var);
821 
822  elem_soln.resize(dof_indices.size());
823 
824  for (std::size_t i=0; i<dof_indices.size(); i++)
825  elem_soln[i] = sys_soln(dof_indices[i]);
826 
828  fe_type,
829  elem,
830  elem_soln,
831  nodal_soln);
832 
833 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
834  // infinite elements should be skipped...
835  if (!elem->infinite())
836 #endif
837  {
838  libmesh_assert_equal_to (nodal_soln.size(), n_vec_dim*elem->n_nodes());
839 
840  for (unsigned int n=0; n<elem->n_nodes(); n++)
841  {
842  for( unsigned int d=0; d < n_vec_dim; d++ )
843  {
844  // For vector-valued elements, all components are in nodal_soln. For each
845  // node, the components are stored in order, i.e. node_0 -> s0_x, s0_y, s0_z
846  parallel_soln.add(nv*(elem->node_id(n)) + (var+d + var_num), nodal_soln[n_vec_dim*n+d]);
847 
848  // Increment the repeat count for this position
849  repeat_count.add(nv*(elem->node_id(n)) + (var+d + var_num), 1);
850  }
851  }
852  }
853  }
854  else // If this variable doesn't exist on this subdomain we have to still increment repeat_count so that we won't divide by 0 later:
855  for (unsigned int n=0; n<elem->n_nodes(); n++)
856  // Only do this if this variable has NO DoFs at this node... it might have some from an ajoining element...
857  if(!elem->node_ptr(n)->n_dofs(sys_num, var))
858  for( unsigned int d=0; d < n_vec_dim; d++ )
859  repeat_count.add(nv*(elem->node_id(n)) + (var+d + var_num), 1);
860 
861  } // end loop over elements
862  } // end loop on variables in this system
863 
864  var_num += nv_sys_split;
865  } // end loop over systems
866 
867  parallel_soln.close();
868  repeat_count.close();
869 
870  // Divide to get the average value at the nodes
871  parallel_soln /= repeat_count;
872 
873  return UniquePtr<NumericVector<Number> >(parallel_soln_ptr.release());
874 }
virtual dof_id_type max_node_id() const =0
static FEFieldType field_type(const FEType &fe_type)
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
const Parallel::Communicator & _communicator
libmesh_assert(j)
static UniquePtr< NumericVector< T > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
std::map< std::string, System * >::const_iterator const_system_iterator
virtual node_iterator local_nodes_end()=0
virtual element_iterator active_local_elements_begin()=0
static unsigned int n_vec_dim(const MeshBase &mesh, const FEType &fe_type)
virtual node_iterator local_nodes_begin()=0
bool verify(const T &r, const Communicator &comm=Communicator_World)
const Parallel::Communicator & comm() const
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
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)
Definition: fe_interface.C:513
virtual dof_id_type n_nodes() const =0
virtual element_iterator active_local_elements_end()=0
uint8_t dof_id_type
Definition: id_types.h:64
std::map< std::string, System * > _systems
void libMesh::EquationSystems::build_solution_vector ( std::vector< Number > &  soln,
const std::string &  system_name,
const std::string &  variable_name = "all_vars" 
) const

Fill the input vector soln with the solution values for the system named name. Note that the input vector soln will only be assembled on processor 0, so this method is only applicable to outputting plot files from processor 0.

Definition at line 581 of file equation_systems.C.

Referenced by libMesh::MeshOutput< MT >::write_equation_systems().

584 {
585  //TODO:[BSK] re-implement this from the method below
586  libmesh_not_implemented();
587 
588  // // Get a reference to the named system
589  // const System & system = this->get_system(system_name);
590 
591  // // Get the number associated with the variable_name we are passed
592  // const unsigned short int variable_num = system.variable_number(variable_name);
593 
594  // // Get the dimension of the current mesh
595  // const unsigned int dim = _mesh.mesh_dimension();
596 
597  // // If we're on processor 0, allocate enough memory to hold the solution.
598  // // Since we're only looking at one variable, there will be one solution value
599  // // for each node in the mesh.
600  // if (_mesh.processor_id() == 0)
601  // soln.resize(_mesh.n_nodes());
602 
603  // // Vector to hold the global solution from all processors
604  // std::vector<Number> sys_soln;
605 
606  // // Update the global solution from all processors
607  // system.update_global_solution (sys_soln, 0);
608 
609  // // Temporary vector to store the solution on an individual element.
610  // std::vector<Number> elem_soln;
611 
612  // // The FE solution interpolated to the nodes
613  // std::vector<Number> nodal_soln;
614 
615  // // The DOF indices for the element
616  // std::vector<dof_id_type> dof_indices;
617 
618  // // Determine the finite/infinite element type used in this system
619  // const FEType & fe_type = system.variable_type(variable_num);
620 
621  // // Define iterators to iterate over all the elements of the mesh
622  // const_active_elem_iterator it (_mesh.elements_begin());
623  // const const_active_elem_iterator end(_mesh.elements_end());
624 
625  // // Loop over elements
626  // for ( ; it != end; ++it)
627  // {
628  // // Convenient shortcut to the element pointer
629  // const Elem * elem = *it;
630 
631  // // Fill the dof_indices vector for this variable
632  // system.get_dof_map().dof_indices(elem,
633  // dof_indices,
634  // variable_num);
635 
636  // // Resize the element solution vector to fit the
637  // // dof_indices for this element.
638  // elem_soln.resize(dof_indices.size());
639 
640  // // Transfer the system solution to the element
641  // // solution by mapping it through the dof_indices vector.
642  // for (std::size_t i=0; i<dof_indices.size(); i++)
643  // elem_soln[i] = sys_soln[dof_indices[i]];
644 
645  // // Using the FE interface, compute the nodal_soln
646  // // for the current elemnt type given the elem_soln
647  // FEInterface::nodal_soln (dim,
648  // fe_type,
649  // elem,
650  // elem_soln,
651  // nodal_soln);
652 
653  // // Sanity check -- make sure that there are the same number
654  // // of entries in the nodal_soln as there are nodes in the
655  // // element!
656  // libmesh_assert_equal_to (nodal_soln.size(), elem->n_nodes());
657 
658  // // Copy the nodal solution over into the correct place in
659  // // the global soln vector which will be returned to the user.
660  // for (unsigned int n=0; n<elem->n_nodes(); n++)
661  // soln[elem->node_id(n)] = nodal_soln[n];
662  // }
663 }
void libMesh::EquationSystems::build_solution_vector ( std::vector< Number > &  soln,
const std::set< std::string > *  system_names = libmesh_nullptr 
) const

Fill the input vector soln with solution values. The entries will be in variable-major format (corresponding to the names from build_variable_names()). If systems_names!=NULL, only include data from the specified systems.

Definition at line 878 of file equation_systems.C.

References build_parallel_solution_vector().

880 {
881  LOG_SCOPE("build_solution_vector()", "EquationSystems");
882 
883  // Call the parallel implementation
884  UniquePtr<NumericVector<Number> > parallel_soln =
885  this->build_parallel_solution_vector(system_names);
886 
887  // Localize the NumericVector into the provided std::vector.
888  parallel_soln->localize_to_one(soln);
889 }
UniquePtr< NumericVector< Number > > build_parallel_solution_vector(const std::set< std::string > *system_names=libmesh_nullptr) const
void libMesh::EquationSystems::build_variable_names ( std::vector< std::string > &  var_names,
const FEType type = libmesh_nullptr,
const std::set< std::string > *  system_names = libmesh_nullptr 
) const

Fill the input vector var_names with the names of the variables for each system. If type is passed, only variables of the specified type will be populated. If systems_names!=NULL, only include names from the specified systems.

Definition at line 478 of file equation_systems.C.

References _systems, end, libMesh::FEInterface::field_type(), get_mesh(), libmesh_nullptr, libMesh::MeshBase::mesh_dimension(), n_vars(), libMesh::FEInterface::n_vec_dim(), and libMesh::TYPE_VECTOR.

Referenced by libMesh::ExodusII_IO::write_discontinuous_exodusII(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), and libMesh::MeshOutput< MT >::write_equation_systems().

481 {
482  unsigned int var_num=0;
483 
484  const_system_iterator pos = _systems.begin();
485  const const_system_iterator end = _systems.end();
486 
487  // Need to size var_names by scalar variables plus all the
488  // vector components for all the vector variables
489  //Could this be replaced by a/some convenience methods?[PB]
490  {
491  unsigned int n_scalar_vars = 0;
492  unsigned int n_vector_vars = 0;
493 
494  for (; pos != end; ++pos)
495  {
496  // Check current system is listed in system_names, and skip pos if not
497  bool use_current_system = (system_names == libmesh_nullptr);
498  if (!use_current_system)
499  use_current_system = system_names->count(pos->first);
500  if (!use_current_system)
501  continue;
502 
503  for (unsigned int vn=0; vn<pos->second->n_vars(); vn++)
504  {
505  if( FEInterface::field_type(pos->second->variable_type(vn)) ==
506  TYPE_VECTOR )
507  n_vector_vars++;
508  else
509  n_scalar_vars++;
510  }
511  }
512 
513  // Here, we're assuming the number of vector components is the same
514  // as the mesh dimension. Will break for mixed dimension meshes.
515  unsigned int dim = this->get_mesh().mesh_dimension();
516  unsigned int nv = n_scalar_vars + dim*n_vector_vars;
517 
518  // We'd better not have more than dim*his->n_vars() (all vector variables)
519  libmesh_assert_less_equal ( nv, dim*this->n_vars() );
520 
521  // Here, we're assuming the number of vector components is the same
522  // as the mesh dimension. Will break for mixed dimension meshes.
523 
524  var_names.resize( nv );
525  }
526 
527  // reset
528  pos = _systems.begin();
529 
530  for (; pos != end; ++pos)
531  {
532  // Check current system is listed in system_names, and skip pos if not
533  bool use_current_system = (system_names == libmesh_nullptr);
534  if (!use_current_system)
535  use_current_system = system_names->count(pos->first);
536  if (!use_current_system)
537  continue;
538 
539  for (unsigned int vn=0; vn<pos->second->n_vars(); vn++)
540  {
541  std::string var_name = pos->second->variable_name(vn);
542  FEType fe_type = pos->second->variable_type(vn);
543 
544  unsigned int n_vec_dim = FEInterface::n_vec_dim( pos->second->get_mesh(), fe_type);
545 
546  // Filter on the type if requested
547  if (type == libmesh_nullptr || (type && *type == fe_type))
548  {
549  if( FEInterface::field_type(fe_type) == TYPE_VECTOR )
550  {
551  switch(n_vec_dim)
552  {
553  case 0:
554  case 1:
555  var_names[var_num++] = var_name;
556  break;
557  case 2:
558  var_names[var_num++] = var_name+"_x";
559  var_names[var_num++] = var_name+"_y";
560  break;
561  case 3:
562  var_names[var_num++] = var_name+"_x";
563  var_names[var_num++] = var_name+"_y";
564  var_names[var_num++] = var_name+"_z";
565  break;
566  default:
567  libmesh_error_msg("Invalid dim in build_variable_names");
568  }
569  }
570  else
571  var_names[var_num++] = var_name;
572  }
573  }
574  }
575  // Now resize again in case we filtered any names
576  var_names.resize(var_num);
577 }
static FEFieldType field_type(const FEType &fe_type)
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
static unsigned int n_vec_dim(const MeshBase &mesh, const FEType &fe_type)
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
unsigned int n_vars() const
const MeshBase & get_mesh() const
std::map< std::string, System * > _systems
void libMesh::EquationSystems::clear ( )
virtual

Returns tha data structure to a pristine state.

Definition at line 73 of file equation_systems.C.

References _systems, libMesh::Parameters::clear(), libmesh_nullptr, parameters, and libMesh::sys.

Referenced by ~EquationSystems().

74 {
75  // Clear any additional parameters
76  parameters.clear ();
77 
78  // clear the systems. We must delete them
79  // since we newed them!
80  while (!_systems.empty())
81  {
82  system_iterator pos = _systems.begin();
83 
84  System * sys = pos->second;
85  delete sys;
86  sys = libmesh_nullptr;
87 
88  _systems.erase (pos);
89  }
90 }
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
std::map< std::string, System * >::iterator system_iterator
virtual void clear()
Definition: parameters.h:321
std::map< std::string, System * > _systems
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
a reference to the Parallel::Communicator object used by this mesh.

Definition at line 87 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_postcheck(), libMesh::__libmesh_petsc_snes_residual(), libMesh::__libmesh_tao_equality_constraints(), libMesh::__libmesh_tao_equality_constraints_jacobian(), libMesh::__libmesh_tao_gradient(), libMesh::__libmesh_tao_hessian(), libMesh::__libmesh_tao_inequality_constraints(), libMesh::__libmesh_tao_inequality_constraints_jacobian(), libMesh::__libmesh_tao_objective(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), _read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::EigenSparseLinearSolver< T >::adjoint_solve(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::Parallel::BinSorter< KeyType, IdxType >::binsort(), libMesh::Parallel::Sort< KeyType, IdxType >::binsort(), libMesh::MeshTools::bounding_box(), libMesh::MeshCommunication::broadcast(), libMesh::SparseMatrix< T >::build(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::Parallel::Histogram< KeyType, IdxType >::build_histogram(), libMesh::PetscNonlinearSolver< T >::build_mat_null_space(), libMesh::BoundaryInfo::build_node_list_from_side_list(), build_parallel_solution_vector(), libMesh::MeshBase::cache_elem_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::DistributedVector< T >::clone(), libMesh::EigenSparseVector< T >::clone(), libMesh::LaspackVector< T >::clone(), libMesh::EpetraVector< T >::clone(), libMesh::PetscVector< T >::clone(), libMesh::EpetraVector< T >::close(), libMesh::Parallel::Sort< KeyType, IdxType >::communicate_bins(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshRefinement::create_parent_error_vector(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::EpetraVector< T >::EpetraVector(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), get_solution(), libMesh::LocationMap< T >::init(), libMesh::PetscDiffSolver::init(), libMesh::TimeSolver::init(), libMesh::TopologyMap::init(), libMesh::TaoOptimizationSolver< T >::init(), libMesh::PetscNonlinearSolver< T >::init(), libMesh::DistributedVector< T >::init(), libMesh::EpetraVector< T >::init(), libMesh::PetscVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_flags(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_p_levels(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshTools::libmesh_assert_valid_unique_ids(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::limit_overrefined_boundary(), libMesh::MeshRefinement::limit_underrefined_boundary(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_new_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_new_nodes_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_unique_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshCommunication::make_p_levels_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::DistributedVector< T >::max(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::DistributedVector< T >::min(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::DistributedMesh::parallel_max_elem_id(), libMesh::DistributedMesh::parallel_max_node_id(), libMesh::ReplicatedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_n_elem(), libMesh::DistributedMesh::parallel_n_nodes(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), libMesh::MetisPartitioner::partition_range(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::SparseMatrix< T >::print(), libMesh::MeshTools::processor_bounding_box(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::MeshTools::subdomain_bounding_box(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::Parallel::sync_node_data_by_element_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::DistributedVector< T >::zero_clone(), libMesh::LaspackVector< T >::zero_clone(), libMesh::EigenSparseVector< T >::zero_clone(), libMesh::EpetraVector< T >::zero_clone(), and libMesh::PetscVector< T >::zero_clone().

88  { return _communicator; }
const Parallel::Communicator & _communicator
bool libMesh::EquationSystems::compare ( const EquationSystems other_es,
const Real  threshold,
const bool  verbose 
) const
virtual
Returns
true when this equation system contains identical data, up to the given threshold. Delegates most of the comparisons to perform to the responsible systems

Definition at line 1176 of file equation_systems.C.

References _systems, libMesh::System::compare(), end, get_system(), n_systems(), and libMesh::out.

Referenced by read().

1179 {
1180  // safety check, whether we handle at least the same number
1181  // of systems
1182  std::vector<bool> os_result;
1183 
1184  if (this->n_systems() != other_es.n_systems())
1185  {
1186  if (verbose)
1187  {
1188  libMesh::out << " Fatal difference. This system handles "
1189  << this->n_systems() << " systems," << std::endl
1190  << " while the other system handles "
1191  << other_es.n_systems()
1192  << " systems." << std::endl
1193  << " Aborting comparison." << std::endl;
1194  }
1195  return false;
1196  }
1197  else
1198  {
1199  // start comparing each system
1200  const_system_iterator pos = _systems.begin();
1201  const const_system_iterator end = _systems.end();
1202 
1203  for (; pos != end; ++pos)
1204  {
1205  const std::string & sys_name = pos->first;
1206  const System & system = *(pos->second);
1207 
1208  // get the other system
1209  const System & other_system = other_es.get_system (sys_name);
1210 
1211  os_result.push_back (system.compare (other_system, threshold, verbose));
1212 
1213  }
1214 
1215  }
1216 
1217 
1218  // sum up the results
1219  if (os_result.size()==0)
1220  return true;
1221  else
1222  {
1223  bool os_identical;
1224  unsigned int n = 0;
1225  do
1226  {
1227  os_identical = os_result[n];
1228  n++;
1229  }
1230  while (os_identical && n<os_result.size());
1231  return os_identical;
1232  }
1233 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
unsigned int n_systems() const
OStreamProxy out(std::cout)
std::map< std::string, System * > _systems
void libMesh::EquationSystems::delete_system ( const std::string &  name)

Remove the system named name from the systems array. This function is now deprecated - write the libmesh-devel mailing list if you need it reimplemented.

Definition at line 434 of file equation_systems.C.

References _systems, and libMesh::Quality::name().

435 {
436  libmesh_deprecated();
437 
438  if (!_systems.count(name))
439  libmesh_error_msg("ERROR: no system named " << name);
440 
441  delete _systems[name];
442 
443  _systems.erase (name);
444 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
std::map< std::string, System * > _systems
void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited
void libMesh::EquationSystems::disable_refine_in_reinit ( )
inline

Calls to reinit() will not try to coarsen or refine the mesh

Definition at line 449 of file equation_systems.h.

References _refine_in_reinit.

450  { this->_refine_in_reinit = false; }
void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

Referenced by libMesh::ReferenceCounter::n_objects().

101 {
102  _enable_print_counter = true;
103  return;
104 }
void libMesh::EquationSystems::enable_refine_in_reinit ( )
inline

Calls to reinit() will also do two-step coarsen-then-refine

Definition at line 443 of file equation_systems.h.

References _refine_in_reinit.

444  { this->_refine_in_reinit = true; }
std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (Counts::iterator it = _counts.begin();
59  it != _counts.end(); ++it)
60  {
61  const std::string name(it->first);
62  const unsigned int creations = it->second.first;
63  const unsigned int destructions = it->second.second;
64 
65  oss << "| " << name << " reference count information:\n"
66  << "| Creations: " << creations << '\n'
67  << "| Destructions: " << destructions << '\n';
68  }
69 
70  oss << " ---------------------------------------------------------------------------- \n";
71 
72  return oss.str();
73 
74 #else
75 
76  return "";
77 
78 #endif
79 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
std::string libMesh::EquationSystems::get_info ( ) const
virtual
Returns
a string containing information about the systems, flags, and parameters.

Definition at line 1237 of file equation_systems.C.

References _systems, end, and n_systems().

Referenced by print_info(), and read().

1238 {
1239  std::ostringstream oss;
1240 
1241  oss << " EquationSystems\n"
1242  << " n_systems()=" << this->n_systems() << '\n';
1243 
1244  // Print the info for the individual systems
1245  const_system_iterator pos = _systems.begin();
1246  const const_system_iterator end = _systems.end();
1247 
1248  for (; pos != end; ++pos)
1249  oss << pos->second->get_info();
1250 
1251 
1252  // // Possibly print the parameters
1253  // if (!this->parameters.empty())
1254  // {
1255  // oss << " n_parameters()=" << this->n_parameters() << '\n';
1256  // oss << " Parameters:\n";
1257 
1258  // for (std::map<std::string, Real>::const_iterator
1259  // param = _parameters.begin(); param != _parameters.end();
1260  // ++param)
1261  // oss << " "
1262  // << "\""
1263  // << param->first
1264  // << "\""
1265  // << "="
1266  // << param->second
1267  // << '\n';
1268  // }
1269 
1270  return oss.str();
1271 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
unsigned int n_systems() const
std::map< std::string, System * > _systems
const MeshBase & libMesh::EquationSystems::get_mesh ( ) const
inline
MeshBase & libMesh::EquationSystems::get_mesh ( )
inline
Returns
a reference to the mesh

Definition at line 536 of file equation_systems.h.

References _mesh.

537 {
538  return _mesh;
539 }
void libMesh::EquationSystems::get_solution ( std::vector< Number > &  soln,
std::vector< std::string > &  names 
) const

Retrieve the solution data for CONSTANT MONOMIALs. If names is populated, only the variables corresponding to those names will be retrieved. This can be used to filter which variables are retrieved.

Definition at line 893 of file equation_systems.C.

References libMesh::ParallelObject::_communicator, _mesh, _systems, libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::Variable::active_on_subdomain(), libMesh::NumericVector< T >::build(), libMesh::NumericVector< T >::close(), libMesh::ParallelObject::comm(), libMesh::CONSTANT, libMesh::System::current_local_solution, libMesh::DofMap::dof_indices(), end, libMesh::System::get_dof_map(), libMesh::DofObject::id(), libMesh::NumericVector< T >::init(), libMesh::libmesh_assert(), libMesh::NumericVector< T >::localize_to_one(), libMesh::MeshBase::max_elem_id(), libMesh::MONOMIAL, libMesh::MeshBase::n_elem(), libMesh::ParallelObject::n_processors(), n_systems(), libMesh::System::n_vars(), libMesh::PARALLEL, libMesh::ParallelObject::processor_id(), libMesh::NumericVector< T >::set(), libMesh::System::solution, libMesh::Elem::subdomain_id(), libMesh::System::update(), libMesh::System::variable(), libMesh::System::variable_name(), libMesh::System::variable_type(), and libMesh::Parallel::verify().

Referenced by libMesh::ExodusII_IO::write_element_data().

895 {
896  // This function must be run on all processors at once
897  parallel_object_only();
898 
899  libmesh_assert (this->n_systems());
900 
901  const dof_id_type ne = _mesh.n_elem();
902 
903  libmesh_assert_equal_to (ne, _mesh.max_elem_id());
904 
905  // If the names vector has entries, we will only populate the soln vector
906  // with names included in that list. Note: The names vector may be
907  // reordered upon exiting this function
908  std::vector<std::string> filter_names = names;
909  bool is_filter_names = !filter_names.empty();
910 
911  soln.clear();
912  names.clear();
913 
914  const FEType type(CONSTANT, MONOMIAL);
915 
916  dof_id_type nv = 0;
917 
918  // Find the total number of variables to output
919  std::vector<std::vector<unsigned> > do_output(_systems.size());
920  {
921  const_system_iterator pos = _systems.begin();
922  const const_system_iterator end = _systems.end();
923  unsigned sys_ctr = 0;
924 
925  for (; pos != end; ++pos, ++sys_ctr)
926  {
927  const System & system = *(pos->second);
928  const unsigned int nv_sys = system.n_vars();
929 
930  do_output[sys_ctr].resize(nv_sys);
931 
932  for (unsigned int var=0; var < nv_sys; ++var)
933  {
934  if (system.variable_type(var) != type ||
935  (is_filter_names && std::find(filter_names.begin(), filter_names.end(), system.variable_name(var)) == filter_names.end()))
936  continue;
937 
938  // Otherwise, this variable should be output
939  nv++;
940  do_output[sys_ctr][var] = 1;
941  }
942  }
943  }
944 
945  // If there are no variables to write out don't do anything...
946  if (!nv)
947  return;
948 
949  // We can handle the case where there are NULLs in the Elem vector
950  // by just having extra zeros in the solution vector.
951  numeric_index_type parallel_soln_global_size = ne*nv;
952 
953  numeric_index_type div = parallel_soln_global_size / this->n_processors();
954  numeric_index_type mod = parallel_soln_global_size % this->n_processors();
955 
956  // Initialize all processors to the average size.
957  numeric_index_type parallel_soln_local_size = div;
958 
959  // The first "mod" processors get an extra entry.
960  if (this->processor_id() < mod)
961  parallel_soln_local_size = div+1;
962 
963  // Create a NumericVector to hold the parallel solution
964  UniquePtr<NumericVector<Number> > parallel_soln_ptr = NumericVector<Number>::build(_communicator);
965  NumericVector<Number> & parallel_soln = *parallel_soln_ptr;
966  parallel_soln.init(parallel_soln_global_size,
967  parallel_soln_local_size,
968  /*fast=*/false,
969  /*ParallelType=*/PARALLEL);
970 
971  dof_id_type var_num = 0;
972 
973  // For each system in this EquationSystems object,
974  // update the global solution and collect the
975  // CONSTANT MONOMIALs. The entries are in variable-major
976  // format.
977  const_system_iterator pos = _systems.begin();
978  const const_system_iterator end = _systems.end();
979  unsigned sys_ctr = 0;
980 
981  for (; pos != end; ++pos, ++sys_ctr)
982  {
983  const System & system = *(pos->second);
984  const unsigned int nv_sys = system.n_vars();
985 
986  // Update the current_local_solution
987  {
988  System & non_const_sys = const_cast<System &>(system);
989  // We used to simply call non_const_sys.solution->close()
990  // here, but that is not allowed when the solution vector is
991  // locked read-only, for example when printing the solution
992  // during during the middle of a solve... So try to be a bit
993  // more careful about calling close() unnecessarily.
994  libmesh_assert(this->comm().verify(non_const_sys.solution->closed()));
995  if (!non_const_sys.solution->closed())
996  non_const_sys.solution->close();
997  non_const_sys.update();
998  }
999 
1000  NumericVector<Number> & sys_soln(*system.current_local_solution);
1001 
1002  // The DOF indices for the finite element
1003  std::vector<dof_id_type> dof_indices;
1004 
1005  // Loop over the variable names and load them in order
1006  for (unsigned int var=0; var < nv_sys; ++var)
1007  {
1008  // Skip this variable if we are not outputting it.
1009  if (!do_output[sys_ctr][var])
1010  continue;
1011 
1012  names.push_back(system.variable_name(var));
1013 
1014  const Variable & variable = system.variable(var);
1015  const DofMap & dof_map = system.get_dof_map();
1016 
1017  MeshBase::element_iterator it = _mesh.active_local_elements_begin();
1018  const MeshBase::element_iterator end_elem = _mesh.active_local_elements_end();
1019 
1020  for ( ; it != end_elem; ++it)
1021  {
1022  const Elem * elem = *it;
1023 
1024  if (variable.active_on_subdomain(elem->subdomain_id()))
1025  {
1026  dof_map.dof_indices (elem, dof_indices, var);
1027 
1028  libmesh_assert_equal_to (1, dof_indices.size());
1029 
1030  parallel_soln.set((ne*var_num)+elem->id(), sys_soln(dof_indices[0]));
1031  }
1032  }
1033 
1034  var_num++;
1035  } // end loop on variables in this system
1036  } // end loop over systems
1037 
1038  parallel_soln.close();
1039 
1040  parallel_soln.localize_to_one(soln);
1041 }
processor_id_type n_processors() const
IterBase * end
const Parallel::Communicator & _communicator
libmesh_assert(j)
static UniquePtr< NumericVector< T > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
std::map< std::string, System * >::const_iterator const_system_iterator
virtual dof_id_type max_elem_id() const =0
dof_id_type numeric_index_type
Definition: id_types.h:92
virtual element_iterator active_local_elements_begin()=0
unsigned int n_systems() const
bool verify(const T &r, const Communicator &comm=Communicator_World)
const Parallel::Communicator & comm() const
virtual dof_id_type n_elem() const =0
virtual element_iterator active_local_elements_end()=0
processor_id_type processor_id() const
uint8_t dof_id_type
Definition: id_types.h:64
std::map< std::string, System * > _systems
template<typename T_sys >
const T_sys & libMesh::EquationSystems::get_system ( const std::string &  name) const
inline
Returns
a constant reference to the system named name. The template argument defines the return type. For example, const SteadySystem & sys = eq.get_system<SteadySystem> ("sys"); is an example of how the method might be used

Definition at line 644 of file equation_systems.h.

References _systems.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), _read_impl(), libMesh::EnsightIO::add_scalar(), add_system(), libMesh::EnsightIO::add_vector(), adjoint_solve(), allgather(), libMesh::ExactSolution::attach_exact_deriv(), libMesh::ExactSolution::attach_exact_value(), libMesh::ExactSolution::attach_reference_solution(), compare(), libMesh::ExactSolution::compute_error(), libMesh::GMVIO::copy_nodal_solution(), libMesh::ExactSolution::error_norm(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactSolution::ExactSolution(), init(), reinit(), sensitivity_solve(), solve(), update(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

645 {
646  const_system_iterator pos = _systems.find(name);
647 
648  // Check for errors
649  if (pos == _systems.end())
650  libmesh_error_msg("ERROR: no system named \"" << name << "\" found!");
651 
652  // Attempt dynamic cast
653  return *cast_ptr<T_sys *>(pos->second);
654 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
std::map< std::string, System * >::const_iterator const_system_iterator
std::map< std::string, System * > _systems
template<typename T_sys >
T_sys & libMesh::EquationSystems::get_system ( const std::string &  name)
inline
Returns
a writeable referene to the system named name. The template argument defines the return type. For example, const SteadySystem & sys = eq.get_system<SteadySystem> ("sys"); is an example of how the method might be used

Definition at line 663 of file equation_systems.h.

References _systems.

664 {
665  system_iterator pos = _systems.find(name);
666 
667  // Check for errors
668  if (pos == _systems.end())
669  libmesh_error_msg("ERROR: no system named " << name << " found!");
670 
671  // Attempt dynamic cast
672  return *cast_ptr<T_sys *>(pos->second);
673 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
std::map< std::string, System * >::iterator system_iterator
std::map< std::string, System * > _systems
template<typename T_sys >
const T_sys & libMesh::EquationSystems::get_system ( const unsigned int  num) const
inline
Returns
a constant reference to system number num. The template argument defines the return type. For example, const SteadySystem & sys = eq.get_system<SteadySystem> (0); is an example of how the method might be used

Definition at line 593 of file equation_systems.h.

References _systems, end, and n_systems().

594 {
595  libmesh_assert_less (num, this->n_systems());
596 
597 
598  const_system_iterator pos = _systems.begin();
599  const const_system_iterator end = _systems.end();
600 
601  for (; pos != end; ++pos)
602  if (pos->second->number() == num)
603  break;
604 
605  // Check for errors
606  if (pos == end)
607  libmesh_error_msg("ERROR: no system number " << num << " found!");
608 
609  // Attempt dynamic cast
610  return *cast_ptr<T_sys *>(pos->second);
611 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
unsigned int n_systems() const
std::map< std::string, System * > _systems
template<typename T_sys >
T_sys & libMesh::EquationSystems::get_system ( const unsigned int  num)
inline
Returns
a writeable referene to the system number num. The template argument defines the return type. For example, const SteadySystem & sys = eq.get_system<SteadySystem> (0); is an example of how the method might be used

Definition at line 618 of file equation_systems.h.

References _systems, end, and n_systems().

619 {
620  libmesh_assert_less (num, this->n_systems());
621 
622  const_system_iterator pos = _systems.begin();
623  const const_system_iterator end = _systems.end();
624 
625  for (; pos != end; ++pos)
626  if (pos->second->number() == num)
627  break;
628 
629  // Check for errors
630  if (pos == end)
631  libmesh_error_msg("ERROR: no system number " << num << " found!");
632 
633  // Attempt dynamic cast
634  return *cast_ptr<T_sys *>(pos->second);
635 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
unsigned int n_systems() const
std::map< std::string, System * > _systems
const System & libMesh::EquationSystems::get_system ( const std::string &  name) const
inline
Returns
a constant reference to the system named name.

Definition at line 682 of file equation_systems.h.

References libMesh::Quality::name().

683 {
684  return this->get_system<System>(name);
685 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
System & libMesh::EquationSystems::get_system ( const std::string &  name)
inline
Returns
a writeable referene to the system named name.

Definition at line 690 of file equation_systems.h.

References libMesh::Quality::name().

691 {
692  return this->get_system<System>(name);
693 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
const System & libMesh::EquationSystems::get_system ( const unsigned int  num) const
inline
Returns
a constant reference to system number num.

Definition at line 698 of file equation_systems.h.

699 {
700  return this->get_system<System>(num);
701 }
System & libMesh::EquationSystems::get_system ( const unsigned int  num)
inline
Returns
a writeable referene to the system number num.

Definition at line 706 of file equation_systems.h.

707 {
708  return this->get_system<System>(num);
709 }
bool libMesh::EquationSystems::has_system ( const std::string &  name) const
inline
Returns
true if the system named name exists within this EquationSystems object.

Definition at line 581 of file equation_systems.h.

References _systems.

Referenced by libMesh::EnsightIO::add_scalar(), libMesh::EnsightIO::add_vector(), libMesh::ExactSolution::compute_error(), and libMesh::ExactSolution::error_norm().

582 {
583  if (_systems.find(name) == _systems.end())
584  return false;
585  return true;
586 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
std::map< std::string, System * > _systems
void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the construction counter. Should be called in the constructor of any derived class that will be reference counted.

Definition at line 160 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

161 {
162  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
163  std::pair<unsigned int, unsigned int> & p = _counts[name];
164 
165  p.first++;
166 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the destruction counter. Should be called in the destructor of any derived class that will be reference counted.

Definition at line 173 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

174 {
175  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
176  std::pair<unsigned int, unsigned int> & p = _counts[name];
177 
178  p.second++;
179 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
void libMesh::EquationSystems::init ( )
virtual

Initialize all the systems

Definition at line 94 of file equation_systems.C.

References _mesh, libMesh::MeshRefinement::clean_refinement_flags(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), get_system(), n_systems(), libMesh::MeshBase::nodes_begin(), and libMesh::MeshBase::nodes_end().

Referenced by _read_impl(), and libMesh::ErrorVector::plot_error().

95 {
96  const unsigned int n_sys = this->n_systems();
97 
98  libmesh_assert_not_equal_to (n_sys, 0);
99 
100  // Tell all the \p DofObject entities how many systems
101  // there are.
102  {
103  MeshBase::node_iterator node_it = _mesh.nodes_begin();
104  const MeshBase::node_iterator node_end = _mesh.nodes_end();
105 
106  for ( ; node_it != node_end; ++node_it)
107  (*node_it)->set_n_systems(n_sys);
108 
109  MeshBase::element_iterator elem_it = _mesh.elements_begin();
110  const MeshBase::element_iterator elem_end = _mesh.elements_end();
111 
112  for ( ; elem_it != elem_end; ++elem_it)
113  (*elem_it)->set_n_systems(n_sys);
114  }
115 
116  for (unsigned int i=0; i != this->n_systems(); ++i)
117  this->get_system(i).init();
118 
119 #ifdef LIBMESH_ENABLE_AMR
120  MeshRefinement mesh_refine(_mesh);
121  mesh_refine.clean_refinement_flags();
122 #endif
123 }
virtual node_iterator nodes_begin()=0
virtual element_iterator elements_begin()=0
virtual element_iterator elements_end()=0
unsigned int n_systems() const
virtual node_iterator nodes_end()=0
const T_sys & get_system(const std::string &name) const
std::size_t libMesh::EquationSystems::n_active_dofs ( ) const

Returns the number of active degrees of freedom for the EquationSystems object.

Definition at line 1323 of file equation_systems.C.

References _systems, and end.

1324 {
1325  std::size_t tot=0;
1326 
1327  const_system_iterator pos = _systems.begin();
1328  const const_system_iterator end = _systems.end();
1329 
1330  for (; pos != end; ++pos)
1331  tot += pos->second->n_active_dofs();
1332 
1333  return tot;
1334 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
std::map< std::string, System * > _systems
std::size_t libMesh::EquationSystems::n_dofs ( ) const
Returns
the total number of degrees of freedom in all systems.

Definition at line 1307 of file equation_systems.C.

References _systems, and end.

1308 {
1309  std::size_t tot=0;
1310 
1311  const_system_iterator pos = _systems.begin();
1312  const const_system_iterator end = _systems.end();
1313 
1314  for (; pos != end; ++pos)
1315  tot += pos->second->n_dofs();
1316 
1317  return tot;
1318 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
std::map< std::string, System * > _systems
static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
the number of processors in the group.

Definition at line 93 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::size().

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::MeshCommunication::broadcast(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::MeshBase::get_info(), get_solution(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::DofMap::local_variable_indices(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::MeshTools::processor_bounding_box(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

94  { return cast_int<processor_id_type>(_communicator.size()); }
unsigned int size() const
Definition: parallel.h:679
const Parallel::Communicator & _communicator
unsigned int libMesh::EquationSystems::n_vars ( ) const
Returns
the total number of variables in all systems.

Definition at line 1292 of file equation_systems.C.

References _systems, and end.

Referenced by build_variable_names().

1293 {
1294  unsigned int tot=0;
1295 
1296  const_system_iterator pos = _systems.begin();
1297  const const_system_iterator end = _systems.end();
1298 
1299  for (; pos != end; ++pos)
1300  tot += pos->second->n_vars();
1301 
1302  return tot;
1303 }
IterBase * end
std::map< std::string, System * >::const_iterator const_system_iterator
std::map< std::string, System * > _systems
void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 88 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

Referenced by libMesh::LibMeshInit::LibMeshInit().

89 {
91 }
static std::string get_info()
void libMesh::EquationSystems::print_info ( std::ostream &  os = libMesh::out) const

Prints information about the equation systems, by default to libMesh::out.

Definition at line 1275 of file equation_systems.C.

References get_info().

Referenced by libMesh::operator<<(), and read().

1276 {
1277  os << this->get_info()
1278  << std::endl;
1279 }
virtual std::string get_info() const
processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
the rank of this processor in the group.

Definition at line 99 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::rank().

Referenced by libMesh::BoundaryInfo::_find_id_maps(), _read_impl(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshRefinement::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::MeshCommunication::broadcast(), build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::ParmetisPartitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DofMap::build_sparsity(), libMesh::DistributedMesh::clear(), libMesh::ExodusII_IO_Helper::close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::DistributedMesh::delete_node(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::DofMap::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), get_solution(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_element_variables(), libMesh::ExodusII_IO_Helper::initialize_global_variables(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::DistributedMesh::insert_elem(), libMesh::DofMap::is_evaluable(), libMesh::SparsityPattern::Build::join(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MetisPartitioner::partition_range(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NumericVector< T >::print_global(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::Parallel::Packing< Node * >::unpack(), libMesh::Parallel::Packing< Elem * >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO_Helper::write_element_values(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), libMesh::ExodusII_IO::write_information_records(), libMesh::ExodusII_IO_Helper::write_information_records(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), and libMesh::ExodusII_IO_Helper::write_timestep().

100  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
unsigned int rank() const
Definition: parallel.h:677
template<typename InValType >
void libMesh::EquationSystems::read ( const std::string &  name,
const XdrMODE  mode,
const unsigned int  read_flags = (READ_HEADER | READ_DATA),
bool  partition_agnostic = true 
)

Read & initialize the systems from disk using the XDR data format. This format allows for machine-independent binary output.

Set which sections of the file to read by bitwise OR'ing the EquationSystems::ReadFlags enumeration together. For example, to read all sections of the file, set read_flags to: (READ_HEADER | READ_DATA | READ_ADDITIONAL_DATA)

Note that the equation system can be defined without initializing the data vectors to any solution values. This can be done by omitting READ_DATA in the read_flags parameter.

If XdrMODE is omitted, it will be inferred as READ for filenames containing .xda or as DECODE for filenames containing .xdr

Parameters
nameName of the file to be read.
read_flagsSingle flag created by bitwise-OR'ing several flags together.
modeControls whether reading is done in binary or ascii mode.
partition_agnosticIf true then the mesh and degrees of freedom will be temporarily renumbered in a partition agnostic way so that files written using "n" mpi processes can be re-read on "m" mpi processes. Note that this renumbering is not compatible with meshes that have two nodes in exactly the same position!

Definition at line 90 of file equation_systems_io.C.

References _mesh, libMesh::MeshRefinement::clean_refinement_flags(), libMesh::Quality::name(), libMesh::out, and TRY_READ_IFEMS.

Referenced by _read_impl(), and read().

94 {
95  // If we have exceptions enabled we can be considerate and try
96  // to read old restart files which contain infinite element
97  // information but do not have the " with infinite elements"
98  // string in the version information.
99 
100  // First try the read the user requested
101  libmesh_try
102  {
103  this->_read_impl<InValType> (name, mode, read_flags, partition_agnostic);
104  }
105 
106  // If that fails, try it again but explicitly request we look for infinite element info
107  libmesh_catch (...)
108  {
109  libMesh::out << "\n*********************************************************************\n"
110  << "READING THE FILE \"" << name << "\" FAILED.\n"
111  << "It is possible this file contains infinite element information,\n"
112  << "but the version string does not contain \" with infinite elements\"\n"
113  << "Let's try this again, but looking for infinite element information...\n"
114  << "*********************************************************************\n"
115  << std::endl;
116 
117  libmesh_try
118  {
119  this->_read_impl<InValType> (name, mode, read_flags | EquationSystems::TRY_READ_IFEMS, partition_agnostic);
120  }
121 
122  // If all that failed, we are out of ideas here...
123  libmesh_catch (...)
124  {
125  libMesh::out << "\n*********************************************************************\n"
126  << "Well, at least we tried!\n"
127  << "Good Luck!!\n"
128  << "*********************************************************************\n"
129  << std::endl;
130  LIBMESH_THROW();
131  }
132  }
133 
134 #ifdef LIBMESH_ENABLE_AMR
135  MeshRefinement mesh_refine(_mesh);
136  mesh_refine.clean_refinement_flags();
137 #endif
138 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
OStreamProxy out(std::cout)
template void libMesh::EquationSystems::read< Real > ( const std::string &  name,
const XdrMODE  mode,
const unsigned int  read_flags = (READ_HEADER|READ_DATA),
bool  partition_agnostic = true 
)
inline

Definition at line 346 of file equation_systems.h.

References libMesh::Quality::name(), read(), READ_DATA, and READ_HEADER.

350  { read<Number>(name, mode, read_flags, partition_agnostic); }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
template<typename InValType >
void libMesh::EquationSystems::read ( const std::string &  name,
const unsigned int  read_flags = (READ_HEADER | READ_DATA),
bool  partition_agnostic = true 
)

Definition at line 72 of file equation_systems_io.C.

References _mesh, libMesh::MeshRefinement::clean_refinement_flags(), libMesh::DECODE, libMesh::READ, and read().

75 {
76  XdrMODE mode = READ;
77  if (name.find(".xdr") != std::string::npos)
78  mode = DECODE;
79  this->read(name, mode, read_flags, partition_agnostic);
80 
81 #ifdef LIBMESH_ENABLE_AMR
82  MeshRefinement mesh_refine(_mesh);
83  mesh_refine.clean_refinement_flags();
84 #endif
85 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
void read(const std::string &name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
template void libMesh::EquationSystems::read< Real > ( const std::string &  name,
const unsigned int  read_flags = (READ_HEADER|READ_DATA),
bool  partition_agnostic = true 
)
inline

Definition at line 357 of file equation_systems.h.

References allgather(), compare(), get_info(), get_mesh(), libMesh::Quality::name(), operator<<, libMesh::out, print_info(), libMesh::Real, write(), and WRITE_DATA.

360  { read<Number>(name, read_flags, partition_agnostic); }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
bool libMesh::EquationSystems::refine_in_reinit_flag ( )
inline
Returns
whether or not calls to reinit() will try to coarsen/refine the mesh

Definition at line 455 of file equation_systems.h.

References _refine_in_reinit.

456  { return this->_refine_in_reinit; }
void libMesh::EquationSystems::reinit ( )
virtual

Reinitialize all the systems

Definition at line 127 of file equation_systems.C.

References _mesh, _refine_in_reinit, libMesh::MeshRefinement::coarsen_elements(), libMesh::MeshBase::contract(), libMesh::DofMap::distribute_dofs(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), libMesh::MeshRefinement::face_level_mismatch_limit(), libMesh::System::get_dof_map(), get_mesh(), get_system(), n_systems(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), libMesh::MeshRefinement::overrefined_boundary_limit(), libMesh::System::prolong_vectors(), libMesh::MeshRefinement::refine_elements(), libMesh::System::reinit_constraints(), libMesh::System::restrict_vectors(), libMesh::DofObject::set_n_systems(), libMesh::sys, and libMesh::MeshRefinement::underrefined_boundary_limit().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), and libMesh::AdjointRefinementEstimator::estimate_error().

128 {
129  parallel_object_only();
130 
131  const unsigned int n_sys = this->n_systems();
132  libmesh_assert_not_equal_to (n_sys, 0);
133 
134  // We may have added new systems since our last
135  // EquationSystems::(re)init call
136  bool _added_new_systems = false;
137  for (unsigned int i=0; i != n_sys; ++i)
138  if (!this->get_system(i).is_initialized())
139  _added_new_systems = true;
140 
141  if (_added_new_systems)
142  {
143  // Our DofObjects will need space for the additional systems
144  MeshBase::node_iterator node_it = _mesh.nodes_begin();
145  const MeshBase::node_iterator node_end = _mesh.nodes_end();
146 
147  for ( ; node_it != node_end; ++node_it)
148  (*node_it)->set_n_systems(n_sys);
149 
150  MeshBase::element_iterator elem_it = _mesh.elements_begin();
151  const MeshBase::element_iterator elem_end = _mesh.elements_end();
152 
153  for ( ; elem_it != elem_end; ++elem_it)
154  (*elem_it)->set_n_systems(n_sys);
155 
156  // And any new systems will need initialization
157  for (unsigned int i=0; i != n_sys; ++i)
158  if (!this->get_system(i).is_initialized())
159  this->get_system(i).init();
160  }
161 
162 
163  // We used to assert that all nodes and elements *already* had
164  // n_systems() properly set; however this is false in the case where
165  // user code has manually added nodes and/or elements to an
166  // already-initialized system.
167 
168  // Make sure all the \p DofObject entities know how many systems
169  // there are.
170  {
171  // All the nodes
172  MeshBase::node_iterator node_it = _mesh.nodes_begin();
173  const MeshBase::node_iterator node_end = _mesh.nodes_end();
174 
175  for ( ; node_it != node_end; ++node_it)
176  {
177  Node * node = *node_it;
178  node->set_n_systems(this->n_systems());
179  }
180 
181  // All the elements
182  MeshBase::element_iterator elem_it = _mesh.elements_begin();
183  const MeshBase::element_iterator elem_end = _mesh.elements_end();
184 
185  for ( ; elem_it != elem_end; ++elem_it)
186  {
187  Elem * elem = *elem_it;
188  elem->set_n_systems(this->n_systems());
189  }
190  }
191 
192  // Localize each system's vectors
193  for (unsigned int i=0; i != this->n_systems(); ++i)
194  this->get_system(i).re_update();
195 
196 #ifdef LIBMESH_ENABLE_AMR
197 
198  bool dof_constraints_created = false;
199  bool mesh_changed = false;
200 
201  // FIXME: For backwards compatibility, assume
202  // refine_and_coarsen_elements or refine_uniformly have already
203  // been called
204  {
205  for (unsigned int i=0; i != this->n_systems(); ++i)
206  {
207  System & sys = this->get_system(i);
208 
209  // Even if the system doesn't have any variables in it we want
210  // consistent behavior; e.g. distribute_dofs should have the
211  // opportunity to count up zero dofs on each processor.
212  //
213  // Who's been adding zero-var systems anyway, outside of my
214  // unit tests? - RHS
215  // if(!sys.n_vars())
216  // continue;
217 
218  sys.get_dof_map().distribute_dofs(_mesh);
219 
220  // Recreate any user or internal constraints
221  sys.reinit_constraints();
222 
223  sys.prolong_vectors();
224  }
225  mesh_changed = true;
226  dof_constraints_created = true;
227  }
228 
229  if (this->_refine_in_reinit)
230  {
231  // Don't override any user refinement settings
232  MeshRefinement mesh_refine(_mesh);
233  mesh_refine.face_level_mismatch_limit() = 0; // unlimited
234  mesh_refine.overrefined_boundary_limit() = -1; // unlimited
235  mesh_refine.underrefined_boundary_limit() = -1; // unlimited
236 
237  // Try to coarsen the mesh, then restrict each system's vectors
238  // if necessary
239  if (mesh_refine.coarsen_elements())
240  {
241  for (unsigned int i=0; i != this->n_systems(); ++i)
242  {
243  System & sys = this->get_system(i);
244  if (!dof_constraints_created)
245  {
246  sys.get_dof_map().distribute_dofs(_mesh);
247  sys.reinit_constraints();
248  }
249  sys.restrict_vectors();
250  }
251  mesh_changed = true;
252  dof_constraints_created = true;
253  }
254 
255  // Once vectors are all restricted, we can delete
256  // children of coarsened elements
257  if (mesh_changed)
258  this->get_mesh().contract();
259 
260  // Try to refine the mesh, then prolong each system's vectors
261  // if necessary
262  if (mesh_refine.refine_elements())
263  {
264  for (unsigned int i=0; i != this->n_systems(); ++i)
265  {
266  System & sys = this->get_system(i);
267  if (!dof_constraints_created)
268  {
269  sys.get_dof_map().distribute_dofs(_mesh);
270  sys.reinit_constraints();
271  }
272  sys.prolong_vectors();
273  }
274  mesh_changed = true;
275  // dof_constraints_created = true;
276  }
277  }
278 
279  // If the mesh has changed, systems will need to create new dof
280  // constraints and update their global solution vectors
281  if (mesh_changed)
282  {
283  for (unsigned int i=0; i != this->n_systems(); ++i)
284  this->get_system(i).reinit();
285  }
286 #endif // #ifdef LIBMESH_ENABLE_AMR
287 }
ImplicitSystem & sys
virtual node_iterator nodes_begin()=0
virtual element_iterator elements_begin()=0
virtual element_iterator elements_end()=0
unsigned int n_systems() const
virtual node_iterator nodes_end()=0
virtual bool contract()=0
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::sensitivity_solve ( const ParameterVector parameters)
virtual

Call sensitivity_solve on all the individual equation systems.

By default this function solves each sensitivity system once, in the order in which in which they were added. For more sophisticated decoupled problems the user may with to override this behavior in a derived class.

Definition at line 458 of file equation_systems.C.

References get_system(), libMesh::libmesh_assert(), and n_systems().

459 {
460  libmesh_assert (this->n_systems());
461 
462  for (unsigned int i=0; i != this->n_systems(); ++i)
463  this->get_system(i).sensitivity_solve(parameters_in);
464 }
libmesh_assert(j)
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::solve ( )
virtual

Call solve on all the individual equation systems.

By default this function solves each equation system once, in the order they were added. For more sophisticated decoupled problems the user may with to override this behavior in a derived class.

Definition at line 448 of file equation_systems.C.

References get_system(), libMesh::libmesh_assert(), and n_systems().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error().

449 {
450  libmesh_assert (this->n_systems());
451 
452  for (unsigned int i=0; i != this->n_systems(); ++i)
453  this->get_system(i).solve();
454 }
libmesh_assert(j)
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::update ( )

Updates local values for all the systems

Definition at line 338 of file equation_systems.C.

References get_system(), and n_systems().

Referenced by _read_impl().

339 {
340  LOG_SCOPE("update()", "EquationSystems");
341 
342  // Localize each system's vectors
343  for (unsigned int i=0; i != this->n_systems(); ++i)
344  this->get_system(i).update();
345 }
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
void libMesh::EquationSystems::write ( const std::string &  name,
const XdrMODE  mode,
const unsigned int  write_flags = (WRITE_DATA),
bool  partition_agnostic = true 
) const

Write the systems to disk using the XDR data format. This format allows for machine-independent binary output.

Set the writing properties using the EquationSystems::WriteFlags enumeration. Set which sections to write out by bitwise OR'ing the enumeration values. Write everything by setting write_flags to: (WRITE_DATA | WRITE_ADDITIONAL_DATA)

Note that the solution data can be omitted by calling this routine with WRITE_DATA omitted in the write_flags argument.

If XdrMODE is omitted, it will be inferred as WRITE for filenames containing .xda or as ENCODE for filenames containing .xdr

Parameters
nameName of the file to be read.
write_flagsSingle flag created by bitwise-OR'ing several flags together.
modeControls whether reading is done in binary or ascii mode.
partition_agnosticIf true then the mesh and degrees of freedom will be temporarily renumbered in a partition agnostic way so that files written using "n" mpi processes can be re-read on "m" mpi processes. Note that this renumbering is not compatible with meshes that have two nodes in exactly the same position!

This program implements the output of an EquationSystems object. This warrants some documentation. The output file essentially consists of 11 sections:

1.) The version header.
2.) The number of individual equation systems (unsigned int)

for each system

3.)  The name of the system (string)
4.)  The type of the system (string)

handled through System::read():

+-------------------------------------------------------------+
|  5.) The number of variables in the system (unsigned int)   |
|                                                             |
|   for each variable in the system                           |
|                                                             |
|    6.) The name of the variable (string)                    |
|                                                             |
|    7.) Combined in an FEType:                               |
|         - The approximation order(s) of the variable (Order |
|           Enum, cast to int/s)                              |
|         - The finite element family/ies of the variable     |
|           (FEFamily Enum, cast to int/s)                    |
|                                                             |
|   end variable loop                                         |
|                                                             |
| 8.) The number of additional vectors (unsigned int),        |
|                                                             |
|    for each additional vector in the equation system object |
|                                                             |
|    9.) the name of the additional vector  (string)          |
+-------------------------------------------------------------+

end system loop


for each system, handled through System::write_{serialized,parallel}_data():

+--------------------------------------------------------------+
| 10.) The global solution vector, re-ordered to be node-major |
|     (More on this later.)                                    |
|                                                              |
|    for each additional vector in the equation system object  |
|                                                              |
|    11.) The global additional vector, re-ordered to be       |
|         node-major (More on this later.)                     |
+--------------------------------------------------------------+

end system loop

Note that the actual IO is handled through the Xdr class (to be renamed later?) which provides a uniform interface to both the XDR (eXternal Data Representation) interface and standard ASCII output. Thus this one section of code will write XDR or ASCII files with no changes.

Definition at line 378 of file equation_systems_io.C.

References _mesh, _systems, libMesh::Xdr::data(), libMesh::get_io_compatibility_version(), get_mesh(), libMesh::MeshTools::Private::globally_renumber_nodes_and_elements(), libMesh::libmesh_assert(), mesh, libMesh::Quality::name(), libMesh::ParallelObject::processor_id(), libMesh::Xdr::set_version(), WRITE_ADDITIONAL_DATA, WRITE_DATA, WRITE_PARALLEL_FILES, and libMesh::Xdr::writing().

Referenced by read(), write(), and libMesh::NameBasedIO::write_equation_systems().

382 {
446  // the EquationSystems::write() method should look constant,
447  // but we need to assign a temporary numbering to the nodes
448  // and elements in the mesh, which requires that we abuse const_cast
449  if(partition_agnostic)
450  {
451  MeshBase & mesh = const_cast<MeshBase &>(this->get_mesh());
453  }
454 
455  // set booleans from write_flags argument
456  const bool write_data = write_flags & EquationSystems::WRITE_DATA;
457  const bool write_additional_data = write_flags & EquationSystems::WRITE_ADDITIONAL_DATA;
458 
459  // always write parallel files if we're instructed to write in
460  // parallel
461  const bool write_parallel_files =
463  // Even if we're on a distributed mesh, we may or may not have a
464  // consistent way of reconstructing the same mesh partitioning
465  // later, but we need the same mesh partitioning if we want to
466  // reread the parallel solution safely, so let's write a serial file
467  // unless specifically requested not to.
468  // ||
469  // // but also write parallel files if we haven't been instructed to
470  // // write in serial and we're on a distributed mesh
471  // (!(write_flags & EquationSystems::WRITE_SERIAL_FILES) &&
472  // !this->get_mesh().is_serial())
473  ;
474 
475  // New scope so that io will close before we try to zip the file
476  {
477  Xdr io((this->processor_id()==0) ? name : "", mode);
478  libmesh_assert (io.writing());
479 
480  LOG_SCOPE("write()", "EquationSystems");
481 
482  const unsigned int proc_id = this->processor_id();
483 
484  unsigned int n_sys = 0;
485  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
486  pos != _systems.end(); ++pos)
487  {
488  if (! pos->second->hide_output()) n_sys++;
489  }
490 
491  // set the version number in the Xdr object
492  io.set_version(LIBMESH_VERSION_ID(LIBMESH_MAJOR_VERSION,
493  LIBMESH_MINOR_VERSION,
494  LIBMESH_MICRO_VERSION));
495 
496  // Only write the header information
497  // if we are processor 0.
498  if (proc_id == 0)
499  {
500  std::string comment;
501  char buf[256];
502 
503  // 1.)
504  // Write the version header
505  std::string version("libMesh-" + libMesh::get_io_compatibility_version());
506  if (write_parallel_files) version += " parallel";
507 
508 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
509  version += " with infinite elements";
510 #endif
511  io.data (version, "# File Format Identifier");
512 
513  // 2.)
514  // Write the number of equation systems
515  io.data (n_sys, "# No. of Equation Systems");
516 
517  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
518  pos != _systems.end(); ++pos)
519  {
520  // Ignore this system if it has been marked as hidden
521  if (pos->second->hide_output()) continue;
522 
523  // 3.)
524  // Write the name of the sys_num-th system
525  {
526  const unsigned int sys_num = pos->second->number();
527  std::string sys_name = pos->first;
528 
529  comment = "# Name, System No. ";
530  std::sprintf(buf, "%u", sys_num);
531  comment += buf;
532 
533  io.data (sys_name, comment.c_str());
534  }
535 
536  // 4.)
537  // Write the type of system handled
538  {
539  const unsigned int sys_num = pos->second->number();
540  std::string sys_type = pos->second->system_type();
541 
542  comment = "# Type, System No. ";
543  std::sprintf(buf, "%u", sys_num);
544  comment += buf;
545 
546  io.data (sys_type, comment.c_str());
547  }
548 
549  // 5.) - 9.)
550  // Let System::write_header() do the job
551  pos->second->write_header (io, version, write_additional_data);
552  }
553  }
554 
555  // Start from the first system, again,
556  // to write vectors to disk, if wanted
557  if (write_data)
558  {
559  // open a parallel buffer if warranted.
560  Xdr local_io (write_parallel_files ? local_file_name(this->processor_id(),name) : "", mode);
561 
562  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
563  pos != _systems.end(); ++pos)
564  {
565  // Ignore this system if it has been marked as hidden
566  if (pos->second->hide_output()) continue;
567 
568  // 10.) + 11.)
569  if (write_parallel_files)
570  pos->second->write_parallel_data (local_io,write_additional_data);
571  else
572  pos->second->write_serialized_data (io,write_additional_data);
573  }
574  }
575  }
576 
577  // the EquationSystems::write() method should look constant,
578  // but we need to undo the temporary numbering of the nodes
579  // and elements in the mesh, which requires that we abuse const_cast
580  if(partition_agnostic)
581  const_cast<MeshBase &>(_mesh).fix_broken_node_and_element_numbering();
582 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
MeshBase & mesh
libmesh_assert(j)
std::string get_io_compatibility_version()
void globally_renumber_nodes_and_elements(MeshBase &)
Definition: mesh_tools.C:1948
const MeshBase & get_mesh() const
processor_id_type processor_id() const
std::map< std::string, System * > _systems
void libMesh::EquationSystems::write ( const std::string &  name,
const unsigned int  write_flags = (WRITE_DATA),
bool  partition_agnostic = true 
) const

Definition at line 366 of file equation_systems_io.C.

References libMesh::ENCODE, libMesh::WRITE, and write().

369 {
370  XdrMODE mode = WRITE;
371  if (name.find(".xdr") != std::string::npos)
372  mode = ENCODE;
373  this->write(name, mode, write_flags, partition_agnostic);
374 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
void write(const std::string &name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const

Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  os,
const EquationSystems es 
)
friend

Same as above, but allows you to also use stream syntax.

Definition at line 1283 of file equation_systems.C.

Referenced by read().

1285 {
1286  es.print_info(os);
1287  return os;
1288 }

Member Data Documentation

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited
bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 134 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

MeshBase& libMesh::EquationSystems::_mesh
protected
Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 128 of file reference_counter.h.

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects. Print the reference count information when the number returns to 0.

Definition at line 123 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

bool libMesh::EquationSystems::_refine_in_reinit
protected

Flag for whether to call coarsen/refine in reinit(). Default value: true

Definition at line 492 of file equation_systems.h.

Referenced by disable_refine_in_reinit(), enable_refine_in_reinit(), EquationSystems(), refine_in_reinit_flag(), and reinit().

std::map<std::string, System *> libMesh::EquationSystems::_systems
protected

The documentation for this class was generated from the following files: