The source file introduction_ex5.C with comments:

Introduction Example 5 - Run-Time Quadrature Rule Selection



This is the fifth example program. It builds on the previous two examples, and extends the use of the \p AutoPtr as a convenient build method to determine the quadrature rule at run time.



C++ include files that we need
        #include <iostream>
        #include <sstream>
        #include <algorithm>
        #include <math.h>
        
Basic include file needed for the mesh functionality.
        #include "libmesh/libmesh.h"
        #include "libmesh/mesh.h"
        #include "libmesh/mesh_generation.h"
        #include "libmesh/exodusII_io.h"
        #include "libmesh/linear_implicit_system.h"
        #include "libmesh/equation_systems.h"
        
Define the Finite Element object.
        #include "libmesh/fe.h"
        
Define the base quadrature class, with which specialized quadrature rules will be built.
        #include "libmesh/quadrature.h"
        
Include the namespace \p QuadratureRules for some handy descriptions.
        #include "libmesh/quadrature_rules.h"
        
Define useful datatypes for finite element matrix and vector components.
        #include "libmesh/sparse_matrix.h"
        #include "libmesh/numeric_vector.h"
        #include "libmesh/dense_matrix.h"
        #include "libmesh/dense_vector.h"
        
Define the DofMap, which handles degree of freedom indexing.
        #include "libmesh/dof_map.h"
        
To impose Dirichlet boundary conditions
        #include "libmesh/dirichlet_boundaries.h"
        #include "libmesh/analytic_function.h"
        
The definition of a geometric element
        #include "libmesh/elem.h"
        
Bring in everything from the libMesh namespace
        using namespace libMesh;
        
        
        
        
        
        
        
Function prototype, as before.
        void assemble_poisson(EquationSystems& es,
                              const std::string& system_name);
        
        
        
Exact solution function prototype, as before.
        Real exact_solution (const Real x,
                             const Real y,
                             const Real z = 0.);
        
Define a wrapper for exact_solution that will be needed below
        void exact_solution_wrapper (DenseVector<Number>& output,
                                     const Point& p,
                                     const Real)
        {
          output(0) = exact_solution(p(0),p(1),p(2));
        }
        
        
The quadrature type the user requests.
        QuadratureType quad_type=INVALID_Q_RULE;
        
        
        
