📄 ex8.php
字号:
<?php $root=""; ?><?php require($root."navigation.php"); ?><html><head> <?php load_style($root); ?></head> <body> <?php make_navigation("ex8",$root)?> <div class="content"><a name="comments"></a> <div class = "comment"><h1>Example 8 - The Wave Equation</h1><br><br>This is the eighth example program. It builds onthe previous example programs. It introduces theNewmarkSystem class. In this example the wave equationis solved using the time integration scheme providedby the NewmarkSystem class.<br><br>This example comes with a cylindrical mesh given in theuniversal file pipe-mesh.unv.The mesh contains HEX8 and PRISM6 elements.<br><br>C++ include files that we need</div><div class ="fragment"><pre> #include <iostream> #include <fstream> #include <algorithm> #include <stdio.h> #include <math.h> </pre></div><div class = "comment">Basic include file needed for the mesh functionality.</div><div class ="fragment"><pre> #include "libmesh.h" #include "mesh.h" #include "gmv_io.h" #include "newmark_system.h" #include "equation_systems.h" </pre></div><div class = "comment">Define the Finite Element object.</div><div class ="fragment"><pre> #include "fe.h" </pre></div><div class = "comment">Define Gauss quadrature rules.</div><div class ="fragment"><pre> #include "quadrature_gauss.h" </pre></div><div class = "comment">Define useful datatypes for finite element</div><div class ="fragment"><pre> #include "dense_matrix.h" #include "dense_vector.h" </pre></div><div class = "comment">Define matrix and vector data types for the global equation system. These are base classes,from which specific implementations, likethe PETSc or LASPACK implementations, are derived.</div><div class ="fragment"><pre> #include "sparse_matrix.h" #include "numeric_vector.h" </pre></div><div class = "comment">Define the DofMap, which handles degree of freedomindexing.</div><div class ="fragment"><pre> #include "dof_map.h" </pre></div><div class = "comment">The definition of a vertex associated with a Mesh.</div><div class ="fragment"><pre> #include "node.h" </pre></div><div class = "comment">The definition of a geometric element</div><div class ="fragment"><pre> #include "elem.h" </pre></div><div class = "comment">Defines the MeshData class, which allows you to storedata about the mesh when reading in files, etc.</div><div class ="fragment"><pre> #include "mesh_data.h" </pre></div><div class = "comment">Function prototype. This is the function that will assemblethe linear system for our problem, governed by the linearwave equation.</div><div class ="fragment"><pre> void assemble_wave(EquationSystems& es, const std::string& system_name); </pre></div><div class = "comment">Function Prototype. This function will be used to apply theinitial conditions.</div><div class ="fragment"><pre> void apply_initial(EquationSystems& es, const std::string& system_name); </pre></div><div class = "comment">Function Prototype. This function imposesDirichlet Boundary conditions via the penaltymethod after the system is assembled.</div><div class ="fragment"><pre> void fill_dirichlet_bc(EquationSystems& es, const std::string& system_name); </pre></div><div class = "comment">The main program</div><div class ="fragment"><pre> int main (int argc, char** argv) {</pre></div><div class = "comment">Initialize Petsc, like in example 2.</div><div class ="fragment"><pre> libMesh::init (argc, argv); </pre></div><div class = "comment">Braces are used to force object scope.</div><div class ="fragment"><pre> {</pre></div><div class = "comment">Check for proper usage.</div><div class ="fragment"><pre> if (argc < 2) { std::cerr << "Usage: " << argv[0] << " [meshfile]" << std::endl; error(); } </pre></div><div class = "comment">Tell the user what we are doing.</div><div class ="fragment"><pre> else { std::cout << "Running " << argv[0]; for (int i=1; i<argc; i++) std::cout << " " << argv[i]; std::cout << std::endl << std::endl; } </pre></div><div class = "comment">LasPack solvers don't work so well for this example(not sure why). Print a warning to the user if PETScis not available, or if they are using LasPack solvers.</div><div class ="fragment"><pre> #ifdef HAVE_PETSC if ((libMesh::on_command_line("--use-laspack")) || (libMesh::on_command_line("--disable-petsc"))) #endif { std::cerr << "WARNING! It appears you are using the\n" << "LasPack solvers. ex8 is known not to converge\n" << "using LasPack, but should work OK with PETSc.\n" << "If possible, download and install the PETSc\n" << "library from www-unix.mcs.anl.gov/petsc/petsc-2/\n" << std::endl; } </pre></div><div class = "comment">Get the name of the mesh filefrom the command line.</div><div class ="fragment"><pre> std::string mesh_file = argv[1]; std::cout << "Mesh file is: " << mesh_file << std::endl; </pre></div><div class = "comment">For now, restrict to dim=3, though thismay easily be changed, see example 4</div><div class ="fragment"><pre> const unsigned int dim = 3; </pre></div><div class = "comment">Create a dim-dimensional mesh.</div><div class ="fragment"><pre> Mesh mesh (dim); MeshData mesh_data(mesh); </pre></div><div class = "comment">Read the meshfile specified in the command line oruse the internal mesh generator to create a uniformgrid on an elongated cube.</div><div class ="fragment"><pre> mesh.read(mesh_file, &mesh_data); </pre></div><div class = "comment">mesh.build_cube (10, 10, 40,-1., 1.,-1., 1.,0., 4.,HEX8);<br><br>Print information about the mesh to the screen.</div><div class ="fragment"><pre> mesh.print_info(); </pre></div><div class = "comment">The node that should be monitored.</div><div class ="fragment"><pre> const unsigned int result_node = 274; </pre></div><div class = "comment">Time stepping issues<br><br>Note that the total current time is stored as a parameterin the \pEquationSystems object.<br><br>the time step size</div><div class ="fragment"><pre> const Real delta_t = .0000625; </pre></div><div class = "comment">The number of time steps.</div><div class ="fragment"><pre> unsigned int n_time_steps = 300; </pre></div><div class = "comment">Create an equation systems object.</div><div class ="fragment"><pre> EquationSystems equation_systems (mesh); </pre></div><div class = "comment">Declare the system and its variables.</div><div class ="fragment"><pre> {</pre></div><div class = "comment">Creates a NewmarkSystem named "Wave"</div><div class ="fragment"><pre> equation_systems.add_system<NewmarkSystem> ("Wave"); </pre></div><div class = "comment">Use a handy reference to this system</div><div class ="fragment"><pre> NewmarkSystem & t_system = equation_systems.get_system<NewmarkSystem> ("Wave"); </pre></div><div class = "comment">Adds the variable "p" to "Wave". "p"will be approximated using first-order approximation.</div><div class ="fragment"><pre> t_system.add_variable("p", FIRST); </pre></div><div class = "comment">Give the system a pointer to the matrix assemblyfunction and the initial condition function definedbelow.</div><div class ="fragment"><pre> t_system.attach_assemble_function (assemble_wave); t_system.attach_init_function (apply_initial); </pre></div><div class = "comment">Set the time step size, and optionally theNewmark parameters, so that \p NewmarkSystem can compute integration constants. Here we simply use pass only the time step and use default values for \p alpha=.25 and \p delta=.5.</div><div class ="fragment"><pre> t_system.set_newmark_parameters(delta_t); </pre></div><div class = "comment">Set the speed of sound and fluid densityas \p EquationSystems parameter,so that \p assemble_wave() can access it.</div><div class ="fragment"><pre> equation_systems.parameters.set<Real>("speed") = 1000.; equation_systems.parameters.set<Real>("fluid density") = 1000.; </pre></div><div class = "comment">Store the current time as an\p EquationSystems parameter, so that\p fill_dirichlet_bc() can access it.</div><div class ="fragment"><pre> equation_systems.parameters.set<Real>("time") = 0.; ⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -