libMesh
Functions
reduced_basis_ex5.C File Reference

Go to the source code of this file.

Functions

void scale_mesh_and_plot (EquationSystems &es, const RBParameters &mu, const std::string &filename)
 
void compute_stresses (EquationSystems &es)
 
int main (int argc, char **argv)
 

Function Documentation

void compute_stresses ( EquationSystems es)

Definition at line 346 of file reduced_basis_ex5.C.

References libMesh::MeshBase::active_local_element_ptr_range(), libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::System::current_solution(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::System::n_vars(), std::pow(), libMesh::DenseMatrix< T >::resize(), libMesh::DenseMatrix< T >::scale(), libMesh::System::solution, libMesh::System::update(), libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

Referenced by scale_mesh_and_plot().

347 {
348  LOG_SCOPE("compute_stresses()", "main");
349 
350  const MeshBase & mesh = es.get_mesh();
351 
352  const unsigned int dim = mesh.mesh_dimension();
353 
354  ElasticityRBConstruction & system = es.get_system<ElasticityRBConstruction>("RBElasticity");
355 
356  unsigned int displacement_vars[3];
357  displacement_vars[0] = system.variable_number ("u");
358  displacement_vars[1] = system.variable_number ("v");
359  displacement_vars[2] = system.variable_number ("w");
360  const unsigned int u_var = system.variable_number ("u");
361 
362  const DofMap & dof_map = system.get_dof_map();
363  FEType fe_type = dof_map.variable_type(u_var);
364  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
365  QGauss qrule (dim, fe_type.default_quadrature_order());
366  fe->attach_quadrature_rule (&qrule);
367 
368  const std::vector<Real> & JxW = fe->get_JxW();
369  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
370 
371  // Also, get a reference to the ExplicitSystem
372  ExplicitSystem & stress_system = es.get_system<ExplicitSystem>("StressSystem");
373  const DofMap & stress_dof_map = stress_system.get_dof_map();
374  unsigned int sigma_vars[3][3];
375  sigma_vars[0][0] = stress_system.variable_number ("sigma_00");
376  sigma_vars[0][1] = stress_system.variable_number ("sigma_01");
377  sigma_vars[0][2] = stress_system.variable_number ("sigma_02");
378  sigma_vars[1][0] = stress_system.variable_number ("sigma_10");
379  sigma_vars[1][1] = stress_system.variable_number ("sigma_11");
380  sigma_vars[1][2] = stress_system.variable_number ("sigma_12");
381  sigma_vars[2][0] = stress_system.variable_number ("sigma_20");
382  sigma_vars[2][1] = stress_system.variable_number ("sigma_21");
383  sigma_vars[2][2] = stress_system.variable_number ("sigma_22");
384  unsigned int vonMises_var = stress_system.variable_number ("vonMises");
385 
386  // Storage for the stress dof indices on each element
387  std::vector<std::vector<dof_id_type>> dof_indices_var(system.n_vars());
388  std::vector<dof_id_type> stress_dof_indices_var;
389 
390  // To store the stress tensor on each element
391  DenseMatrix<Number> elem_sigma;
392 
393  for (const auto & elem : mesh.active_local_element_ptr_range())
394  {
395  for (unsigned int var=0; var<3; var++)
396  dof_map.dof_indices (elem, dof_indices_var[var], displacement_vars[var]);
397 
398  fe->reinit (elem);
399 
400  elem_sigma.resize(3, 3);
401 
402  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
403  for (unsigned int C_i=0; C_i<3; C_i++)
404  for (unsigned int C_j=0; C_j<3; C_j++)
405  for (unsigned int C_k=0; C_k<3; C_k++)
406  {
407  const unsigned int n_var_dofs = dof_indices_var[C_k].size();
408 
409  // Get the gradient at this quadrature point
410  Gradient displacement_gradient;
411  for (unsigned int l=0; l<n_var_dofs; l++)
412  displacement_gradient.add_scaled(dphi[l][qp], system.current_solution(dof_indices_var[C_k][l]));
413 
414  for (unsigned int C_l=0; C_l<3; C_l++)
415  elem_sigma(C_i,C_j) += JxW[qp] * (elasticity_tensor(C_i, C_j, C_k, C_l) * displacement_gradient(C_l));
416  }
417 
418  // Get the average stresses by dividing by the element volume
419  elem_sigma.scale(1./elem->volume());
420 
421  // load elem_sigma data into stress_system
422  for (unsigned int i=0; i<3; i++)
423  for (unsigned int j=0; j<3; j++)
424  {
425  stress_dof_map.dof_indices (elem, stress_dof_indices_var, sigma_vars[i][j]);
426 
427  // We are using CONSTANT MONOMIAL basis functions, hence we only need to get
428  // one dof index per variable
429  dof_id_type dof_index = stress_dof_indices_var[0];
430 
431  if ((stress_system.solution->first_local_index() <= dof_index) &&
432  (dof_index < stress_system.solution->last_local_index()))
433  {
434  stress_system.solution->set(dof_index, elem_sigma(i,j));
435  }
436 
437  }
438 
439  // Also, the von Mises stress
440  Number vonMises_value = std::sqrt(0.5*(pow(elem_sigma(0,0) - elem_sigma(1,1),2.) +
441  pow(elem_sigma(1,1) - elem_sigma(2,2),2.) +
442  pow(elem_sigma(2,2) - elem_sigma(0,0),2.) +
443  6.*(pow(elem_sigma(0,1),2.) + pow(elem_sigma(1,2),2.) + pow(elem_sigma(2,0),2.))
444  ));
445  stress_dof_map.dof_indices (elem, stress_dof_indices_var, vonMises_var);
446  dof_id_type dof_index = stress_dof_indices_var[0];
447  if ((stress_system.solution->first_local_index() <= dof_index) &&
448  (dof_index < stress_system.solution->last_local_index()))
449  {
450  stress_system.solution->set(dof_index, vonMises_value);
451  }
452 
453  }
454 
455  // Should call close and update when we set vector entries directly
456  stress_system.solution->close();
457  stress_system.update();
458 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
void add_scaled(const TypeVector< T2 > &, const T)
Add a scaled value to this vector without creating a temporary.
Definition: type_vector.h:624
unsigned int dim
Real elasticity_tensor(unsigned int i, unsigned int j, unsigned int k, unsigned int l)
Definition: assembly.C:38
MeshBase & mesh
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
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1263
const DofMap & get_dof_map() const
Definition: system.h:2030
Order default_quadrature_order() const
Definition: fe_type.h:332
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
double pow(double a, int b)
void scale(const T factor)
Multiplies every element in the matrix by factor.
Definition: dense_matrix.h:851
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:425
Number current_solution(const dof_id_type global_dof_number) const
Definition: system.C:192
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
unsigned int n_vars() const
Definition: system.h:2086
The ExplicitSystem provides only "right hand side" storage, which should be sufficient for solving mo...
uint8_t dof_id_type
Definition: id_types.h:64
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 82 of file reduced_basis_ex5.C.

References libMesh::BoundaryInfo::add_node(), libMesh::EquationSystems::add_system(), libMesh::System::add_variable(), libMesh::MeshTools::Generation::build_cube(), libMesh::CG, libMesh::LibMeshInit::comm(), libMesh::ParallelObject::comm(), libMesh::CONSTANT, libMesh::default_solver_package(), dim, libMesh::EIGEN_SOLVERS, libMesh::MeshBase::element_ptr_range(), libMesh::MeshBase::get_boundary_info(), libMesh::LinearImplicitSystem::get_linear_solver(), libMesh::RBConstruction::get_rb_evaluation(), libMesh::BoundaryInfo::has_boundary_id(), libMesh::HEX8, libMesh::TriangleWrapper::init(), libMesh::EquationSystems::init(), libMesh::RBConstruction::initialize_rb_construction(), libMesh::INVALID_SOLVER_PACKAGE, libMesh::RBEvaluation::legacy_write_offline_data_to_files(), libMesh::RBConstruction::load_rb_solution(), mesh, libMesh::MONOMIAL, libMesh::RBConstruction::print_info(), libMesh::EquationSystems::print_info(), libMesh::MeshBase::print_info(), libMesh::RBConstruction::process_parameters_file(), libMesh::RBDataDeserialization::RBEvaluationDeserialization::read_from_file(), libMesh::Real, scale_mesh_and_plot(), libMesh::RBConstruction::set_rb_evaluation(), libMesh::LinearSolver< T >::set_solver_type(), libMesh::RBParameters::set_value(), side, libMesh::RBConstruction::train_reduced_basis(), libMesh::TRILINOS_SOLVERS, libMesh::RBEvaluation::write_out_basis_functions(), and libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file().

83 {
84  // Initialize libMesh.
85  LibMeshInit init (argc, argv);
86 
87  // This example requires a linear solver package.
88  libmesh_example_requires(libMesh::default_solver_package() != INVALID_SOLVER_PACKAGE,
89  "--enable-petsc, --enable-trilinos, or --enable-eigen");
90 
91 #if !defined(LIBMESH_HAVE_XDR)
92  // We need XDR support to write out reduced bases
93  libmesh_example_requires(false, "--enable-xdr");
94 #elif defined(LIBMESH_DEFAULT_SINGLE_PRECISION)
95  // XDR binary support requires double precision
96  libmesh_example_requires(false, "--disable-singleprecision");
97 #elif defined(LIBMESH_DEFAULT_TRIPLE_PRECISION)
98  // I have no idea why long double isn't working here... [RHS]
99  libmesh_example_requires(false, "double precision");
100 #endif
101 
102  // Eigen can take forever to solve the offline mode portion of this
103  // example
104  libmesh_example_requires(libMesh::default_solver_package() != EIGEN_SOLVERS, "--enable-petsc or --enable-laspack");
105 
106  // Trilinos reports "true residual is too large" in the offline mode
107  // portion of this example
108  libmesh_example_requires(libMesh::default_solver_package() != TRILINOS_SOLVERS, "--enable-petsc or --enable-laspack");
109 
110  // This example only works if libMesh was compiled for 3D
111  const unsigned int dim = 3;
112  libmesh_example_requires(dim == LIBMESH_DIM, "3D support");
113 
114  std::string parameters_filename = "reduced_basis_ex5.in";
115  GetPot infile(parameters_filename);
116 
117  unsigned int n_elem_x = infile("n_elem_x", 0);
118  unsigned int n_elem_y = infile("n_elem_y", 0);
119  unsigned int n_elem_z = infile("n_elem_z", 0);
120  Real x_size = infile("x_size", 0.);
121  Real y_size = infile("y_size", 0.);
122  Real z_size = infile("z_size", 0.);
123 
124  bool store_basis_functions = infile("store_basis_functions", true);
125 
126  // Read the "online_mode" flag from the command line
127  GetPot command_line(argc, argv);
128  int online_mode = 0;
129  if (command_line.search(1, "-online_mode"))
130  online_mode = command_line.next(online_mode);
131 
132 
133  Mesh mesh (init.comm(), dim);
135  n_elem_x,
136  n_elem_y,
137  n_elem_z,
138  0., x_size,
139  0., y_size,
140  0., z_size,
141  HEX8);
142 
143  // Let's add a node boundary condition so that we can impose a point
144  // load.
145  // Each processor should know about each boundary condition it can
146  // see, so we loop over all elements, not just local elements.
147  for (const auto & elem : mesh.element_ptr_range())
148  {
149  unsigned int side_max_x = 0, side_max_y = 0, side_max_z = 0;
150  bool found_side_max_x = false, found_side_max_y = false, found_side_max_z = false;
151  for (auto side : elem->side_index_range())
152  {
153  if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MAX_X))
154  {
155  side_max_x = side;
156  found_side_max_x = true;
157  }
158 
159  if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MAX_Y))
160  {
161  side_max_y = side;
162  found_side_max_y = true;
163  }
164 
165  if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MAX_Z))
166  {
167  side_max_z = side;
168  found_side_max_z = true;
169  }
170  }
171 
172  // If elem has sides on boundaries
173  // BOUNDARY_ID_MAX_X, BOUNDARY_ID_MAX_Y, BOUNDARY_ID_MAX_Z
174  // then let's set a node boundary condition
175  if (found_side_max_x && found_side_max_y && found_side_max_z)
176  {
177  for (auto n : elem->node_index_range())
178  {
179  if (elem->is_node_on_side(n, side_max_x) &&
180  elem->is_node_on_side(n, side_max_y) &&
181  elem->is_node_on_side(n, side_max_z))
182  {
183  mesh.get_boundary_info().add_node(elem->node_ptr(n), NODE_BOUNDARY_ID);
184  }
185  }
186  }
187  }
188 
189  mesh.print_info();
190 
191  // Create an equation systems object.
192  EquationSystems equation_systems(mesh);
193 
194  // We override RBConstruction with ElasticityRBConstruction in order to
195  // specialize a few functions for this particular problem.
196  ElasticityRBConstruction & rb_con =
197  equation_systems.add_system<ElasticityRBConstruction>("RBElasticity");
198 
199  // Also, initialize an ExplicitSystem to store stresses
200  ExplicitSystem & stress_system =
201  equation_systems.add_system<ExplicitSystem> ("StressSystem");
202  stress_system.add_variable("sigma_00", CONSTANT, MONOMIAL);
203  stress_system.add_variable("sigma_01", CONSTANT, MONOMIAL);
204  stress_system.add_variable("sigma_02", CONSTANT, MONOMIAL);
205  stress_system.add_variable("sigma_10", CONSTANT, MONOMIAL);
206  stress_system.add_variable("sigma_11", CONSTANT, MONOMIAL);
207  stress_system.add_variable("sigma_12", CONSTANT, MONOMIAL);
208  stress_system.add_variable("sigma_20", CONSTANT, MONOMIAL);
209  stress_system.add_variable("sigma_21", CONSTANT, MONOMIAL);
210  stress_system.add_variable("sigma_22", CONSTANT, MONOMIAL);
211  stress_system.add_variable("vonMises", CONSTANT, MONOMIAL);
212 
213  // Initialize the data structures for the equation system.
214  equation_systems.init ();
215  equation_systems.print_info();
216 
217  // Build a new RBEvaluation object which will be used to perform
218  // Reduced Basis calculations. This is required in both the
219  // "Offline" and "Online" stages.
220  ElasticityRBEvaluation rb_eval(mesh.comm());
221 
222  // We need to give the RBConstruction object a pointer to
223  // our RBEvaluation object
224  rb_con.set_rb_evaluation(rb_eval);
225 
226  if (!online_mode) // Perform the Offline stage of the RB method
227  {
228  // Use CG solver. This echoes the Petsc command line options set
229  // in run.sh, but also works for non-PETSc linear solvers.
231 
232  // Read in the data that defines this problem from the specified text file
233  rb_con.process_parameters_file(parameters_filename);
234 
235  // Print out info that describes the current setup of rb_con
236  rb_con.print_info();
237 
238  // Prepare rb_con for the Construction stage of the RB method.
239  // This sets up the necessary data structures and performs
240  // initial assembly of the "truth" affine expansion of the PDE.
242 
243  // Compute the reduced basis space by computing "snapshots", i.e.
244  // "truth" solves, at well-chosen parameter values and employing
245  // these snapshots as basis functions.
246  rb_con.train_reduced_basis();
247 
248  // Write out the data that will subsequently be required for the Evaluation stage
249 #if defined(LIBMESH_HAVE_CAPNPROTO)
251  rb_eval_writer.write_to_file("rb_eval.bin");
252 #else
254 #endif
255 
256  // If requested, write out the RB basis functions for visualization purposes
257  if (store_basis_functions)
258  {
259  // Write out the basis functions
261  }
262  }
263  else // Perform the Online stage of the RB method
264  {
265  // Read in the reduced basis data
266 #if defined(LIBMESH_HAVE_CAPNPROTO)
268  rb_eval_reader.read_from_file("rb_eval.bin", /*read_error_bound_data*/ true);
269 #else
270  rb_eval.legacy_read_offline_data_from_files();
271 #endif
272 
273  // Initialize online parameters
274  Real online_x_scaling = infile("online_x_scaling", 0.);
275  Real online_load_Fx = infile("online_load_Fx", 0.);
276  Real online_load_Fy = infile("online_load_Fy", 0.);
277  Real online_load_Fz = infile("online_load_Fz", 0.);
278  Real online_point_load_Fx = infile("online_point_load_Fx", 0.);
279  Real online_point_load_Fy = infile("online_point_load_Fy", 0.);
280  Real online_point_load_Fz = infile("online_point_load_Fz", 0.);
281  RBParameters online_mu;
282  online_mu.set_value("x_scaling", online_x_scaling);
283  online_mu.set_value("load_Fx", online_load_Fx);
284  online_mu.set_value("load_Fy", online_load_Fy);
285  online_mu.set_value("load_Fz", online_load_Fz);
286  online_mu.set_value("point_load_Fx", online_point_load_Fx);
287  online_mu.set_value("point_load_Fy", online_point_load_Fy);
288  online_mu.set_value("point_load_Fz", online_point_load_Fz);
289  rb_eval.set_parameters(online_mu);
290  rb_eval.print_parameters();
291 
292  // Now do the Online solve using the precomputed reduced basis
293  rb_eval.rb_solve(rb_eval.get_n_basis_functions());
294 
295  if (store_basis_functions)
296  {
297  // Read in the basis functions
298  rb_eval.read_in_basis_functions(rb_con);
299 
300  // Plot the solution
301  rb_con.load_rb_solution();
302 
303  const RBParameters & rb_eval_params = rb_eval.get_parameters();
304  scale_mesh_and_plot(equation_systems, rb_eval_params, "RB_sol.e");
305  }
306  }
307 
308  return 0;
309 }
virtual void initialize_rb_construction(bool skip_matrix_assembly=false, bool skip_vector_assembly=false)
Allocate all the data structures necessary for the construction stage of the RB method.
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:117
This is the EquationSystems class.
virtual void write_out_basis_functions(System &sys, const std::string &directory_name="offline_data", const bool write_binary_basis_functions=true)
Write out all the basis functions to file.
virtual Real train_reduced_basis(const bool resize_rb_eval_data=true)
Train the reduced basis.
This class serializes an RBEvaluation object using the Cap&#39;n Proto library.
unsigned int dim
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
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
SolverPackage default_solver_package()
Definition: libmesh.C:995
void add_node(const Node *node, const boundary_id_type id)
Add Node node with boundary id id to the boundary information data structures.
virtual LinearSolver< Number > * get_linear_solver() const libmesh_override
virtual SimpleRange< element_iterator > element_ptr_range()=0
virtual void legacy_write_offline_data_to_files(const std::string &directory_name="offline_data", const bool write_binary_data=true)
Write out all the data to text files in order to segregate the Offline stage from the Online stage...
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
virtual void print_info()
Print out info that describes the current setup of this RBConstruction.
void set_solver_type(const SolverType st)
Sets the type of solver to use.
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:42
void scale_mesh_and_plot(EquationSystems &es, const RBParameters &mu, const std::string &filename)
void write_to_file(const std::string &path)
Write the Cap&#39;n&#39;Proto buffer to disk.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void set_value(const std::string &param_name, Real value)
Set the value of the specified parameter.
Definition: rb_parameters.C:57
const Parallel::Communicator & comm() const
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
virtual void process_parameters_file(const std::string &parameters_filename)
Read in from the file specified by parameters_filename and set the this system&#39;s member variables acc...
This class de-serializes an RBEvaluation object using the Cap&#39;n Proto library.
RBEvaluation & get_rb_evaluation()
Get a reference to the RBEvaluation object.
virtual void load_rb_solution()
Load the RB solution from the most recent solve with rb_eval into this system&#39;s solution vector...
void set_rb_evaluation(RBEvaluation &rb_eval_in)
Set the RBEvaluation object.
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
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.
The ExplicitSystem provides only "right hand side" storage, which should be sufficient for solving mo...
void scale_mesh_and_plot ( EquationSystems es,
const RBParameters mu,
const std::string &  filename 
)