Begin the main program.
        int main (int argc, char** argv)
        {
Initialize libMesh and any dependent libaries, like in example 2.
          LibMeshInit init (argc, argv);
        
Check for proper usage. The quadrature rule must be given at run time.
          if (argc < 3)
            {
              if (libMesh::processor_id() == 0)
                {
                  std::cerr << "Usage: " << argv[0] << " -q n"
                            << std::endl;
                  std::cerr << "  where n stands for:" << std::endl;
        
        
Note that only some of all quadrature rules are valid choices. For example, the Jacobi quadrature is actually a "helper" for higher-order rules, included in QGauss.
                  for (unsigned int n=0; n<QuadratureRules::num_valid_elem_rules; n++)
                    std::cerr << "  " << QuadratureRules::valid_elem_rules[n] << "    "
                              << QuadratureRules::name(QuadratureRules::valid_elem_rules[n])
                              << std::endl;
        
                  std::cerr << std::endl;
                }
        
              libmesh_error();
            }
        
        
Tell the user what we are doing.
          else
            {
              std::cout << "Running " << argv[0];
        
              for (int i=1; i<argc; i++)
                std::cout << " " << argv[i];
        
              std::cout << std::endl << std::endl;
            }
        
        
Set the quadrature rule type that the user wants from argv[2]
          quad_type = static_cast<QuadratureType>(std::atoi(argv[2]));
        
Skip this 3D example if libMesh was compiled as 1D-only.
          libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
        
The following is identical to example 4, and therefore not commented. Differences are mentioned when present.
          Mesh mesh(init.comm());
        
We will use a linear approximation space in this example, hence 8-noded hexahedral elements are sufficient. This is different than example 4 where we used 27-noded hexahedral elements to support a second-order approximation space.
          MeshTools::Generation::build_cube (mesh,
                                             16, 16, 16,
                                             -1., 1.,
                                             -1., 1.,
                                             -1., 1.,
                                             HEX8);
        
          mesh.print_info();
        
          EquationSystems equation_systems (mesh);
        
          equation_systems.add_system<LinearImplicitSystem> ("Poisson");
        
          unsigned int u_var = equation_systems.get_system("Poisson").add_variable("u", FIRST);
        
          equation_systems.get_system("Poisson").attach_assemble_function (assemble_poisson);
        
Construct a Dirichlet boundary condition object

Indicate which boundary IDs we impose the BC on We either build a line, a square or a cube, and here we indicate the boundaries IDs in each case
          std::set<boundary_id_type> boundary_ids;
the dim==1 mesh has two boundaries with IDs 0 and 1
          boundary_ids.insert(0);
          boundary_ids.insert(1);
          boundary_ids.insert(2);
          boundary_ids.insert(3);
          boundary_ids.insert(4);
          boundary_ids.insert(5);
        
Create a vector storing the variable numbers which the BC applies to
          std::vector<unsigned int> variables(1);
          variables[0] = u_var;
        
Create an AnalyticFunction object that we use to project the BC This function just calls the function exact_solution via exact_solution_wrapper
          AnalyticFunction<> exact_solution_object(exact_solution_wrapper);
        
          DirichletBoundary dirichlet_bc(boundary_ids,
                                         variables,
                                         &exact_solution_object);
        
We must add the Dirichlet boundary condition _before_ we call equation_systems.init()
          equation_systems.get_system("Poisson").get_dof_map().add_dirichlet_boundary(dirichlet_bc);
        
          equation_systems.init();
        
          equation_systems.print_info();
        
          equation_systems.get_system("Poisson").solve();
        
"Personalize" the output, with the number of the quadrature rule appended.
          std::ostringstream f_name;
          f_name << "out_" << quad_type << ".e";
        
        #ifdef LIBMESH_HAVE_EXODUS_API
          ExodusII_IO(mesh).write_equation_systems (f_name.str(),
                                              equation_systems);
        #endif // #ifdef LIBMESH_HAVE_EXODUS_API
        
All done.
          return 0;
        }
        
        
        
        
        void assemble_poisson(EquationSystems& es,
                              const std::string& system_name)
        {
          libmesh_assert_equal_to (system_name, "Poisson");
        
          const MeshBase& mesh = es.get_mesh();
        
          const unsigned int dim = mesh.mesh_dimension();
        
          LinearImplicitSystem& system = es.get_system<LinearImplicitSystem>("Poisson");
        
          const DofMap& dof_map = system.get_dof_map();
        
          FEType fe_type = dof_map.variable_type(0);
        
        
Build a Finite Element object of the specified type. Since the \p FEBase::build() member dynamically creates memory we will store the object as an \p AutoPtr. Below, the functionality of \p AutoPtr's is described more detailed in the context of building quadrature rules.
          AutoPtr<FEBase> fe (FEBase::build(dim, fe_type));
        
        
Now this deviates from example 4. we create a 5th order quadrature rule of user-specified type for numerical integration. Note that not all quadrature rules support this order.
          AutoPtr<QBase> qrule(QBase::build(quad_type, dim, THIRD));
        
        
        
Tell the finte element object to use our quadrature rule. Note that a \p AutoPtr returns a QBase* pointer to the object it handles with \p get(). However, using \p get(), the \p AutoPtr \p qrule is still in charge of this pointer. I.e., when \p qrule goes out of scope, it will safely delete the \p QBase object it points to. This behavior may be overridden using \p AutoPtr::release(), but is currently not recommended.
          fe->attach_quadrature_rule (qrule.get());
        
        
Declare a special finite element object for boundary integration.
          AutoPtr<FEBase> fe_face (FEBase::build(dim, fe_type));
        
        
As already seen in example 3, boundary integration requires a quadrature rule. Here, however, we use the more convenient way of building this rule at run-time using \p quad_type. Note that one could also have initialized the face quadrature rules with the type directly determined from \p qrule, namely through: \verbatim AutoPtr qface (QBase::build(qrule->type(), dim-1, THIRD)); \endverbatim And again: using the \p AutoPtr relaxes the need to delete the object afterwards, they clean up themselves.
          AutoPtr<QBase>  qface (QBase::build(quad_type,
                                              dim-1,
                                              THIRD));
        
        
Tell the finte element object to use our quadrature rule. Note that a \p AutoPtr returns a \p QBase* pointer to the object it handles with \p get(). However, using \p get(), the \p AutoPtr \p qface is still in charge of this pointer. I.e., when \p qface goes out of scope, it will safely delete the \p QBase object it points to. This behavior may be overridden using \p AutoPtr::release(), but is not recommended.
          fe_face->attach_quadrature_rule (qface.get());
        
        
        
This is again identical to example 4, and not commented.
          const std::vector<Real>& JxW = fe->get_JxW();
        
          const std::vector<Point>& q_point = fe->get_xyz();
        
          const std::vector<std::vector<Real> >& phi = fe->get_phi();
        
          const std::vector<std::vector<RealGradient> >& dphi = fe->get_dphi();
        
          DenseMatrix<Number> Ke;
          DenseVector<Number> Fe;
        
          std::vector<dof_id_type> dof_indices;
        
        
        
        
        
Now we will loop over all the elements in the mesh. See example 3 for details.
          MeshBase::const_element_iterator       el     = mesh.active_local_elements_begin();
          const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
        
          for ( ; el != end_el; ++el)
            {
              const Elem* elem = *el;
        
              dof_map.dof_indices (elem, dof_indices);
        
              fe->reinit (elem);
        
              Ke.resize (dof_indices.size(),
                         dof_indices.size());
        
              Fe.resize (dof_indices.size());
        
        
        
        
Now loop over the quadrature points. This handles the numeric integration. Note the slightly different access to the QBase members!
              for (unsigned int qp=0; qp<qrule->n_points(); qp++)
                {
Add the matrix contribution
                  for (unsigned int i=0; i<phi.size(); i++)
                    for (unsigned int j=0; j<phi.size(); j++)
                      Ke(i,j) += JxW[qp]*(dphi[i][qp]*dphi[j][qp]);
        
        
fxy is the forcing function for the Poisson equation. In this case we set fxy to be a finite difference Laplacian approximation to the (known) exact solution.

We will use the second-order accurate FD Laplacian approximation, which in 2D on a structured grid is

u_xx + u_yy = (u(i-1,j) + u(i+1,j) + u(i,j-1) + u(i,j+1) + -4*u(i,j))/h^2

Since the value of the forcing function depends only on the location of the quadrature point (q_point[qp]) we will compute it here, outside of the i-loop
                  const Real x = q_point[qp](0);
                  const Real y = q_point[qp](1);
                  const Real z = q_point[qp](2);
                  const Real eps = 1.e-3;
        
                  const Real uxx = (exact_solution(x-eps,y,z) +
                                    exact_solution(x+eps,y,z) +
                                    -2.*exact_solution(x,y,z))/eps/eps;
        
                  const Real uyy = (exact_solution(x,y-eps,z) +
                                    exact_solution(x,y+eps,z) +
                                    -2.*exact_solution(x,y,z))/eps/eps;
        
                  const Real uzz = (exact_solution(x,y,z-eps) +
                                    exact_solution(x,y,z+eps) +
                                    -2.*exact_solution(x,y,z))/eps/eps;
        
                  const Real fxy = - (uxx + uyy + ((dim==2) ? 0. : uzz));
        
        
Add the RHS contribution
                  for (unsigned int i=0; i<phi.size(); i++)
                    Fe(i) += JxW[qp]*fxy*phi[i][qp];
                }
        
If this assembly program were to be used on an adaptive mesh, we would have to apply any hanging node constraint equations Call heterogenously_constrain_element_matrix_and_vector to impose non-homogeneous Dirichlet BCs
              dof_map.heterogenously_constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
        
The element matrix and right-hand-side are now built for this element. Add them to the global matrix and right-hand-side vector. The \p SparseMatrix::add_matrix() and \p NumericVector::add_vector() members do this for us.
              system.matrix->add_matrix (Ke, dof_indices);
              system.rhs->add_vector    (Fe, dof_indices);
        
            } // end of element loop
        
        
        
        
All done!
          return;
        }



