libMesh
miscellaneous_ex2.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 
21 // <h1>Miscellaneous Example 2 - Complex Numbers and the "FrequencySystem"</h1>
22 // \author Steffen Petersen
23 // \date 2003
24 //
25 // This example program introduces complex numbers and the
26 // FrequencySystem class to solve a simple Helmholtz equation,
27 // Laplacian(p) + (omega/c)^2*p = 0,
28 // for multiple frequencies rather efficiently.
29 //
30 // The FrequencySystem class offers two solution styles:
31 // 1.) Solve large systems once, and
32 // 2.) Solve moderately-sized systems for multiple frequencies.
33 // The latter approach is implemented here.
34 //
35 // This example uses an L-shaped domain and nodal boundary data given
36 // in the files lshape.unv and lshape_data.unv. For this example, the
37 // library has to be compiled with complex numbers enabled.
38 
39 // C++ include files that we need
40 #include <iostream>
41 #include <algorithm>
42 #include <stdio.h>
43 
44 // Basic include files needed for overall functionality.
45 #include "libmesh/libmesh.h"
46 #include "libmesh/libmesh_logging.h"
47 #include "libmesh/mesh.h"
48 #include "libmesh/mesh_generation.h"
49 #include "libmesh/exodusII_io.h"
50 #include "libmesh/unv_io.h"
51 #include "libmesh/equation_systems.h"
52 #include "libmesh/elem.h"
53 
54 // Include FrequencySystem. This class offers added functionality for
55 // the solution of frequency-dependent systems.
56 #include "libmesh/frequency_system.h"
57 
58 // Define the FE object.
59 #include "libmesh/fe.h"
60 
61 // Define the QGauss quadrature rule objects.
62 #include "libmesh/quadrature_gauss.h"
63 
64 // Useful datatypes for finite element assembly.
65 #include "libmesh/dense_matrix.h"
66 #include "libmesh/dense_vector.h"
67 
68 // Sparse matrix and vector data types for parallel linear algebra.
69 #include "libmesh/sparse_matrix.h"
70 #include "libmesh/numeric_vector.h"
71 
72 // Define the DofMap, which handles degree of freedom indexing.
73 #include "libmesh/dof_map.h"
74 
75 // Bring in everything from the libMesh namespace
76 using namespace libMesh;
77 
78 // This problem is only defined on complex-valued fields, for
79 // which libMesh must be configured with Number == Complex.
80 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
81 
82 // Function prototype. This is the function that will assemble
83 // the mass, damping and stiffness matrices. It will not
84 // form an overall system matrix ready for solution.
86  const std::string & system_name);
87 
88 // Function prototype. This is the function that will combine
89 // the previously-assembled mass, damping and stiffness matrices
90 // to the overall matrix, which then renders ready for solution.
92  const std::string & system_name);
93 #endif
94 
95 // This example only works correctly if libmesh has been configured
96 // with --enable-complex. If you wish to use PETSc, you must also
97 // build PETSc with complex number support by configuring with
98 // --with-scalar-type=complex --with-clanguage=cxx, and using the same
99 // C++ compiler to build both PETSc and libmesh. This example also
100 // works with the Eigen sparse linear solvers which are provided in
101 // libmesh's contrib directory.
102 int main (int argc, char ** argv)
103 {
104  // The libMeshInit object initializes MPI, PETSc, etc, and must be
105  // constructed before all other objects.
106  LibMeshInit init (argc, argv);
107 
108  // This example is designed for complex numbers.
109 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
110  libmesh_example_requires(false, "--enable-complex");
111 #else
112 
113  // This example requires at least 2D support.
114  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");
115 
116  // Check for proper usage. frequency has two terms:
117  if ( argc <4 )
118  libmesh_error_msg("Usage: " << argv[0] << " -f [real_frequency imag_frequency]");
119 
120  if (init.comm().size() > 1)
121  {
122  if (init.comm().rank() == 0)
123  {
124  libMesh::err << "TODO: This example should be able to run in parallel."
125  << std::endl;
126  }
127  return 0;
128  }
129 
130  // Tell the user what we are doing.
131  else
132  {
133  libMesh::out << "Running " << argv[0];
134 
135  for (int i=1; i<argc; i++)
136  libMesh::out << " " << argv[i];
137 
138  libMesh::out << std::endl << std::endl;
139  }
140 
141  // Get the frequency from argv[2] as a float, solve for 1/3rd, 2/3rd
142  // and 1/1th of the given frequency.
143 
144  const Number frequency_in (atof(argv[2]), atof(argv[3]));
145 
146  // this must work as well, but is deprecated.
147  // const Real frequency_in = atof(argv[2]);
148 
149  const unsigned int n_frequencies = 3;
150 
151  // Create a mesh, with dimension to be overridden later, distributed
152  // across the default MPI communicator.
153  Mesh mesh(init.comm());
154 
155  // Read the mesh file. Here the file lshape.unv contains
156  // an L-shaped domain in .unv format.
157  UNVIO unvio(mesh);
158  unvio.read("lshape.unv");
159 
160  // Manually prepare the mesh for use, this is not done automatically
161  // by the UNVIO reader, so we have to do it here. Note that calling
162  // this function renumbers the nodes and elements of the Mesh, but
163  // the original numbering is still stored in the UNVIO object for
164  // correctly mapping dataset values (see below) to the correct nodes.
166 
167  // Read the dataset accompanying this problem. The load on the
168  // boundary of the domain is stored in the .unv data file
169  // lshape_data.unv. The data are given as complex-valued normal
170  // velocities.
171  unvio.read_dataset("lshape_data.unv");
172 
173  // Print information about the mesh to the screen.
174  mesh.print_info();
175 
176  // Create an EquationSystems object.
177  EquationSystems equation_systems (mesh);
178 
179  // Create a FrequencySystem named "Helmholtz" and store a
180  // reference to it.
181  FrequencySystem & f_system =
182  equation_systems.add_system<FrequencySystem> ("Helmholtz");
183 
184  // Add the variable "p" to the "Helmholtz" system. "p"
185  // will be approximated using second-order approximation.
186  f_system.add_variable("p", SECOND);
187 
188  // The FrequencySystem requires two user-provided functions: one for
189  // assembling the different operators and another that specifies how
190  // to combine them before the solve.
193 
194  // To enable the fast solution scheme, additional sparse matrices
195  // and one vector have to be added. The FrequencySystem object
196  // takes care of sizing the additional objects. The user should
197  // still set the sparsity structure of the f_system.matrix, so that
198  // the fast matrix addition method can be used. The procedure for
199  // this is shown in detail in the assembly function.
200  f_system.add_matrix ("stiffness");
201  f_system.add_matrix ("damping");
202  f_system.add_matrix ("mass");
203  f_system.add_vector ("rhs");
204 
205  // Communicate the frequencies to the system. Note that the
206  // frequency system stores the frequencies as parameters in the
207  // EquationSystems object, so that our assemble and solve functions
208  // may directly access them. Must be called before the
209  // EquationSystems object is initialized. Will solve for 1/3,
210  // 2/3, and 1 times the given frequency.
211  f_system.set_frequencies_by_steps (frequency_in/static_cast<Number>(n_frequencies),
212  frequency_in,
213  n_frequencies);
214 
215  // Set the wave velocity and fluid density parameter values,
216  // otherwise default values will be used.
217  equation_systems.parameters.set<Real> ("wave speed") = 1.;
218  equation_systems.parameters.set<Real> ("rho") = 1.;
219 
220  // Initialize the EquationSystems object.
221  equation_systems.init ();
222 
223  // Set values in the "rhs" vector based on the entries stored in the
224  // UNVIO object from the dataset we read in. These values only need
225  // to be set once, as they are the same for every frequency. We can
226  // only do this once equation_systems.init() has been called...
227  {
228  NumericVector<Number> & freq_indep_rhs = f_system.get_vector("rhs");
229 
231  const MeshBase::const_node_iterator node_end = mesh.nodes_end();
232 
233  for ( ; node_it != node_end; ++node_it)
234  {
235  // the current node pointer
236  Node * node = *node_it;
237 
238  // Get the data read in from the dataset for the current Node, if any.
239  const std::vector<Number> * nodal_data = unvio.get_data(node);
240 
241  // Set the rhs value based on values read in from the dataset.
242  if (nodal_data)
243  {
244  unsigned int dn = node->dof_number(/*system=*/0,
245  /*variable=*/0,
246  /*component=*/0);
247  freq_indep_rhs.set(dn, (*nodal_data)[0]);
248  }
249  }
250  }
251 
252  // Print information about the system to the screen.
253  equation_systems.print_info ();
254 
255  for (unsigned int n=0; n < n_frequencies; n++)
256  {
257  // Solve the "Helmholtz" system for the n-th frequency.
258  // Since we attached an assemble() function to the system,
259  // the mass, damping, and stiffness contributions will only
260  // be assembled once. Then, the system is solved for the
261  // given frequencies. Note that solve() may also solve
262  // the system only for specific frequencies.
263  f_system.solve (n, n);
264 
265  // After solving the system, write the solution to an ExodusII
266  // file for every frequency.
267 #ifdef LIBMESH_HAVE_EXODUS_API
268  std::ostringstream file_name;
269 
270  file_name << "out_"
271  << std::setw(4)
272  << std::setfill('0')
273  << std::right
274  << n
275  << ".e";
276 
277  ExodusII_IO(mesh).write_equation_systems (file_name.str(),
278  equation_systems);
279 #endif
280  }
281 
282  // Alternatively, the whole EquationSystems object can be
283  // written to disk. By default, the additional vectors are also
284  // saved.
285  equation_systems.write ("eqn_sys.dat", WRITE);
286 
287  // All done.
288  return 0;
289 
290 #endif
291 }
292 
293 
294 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
295 // Here we define the matrix assembly routine for
296 // the Helmholtz system. This function will be
297 // called to form the stiffness matrix and right-hand side.
299  const std::string & system_name)
300 {
301 
302  // It is a good idea to make sure we are assembling
303  // the proper system.
304  libmesh_assert_equal_to (system_name, "Helmholtz");
305 
306  // Get a constant reference to the mesh object.
307  const MeshBase & mesh = es.get_mesh();
308 
309  // The maximum dimension of the elements stored in the mesh.
310  const unsigned int dim = mesh.mesh_dimension();
311 
312  // Get a reference to our system, as before
313  FrequencySystem & f_system =
314  es.get_system<FrequencySystem> (system_name);
315 
316  // A const reference to the DofMap object for this system. The DofMap
317  // object handles the index translation from node and element numbers
318  // to degree of freedom numbers.
319  const DofMap & dof_map = f_system.get_dof_map();
320 
321  // Get a constant reference to the finite element type
322  // for the first (and only) variable in the system.
323  const FEType & fe_type = dof_map.variable_type(0);
324 
325  // The fluid density is used by the admittance boundary condition.
326  const Real rho = es.parameters.get<Real>("rho");
327 
328  // In this example, we assemble the element matrices into the
329  // additional matrices "stiffness_mass", "damping", and the element
330  // vectors into "rhs". We get writable references to these objects
331  // now.
332  SparseMatrix<Number> & stiffness = f_system.get_matrix("stiffness");
333  SparseMatrix<Number> & damping = f_system.get_matrix("damping");
334  SparseMatrix<Number> & mass = f_system.get_matrix("mass");
335 
336  // Build a finite element object of the specified type. Since the
337  // FEBase::build() member dynamically creates memory we will
338  // store the object as a UniquePtr<FEBase>. This can be thought
339  // of as a pointer that will clean up after itself.
340  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
341 
342  // A 5th-order Gauss quadrature rule for numerical integration.
343  QGauss qrule (dim, FIFTH);
344 
345  // Tell the finite element object to use our quadrature rule.
346  fe->attach_quadrature_rule (&qrule);
347 
348  // The element Jacobian times the quadrature weight at each integration point.
349  const std::vector<Real> & JxW = fe->get_JxW();
350 
351  // The element shape functions evaluated at the quadrature points.
352  const std::vector<std::vector<Real>> & phi = fe->get_phi();
353 
354  // The element shape function gradients evaluated at the quadrature
355  // points.
356  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
357 
358  // Here we do not assemble directly in the System matrix, but to the
359  // additional matrices "stiffness_mass" and "damping". The same
360  // approach is followed for the right-hand-side vector Fe, which we
361  // will later on store in the additional vector "rhs".
362  // The zero_matrix is used to explicitly induce the same sparsity
363  // structure in the overall matrix. The mass and stiffness matrices
364  // are real-valued. Therefore, it would be possible to store them
365  // as a single complex-valued matrix to save on memory, but we do
366  // not follow this approach here.
367  DenseMatrix<Number> Ke, Ce, Me, zero_matrix;
369 
370  // This vector will hold the degree of freedom indices for
371  // the element. These define where in the global system
372  // the element degrees of freedom get mapped.
373  std::vector<dof_id_type> dof_indices;
374 
375  // Now we will loop over all the elements in the mesh, and compute
376  // the element matrix and right-hand-side contributions.
377  for (const auto & elem : mesh.active_local_element_ptr_range())
378  {
379  // Start logging the element initialization.
380  START_LOG("elem init", "assemble_helmholtz");
381 
382  // Get the degree of freedom indices for the
383  // current element. These define where in the global
384  // matrix and right-hand-side this element will
385  // contribute to.
386  dof_map.dof_indices (elem, dof_indices);
387 
388  // Compute the element-specific data for the current
389  // element. This involves computing the location of the
390  // quadrature points (q_point) and the shape functions
391  // (phi, dphi) for the current element.
392  fe->reinit (elem);
393 
394  // Zero & resize the element matrix and right-hand side before
395  // summing them, with different element types in the mesh this
396  // is quite necessary.
397  {
398  const unsigned int n_dof_indices = dof_indices.size();
399 
400  Ke.resize (n_dof_indices, n_dof_indices);
401  Ce.resize (n_dof_indices, n_dof_indices);
402  Me.resize (n_dof_indices, n_dof_indices);
403  zero_matrix.resize (n_dof_indices, n_dof_indices);
404  Fe.resize (n_dof_indices);
405  }
406 
407  // Stop logging the element initialization.
408  STOP_LOG("elem init", "assemble_helmholtz");
409 
410  // Now loop over the quadrature points. This handles
411  // the numeric integration.
412  START_LOG("stiffness & mass", "assemble_helmholtz");
413 
414  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
415  {
416  // Now we will build the element matrix. This involves
417  // a double loop to integrate the test functions (i) against
418  // the trial functions (j).
419  for (std::size_t i=0; i<phi.size(); i++)
420  for (std::size_t j=0; j<phi.size(); j++)
421  {
422  Ke(i,j) += JxW[qp]*(dphi[i][qp]*dphi[j][qp]);
423  Me(i,j) += JxW[qp]*(phi[i][qp]*phi[j][qp]);
424  }
425  }
426 
427  STOP_LOG("stiffness & mass", "assemble_helmholtz");
428 
429  // Now compute the contribution to the element matrix
430  // (due to mixed boundary conditions) if the current
431  // element lies on the boundary.
432  //
433  // The following loops over the sides of the element.
434  // If the element has no neighbor on a side then that
435  // side MUST live on a boundary of the domain.
436  for (auto side : elem->side_index_range())
437  if (elem->neighbor_ptr(side) == libmesh_nullptr)
438  {
439  LOG_SCOPE("damping", "assemble_helmholtz");
440 
441  // Declare a special finite element object for
442  // boundary integration.
443  UniquePtr<FEBase> fe_face (FEBase::build(dim, fe_type));
444 
445  // Boundary integration requires one quadrature rule,
446  // with dimensionality one less than the dimensionality
447  // of the element.
448  QGauss qface(dim-1, SECOND);
449 
450  // Tell the finite element object to use our
451  // quadrature rule.
452  fe_face->attach_quadrature_rule (&qface);
453 
454  // The value of the shape functions at the quadrature
455  // points.
456  const std::vector<std::vector<Real>> & phi_face =
457  fe_face->get_phi();
458 
459  // The Jacobian times the quadrature weight at the quadrature
460  // points on the face.
461  const std::vector<Real> & JxW_face = fe_face->get_JxW();
462 
463  // Compute the shape function values on the element
464  // face.
465  fe_face->reinit(elem, side);
466 
467  // For the Robin BCs, consider a normal admittance an=1
468  // at some parts of the boundary
469  const Real an_value = 1.;
470 
471  // Loop over the face quadrature points for integration.
472  for (unsigned int qp=0; qp<qface.n_points(); qp++)
473  {
474  // Element matrix contribution due to prescribed
475  // admittance boundary conditions.
476  for (std::size_t i=0; i<phi_face.size(); i++)
477  for (std::size_t j=0; j<phi_face.size(); j++)
478  Ce(i,j) += rho * an_value * JxW_face[qp] * phi_face[i][qp] * phi_face[j][qp];
479  }
480  }
481 
482  // If this assembly program were to be used on an adaptive mesh,
483  // we would have to apply any hanging node constraint equations
484  // by uncommenting the following lines:
485  // std::vector<unsigned int> dof_indicesC = dof_indices;
486  // std::vector<unsigned int> dof_indicesM = dof_indices;
487  // dof_map.constrain_element_matrix (Ke, dof_indices);
488  // dof_map.constrain_element_matrix (Ce, dof_indicesC);
489  // dof_map.constrain_element_matrix (Me, dof_indicesM);
490 
491  // Finally, simply add the contributions to the additional
492  // matrices and vector.
493  stiffness.add_matrix (Ke, dof_indices);
494  damping.add_matrix (Ce, dof_indices);
495  mass.add_matrix (Me, dof_indices);
496 
497  // For the overall matrix, explicitly zero the entries where
498  // we added values in the other ones, so that we have
499  // identical sparsity footprints.
500  f_system.matrix->add_matrix(zero_matrix, dof_indices);
501  } // end loop over elements
502 }
503 
504 
505 // We now define the function which will combine the previously-assembled
506 // mass, stiffness, and damping matrices into a single system matrix.
508  const std::string & system_name)
509 {
510  START_LOG("init phase", "add_M_C_K_helmholtz");
511 
512  // Verify that we are assembling the system we think we are.
513  libmesh_assert_equal_to (system_name, "Helmholtz");
514 
515  // Get a reference to the FrequencySystem.
516  FrequencySystem & f_system =
517  es.get_system<FrequencySystem> (system_name);
518 
519  // Get the frequency, fluid density, and sound speed of the current solve.
520  const Number frequency = es.parameters.get<Number> ("current frequency");
521  const Real rho = es.parameters.get<Real> ("rho");
522  const Real speed = es.parameters.get<Real> ("wave speed");
523 
524  // Compute angular frequency omega and wave number k
525  const Number omega = 2.0*libMesh::pi*frequency;
526  const Number k = omega / speed;
527 
528  // Get writable references to the system matrix and vector, where the
529  // frequency-dependent system is to be collected.
530  SparseMatrix<Number> & matrix = *f_system.matrix;
531  NumericVector<Number> & rhs = *f_system.rhs;
532 
533  // Get writable references to the frequency-independent matrices
534  // and rhs, though we only need to extract values. This write access
535  // is necessary, since solver packages have to close the data structure
536  // before they can extract values for computation.
537  SparseMatrix<Number> & stiffness = f_system.get_matrix("stiffness");
538  SparseMatrix<Number> & damping = f_system.get_matrix("damping");
539  SparseMatrix<Number> & mass = f_system.get_matrix("mass");
540  NumericVector<Number> & freq_indep_rhs = f_system.get_vector("rhs");
541 
542  Number unit_i (0,1);
543 
544  // Compute the scale values for the different operators.
545  const Number scale_stiffness (1., 0.);
546  const Number scale_damping=unit_i*omega; // I is imaginary unit (from complex.h)
547  const Number scale_mass=-k*k;
548  const Number scale_rhs=-unit_i*rho*omega;
549 
550  // Now simply add the matrices together, store the result
551  // in matrix and rhs. Clear them first.
552  matrix.close ();
553  matrix.zero ();
554  rhs.close ();
555  rhs.zero ();
556 
557  // The matrices from which values are added to another matrix
558  // have to be closed. The add() method does take care of
559  // that, but let us do it explicitly.
560  stiffness.close ();
561  damping.close ();
562  mass.close ();
563 
564  STOP_LOG("init phase", "add_M_C_K_helmholtz");
565 
566  START_LOG("global matrix & vector additions", "add_M_C_K_helmholtz");
567 
568  // Add the stiffness and mass with the proper frequency to the
569  // overall system. For this to work properly, the matrix has
570  // to be not only initialized, but filled with the identical
571  // sparsity structure as the matrix added to it, otherwise
572  // solver packages like PETSc crash.
573  //
574  // Note that we have to add the mass and stiffness contributions one
575  // at a time, otherwise the real part of matrix would be fine, but
576  // the imaginary part would be cluttered with unwanted products.
577  matrix.add (scale_stiffness, stiffness);
578  matrix.add (scale_mass, mass);
579  matrix.add (scale_damping, damping);
580  rhs.add (scale_rhs, freq_indep_rhs);
581 
582  STOP_LOG("global matrix & vector additions", "add_M_C_K_helmholtz");
583 
584  // The linear system involving "matrix" and "rhs" is now ready to be solved.
585 }
586 
587 #endif // LIBMESH_USE_COMPLEX_NUMBERS
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
OStreamProxy err
This is the EquationSystems class.
A Node is like a Point, but with more information.
Definition: node.h:52
int main(int argc, char **argv)
void write(const std::string &name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
Write the systems to disk using the XDR data format.
virtual void solve() libmesh_override
Solves the system for all frequencies.
unsigned int size() const
Definition: parallel.h:726
FrequencySystem provides a specific system class for frequency-dependent (linear) systems...
unsigned int dim
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
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
SparseMatrix< Number > & add_matrix(const std::string &mat_name)
Adds the additional matrix mat_name to this system.
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
NumericVector< Number > * rhs
The system matrix.
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
Definition: dof_object.h:810
const T & get(const std::string &) const
Definition: parameters.h:430
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.
virtual void add(const numeric_index_type i, const numeric_index_type j, const T value)=0
Add value to the element (i,j).
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual node_iterator nodes_begin()=0
Iterate over all the nodes in the Mesh.
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
void attach_solve_function(void fptr(EquationSystems &es, const std::string &name))
Register a required user function to use in assembling/solving the system.
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.
NumericVector< Number > & add_vector(const std::string &vec_name, const bool projections=true, const ParallelType type=PARALLEL)
Adds the additional vector vec_name to this system.
Definition: system.C:681
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual void zero()=0
Set all entries to 0.
void print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
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 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.
void prepare_for_use(const bool skip_renumber_nodes_and_elements=false, const bool skip_find_neighbors=false)
Prepare a newly created (or read) mesh for use.
Definition: mesh_base.C:174
void set_frequencies_by_steps(const Number base_freq, const Number freq_step=0., const unsigned int n_freq=1, const bool allocate_solution_duplicates=true)
Set the frequency range for which the system should be solved.
const Parallel::Communicator & comm() const
Definition: libmesh.h:81
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 add_M_C_K_helmholtz(EquationSystems &es, const std::string &system_name)
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 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
virtual node_iterator nodes_end()=0
T & set(const std::string &)
Definition: parameters.h:469
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
The definition of the const_node_iterator struct.
Definition: mesh_base.h:1548
This class implements specific orders of Gauss quadrature.
Parameters parameters
Data structure holding arbitrary parameters.
const MeshBase & get_mesh() const
unsigned int rank() const
Definition: parallel.h:724
void assemble_helmholtz(EquationSystems &es, const std::string &system_name)
const T_sys & get_system(const std::string &name) const
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
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
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
The UNVIO class implements the Ideas UNV universal file format.
Definition: unv_io.h:52
const Real pi
.
Definition: libmesh.h:172
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.