equation_systems.h

Go to the documentation of this file.

// The libMesh Finite Element Library.
// Copyright (C) 2002-2018 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner

// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.

// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.

// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA



#ifndef LIBMESH_EQUATION_SYSTEMS_H
#define LIBMESH_EQUATION_SYSTEMS_H

// Local Includes
#include "libmesh/libmesh_common.h"
#include "libmesh/parameters.h"
#include "libmesh/system.h"
#include "libmesh/parallel_object.h"

#ifdef LIBMESH_FORWARD_DECLARE_ENUMS
namespace libMesh
{
enum XdrMODE : int;
}
#else
#include "libmesh/enum_xdr_mode.h"
#endif

// HP aCC needs these for some reason
#ifdef __HP_aCC
# include "libmesh/frequency_system.h"
# include "libmesh/transient_system.h"
# include "libmesh/newmark_system.h"
# include "libmesh/steady_system.h"
#endif

// C++ includes
#include <cstddef>
#include <map>
#include <set>
#include <string>
#include <vector>

namespace libMesh
{

// Forward Declarations
class Elem;
class MeshBase;

class EquationSystems : public ReferenceCountedObject<EquationSystems>,
                        public ParallelObject

{
public:

  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 };

  EquationSystems (MeshBase & mesh);

  virtual ~EquationSystems ();

  virtual void clear ();

  virtual void init ();

  virtual void reinit ();

  virtual void enable_default_ghosting (bool enable);

  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 System & get_system (const std::string & name) const;

  System & get_system (const std::string & name);

  const System & get_system (const unsigned int num) const;

  System & get_system (const unsigned int num);

  virtual System & add_system (const std::string & system_type,
                               const std::string & name);

  template <typename T_sys>
  T_sys & add_system (const std::string & name);

#ifdef LIBMESH_ENABLE_DEPRECATED
  void delete_system (const std::string & name);
#endif

  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=nullptr,
                             const std::set<std::string> * system_names=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=nullptr) const;

  std::unique_ptr<NumericVector<Number>>
  build_parallel_solution_vector(const std::set<std::string> * system_names=nullptr) const;

  void get_vars_active_subdomains(const std::vector<std::string> & names,
                                  std::vector<std::set<subdomain_id_type>> & vars_active_subdomains) const;

  void get_solution (std::vector<Number> & soln,
                     std::vector<std::string> & names) const;

  void build_elemental_solution_vector (std::vector<Number> & soln,
                                        std::vector<std::string> & names) const;

  std::unique_ptr<NumericVector<Number>>
  build_parallel_elemental_solution_vector (std::vector<std::string> & names) const;

  void build_discontinuous_solution_vector (std::vector<Number> & soln,
                                            const std::set<std::string> * system_names=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)
  { read<Number>(name, mode, read_flags, partition_agnostic); }

  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)
  { read<Number>(name, read_flags, partition_agnostic); }


  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;

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

  const MeshBase & get_mesh() const;

  MeshBase & get_mesh();

  void allgather ();

  void enable_refine_in_reinit() { this->_refine_in_reinit = true; }

  void disable_refine_in_reinit() { this->_refine_in_reinit = false; }

  bool refine_in_reinit_flag() { return this->_refine_in_reinit; }

  bool reinit_solutions ();

  virtual void reinit_systems ();

  Parameters parameters;


protected:


  MeshBase & _mesh;

  std::map<std::string, System *> _systems;

  typedef std::map<std::string, System *>::iterator       system_iterator;

  typedef std::map<std::string, System *>::const_iterator const_system_iterator;

  bool _refine_in_reinit;

  bool _enable_default_ghosting;

private:

  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();

  void _remove_default_ghosting(unsigned int sys_num);
};



// ------------------------------------------------------------
// EquationSystems inline methods
inline
const MeshBase & EquationSystems::get_mesh () const
{
  return _mesh;
}



inline
MeshBase & EquationSystems::get_mesh ()
{
  return _mesh;
}


inline
unsigned int EquationSystems::n_systems () const
{
  return cast_int<unsigned int>(_systems.size());
}




template <typename T_sys>
inline
T_sys & EquationSystems::add_system (const std::string & name)
{
  T_sys * ptr = nullptr;

  if (!_systems.count(name))
    {
      const unsigned int sys_num = this->n_systems();
      ptr = new T_sys(*this, name, sys_num);

      _systems.insert (std::make_pair(name, ptr));

      if (!_enable_default_ghosting)
        this->_remove_default_ghosting(sys_num);

      // Tell all the \p DofObject entities to add a system.
      this->_add_system_to_nodes_and_elems();
    }
  else
    {
      // We now allow redundant add_system calls, to make it
      // easier to load data from files for user-derived system
      // subclasses
      ptr = &(this->get_system<T_sys>(name));
    }

  // Return a dynamically casted reference to the newly added System.
  return *ptr;
}



inline
bool EquationSystems::has_system (const std::string & name) const
{
  if (_systems.find(name) == _systems.end())
    return false;
  return true;
}




template <typename T_sys>
inline
const T_sys & EquationSystems::get_system (const unsigned int num) const
{
  libmesh_assert_less (num, this->n_systems());


  const_system_iterator       pos = _systems.begin();
  const const_system_iterator end = _systems.end();

  for (; pos != end; ++pos)
    if (pos->second->number() == num)
      break;

  // Check for errors
  if (pos == end)
    libmesh_error_msg("ERROR: no system number " << num << " found!");

  // Attempt dynamic cast
  return *cast_ptr<T_sys *>(pos->second);
}




template <typename T_sys>
inline
T_sys & EquationSystems::get_system (const unsigned int num)
{
  libmesh_assert_less (num, this->n_systems());

  const_system_iterator       pos = _systems.begin();
  const const_system_iterator end = _systems.end();

  for (; pos != end; ++pos)
    if (pos->second->number() == num)
      break;

  // Check for errors
  if (pos == end)
    libmesh_error_msg("ERROR: no system number " << num << " found!");

  // Attempt dynamic cast
  return *cast_ptr<T_sys *>(pos->second);
}






template <typename T_sys>
inline
const T_sys & EquationSystems::get_system (const std::string & name) const
{
  const_system_iterator pos = _systems.find(name);

  // Check for errors
  if (pos == _systems.end())
    libmesh_error_msg("ERROR: no system named \"" << name << "\" found!");

  // Attempt dynamic cast
  return *cast_ptr<T_sys *>(pos->second);
}






template <typename T_sys>
inline
T_sys & EquationSystems::get_system (const std::string & name)
{
  system_iterator pos = _systems.find(name);

  // Check for errors
  if (pos == _systems.end())
    libmesh_error_msg("ERROR: no system named " << name << " found!");

  // Attempt dynamic cast
  return *cast_ptr<T_sys *>(pos->second);
}







inline
const System & EquationSystems::get_system (const std::string & name) const
{
  return this->get_system<System>(name);
}



inline
System & EquationSystems::get_system (const std::string & name)
{
  return this->get_system<System>(name);
}



inline
const System & EquationSystems::get_system (const unsigned int num) const
{
  return this->get_system<System>(num);
}



inline
System & EquationSystems::get_system (const unsigned int num)
{
  return this->get_system<System>(num);
}


} // namespace libMesh


#endif // LIBMESH_EQUATION_SYSTEMS_H