The source file introduction_ex5.C without comments:

 
  
  
  #include <iostream>
  #include <sstream>
  #include <algorithm>
  #include <math.h>
  
  #include "libmesh/libmesh.h"
  #include "libmesh/mesh.h"
  #include "libmesh/mesh_generation.h"
  #include "libmesh/exodusII_io.h"
  #include "libmesh/linear_implicit_system.h"
  #include "libmesh/equation_systems.h"
  
  #include "libmesh/fe.h"
  
  #include "libmesh/quadrature.h"
  
  #include "libmesh/quadrature_rules.h"
  
  #include "libmesh/sparse_matrix.h"
  #include "libmesh/numeric_vector.h"
  #include "libmesh/dense_matrix.h"
  #include "libmesh/dense_vector.h"
  
  #include "libmesh/dof_map.h"
  
  #include "libmesh/dirichlet_boundaries.h"
  #include "libmesh/analytic_function.h"
  
  #include "libmesh/elem.h"
  
  using namespace libMesh;
  
  
  
  
  
  
  
  void assemble_poisson(EquationSystems& es,
                        const std::string& system_name);
  
  
  
  Real exact_solution (const Real x,
                       const Real y,
                       const Real z = 0.);
  
  void exact_solution_wrapper (DenseVector<Number>& output,
                               const Point& p,
                               const Real)
  {
    output(0) = exact_solution(p(0),p(1),p(2));
  }
  
  
  QuadratureType quad_type=INVALID_Q_RULE;
  
  
  
  int main (int argc, char** argv)
  {
    LibMeshInit init (argc, argv);
  
    if (argc < 3)
      {
        if (libMesh::processor_id() == 0)
          {
            std::cerr << "Usage: " << argv[0] << " -q n"
                      << std::endl;
            std::cerr << "  where n stands for:" << std::endl;
  
  
            for (unsigned int n=0; n<QuadratureRules::num_valid_elem_rules; n++)
              std::cerr << "  " << QuadratureRules::valid_elem_rules[n] << "    "
                        << QuadratureRules::name(QuadratureRules::valid_elem_rules[n])
                        << std::endl;
  
            std::cerr << std::endl;
          }
  
        libmesh_error();
      }
  
  
    else
      {
        std::cout << "Running " << argv[0];
  
        for (int i=1; i<argc; i++)
          std::cout << " " << argv[i];
  
        std::cout << std::endl << std::endl;
      }
  
  
    quad_type = static_cast<QuadratureType>(std::atoi(argv[2]));
  
    libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
  
    Mesh mesh(init.comm());
  
    MeshTools::Generation::build_cube (mesh,
                                       16, 16, 16,
                                       -1., 1.,
                                       -1., 1.,
                                       -1., 1.,
                                       HEX8);
  
    mesh.print_info();
  
    EquationSystems equation_systems (mesh);
  
    equation_systems.add_system<LinearImplicitSystem> ("Poisson");
  
    unsigned int u_var = equation_systems.get_system("Poisson").add_variable("u", FIRST);
  
    equation_systems.get_system("Poisson").attach_assemble_function (assemble_poisson);
  
  
    std::set<boundary_id_type> boundary_ids;
    boundary_ids.insert(0);
    boundary_ids.insert(1);
    boundary_ids.insert(2);
    boundary_ids.insert(3);
    boundary_ids.insert(4);
    boundary_ids.insert(5);
  
    std::vector<unsigned int> variables(1);
    variables[0] = u_var;
  
    AnalyticFunction<> exact_solution_object(exact_solution_wrapper);
  
    DirichletBoundary dirichlet_bc(boundary_ids,
                                   variables,
                                   &exact_solution_object);
  
    equation_systems.get_system("Poisson").get_dof_map().add_dirichlet_boundary(dirichlet_bc);
  
    equation_systems.init();
  
    equation_systems.print_info();
  
    equation_systems.get_system("Poisson").solve();
  
    std::ostringstream f_name;
    f_name << "out_" << quad_type << ".e";
  
  #ifdef LIBMESH_HAVE_EXODUS_API
    ExodusII_IO(mesh).write_equation_systems (f_name.str(),
                                        equation_systems);
  #endif // #ifdef LIBMESH_HAVE_EXODUS_API
  
    return 0;
  }
  
  
  
  
  void assemble_poisson(EquationSystems& es,
                        const std::string& system_name)
  {
    libmesh_assert_equal_to (system_name, "Poisson");
  
    const MeshBase& mesh = es.get_mesh();
  
    const unsigned int dim = mesh.mesh_dimension();
  
    LinearImplicitSystem& system = es.get_system<LinearImplicitSystem>("Poisson");
  
    const DofMap& dof_map = system.get_dof_map();
  
    FEType fe_type = dof_map.variable_type(0);
  
  
    AutoPtr<FEBase> fe (FEBase::build(dim, fe_type));
  
  
    AutoPtr<QBase> qrule(QBase::build(quad_type, dim, THIRD));
  
  
  
    fe->attach_quadrature_rule (qrule.get());
  
  
    AutoPtr<FEBase> fe_face (FEBase::build(dim, fe_type));
  
  
    AutoPtr<QBase>  qface (QBase::build(quad_type,
                                        dim-1,
                                        THIRD));
  
  
    fe_face->attach_quadrature_rule (qface.get());
  
  
  
    const std::vector<Real>& JxW = fe->get_JxW();
  
    const std::vector<Point>& q_point = fe->get_xyz();
  
    const std::vector<std::vector<Real> >& phi = fe->get_phi();
  
    const std::vector<std::vector<RealGradient> >& dphi = fe->get_dphi();
  
    DenseMatrix<Number> Ke;
    DenseVector<Number> Fe;
  
    std::vector<dof_id_type> dof_indices;
  
  
  
  
  
    MeshBase::const_element_iterator       el     = mesh.active_local_elements_begin();
    const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
  
    for ( ; el != end_el; ++el)
      {
        const Elem* elem = *el;
  
        dof_map.dof_indices (elem, dof_indices);
  
        fe->reinit (elem);
  
        Ke.resize (dof_indices.size(),
                   dof_indices.size());
  
        Fe.resize (dof_indices.size());
  
  
  
  
        for (unsigned int qp=0; qp<qrule->n_points(); qp++)
          {
            for (unsigned int i=0; i<phi.size(); i++)
              for (unsigned int j=0; j<phi.size(); j++)
                Ke(i,j) += JxW[qp]*(dphi[i][qp]*dphi[j][qp]);
  
  
            const Real x = q_point[qp](0);
            const Real y = q_point[qp](1);
            const Real z = q_point[qp](2);
            const Real eps = 1.e-3;
  
            const Real uxx = (exact_solution(x-eps,y,z) +
                              exact_solution(x+eps,y,z) +
                              -2.*exact_solution(x,y,z))/eps/eps;
  
            const Real uyy = (exact_solution(x,y-eps,z) +
                              exact_solution(x,y+eps,z) +
                              -2.*exact_solution(x,y,z))/eps/eps;
  
            const Real uzz = (exact_solution(x,y,z-eps) +
                              exact_solution(x,y,z+eps) +
                              -2.*exact_solution(x,y,z))/eps/eps;
  
            const Real fxy = - (uxx + uyy + ((dim==2) ? 0. : uzz));
  
  
            for (unsigned int i=0; i<phi.size(); i++)
              Fe(i) += JxW[qp]*fxy*phi[i][qp];
          }
  
        dof_map.heterogenously_constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
  
        system.matrix->add_matrix (Ke, dof_indices);
        system.rhs->add_vector    (Fe, dof_indices);
  
      } // end of element loop
  
  
  
  
    return;
  }



