libMesh
miscellaneous_ex3.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // <h1>Miscellaneous Example 3 - 2D Laplace-Young Problem Using Nonlinear Solvers</h1>
21 // \author Derek Gaston
22 // \date 2008
23 //
24 // This example shows how to use the NonlinearImplicitSystem class to
25 // solve nonlinear problems in libMesh. The NonlinearImplicitSystem
26 // class employs the user's ComputeResidual and ComputeJacobian
27 // objects during the solve. In this particular example, the
28 // LaplaceYoung object provides this functionality.
29 //
30 // You can turn on preconditioning of the matrix-free system using the
31 // jacobian by passing "-pre" on the command line. Currently, this
32 // feature only works with Petsc, so this isn't used by "make run".
33 //
34 // This example also runs with the experimental Trilinos NOX solvers
35 // by specifying the --use-trilinos command line argument.
36 
37 
38 // C++ include files that we need
39 #include <iostream>
40 #include <algorithm>
41 #include <cmath>
42 
43 // Various include files needed for the mesh & solver functionality.
44 #include "libmesh/libmesh.h"
45 #include "libmesh/mesh.h"
46 #include "libmesh/mesh_refinement.h"
47 #include "libmesh/exodusII_io.h"
48 #include "libmesh/equation_systems.h"
49 #include "libmesh/fe.h"
50 #include "libmesh/quadrature_gauss.h"
51 #include "libmesh/dof_map.h"
52 #include "libmesh/sparse_matrix.h"
53 #include "libmesh/numeric_vector.h"
54 #include "libmesh/dense_matrix.h"
55 #include "libmesh/dense_vector.h"
56 #include "libmesh/elem.h"
57 #include "libmesh/string_to_enum.h"
58 #include "libmesh/getpot.h"
59 
60 // The nonlinear solver and system we will be using
61 #include "libmesh/nonlinear_solver.h"
62 #include "libmesh/nonlinear_implicit_system.h"
63 
64 // Necessary for programmatically setting petsc options
65 #ifdef LIBMESH_HAVE_PETSC
66 #include <petsc.h>
67 #include "libmesh/petsc_macro.h"
68 #endif
69 
70 // Bring in everything from the libMesh namespace
71 using namespace libMesh;
72 
77 class LaplaceYoung :
81 {
82 public:
84  _kappa(1.),
85  _sigma(0.2),
86  _gamma(1.0)
87  {}
88 
92  virtual void residual (const NumericVector<Number> & X,
95 
99  virtual void jacobian (const NumericVector<Number> & X,
102 
116  virtual void postcheck (const NumericVector<Number> & old_soln,
117  NumericVector<Number> & search_direction,
118  NumericVector<Number> & new_soln,
119  bool & changed_search_direction,
120  bool & changed_new_soln,
122 
123 private:
126 
127  // Damping factor used for the solve postcheck
129 };
130 
131 
132 
133 // Begin the main program.
134 int main (int argc, char ** argv)
135 {
136  // Initialize libMesh and any dependent libraries, like in example 2.
137  LibMeshInit init (argc, argv);
138 
139  // This example requires a NonlinearSolver.
140 #if !defined(LIBMESH_HAVE_PETSC) && (!defined(LIBMESH_TRILINOS_HAVE_NOX) || !defined(LIBMESH_TRILINOS_HAVE_EPETRA))
141  libmesh_example_requires(false, "--enable-petsc or --enable-trilinos");
142 #endif
143 
144  if (libMesh::on_command_line ("--use-eigen"))
145  {
146  libMesh::err << "This example requires a NonlinearSolver, and therefore does not "
147  << "support --use-eigen on the command line."
148  << std::endl;
149  return 0;
150  }
151 
152 #ifndef LIBMESH_ENABLE_AMR
153  libmesh_example_requires(false, "--enable-amr");
154 #else
155 
156  // Create a GetPot object to parse the command line
157  GetPot command_line (argc, argv);
158 
159  // Check for proper calling arguments.
160  if (argc < 3)
161  {
162  // This handy function will print the file name, line number,
163  // specified message, and then throw an exception.
164  libmesh_error_msg("Usage:\n" << "\t " << argv[0] << " -r 2");
165  }
166 
167  // Brief message to the user regarding the program name
168  // and command line arguments.
169  else
170  {
171  libMesh::out << "Running " << argv[0];
172 
173  for (int i=1; i<argc; i++)
174  libMesh::out << " " << argv[i];
175 
176  libMesh::out << std::endl << std::endl;
177  }
178 
179 
180  // Read number of refinements
181  int nr = 2;
182  if (command_line.search(1, "-r"))
183  nr = command_line.next(nr);
184 
185  // Read FE order from command line
186  std::string order = "FIRST";
187  if (command_line.search(2, "-Order", "-o"))
188  order = command_line.next(order);
189 
190  // Read FE Family from command line
191  std::string family = "LAGRANGE";
192  if (command_line.search(2, "-FEFamily", "-f"))
193  family = command_line.next(family);
194 
195  // Cannot use discontinuous basis.
196  if ((family == "MONOMIAL") || (family == "XYZ"))
197  libmesh_error_msg("This example requires a C^0 (or higher) FE basis.");
198 
199  if (command_line.search(1, "-pre"))
200  {
201 #ifdef LIBMESH_HAVE_PETSC
202  //Use the jacobian for preconditioning.
203 # if PETSC_VERSION_LESS_THAN(3,7,0)
204  PetscOptionsSetValue("-snes_mf_operator", PETSC_NULL);
205 # else
206  PetscOptionsSetValue(PETSC_NULL, "-snes_mf_operator", PETSC_NULL);
207 # endif
208 #else
209  libMesh::err << "Must be using PETSc to use jacobian based preconditioning" << std::endl;
210 
211  //returning zero so that "make run" won't fail if we ever enable this capability there.
212  return 0;
213 #endif //LIBMESH_HAVE_PETSC
214  }
215 
216  // Skip this 2D example if libMesh was compiled as 1D-only.
217  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");
218 
219  // Create a mesh, with dimension to be overridden by the file,
220  // distributed across the default MPI communicator.
221  Mesh mesh(init.comm());
222 
223  mesh.read ("lshaped.xda");
224 
225  if (order != "FIRST")
227 
229 
230  // Print information about the mesh to the screen.
231  mesh.print_info();
232 
233  // Create an equation systems object.
234  EquationSystems equation_systems (mesh);
235 
236  // Declare the system and its variables.
237 
238  // Creates a system named "Laplace-Young"
239  NonlinearImplicitSystem & system =
240  equation_systems.add_system<NonlinearImplicitSystem> ("Laplace-Young");
241 
242  // Here we specify the tolerance for the nonlinear solver and
243  // the maximum of nonlinear iterations.
244  equation_systems.parameters.set<Real> ("nonlinear solver tolerance") = 1.e-12;
245  equation_systems.parameters.set<unsigned int> ("nonlinear solver maximum iterations") = 50;
246 
247 
248  // Adds the variable "u" to "Laplace-Young". "u"
249  // will be approximated using second-order approximation.
250  system.add_variable("u",
251  Utility::string_to_enum<Order> (order),
252  Utility::string_to_enum<FEFamily>(family));
253 
254  // Construct object which provides the residual and jacobian
255  // computations and tell the solver to use it.
256  LaplaceYoung laplace_young;
257  system.nonlinear_solver->residual_object = &laplace_young;
258  system.nonlinear_solver->jacobian_object = &laplace_young;
259  system.nonlinear_solver->postcheck_object = &laplace_young;
260 
261  // Initialize the data structures for the equation system.
262  equation_systems.init();
263 
264  // Prints information about the system to the screen.
265  equation_systems.print_info();
266 
267  // Solve the system "Laplace-Young", print the number of iterations
268  // and final residual
269  equation_systems.get_system("Laplace-Young").solve();
270 
271  // Print out final convergence information. This duplicates some
272  // output from during the solve itself, but demonstrates another way
273  // to get this information after the solve is complete.
274  libMesh::out << "Laplace-Young system solved at nonlinear iteration "
275  << system.n_nonlinear_iterations()
276  << " , final nonlinear residual norm: "
277  << system.final_nonlinear_residual()
278  << std::endl;
279 
280 #ifdef LIBMESH_HAVE_EXODUS_API
281  // After solving the system write the solution
283  equation_systems);
284 #endif // #ifdef LIBMESH_HAVE_EXODUS_API
285 #endif // #ifndef LIBMESH_ENABLE_AMR
286 
287  // All done.
288  return 0;
289 }
290 
291 
292 
293 // Residual assembly function for the Laplace-Young system
295  NumericVector<Number> & residual,
297 {
299 
300  // Get a constant reference to the mesh object.
301  const MeshBase & mesh = es.get_mesh();
302 
303  // The dimension that we are running
304  const unsigned int dim = mesh.mesh_dimension();
305  libmesh_assert_equal_to (dim, 2);
306 
307  // Get a reference to the NonlinearImplicitSystem we are solving
308  NonlinearImplicitSystem & system =
309  es.get_system<NonlinearImplicitSystem>("Laplace-Young");
310 
311  // A reference to the DofMap object for this system. The DofMap
312  // object handles the index translation from node and element numbers
313  // to degree of freedom numbers. We will talk more about the DofMap
314  // in future examples.
315  const DofMap & dof_map = system.get_dof_map();
316 
317  // Get a constant reference to the Finite Element type
318  // for the first (and only) variable in the system.
319  FEType fe_type = dof_map.variable_type(0);
320 
321  // Build a Finite Element object of the specified type. Since the
322  // FEBase::build() member dynamically creates memory we will
323  // store the object as a UniquePtr<FEBase>. This can be thought
324  // of as a pointer that will clean up after itself.
325  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
326 
327  // A 5th order Gauss quadrature rule for numerical integration.
328  QGauss qrule (dim, FIFTH);
329 
330  // Tell the finite element object to use our quadrature rule.
331  fe->attach_quadrature_rule (&qrule);
332 
333  // Declare a special finite element object for
334  // boundary integration.
335  UniquePtr<FEBase> fe_face (FEBase::build(dim, fe_type));
336 
337  // Boundary integration requires one quadrature rule,
338  // with dimensionality one less than the dimensionality
339  // of the element.
340  QGauss qface(dim-1, FIFTH);
341 
342  // Tell the finite element object to use our
343  // quadrature rule.
344  fe_face->attach_quadrature_rule (&qface);
345 
346  // Here we define some references to cell-specific data that
347  // will be used to assemble the linear system.
348  // We begin with the element Jacobian * quadrature weight at each
349  // integration point.
350  const std::vector<Real> & JxW = fe->get_JxW();
351 
352  // The element shape functions evaluated at the quadrature points.
353  const std::vector<std::vector<Real>> & phi = fe->get_phi();
354 
355  // The element shape function gradients evaluated at the quadrature
356  // points.
357  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
358 
359  // Define data structures to contain the residual contributions
361 
362  // This vector will hold the degree of freedom indices for
363  // the element. These define where in the global system
364  // the element degrees of freedom get mapped.
365  std::vector<dof_id_type> dof_indices;
366 
367  // Now we will loop over all the active elements in the mesh which
368  // are local to this processor.
369  // We will compute the element residual.
370  residual.zero();
371 
372  for (const auto & elem : mesh.active_local_element_ptr_range())
373  {
374  // Get the degree of freedom indices for the
375  // current element. These define where in the global
376  // matrix and right-hand-side this element will
377  // contribute to.
378  dof_map.dof_indices (elem, dof_indices);
379 
380  // Compute the element-specific data for the current
381  // element. This involves computing the location of the
382  // quadrature points (q_point) and the shape functions
383  // (phi, dphi) for the current element.
384  fe->reinit (elem);
385 
386  // We use the resize member here because
387  // the number of degrees of freedom might have changed from
388  // the last element. Note that this will be the case if the
389  // element type is different (i.e. the last element was a
390  // triangle, now we are on a quadrilateral).
391  Re.resize (dof_indices.size());
392 
393  // Now we will build the residual. This involves
394  // the construction of the matrix K and multiplication of it
395  // with the current solution x. We rearrange this into two loops:
396  // In the first, we calculate only the contribution of
397  // K_ij*x_j which is independent of the row i. In the second loops,
398  // we multiply with the row-dependent part and add it to the element
399  // residual.
400 
401  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
402  {
403  Number u = 0;
404  Gradient grad_u;
405 
406  for (std::size_t j=0; j<phi.size(); j++)
407  {
408  u += phi[j][qp]*soln(dof_indices[j]);
409  grad_u += dphi[j][qp]*soln(dof_indices[j]);
410  }
411 
412  const Number K = 1./std::sqrt(1. + grad_u*grad_u);
413 
414  for (std::size_t i=0; i<phi.size(); i++)
415  Re(i) += JxW[qp]*(
416  K*(dphi[i][qp]*grad_u) +
417  _kappa*phi[i][qp]*u
418  );
419  }
420 
421  // At this point the interior element integration has
422  // been completed. However, we have not yet addressed
423  // boundary conditions.
424 
425  // The following loops over the sides of the element.
426  // If the element has no neighbor on a side then that
427  // side MUST live on a boundary of the domain.
428  for (auto side : elem->side_index_range())
429  if (elem->neighbor_ptr(side) == libmesh_nullptr)
430  {
431  // The value of the shape functions at the quadrature
432  // points.
433  const std::vector<std::vector<Real>> & phi_face = fe_face->get_phi();
434 
435  // The Jacobian * Quadrature Weight at the quadrature
436  // points on the face.
437  const std::vector<Real> & JxW_face = fe_face->get_JxW();
438 
439  // Compute the shape function values on the element face.
440  fe_face->reinit(elem, side);
441 
442  // Loop over the face quadrature points for integration.
443  for (unsigned int qp=0; qp<qface.n_points(); qp++)
444  {
445  // This is the right-hand-side contribution (f),
446  // which has to be subtracted from the current residual
447  for (std::size_t i=0; i<phi_face.size(); i++)
448  Re(i) -= JxW_face[qp]*_sigma*phi_face[i][qp];
449  }
450  }
451 
452  dof_map.constrain_element_vector (Re, dof_indices);
453  residual.add_vector (Re, dof_indices);
454  }
455 
456  // That's it.
457 }
458 
459 
460 
461 // Jacobian assembly function for the Laplace-Young system
463  SparseMatrix<Number> & jacobian,
465 {
466  // Get a reference to the equation system.
468 
469  // Get a constant reference to the mesh object.
470  const MeshBase & mesh = es.get_mesh();
471 
472  // The dimension that we are running
473  const unsigned int dim = mesh.mesh_dimension();
474 
475  // Get a reference to the NonlinearImplicitSystem we are solving
476  NonlinearImplicitSystem & system =
477  es.get_system<NonlinearImplicitSystem>("Laplace-Young");
478 
479  // A reference to the DofMap object for this system. The DofMap
480  // object handles the index translation from node and element numbers
481  // to degree of freedom numbers. We will talk more about the DofMap
482  // in future examples.
483  const DofMap & dof_map = system.get_dof_map();
484 
485  // Get a constant reference to the Finite Element type
486  // for the first (and only) variable in the system.
487  FEType fe_type = dof_map.variable_type(0);
488 
489  // Build a Finite Element object of the specified type. Since the
490  // FEBase::build() member dynamically creates memory we will
491  // store the object as a UniquePtr<FEBase>. This can be thought
492  // of as a pointer that will clean up after itself.
493  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
494 
495  // A 5th order Gauss quadrature rule for numerical integration.
496  QGauss qrule (dim, FIFTH);
497 
498  // Tell the finite element object to use our quadrature rule.
499  fe->attach_quadrature_rule (&qrule);
500 
501  // Here we define some references to cell-specific data that
502  // will be used to assemble the linear system.
503  // We begin with the element Jacobian * quadrature weight at each
504  // integration point.
505  const std::vector<Real> & JxW = fe->get_JxW();
506 
507  // The element shape functions evaluated at the quadrature points.
508  const std::vector<std::vector<Real>> & phi = fe->get_phi();
509 
510  // The element shape function gradients evaluated at the quadrature
511  // points.
512  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
513 
514  // Define data structures to contain the Jacobian element matrix.
515  // Following basic finite element terminology we will denote these
516  // "Ke". More detail is in example 3.
518 
519  // This vector will hold the degree of freedom indices for
520  // the element. These define where in the global system
521  // the element degrees of freedom get mapped.
522  std::vector<dof_id_type> dof_indices;
523 
524  // Now we will loop over all the active elements in the mesh which
525  // are local to this processor.
526  // We will compute the element Jacobian contribution.
527  for (const auto & elem : mesh.active_local_element_ptr_range())
528  {
529  // Get the degree of freedom indices for the
530  // current element. These define where in the global
531  // matrix and right-hand-side this element will
532  // contribute to.
533  dof_map.dof_indices (elem, dof_indices);
534 
535  // Compute the element-specific data for the current
536  // element. This involves computing the location of the
537  // quadrature points (q_point) and the shape functions
538  // (phi, dphi) for the current element.
539  fe->reinit (elem);
540 
541  // Zero the element Jacobian before
542  // summing them. We use the resize member here because
543  // the number of degrees of freedom might have changed from
544  // the last element. Note that this will be the case if the
545  // element type is different (i.e. the last element was a
546  // triangle, now we are on a quadrilateral).
547  Ke.resize (dof_indices.size(),
548  dof_indices.size());
549 
550  // Now we will build the element Jacobian. This involves
551  // a double loop to integrate the test functions (i) against
552  // the trial functions (j). Note that the Jacobian depends
553  // on the current solution x, which we access using the soln
554  // vector.
555  //
556  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
557  {
558  Gradient grad_u;
559 
560  for (std::size_t i=0; i<phi.size(); i++)
561  grad_u += dphi[i][qp]*soln(dof_indices[i]);
562 
563  const Number
564  sa = 1. + grad_u*grad_u,
565  K = 1. / std::sqrt(sa),
566  dK = -K / sa;
567 
568  for (std::size_t i=0; i<phi.size(); i++)
569  for (std::size_t j=0; j<phi.size(); j++)
570  Ke(i,j) += JxW[qp]*(
571  K * (dphi[i][qp]*dphi[j][qp]) +
572  dK * (grad_u*dphi[j][qp]) * (grad_u*dphi[i][qp]) +
573  _kappa * phi[i][qp] * phi[j][qp]
574  );
575  }
576 
577  dof_map.constrain_element_matrix (Ke, dof_indices);
578 
579  // Add the element matrix to the system Jacobian.
580  jacobian.add_matrix (Ke, dof_indices);
581  }
582 
583  // That's it.
584 }
585 
586 
587 
588 // Jacobian assembly function for the Laplace-Young system
590  NumericVector<Number> & search_direction,
591  NumericVector<Number> & new_soln,
592  bool & /*changed_search_direction*/,
593  bool & changed_new_soln,
594  NonlinearImplicitSystem & /*S*/)
595 {
596  // Back up along the search direction by some amount. Since Newton
597  // already works well for this problem, the only affect of this is
598  // to degrade the rate of convergence.
599  //
600  // The minus sign is due to the sign of the "search_direction"
601  // vector which comes in from the Newton solve. The RHS of the
602  // nonlinear system, i.e. the residual, is generally not multiplied
603  // by -1, so the solution vector, i.e. the search_direction, has a
604  // factor -1.
605  if (_gamma != 1.0)
606  {
607  new_soln = old_soln;
608  new_soln.add(-_gamma, search_direction);
609  changed_new_soln = true;
610  }
611  else
612  changed_new_soln = false;
613 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
OStreamProxy err
This is the EquationSystems class.
int main(int argc, char **argv)
unsigned int dim
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
ImplicitSystem & sys
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 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
unsigned short int side
Definition: xdr_io.C:49
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
virtual void postcheck(const NumericVector< Number > &old_soln, NumericVector< Number > &search_direction, NumericVector< Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln, NonlinearImplicitSystem &S)
Function which performs a postcheck on the solution.
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
The libMesh namespace provides an interface to certain functionality in the library.
This is the MeshBase class.
Definition: mesh_base.h:68
virtual void zero()=0
Set all entries to zero.
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
Abstract base class to be used to calculate the residual of a nonlinear system.
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.
UniquePtr< NonlinearSolver< Number > > nonlinear_solver
The NonlinearSolver defines the default interface used to solve the nonlinear_implicit system...
const DofMap & get_dof_map() const
Definition: system.h:2030
This is the MeshRefinement class.
Abstract base class to be used to calculate the Jacobian of a nonlinear system.
virtual void all_second_order(const bool full_ordered=true)=0
Converts a (conforming, non-refined) mesh with linear elements into a mesh with second-order elements...
void print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
This class provides a specific system class.
virtual System & add_system(const std::string &system_type, const std::string &name)
Add the system of type system_type named name to the systems array.
const Parallel::Communicator & comm() const
Definition: libmesh.h:81
A class which provides the residual and jacobian assembly functions for the Laplace-Young system of e...
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
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &J, NonlinearImplicitSystem &S)
Function which computes the jacobian.
T & set(const std::string &)
Definition: parameters.h:469
bool on_command_line(const std::string &arg)
Definition: libmesh.C:921
OStreamProxy out
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 EquationSystems & get_equation_systems() const
Definition: system.h:712
Parameters parameters
Data structure holding arbitrary parameters.
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
virtual void init()
Initialize all the systems.
virtual void add(const numeric_index_type i, const T value)=0
Adds value to each entry of the vector.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
virtual void residual(const NumericVector< Number > &X, NumericVector< Number > &R, NonlinearImplicitSystem &S)
Function which computes the residual.
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.
Abstract base class to be used for applying user modifications to the solution vector and/or Newton u...
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448
unsigned int n_points() const
Definition: quadrature.h:113
void uniformly_refine(unsigned int n=1)
Uniformly refines the mesh n times.
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.