libMesh
Functions
amr.C File Reference

Go to the source code of this file.

Functions

void assemble (EquationSystems &es, const std::string &system_name)
 
int main (int argc, char **argv)
 
void assemble (EquationSystems &es, const std::string &libmesh_dbg_var(system_name))
 

Function Documentation

void assemble ( EquationSystems es,
const std::string &  system_name 
)

Definition at line 232 of file miscellaneous_ex4.C.

References libMesh::MeshBase::active_local_element_ptr_range(), libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::NumericVector< T >::close(), libMesh::SparseMatrix< T >::close(), libMesh::DofMap::constrain_element_dyad_matrix(), libMesh::DofMap::constrain_element_matrix_and_vector(), dim, libMesh::DofMap::dof_indices(), libMesh::System::get_dof_map(), libMesh::ImplicitSystem::get_matrix(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::System::get_vector(), libMesh::libmesh_ignore(), libmesh_nullptr, libMesh::ImplicitSystem::matrix, mesh, libMesh::MeshBase::mesh_dimension(), libMesh::Real, libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, and libMesh::System::variable_type().

Referenced by main().

234 {
235  // Ignore unused parameter warnings when !LIBMESH_ENABLE_AMR.
236  libmesh_ignore(es);
237  libmesh_ignore(system_name);
238 
239 #ifdef LIBMESH_ENABLE_AMR
240  // It is a good idea to make sure we are assembling
241  // the proper system.
242  libmesh_assert_equal_to (system_name, "System");
243 
244  // Get a constant reference to the mesh object.
245  const MeshBase & mesh = es.get_mesh();
246 
247  // The dimension that we are running
248  const unsigned int dim = mesh.mesh_dimension();
249 
250  // Get a reference to the Convection-Diffusion system object.
251  LinearImplicitSystem & system =
252  es.get_system<LinearImplicitSystem> ("System");
253 
254  // Get the Finite Element type for the first (and only)
255  // variable in the system.
256  FEType fe_type = system.variable_type(0);
257 
258  // Build a Finite Element object of the specified type. Since the
259  // FEBase::build() member dynamically creates memory we will
260  // store the object as a UniquePtr<FEBase>. This can be thought
261  // of as a pointer that will clean up after itself.
262  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
263  UniquePtr<FEBase> fe_face (FEBase::build(dim, fe_type));
264 
265  // A Gauss quadrature rule for numerical integration.
266  // Let the FEType object decide what order rule is appropriate.
267  QGauss qrule (dim, fe_type.default_quadrature_order());
268  QGauss qface (dim-1, fe_type.default_quadrature_order());
269 
270  // Tell the finite element object to use our quadrature rule.
271  fe->attach_quadrature_rule (&qrule);
272  fe_face->attach_quadrature_rule (&qface);
273 
274  // Here we define some references to cell-specific data that
275  // will be used to assemble the linear system. We will start
276  // with the element Jacobian * quadrature weight at each integration point.
277  const std::vector<Real> & JxW = fe->get_JxW();
278  const std::vector<Real> & JxW_face = fe_face->get_JxW();
279 
280  // The element shape functions evaluated at the quadrature points.
281  const std::vector<std::vector<Real>> & phi = fe->get_phi();
282  const std::vector<std::vector<Real>> & psi = fe_face->get_phi();
283 
284  // The element shape function gradients evaluated at the quadrature
285  // points.
286  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
287 
288  // The XY locations of the quadrature points used for face integration
289  //const std::vector<Point>& qface_points = fe_face->get_xyz();
290 
291  // A reference to the DofMap object for this system. The DofMap
292  // object handles the index translation from node and element numbers
293  // to degree of freedom numbers. We will talk more about the DofMap
294  // in future examples.
295  const DofMap & dof_map = system.get_dof_map();
296 
297  // Define data structures to contain the element matrix
298  // and right-hand-side vector contribution. Following
299  // basic finite element terminology we will denote these
300  // "Ke" and "Fe".
303 
304  // Analogous data structures for thw two vectors v and w that form
305  // the tensor shell matrix.
308 
309  // This vector will hold the degree of freedom indices for
310  // the element. These define where in the global system
311  // the element degrees of freedom get mapped.
312  std::vector<dof_id_type> dof_indices;
313 
314  // Now we will loop over all the elements in the mesh that
315  // live on the local processor. We will compute the element
316  // matrix and right-hand-side contribution. Since the mesh
317  // will be refined we want to only consider the ACTIVE elements,
318  // hence we use a variant of the active_elem_iterator.
319  for (const auto & elem : mesh.active_local_element_ptr_range())
320  {
321  // Get the degree of freedom indices for the
322  // current element. These define where in the global
323  // matrix and right-hand-side this element will
324  // contribute to.
325  dof_map.dof_indices (elem, dof_indices);
326 
327  // Compute the element-specific data for the current
328  // element. This involves computing the location of the
329  // quadrature points (q_point) and the shape functions
330  // (phi, dphi) for the current element.
331  fe->reinit (elem);
332 
333  // Zero the element matrix and right-hand side before
334  // summing them. We use the resize member here because
335  // the number of degrees of freedom might have changed from
336  // the last element. Note that this will be the case if the
337  // element type is different (i.e. the last element was a
338  // triangle, now we are on a quadrilateral).
339  Ke.resize (dof_indices.size(),
340  dof_indices.size());
341 
342  Fe.resize (dof_indices.size());
343  Ve.resize (dof_indices.size());
344  We.resize (dof_indices.size());
345 
346  // Now we will build the element matrix and right-hand-side.
347  // Constructing the RHS requires the solution and its
348  // gradient from the previous timestep. This myst be
349  // calculated at each quadrature point by summing the
350  // solution degree-of-freedom values by the appropriate
351  // weight functions.
352  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
353  {
354  // Now compute the element matrix and RHS contributions.
355  for (std::size_t i=0; i<phi.size(); i++)
356  {
357  // The RHS contribution
358  Fe(i) += JxW[qp]*phi[i][qp];
359 
360  for (std::size_t j=0; j<phi.size(); j++)
361  {
362  // The matrix contribution
363  Ke(i,j) += JxW[qp]*(
364  // Stiffness matrix
365  (dphi[i][qp]*dphi[j][qp])
366  );
367  }
368 
369  // V and W are the same for this example.
370  Ve(i) += JxW[qp]*phi[i][qp];
371  We(i) += JxW[qp]*phi[i][qp];
372  }
373  }
374 
375  // At this point the interior element integration has
376  // been completed. However, we have not yet addressed
377  // boundary conditions. For this example we will only
378  // consider simple Dirichlet boundary conditions imposed
379  // via the penalty method.
380  //
381  // The following loops over the sides of the element.
382  // If the element has no neighbor on a side then that
383  // side MUST live on a boundary of the domain.
384  {
385  // The penalty value.
386  const Real penalty = 1.e10;
387 
388  // The following loops over the sides of the element.
389  // If the element has no neighbor on a side then that
390  // side MUST live on a boundary of the domain.
391  for (auto s : elem->side_index_range())
392  if (elem->neighbor_ptr(s) == libmesh_nullptr)
393  {
394  fe_face->reinit(elem, s);
395 
396  for (unsigned int qp=0; qp<qface.n_points(); qp++)
397  {
398  // Matrix contribution
399  for (std::size_t i=0; i<psi.size(); i++)
400  for (std::size_t j=0; j<psi.size(); j++)
401  Ke(i,j) += penalty*JxW_face[qp]*psi[i][qp]*psi[j][qp];
402  }
403  }
404  }
405 
406 
407  // We have now built the element matrix and RHS vector in terms
408  // of the element degrees of freedom. However, it is possible
409  // that some of the element DOFs are constrained to enforce
410  // solution continuity, i.e. they are not really "free". We need
411  // to constrain those DOFs in terms of non-constrained DOFs to
412  // ensure a continuous solution. The
413  // DofMap::constrain_element_matrix_and_vector() method does
414  // just that.
415 
416  // However, constraining both the sparse matrix (and right hand
417  // side) plus the rank 1 matrix is tricky. The dof_indices
418  // vector has to be backed up for that because the constraining
419  // functions modify it.
420 
421  std::vector<dof_id_type> dof_indices_backup(dof_indices);
422  dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
423  dof_indices = dof_indices_backup;
424  dof_map.constrain_element_dyad_matrix(Ve, We, dof_indices);
425 
426  // The element matrix and right-hand-side are now built
427  // for this element. Add them to the global matrix and
428  // right-hand-side vector. The SparseMatrix::add_matrix()
429  // and NumericVector::add_vector() members do this for us.
430  system.matrix->add_matrix (Ke, dof_indices);
431  system.get_matrix("Preconditioner").add_matrix (Ke, dof_indices);
432  system.rhs->add_vector (Fe, dof_indices);
433  system.get_vector("v").add_vector(Ve, dof_indices);
434  system.get_vector("w").add_vector(We, dof_indices);
435  }
436  // Finished computing the system matrix and right-hand side.
437 
438  // Matrices and vectors must be closed manually. This is necessary
439  // because the matrix is not directly used as the system matrix (in
440  // which case the solver closes it) but as a part of a shell matrix.
441  system.matrix->close();
442  system.get_matrix("Preconditioner").close();
443  system.rhs->close();
444  system.get_vector("v").close();
445  system.get_vector("w").close();
446 
447 #endif // #ifdef LIBMESH_ENABLE_AMR
448 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
unsigned int dim
This class provides a specific system class.
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
NumericVector< Number > * rhs
The system matrix.
const FEType & variable_type(const unsigned int i) const
Definition: system.h:2164
This is the MeshBase class.
Definition: mesh_base.h:68
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
const DofMap & get_dof_map() const
Definition: system.h:2030
void constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix and vector.
Definition: dof_map.h:1806
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:794
virtual void close()=0
Calls the NumericVector&#39;s internal assembly routines, ensuring that the values are consistent across ...
const SparseMatrix< Number > & get_matrix(const std::string &mat_name) const
void libmesh_ignore(const T &)
virtual void close()=0
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
SparseMatrix< Number > * matrix
The system matrix.
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
This class implements specific orders of Gauss quadrature.
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
void constrain_element_dyad_matrix(DenseVector< Number > &v, DenseVector< Number > &w, std::vector< dof_id_type > &row_dofs, bool asymmetric_constraint_rows=true) const
Constrains a dyadic element matrix B = v w&#39;.
Definition: dof_map.h:1811
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
void assemble ( EquationSystems es,
const std::string &  libmesh_dbg_varsystem_name 
)

Definition at line 131 of file amr.C.

References libMesh::MeshBase::active_local_element_ptr_range(), libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::constrain_element_matrix_and_vector(), dim, libMesh::DofMap::dof_indices(), libMesh::FIFTH, libMesh::FIRST, libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libmesh_nullptr, libMesh::ImplicitSystem::matrix, mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::out, libMesh::Real, libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, and side.

133 {
134  libmesh_assert_equal_to (system_name, "primary");
135 
136  const MeshBase & mesh = es.get_mesh();
137  const unsigned int dim = mesh.mesh_dimension();
138 
139  // Also use a 3x3x3 quadrature rule (3D). Then tell the FE
140  // about the geometry of the problem and the quadrature rule
141  FEType fe_type (FIRST);
142 
143  UniquePtr<FEBase> fe(FEBase::build(dim, fe_type));
144  QGauss qrule(dim, FIFTH);
145 
146  fe->attach_quadrature_rule (&qrule);
147 
148  UniquePtr<FEBase> fe_face(FEBase::build(dim, fe_type));
149  QGauss qface(dim-1, FIFTH);
150 
151  fe_face->attach_quadrature_rule(&qface);
152 
153  LinearImplicitSystem & system =
154  es.get_system<LinearImplicitSystem>("primary");
155 
156 
157  // These are references to cell-specific data
158  const std::vector<Real> & JxW_face = fe_face->get_JxW();
159  const std::vector<Real> & JxW = fe->get_JxW();
160  const std::vector<Point> & q_point = fe->get_xyz();
161  const std::vector<std::vector<Real>> & phi = fe->get_phi();
162  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
163 
164  std::vector<dof_id_type> dof_indices_U;
165  std::vector<dof_id_type> dof_indices_V;
166  const DofMap & dof_map = system.get_dof_map();
167 
172 
173  Real vol=0., area=0.;
174 
175  for (const auto & elem : mesh.active_local_element_ptr_range())
176  {
177  // recompute the element-specific data for the current element
178  fe->reinit (elem);
179 
180 
181  //fe->print_info();
182 
183  dof_map.dof_indices(elem, dof_indices_U, 0);
184  dof_map.dof_indices(elem, dof_indices_V, 1);
185 
186  // zero the element matrix and vector
187  Kuu.resize (phi.size(),
188  phi.size());
189 
190  Kvv.resize (phi.size(),
191  phi.size());
192 
193  Fu.resize (phi.size());
194  Fv.resize (phi.size());
195 
196  // standard stuff... like in code 1.
197  for (unsigned int gp=0; gp<qrule.n_points(); gp++)
198  {
199  for (std::size_t i=0; i<phi.size(); ++i)
200  {
201  // this is tricky. ig is the _global_ dof index corresponding
202  // to the _global_ vertex number elem->node_id(i). Note that
203  // in general these numbers will not be the same (except for
204  // the case of one unknown per node on one subdomain) so
205  // we need to go through the dof_map
206 
207  const Real f = q_point[gp]*q_point[gp];
208  // const Real f = (q_point[gp](0) +
209  // q_point[gp](1) +
210  // q_point[gp](2));
211 
212  // add jac*weight*f*phi to the RHS in position ig
213 
214  Fu(i) += JxW[gp]*f*phi[i][gp];
215  Fv(i) += JxW[gp]*f*phi[i][gp];
216 
217  for (std::size_t j=0; j<phi.size(); ++j)
218  {
219 
220  Kuu(i,j) += JxW[gp]*((phi[i][gp])*(phi[j][gp]));
221 
222  Kvv(i,j) += JxW[gp]*((phi[i][gp])*(phi[j][gp]) +
223  1.*((dphi[i][gp])*(dphi[j][gp])));
224  };
225  };
226  vol += JxW[gp];
227  };
228 
229 
230  // You can't compute "area" (perimeter) if you are in 2D
231  if (dim == 3)
232  {
233  for (auto side : elem->side_index_range())
234  if (elem->neighbor_ptr(side) == libmesh_nullptr)
235  {
236  fe_face->reinit (elem, side);
237 
238  //fe_face->print_info();
239 
240  for (std::size_t gp=0; gp<JxW_face.size(); gp++)
241  area += JxW_face[gp];
242  }
243  }
244 
245  // Constrain the DOF indices.
246  dof_map.constrain_element_matrix_and_vector(Kuu, Fu, dof_indices_U);
247  dof_map.constrain_element_matrix_and_vector(Kvv, Fv, dof_indices_V);
248 
249 
250  system.rhs->add_vector(Fu, dof_indices_U);
251  system.rhs->add_vector(Fv, dof_indices_V);
252 
253  system.matrix->add_matrix(Kuu, dof_indices_U);
254  system.matrix->add_matrix(Kvv, dof_indices_V);
255  }
256 
257  libMesh::out << "Vol=" << vol << std::endl;
258 
259  if (dim == 3)
260  libMesh::out << "Area=" << area << std::endl;
261 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
unsigned int dim
This class provides a specific system class.
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
unsigned short int side
Definition: xdr_io.C:49
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
NumericVector< Number > * rhs
The system matrix.
This is the MeshBase class.
Definition: mesh_base.h:68
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
const DofMap & get_dof_map() const
Definition: system.h:2030
void constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix and vector.
Definition: dof_map.h:1806
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OStreamProxy out
SparseMatrix< Number > * matrix
The system matrix.
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
This class implements specific orders of Gauss quadrature.
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
int main ( int  argc,
char **  argv 
)

Definition at line 46 of file amr.C.

References libMesh::DofMap::_dof_coupling, libMesh::EquationSystems::add_system(), libMesh::System::add_variable(), assemble(), libMesh::System::attach_assemble_function(), libMesh::LibMeshInit::comm(), dim, libMesh::MeshBase::elem_ref(), libMesh::FIRST, libMesh::System::get_dof_map(), libMesh::TriangleWrapper::init(), libMesh::EquationSystems::init(), mesh, libMesh::DofMap::print_dof_constraints(), libMesh::EquationSystems::print_info(), libMesh::MeshBase::print_info(), libMesh::MeshBase::read(), libMesh::Elem::REFINE, libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::EquationSystems::reinit(), libMesh::CouplingMatrix::resize(), libMesh::Elem::set_refinement_flag(), libMesh::LinearImplicitSystem::solve(), libMesh::MeshRefinement::uniformly_refine(), libMesh::GMVIO::write(), and libMesh::MeshOutput< MT >::write_equation_systems().

47 {
48  LibMeshInit init(argc, argv);
49 
50  if (argc < 4)
51  libmesh_error_msg("Usage: ./prog -d DIM filename");
52 
53  // Variables to get us started
54  const unsigned int dim = atoi(argv[2]);
55 
56  std::string meshname (argv[3]);
57 
58  // declare a mesh...
59  Mesh mesh(init.comm(), dim);
60 
61  // Read a mesh
62  mesh.read(meshname);
63 
64  GMVIO(mesh).write ("out_0.gmv");
65 
66  mesh.elem_ref(0).set_refinement_flag (Elem::REFINE);
67 
68  MeshRefinement mesh_refinement (mesh);
69 
70  mesh_refinement.refine_and_coarsen_elements ();
71  mesh_refinement.uniformly_refine (2);
72 
73  mesh.print_info();
74 
75 
76  // Set up the equation system(s)
77  EquationSystems es (mesh);
78 
79  LinearImplicitSystem & primary =
80  es.add_system<LinearImplicitSystem>("primary");
81 
82  primary.add_variable ("U", FIRST);
83  primary.add_variable ("V", FIRST);
84 
85  primary.get_dof_map()._dof_coupling->resize(2);
86  (*primary.get_dof_map()._dof_coupling)(0,0) = 1;
87  (*primary.get_dof_map()._dof_coupling)(1,1) = 1;
88 
90 
91  es.init ();
92 
93  es.print_info ();
95 
96  // call the solver.
97  primary.solve ();
98 
99  GMVIO(mesh).write_equation_systems ("out_1.gmv",
100  es);
101 
102 
103 
104  // Refine uniformly
105  mesh_refinement.uniformly_refine (1);
106  es.reinit ();
107 
108  // Write out the projected solution
109  GMVIO(mesh).write_equation_systems ("out_2.gmv",
110  es);
111 
112  // Solve again. Output the refined solution
113  primary.solve ();
114  GMVIO(mesh).write_equation_systems ("out_3.gmv",
115  es);
116 
117  return 0;
118 }
This is the EquationSystems class.
unsigned int dim
This class provides a specific system class.
unsigned int add_variable(const std::string &var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=libmesh_nullptr)
Adds the variable var to the list of variables for this system.
Definition: system.C:1101
MeshBase & mesh
void set_refinement_flag(const RefinementState rflag)
Sets the value of the refinement flag for the element.
Definition: elem.h:2513
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
This class implements writing meshes in the GMV format.
Definition: gmv_io.h:46
void assemble(EquationSystems &es, const std::string &system_name)
void print_dof_constraints(std::ostream &os=libMesh::out, bool print_nonlocal=false) const
Prints (from processor 0) all DoF and Node constraints.
const DofMap & get_dof_map() const
Definition: system.h:2030
This is the MeshRefinement class.
void attach_assemble_function(void fptr(EquationSystems &es, const std::string &name))
Register a user function to use in assembling the system matrix and RHS.
Definition: system.C:1795
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
CouplingMatrix * _dof_coupling
Degree of freedom coupling.
Definition: dof_map.h:1240
virtual void write(const std::string &) libmesh_override
This method implements writing a mesh to a specified file.
Definition: gmv_io.C:267
virtual void write_equation_systems(const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
This method implements writing a mesh with data to a specified file where the data is taken from the ...
Definition: mesh_output.C:31
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:490
void resize(const unsigned int n)
Resizes the matrix and initializes all entries to be 0.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
virtual void read(const std::string &name, void *mesh_data=libmesh_nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false)=0
Interfaces for reading/writing a mesh to/from a file.
virtual void solve() libmesh_override
Assembles & solves the linear system A*x=b.
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448