The console output of the program:

make[4]: Entering directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/introduction/introduction_ex5'
***************************************************************
* Running Example introduction_ex5:
*  mpirun -np 4 example-devel -q 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
 
Running /net/spark/workspace/roystgnr/libmesh/git/devel/examples/introduction/introduction_ex5/.libs/lt-example-devel -q 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc

 Mesh Information:
  mesh_dimension()=3
  spatial_dimension()=3
  n_nodes()=4913
    n_local_nodes()=1388
  n_elem()=4096
    n_local_elem()=1024
    n_active_elem()=4096
  n_subdomains()=1
  n_partitions()=4
  n_processors()=4
  n_threads()=1
  processor_id()=0

 EquationSystems
  n_systems()=1
   System #0, "Poisson"
    Type "LinearImplicit"
    Variables="u" 
    Finite Element Types="LAGRANGE", "JACOBI_20_00" 
    Infinite Element Mapping="CARTESIAN" 
    Approximation Orders="FIRST", "THIRD" 
    n_dofs()=4913
    n_local_dofs()=1388
    n_constrained_dofs()=1538
    n_local_constrained_dofs()=418
    n_vectors()=1
    n_matrices()=1
    DofMap Sparsity
      Average  On-Processor Bandwidth <= 23.0124
      Average Off-Processor Bandwidth <= 2.04234
      Maximum  On-Processor Bandwidth <= 36
      Maximum Off-Processor Bandwidth <= 21
    DofMap Constraints
      Number of DoF Constraints = 1538
      Number of Heterogenous Constraints= 1474
      Average DoF Constraint Length= 0
      Number of Node Constraints = 0


 -------------------------------------------------------------------------------------------------------------------