Definition at line 311 of file reduced_basis_ex5.C.

References compute_stresses(), libMesh::EquationSystems::get_mesh(), libMesh::RBParameters::get_value(), mesh, libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), and libMesh::MeshOutput< MT >::write_equation_systems().

Referenced by main().

314 {
315  // Loop over the mesh nodes and move them!
316  MeshBase & mesh = es.get_mesh();
317 
318  MeshBase::node_iterator node_it = mesh.nodes_begin();
319  const MeshBase::node_iterator node_end = mesh.nodes_end();
320 
321  for ( ; node_it != node_end; node_it++)
322  {
323  Node * node = *node_it;
324 
325  (*node)(0) *= mu.get_value("x_scaling");
326  }
327 
328  // Post-process the solution to compute the stresses
329  compute_stresses(es);
330 
331 #ifdef LIBMESH_HAVE_EXODUS_API
332  ExodusII_IO (mesh).write_equation_systems (filename, es);
333 #endif
334 
335  // Loop over the mesh nodes and move them!
336  node_it = mesh.nodes_begin();
337 
338  for ( ; node_it != node_end; node_it++)
339  {
340  Node * node = *node_it;
341 
342  (*node)(0) /= mu.get_value("x_scaling");
343  }
344 }
A Node is like a Point, but with more information.
Definition: node.h:52
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
Real get_value(const std::string &param_name) const
Get the value of the specific parameter.
Definition: rb_parameters.C:45
MeshBase & mesh
This is the MeshBase class.
Definition: mesh_base.h:68
virtual node_iterator nodes_begin()=0
Iterate over all the nodes in the Mesh.
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
The definition of the node_iterator struct.
Definition: mesh_base.h:1528
virtual node_iterator nodes_end()=0
const MeshBase & get_mesh() const
void compute_stresses(EquationSystems &es)