libMesh
Functions
subdomains_ex2.C File Reference

Go to the source code of this file.

Functions

void assemble_poisson (EquationSystems &es, const std::string &system_name)
 
Real exact_solution (const Real x, const Real y=0., const Real z=0.)
 This is the exact solution that we are trying to obtain. More...
 
int main (int argc, char **argv)
 
void assemble_poisson (EquationSystems &es, const std::string &libmesh_dbg_var(system_name))
 

Function Documentation

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

Referenced by main().

void assemble_poisson ( EquationSystems es,
const std::string &  libmesh_dbg_varsystem_name 
)

Definition at line 287 of file subdomains_ex2.C.

References libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), dim, exact_solution(), libMesh::FIFTH, libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libmesh_nullptr, libMesh::MeshBase::local_elements_begin(), libMesh::MeshBase::local_elements_end(), libMesh::ImplicitSystem::matrix, std::max(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::Elem::neighbor_ptr(), libMesh::pi, libMesh::PerfLog::pop(), libMesh::PerfLog::push(), libMesh::Real, libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, side, libMesh::Elem::side_index_range(), libMesh::Elem::subdomain_id(), value, and libMesh::x.

289 {
290  // It is a good idea to make sure we are assembling
291  // the proper system.
292  libmesh_assert_equal_to (system_name, "Poisson");
293 
294  // Declare a performance log. Give it a descriptive
295  // string to identify what part of the code we are
296  // logging, since there may be many PerfLogs in an
297  // application.
298  PerfLog perf_log ("Matrix Assembly");
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 
306  // Get a reference to the LinearImplicitSystem we are solving
307  LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>("Poisson");
308 
309  // A reference to the DofMap object for this system. The DofMap
310  // object handles the index translation from node and element numbers
311  // to degree of freedom numbers. We will talk more about the DofMap
312  // in future examples.
313  const DofMap & dof_map = system.get_dof_map();
314 
315  // Get a constant reference to the Finite Element type
316  // for the first (and only) variable in the system.
317  FEType fe_type = dof_map.variable_type(0);
318 
319  // Build a Finite Element object of the specified type. Since the
320  // FEBase::build() member dynamically creates memory we will
321  // store the object as a UniquePtr<FEBase>. This can be thought
322  // of as a pointer that will clean up after itself.
323  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
324 
325  // A 5th order Gauss quadrature rule for numerical integration.
326  QGauss qrule (dim, FIFTH);
327 
328  // Tell the finite element object to use our quadrature rule.
329  fe->attach_quadrature_rule (&qrule);
330 
331  // Declare a special finite element object for
332  // boundary integration.
333  UniquePtr<FEBase> fe_face (FEBase::build(dim, fe_type));
334 
335  // Boundary integration requires one quadrature rule,
336  // with dimensionality one less than the dimensionality
337  // of the element.
338  QGauss qface(dim-1, FIFTH);
339 
340  // Tell the finite element object to use our
341  // quadrature rule.
342  fe_face->attach_quadrature_rule (&qface);
343 
344  // Here we define some references to cell-specific data that
345  // will be used to assemble the linear system.
346  // We begin with the element Jacobian * quadrature weight at each
347  // integration point.
348  const std::vector<Real> & JxW = fe->get_JxW();
349 
350  // The physical XY locations of the quadrature points on the element.
351  // These might be useful for evaluating spatially varying material
352  // properties at the quadrature points.
353  const std::vector<Point> & q_point = fe->get_xyz();
354 
355  // The element shape functions evaluated at the quadrature points.
356  const std::vector<std::vector<Real>> & phi = fe->get_phi();
357 
358  // The element shape function gradients evaluated at the quadrature
359  // points.
360  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
361 
362  // Define data structures to contain the element matrix
363  // and right-hand-side vector contribution. Following
364  // basic finite element terminology we will denote these
365  // "Ke" and "Fe". More detail is in example 3.
368 
369  // This vector will hold the degree of freedom indices for
370  // the element. These define where in the global system
371  // the element degrees of freedom get mapped.
372  std::vector<dof_id_type> dof_indices, dof_indices2;
373 
374  // Now we will loop over all the elements in the mesh.
375  // We will compute the element matrix and right-hand-side
376  // contribution. See example 3 for a discussion of the
377  // element iterators. Here we use the const_local_elem_iterator
378  // to indicate we only want to loop over elements that are assigned
379  // to the local processor. This allows each processor to compute
380  // its components of the global matrix.
381  //
382  // "PARALLEL CHANGE"
385 
386  for ( ; el != end_el; ++el)
387  {
388  // Start logging the shape function initialization.
389  // This is done through a simple function call with
390  // the name of the event to log.
391  perf_log.push("elem init");
392 
393  // Store a pointer to the element we are currently
394  // working on. This allows for nicer syntax later.
395  const Elem * elem = *el;
396 
397  // Get the degree of freedom indices for the
398  // current element. These define where in the global
399  // matrix and right-hand-side this element will
400  // contribute to.
401  dof_map.dof_indices (elem, dof_indices, 0);
402  dof_map.dof_indices (elem, dof_indices2, 1);
403 
404  // libMesh::out << "dof_indices.size()="
405  // << dof_indices.size()
406  // << ", dof_indices2.size()="
407  // << dof_indices2.size()
408  // << std::endl;
409 
410  // Compute the element-specific data for the current
411  // element. This involves computing the location of the
412  // quadrature points (q_point) and the shape functions
413  // (phi, dphi) for the current element.
414  fe->reinit (elem);
415 
416  // Zero the element matrix and right-hand side before
417  // summing them. We use the resize member here because
418  // the number of degrees of freedom might have changed from
419  // the last element. Note that this will be the case if the
420  // element type is different (i.e. the last element was a
421  // triangle, now we are on a quadrilateral).
422  Ke.resize (std::max(dof_indices.size(), dof_indices2.size()),
423  std::max(dof_indices.size(), dof_indices2.size()));
424 
425  Fe.resize (std::max(dof_indices.size(), dof_indices2.size()));
426 
427  // Stop logging the shape function initialization.
428  // If you forget to stop logging an event the PerfLog
429  // object will probably catch the error and abort.
430  perf_log.pop("elem init");
431 
432  // Now we will build the element matrix. This involves
433  // a double loop to integrate the test functions (i) against
434  // the trial functions (j).
435  //
436  // We have split the numeric integration into two loops
437  // so that we can log the matrix and right-hand-side
438  // computation separately.
439  //
440  // Now start logging the element matrix computation
441  perf_log.push ("Ke");
442 
443  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
444  for (std::size_t i=0; i<phi.size(); i++)
445  for (std::size_t j=0; j<phi.size(); j++)
446  Ke(i,j) += JxW[qp]*(dphi[i][qp]*dphi[j][qp]);
447 
448 
449  // Stop logging the matrix computation
450  perf_log.pop ("Ke");
451 
452  // Now we build the element right-hand-side contribution.
453  // This involves a single loop in which we integrate the
454  // "forcing function" in the PDE against the test functions.
455  //
456  // Start logging the right-hand-side computation
457  perf_log.push ("Fe");
458 
459  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
460  {
461  // fxy is the forcing function for the Poisson equation.
462  // In this case we set fxy to be a finite difference
463  // Laplacian approximation to the (known) exact solution.
464  //
465  // We will use the second-order accurate FD Laplacian
466  // approximation, which in 2D on a structured grid is
467  //
468  // u_xx + u_yy = (u(i-1,j) + u(i+1,j) +
469  // u(i,j-1) + u(i,j+1) +
470  // -4*u(i,j))/h^2
471  //
472  // Since the value of the forcing function depends only
473  // on the location of the quadrature point (q_point[qp])
474  // we will compute it here, outside of the i-loop
475  const Real x = q_point[qp](0);
476 #if LIBMESH_DIM > 1
477  const Real y = q_point[qp](1);
478 #else
479  const Real y = 0;
480 #endif
481 #if LIBMESH_DIM > 2
482  const Real z = q_point[qp](2);
483 #else
484  const Real z = 0;
485 #endif
486  const Real eps = 1.e-3;
487 
488  const Real uxx = (exact_solution(x-eps, y, z) +
489  exact_solution(x+eps, y, z) +
490  -2.*exact_solution(x, y, z))/eps/eps;
491 
492  const Real uyy = (exact_solution(x, y-eps, z) +
493  exact_solution(x, y+eps, z) +
494  -2.*exact_solution(x, y, z))/eps/eps;
495 
496  const Real uzz = (exact_solution(x, y, z-eps) +
497  exact_solution(x, y, z+eps) +
498  -2.*exact_solution(x, y, z))/eps/eps;
499 
500  Real fxy;
501  if (dim==1)
502  {
503  // In 1D, compute the rhs by differentiating the
504  // exact solution twice.
505  const Real pi = libMesh::pi;
506  fxy = (0.25*pi*pi)*sin(.5*pi*x);
507  }
508  else
509  {
510  fxy = - (uxx + uyy + ((dim==2) ? 0. : uzz));
511  }
512 
513  // Add the RHS contribution
514  for (std::size_t i=0; i<phi.size(); i++)
515  Fe(i) += JxW[qp]*fxy*phi[i][qp];
516  }
517 
518  // Stop logging the right-hand-side computation
519  perf_log.pop ("Fe");
520 
521  // At this point the interior element integration has
522  // been completed. However, we have not yet addressed
523  // boundary conditions. For this example we will only
524  // consider simple Dirichlet boundary conditions imposed
525  // via the penalty method. This is discussed at length in
526  // example 3.
527  {
528 
529  // Start logging the boundary condition computation
530  perf_log.push ("BCs");
531 
532  // The following loops over the sides of the element.
533  // If the element has no neighbor on a side then that
534  // side MUST live on a boundary of the domain.
535  for (auto side : elem->side_index_range())
536  if ((elem->neighbor_ptr(side) == libmesh_nullptr) ||
537  (elem->neighbor_ptr(side)->subdomain_id() != elem->subdomain_id()))
538  {
539 
540  // The penalty value. \frac{1}{\epsilon}
541  // in the discussion above.
542  const Real penalty = 1.e10;
543 
544  // The value of the shape functions at the quadrature
545  // points.
546  const std::vector<std::vector<Real>> & phi_face = fe_face->get_phi();
547 
548  // The Jacobian * Quadrature Weight at the quadrature
549  // points on the face.
550  const std::vector<Real> & JxW_face = fe_face->get_JxW();
551 
552  // The XYZ locations (in physical space) of the
553  // quadrature points on the face. This is where
554  // we will interpolate the boundary value function.
555  const std::vector<Point> & qface_point = fe_face->get_xyz();
556 
557  // Compute the shape function values on the element
558  // face.
559  fe_face->reinit(elem, side);
560 
561  // Loop over the face quadrature points for integration.
562  for (unsigned int qp=0; qp<qface.n_points(); qp++)
563  {
564  // The location on the boundary of the current
565  // face quadrature point.
566  const Real xf = qface_point[qp](0);
567 #if LIBMESH_DIM > 1
568  const Real yf = qface_point[qp](1);
569 #else
570  const Real yf = 0.;
571 #endif
572 #if LIBMESH_DIM > 2
573  const Real zf = qface_point[qp](2);
574 #else
575  const Real zf = 0.;
576 #endif
577 
578 
579  // The boundary value.
580  const Real value = exact_solution(xf, yf, zf);
581 
582  // Matrix contribution of the L2 projection.
583  for (std::size_t i=0; i<phi_face.size(); i++)
584  for (std::size_t j=0; j<phi_face.size(); j++)
585  Ke(i,j) += JxW_face[qp]*penalty*phi_face[i][qp]*phi_face[j][qp];
586 
587  // Right-hand-side contribution of the L2
588  // projection.
589  for (std::size_t i=0; i<phi_face.size(); i++)
590  Fe(i) += JxW_face[qp]*penalty*value*phi_face[i][qp];
591  }
592  }
593 
594 
595  // Stop logging the boundary condition computation
596  perf_log.pop ("BCs");
597  }
598 
599 
600  // The element matrix and right-hand-side are now built
601  // for this element. Add them to the global matrix and
602  // right-hand-side vector. The PetscMatrix::add_matrix()
603  // and PetscVector::add_vector() members do this for us.
604  // Start logging the insertion of the local (element)
605  // matrix and vector into the global matrix and vector
606  perf_log.push ("matrix insertion");
607 
608  if (dof_indices.size())
609  {
610  system.matrix->add_matrix (Ke, dof_indices);
611  system.rhs->add_vector (Fe, dof_indices);
612  }
613 
614  if (dof_indices2.size())
615  {
616  system.matrix->add_matrix (Ke, dof_indices2);
617  system.rhs->add_vector (Fe, dof_indices2);
618  }
619 
620  // Start logging the insertion of the local (element)
621  // matrix and vector into the global matrix and vector
622  perf_log.pop ("matrix insertion");
623  }
624 
625  // That's it. We don't need to do anything else to the
626  // PerfLog. When it goes out of scope (at this function return)
627  // it will print its log to the screen. Pretty easy, huh?
628 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
The definition of the const_element_iterator struct.
Definition: mesh_base.h:1494
unsigned int dim
IntRange< unsigned short > side_index_range() const
Definition: elem.h:2083
This class provides a specific system class.
virtual element_iterator local_elements_begin()=0
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
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
NumericVector< Number > * rhs
The system matrix.
long double max(long double a, double b)
This is the MeshBase class.
Definition: mesh_base.h:68
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const Elem * neighbor_ptr(unsigned int i) const
Definition: elem.h:1967
The PerfLog class allows monitoring of specific events.
Definition: perf_log.h:125
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.
PetscErrorCode Vec x
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual element_iterator local_elements_end()=0
subdomain_id_type subdomain_id() const
Definition: elem.h:1951
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
static const bool value
Definition: xdr_io.C:108
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
Real exact_solution(const Real x, const Real y=0., const Real z=0.)
This is the exact solution that we are trying to obtain.
const Real pi
.
Definition: libmesh.h:172
Real exact_solution ( const Real  x,
const Real  y,
const Real  t 
)

This is the exact solution that we are trying to obtain.

We will solve

  • (u_xx + u_yy) = f

and take a finite difference approximation using this function to get f. This is the well-known "method of manufactured solutions".

Definition at line 43 of file exact_solution.C.

Referenced by assemble_poisson().

46 {
47  static const Real pi = acos(-1.);
48 
49  return cos(.5*pi*x)*sin(.5*pi*y)*cos(.5*pi*z);
50 }
PetscErrorCode Vec x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Real pi
.
Definition: libmesh.h:172
int main ( int  argc,
char **  argv 
)

Definition at line 91 of file subdomains_ex2.C.

References libMesh::EquationSystems::add_system(), libMesh::System::add_variable(), assemble_poisson(), libMesh::System::attach_assemble_function(), libMesh::MeshTools::Generation::build_cube(), libMesh::LinearImplicitSystem::clear(), libMesh::LibMeshInit::comm(), libMesh::default_solver_package(), dim, libMesh::EDGE2, libMesh::EDGE3, libMesh::MeshBase::element_ptr_range(), libMesh::EquationSystems::get_system(), libMesh::GnuPlotIO::GRID_ON, libMesh::HEX27, libMesh::HEX8, libMesh::TriangleWrapper::init(), libMesh::EquationSystems::init(), libMesh::INVALID_SOLVER_PACKAGE, mesh, libMesh::out, libMesh::EquationSystems::print_info(), libMesh::MeshBase::print_info(), libMesh::ParallelObject::processor_id(), libMesh::QUAD4, libMesh::QUAD9, libMesh::Real, and libMesh::MeshOutput< MT >::write_equation_systems().

92 {
93  // Initialize libMesh and any dependent libraries, like in example 2.
94  LibMeshInit init (argc, argv);
95 
96  // This example requires a linear solver package.
97  libmesh_example_requires(libMesh::default_solver_package() != INVALID_SOLVER_PACKAGE,
98  "--enable-petsc, --enable-trilinos, or --enable-eigen");
99 
100  // Declare a performance log for the main program
101  // PerfLog perf_main("Main Program");
102 
103  // Create a GetPot object to parse the command line
104  GetPot command_line (argc, argv);
105 
106  // Check for proper calling arguments.
107  if (argc < 3)
108  {
109  // This handy function will print the file name, line number,
110  // specified message, and then throw an exception.
111  libmesh_error_msg("Usage:\n" << "\t " << argv[0] << " -d 2(3)" << " -n 15");
112  }
113 
114  // Brief message to the user regarding the program name
115  // and command line arguments.
116  else
117  {
118  libMesh::out << "Running " << argv[0];
119 
120  for (int i=1; i<argc; i++)
121  libMesh::out << " " << argv[i];
122 
123  libMesh::out << std::endl << std::endl;
124  }
125 
126 
127  // Read problem dimension from command line. Use int
128  // instead of unsigned since the GetPot overload is ambiguous
129  // otherwise.
130  int dim = 2;
131  if (command_line.search(1, "-d"))
132  dim = command_line.next(dim);
133 
134  // Skip higher-dimensional examples on a lower-dimensional libMesh build
135  libmesh_example_requires(dim <= LIBMESH_DIM, "2D/3D support");
136 
137  // Create a mesh with user-defined dimension on the default MPI
138  // communicator.
139  Mesh mesh (init.comm(), dim);
140 
141  // Read number of elements from command line
142  int ps = 15;
143  if (command_line.search(1, "-n"))
144  ps = command_line.next(ps);
145 
146  // Read FE order from command line
147  std::string order = "SECOND";
148  if (command_line.search(2, "-Order", "-o"))
149  order = command_line.next(order);
150 
151  // Read FE Family from command line
152  std::string family = "LAGRANGE";
153  if (command_line.search(2, "-FEFamily", "-f"))
154  family = command_line.next(family);
155 
156  // Cannot use discontinuous basis.
157  if ((family == "MONOMIAL") || (family == "XYZ"))
158  {
159  if (mesh.processor_id() == 0)
160  libmesh_error_msg("This example requires a C^0 (or higher) FE basis.");
161  }
162 
163  // Use the MeshTools::Generation mesh generator to create a uniform
164  // grid on the square [-1,1]^D. We instruct the mesh generator
165  // to build a mesh of 8x8 Quad9 elements in 2D, or Hex27
166  // elements in 3D. Building these higher-order elements allows
167  // us to use higher-order approximation, as in example 3.
168 
169  Real halfwidth = dim > 1 ? 1. : 0.;
170  Real halfheight = dim > 2 ? 1. : 0.;
171 
172  if ((family == "LAGRANGE") && (order == "FIRST"))
173  {
174  // No reason to use high-order geometric elements if we are
175  // solving with low-order finite elements.
177  ps,
178  (dim>1) ? ps : 0,
179  (dim>2) ? ps : 0,
180  -1., 1.,
181  -halfwidth, halfwidth,
182  -halfheight, halfheight,
183  (dim==1) ? EDGE2 :
184  ((dim == 2) ? QUAD4 : HEX8));
185  }
186 
187  else
188  {
190  ps,
191  (dim>1) ? ps : 0,
192  (dim>2) ? ps : 0,
193  -1., 1.,
194  -halfwidth, halfwidth,
195  -halfheight, halfheight,
196  (dim==1) ? EDGE3 :
197  ((dim == 2) ? QUAD9 : HEX27));
198  }
199 
200  for (auto & elem : mesh.element_ptr_range())
201  {
202  const Point cent = elem->centroid();
203  if (dim > 1)
204  {
205  if ((cent(0) > 0) == (cent(1) > 0))
206  elem->subdomain_id() = 1;
207  }
208  else if (cent(0) > 0)
209  elem->subdomain_id() = 1;
210  }
211 
212  // Print information about the mesh to the screen.
213  mesh.print_info();
214 
215  // Create an equation systems object.
216  EquationSystems equation_systems (mesh);
217 
218  // Declare the system and its variables.
219  // Create a system named "Poisson"
220  LinearImplicitSystem & system =
221  equation_systems.add_system<LinearImplicitSystem> ("Poisson");
222 
223 
224  std::set<subdomain_id_type> active_subdomains;
225 
226 
227  // Add the variable "u" to "Poisson". "u"
228  // will be approximated using second-order approximation.
229  active_subdomains.clear(); active_subdomains.insert(0);
230  system.add_variable("u",
231  Utility::string_to_enum<Order> (order),
232  Utility::string_to_enum<FEFamily>(family),
233  &active_subdomains);
234 
235  // Add the variable "v" to "Poisson". "v"
236  // will be approximated using second-order approximation.
237  active_subdomains.clear(); active_subdomains.insert(1);
238  system.add_variable("v",
239  Utility::string_to_enum<Order> (order),
240  Utility::string_to_enum<FEFamily>(family),
241  &active_subdomains);
242 
243  // Give the system a pointer to the matrix assembly
244  // function.
246 
247  // Initialize the data structures for the equation system.
248  equation_systems.init();
249 
250  // Print information about the system to the screen.
251  equation_systems.print_info();
252  mesh.print_info();
253 
254  // Solve the system "Poisson", just like example 2.
255  equation_systems.get_system("Poisson").solve();
256 
257  // After solving the system write the solution
258  // to a GMV-formatted plot file.
259  if (dim == 1)
260  {
261  GnuPlotIO plot(mesh, "Subdomains Example 2, 1D", GnuPlotIO::GRID_ON);
262  plot.write_equation_systems("gnuplot_script", equation_systems);
263  }
264  else
265  {
266 #ifdef LIBMESH_HAVE_EXODUS_API
267  ExodusII_IO (mesh).write_equation_systems ((dim == 3) ?
268  "out_3.e" : "out_2.e", equation_systems);
269 #endif // #ifdef LIBMESH_HAVE_EXODUS_API
270  }
271 
272  // All done.
273  return 0;
274 }
This is the EquationSystems class.
unsigned int dim
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
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
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
virtual void clear() libmesh_override
Clear all the data structures associated with the system.
SolverPackage default_solver_package()
Definition: libmesh.C:995
This class implements writing meshes using GNUplot, designed for use only with 1D meshes...
Definition: gnuplot_io.h:43
virtual SimpleRange< element_iterator > element_ptr_range()=0
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.
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
OStreamProxy out
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.
processor_id_type processor_id() const
void assemble_poisson(EquationSystems &es, const std::string &system_name)