| Processor id:   0                                                                                                 |
| Num Processors: 4                                                                                                 |
| Time:           Fri Apr 19 11:44:55 2013                                                                          |
| OS:             Linux                                                                                             |
| HostName:       spark.ices.utexas.edu                                                                             |
| OS Release:     2.6.32-279.22.1.el6.x86_64                                                                        |
| OS Version:     #1 SMP Tue Feb 5 14:33:39 CST 2013                                                                |
| Machine:        x86_64                                                                                            |
| Username:       roystgnr                                                                                          |
| Configuration:  ../configure  '--enable-everything'                                                               |
|  'METHODS=devel'                                                                                                  |
|  '--prefix=/h2/roystgnr/libmesh-test'                                                                             |
|  'CXX=distcc /usr/bin/g++'                                                                                        |
|  'CC=distcc /usr/bin/gcc'                                                                                         |
|  'FC=distcc /usr/bin/gfortran'                                                                                    |
|  'F77=distcc /usr/bin/gfortran'                                                                                   |
|  'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2'                                                          |
|  'PETSC_ARCH=gcc-system-mkl-gf-10.3.12.361-mpich2-1.4.1p1-cxx-opt'                                                |
|  'SLEPC_DIR=/opt/apps/ossw/libraries/slepc/slepc-3.3-p2-petsc-3.3-p2-cxx-opt'                                     |
|  'TRILINOS_DIR=/opt/apps/ossw/libraries/trilinos/trilinos-10.12.2/sl6/gcc-system/mpich2-1.4.1p1/mkl-gf-10.3.12.361'|
|  'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/gcc-system'                                                 |
|  'HDF5_DIR=/opt/apps/ossw/libraries/hdf5/hdf5-1.8.9/sl6/gcc-system'                                               |
 -------------------------------------------------------------------------------------------------------------------
 --------------------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=0.71082, Active time=0.607454                                                      |
 --------------------------------------------------------------------------------------------------------------------
| Event                                  nCalls    Total Time  Avg Time    Total Time  Avg Time    % of Active Time  |
|                                                  w/o Sub     w/o Sub     With Sub    With Sub    w/o S    With S   |
|--------------------------------------------------------------------------------------------------------------------|
|                                                                                                                    |
|                                                                                                                    |
| DofMap                                                                                                             |
|   add_neighbors_to_send_list()         1         0.0098      0.009805    0.0149      0.014925    1.61     2.46     |
|   build_constraint_matrix_and_vector() 1024      0.0018      0.000002    0.0018      0.000002    0.30     0.30     |
|   build_sparsity()                     1         0.0044      0.004387    0.0174      0.017384    0.72     2.86     |
|   create_dof_constraints()             1         0.0201      0.020071    0.0530      0.052972    3.30     8.72     |
|   distribute_dofs()                    1         0.0127      0.012660    0.0363      0.036348    2.08     5.98     |
|   dof_indices()                        7738      0.0531      0.000007    0.0531      0.000007    8.75     8.75     |
|   hetero_cnstrn_elem_mat_vec()         1024      0.0068      0.000007    0.0068      0.000007    1.12     1.12     |
|   prepare_send_list()                  1         0.0001      0.000095    0.0001      0.000095    0.02     0.02     |
|   reinit()                             1         0.0232      0.023249    0.0232      0.023249    3.83     3.83     |
|                                                                                                                    |
| EquationSystems                                                                                                    |
|   build_solution_vector()              1         0.0018      0.001841    0.0132      0.013197    0.30     2.17     |
|                                                                                                                    |
| ExodusII_IO                                                                                                        |
|   write_nodal_data()                   1         0.0522      0.052181    0.0522      0.052181    8.59     8.59     |
|                                                                                                                    |
| FE                                                                                                                 |
|   compute_shape_functions()            1024      0.0067      0.000007    0.0067      0.000007    1.10     1.10     |
|   init_shape_functions()               1         0.0000      0.000031    0.0000      0.000031    0.01     0.01     |
|                                                                                                                    |
| FEMap                                                                                                              |
|   compute_affine_map()                 1024      0.0036      0.000003    0.0036      0.000003    0.58     0.58     |
|   init_reference_to_physical_map()     1         0.0000      0.000038    0.0000      0.000038    0.01     0.01     |
|                                                                                                                    |
| Mesh                                                                                                               |
|   find_neighbors()                     1         0.0216      0.021566    0.0217      0.021664    3.55     3.57     |
|   renumber_nodes_and_elem()            2         0.0018      0.000880    0.0018      0.000880    0.29     0.29     |
|                                                                                                                    |
| MeshCommunication                                                                                                  |
|   compute_hilbert_indices()            2         0.0347      0.017365    0.0347      0.017365    5.72     5.72     |
|   find_global_indices()                2         0.0051      0.002533    0.0413      0.020663    0.83     6.80     |
|   parallel_sort()                      2         0.0010      0.000488    0.0011      0.000536    0.16     0.18     |
|                                                                                                                    |
| MeshOutput                                                                                                         |
|   write_equation_systems()             1         0.0001      0.000059    0.0655      0.065508    0.01     10.78    |
|                                                                                                                    |
| MeshTools::Generation                                                                                              |
|   build_cube()                         1         0.0046      0.004573    0.0046      0.004573    0.75     0.75     |
|                                                                                                                    |
| MetisPartitioner                                                                                                   |
|   partition()                          1         0.1961      0.196062    0.2165      0.216493    32.28    35.64    |
|                                                                                                                    |
| Parallel                                                                                                           |
|   allgather()                          9         0.0001      0.000017    0.0002      0.000021    0.02     0.03     |
|   max(bool)                            1         0.0000      0.000004    0.0000      0.000004    0.00     0.00     |
|   max(scalar)                          105       0.0006      0.000006    0.0006      0.000006    0.10     0.10     |
|   max(vector)                          24        0.0002      0.000006    0.0005      0.000019    0.03     0.08     |
|   min(bool)                            121       0.0005      0.000004    0.0005      0.000004    0.08     0.08     |
|   min(scalar)                          99        0.0257      0.000260    0.0257      0.000260    4.23     4.23     |
|   min(vector)                          24        0.0002      0.000009    0.0005      0.000022    0.04     0.09     |
|   probe()                              36        0.0004      0.000011    0.0004      0.000011    0.07     0.07     |
|   receive()                            36        0.0002      0.000005    0.0006      0.000016    0.03     0.10     |
|   send()                               36        0.0001      0.000003    0.0001      0.000003    0.02     0.02     |
|   send_receive()                       40        0.0003      0.000007    0.0011      0.000026    0.04     0.17     |
|   sum()                                20        0.0004      0.000019    0.0063      0.000316    0.06     1.04     |
|                                                                                                                    |
| Parallel::Request                                                                                                  |
|   wait()                               36        0.0001      0.000001    0.0001      0.000001    0.01     0.01     |
|                                                                                                                    |
| Partitioner                                                                                                        |
|   set_node_processor_ids()             1         0.0038      0.003842    0.0040      0.004024    0.63     0.66     |
|   set_parent_processor_ids()           1         0.0012      0.001172    0.0012      0.001172    0.19     0.19     |
|                                                                                                                    |
| PetscLinearSolver                                                                                                  |
|   solve()                              1         0.0931      0.093097    0.0931      0.093097    15.33    15.33    |
|                                                                                                                    |
| System                                                                                                             |
|   assemble()                           1         0.0195      0.019491    0.0445      0.044461    3.21     7.32     |
 --------------------------------------------------------------------------------------------------------------------
| Totals:                                12447     0.6075                                          100.00            |
 --------------------------------------------------------------------------------------------------------------------

 
***************************************************************
* Done Running Example introduction_ex5:
*  mpirun -np 4 example-devel -q 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
make[4]: Leaving directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/introduction/introduction_ex5'

Site Created By: libMesh Developers
Last modified: April 23 2013 04:19:30 UTC

Hosted By:
SourceForge.net Logo