The source file exact_solution.C with comments:
#include <math.h>
Mesh library includes
#include "libmesh/libmesh_common.h"
Bring in everything from the libMesh namespace
using namespace libMesh;
/**
* This is the exact solution that
* we are trying to obtain. We will solve
*
* - (u_xx + u_yy) = f
*
* and take a finite difference approximation using this
* function to get f. This is the well-known "method of
* manufactured solutions".
*/
Real exact_solution (const Real x,
const Real y,
const Real z = 0.)
{
static const Real pi = acos(-1.);
return cos(.5*pi*x)*sin(.5*pi*y)*cos(.5*pi*z);
}
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;
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
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
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 exact_solution.C without comments:
#include <math.h>
#include "libmesh/libmesh_common.h"
using namespace libMesh;
/**
* This is the exact solution that
* we are trying to obtain. We will solve
*
* - (u_xx + u_yy) = f
*
* and take a finite difference approximation using this
* function to get f. This is the well-known "method of
* manufactured solutions".
*/
Real exact_solution (const Real x,
const Real y,
const Real z = 0.)
{
static const Real pi = acos(-1.);
return cos(.5*pi*x)*sin(.5*pi*y)*cos(.5*pi*z);
}
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;
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:
***************************************************************
* Running Example introduction_ex5:
* mpirun -np 2 example-devel -q 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
Running /workspace/libmesh/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 -log_summary
Mesh Information:
mesh_dimension()=3
spatial_dimension()=3
n_nodes()=4913
n_local_nodes()=2618
n_elem()=4096
n_local_elem()=2048
n_active_elem()=4096
n_subdomains()=1
n_partitions()=2
n_processors()=2
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()=2618
n_constrained_dofs()=1538
n_local_constrained_dofs()=803
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 23.4777
Average Off-Processor Bandwidth <= 1.0173
Maximum On-Processor Bandwidth <= 33
Maximum Off-Processor Bandwidth <= 12
DofMap Constraints
Number of DoF Constraints = 1538
Number of Heterogenous Constraints= 1474
Average DoF Constraint Length= 0
Number of Node Constraints = 0
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/introduction/introduction_ex5/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:38:06 2013
Using Petsc Release Version 3.3.0, Patch 2, Fri Jul 13 15:42:00 CDT 2012
Max Max/Min Avg Total
Time (sec): 1.610e+00 1.00000 1.610e+00
Objects: 6.200e+01 1.00000 6.200e+01
Flops: 4.897e+07 1.27983 4.362e+07 8.724e+07
Flops/sec: 3.042e+07 1.27983 2.710e+07 5.419e+07
MPI Messages: 3.250e+01 1.00000 3.250e+01 6.500e+01
MPI Message Lengths: 1.524e+05 1.00000 4.688e+03 3.047e+05
MPI Reductions: 9.500e+01 1.00000
Flop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)
e.g., VecAXPY() for real vectors of length N --> 2N flops
and VecAXPY() for complex vectors of length N --> 8N flops
Summary of Stages: ----- Time ------ ----- Flops ----- --- Messages --- -- Message Lengths -- -- Reductions --
Avg %Total Avg %Total counts %Total Avg %Total counts %Total
0: Main Stage: 1.6097e+00 100.0% 8.7239e+07 100.0% 6.500e+01 100.0% 4.688e+03 100.0% 9.400e+01 98.9%
------------------------------------------------------------------------------------------------------------------------
See the 'Profiling' chapter of the users' manual for details on interpreting output.
Phase summary info:
Count: number of times phase was executed
Time and Flops: Max - maximum over all processors
Ratio - ratio of maximum to minimum over all processors
Mess: number of messages sent
Avg. len: average message length
Reduct: number of global reductions
Global: entire computation
Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().
%T - percent time in this phase %f - percent flops in this phase
%M - percent messages in this phase %L - percent message lengths in this phase
%R - percent reductions in this phase
Total Mflop/s: 10e-6 * (sum of flops over all processors)/(max time over all processors)
------------------------------------------------------------------------------------------------------------------------
Event Count Time (sec) Flops --- Global --- --- Stage --- Total
Max Ratio Max Ratio Max Ratio Mess Avg len Reduct %T %f %M %L %R %T %f %M %L %R Mflop/s
------------------------------------------------------------------------------------------------------------------------
--- Event Stage 0: Main Stage
VecMDot 24 1.0 1.9994e-03 4.2 1.57e+06 1.1 0.0e+00 0.0e+00 2.4e+01 0 3 0 0 25 0 3 0 0 26 1474
VecNorm 26 1.0 2.1601e-04 1.7 1.36e+05 1.1 0.0e+00 0.0e+00 2.6e+01 0 0 0 0 27 0 0 0 0 28 1183
VecScale 25 1.0 4.0531e-05 1.0 6.54e+04 1.1 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 3030
VecCopy 2 1.0 6.9141e-06 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
VecSet 32 1.0 4.5538e-05 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
VecAXPY 2 1.0 2.1935e-05 1.4 1.05e+04 1.1 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 896
VecMAXPY 25 1.0 4.5085e-04 1.1 1.70e+06 1.1 0.0e+00 0.0e+00 0.0e+00 0 4 0 0 0 0 4 0 0 0 7061
VecAssemblyBegin 3 1.0 7.2956e-05 1.0 0.00e+00 0.0 2.0e+00 6.9e+03 9.0e+00 0 0 3 5 9 0 0 3 5 10 0
VecAssemblyEnd 3 1.0 2.1696e-05 1.2 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
VecScatterBegin 26 1.0 1.2279e-04 1.2 0.00e+00 0.0 5.2e+01 2.6e+03 0.0e+00 0 0 80 45 0 0 0 80 45 0 0
VecScatterEnd 26 1.0 1.4198e-02313.4 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
VecNormalize 25 1.0 2.7084e-04 1.1 1.96e+05 1.1 0.0e+00 0.0e+00 2.5e+01 0 0 0 0 26 0 0 0 0 27 1360
MatMult 25 1.0 1.5721e-02 9.0 3.08e+06 1.1 5.0e+01 2.6e+03 0.0e+00 1 7 77 42 0 1 7 77 42 0 366
MatSolve 26 1.0 8.3320e-03 1.2 2.03e+07 1.2 0.0e+00 0.0e+00 0.0e+00 0 42 0 0 0 0 42 0 0 0 4434
MatLUFactorNum 1 1.0 1.4023e-02 1.3 2.21e+07 1.4 0.0e+00 0.0e+00 0.0e+00 1 44 0 0 0 1 44 0 0 0 2710
MatILUFactorSym 1 1.0 5.3577e-02 1.3 0.00e+00 0.0 0.0e+00 0.0e+00 3.0e+00 3 0 0 0 3 3 0 0 0 3 0
MatAssemblyBegin 2 1.0 1.0502e-03 4.1 0.00e+00 0.0 3.0e+00 4.9e+04 4.0e+00 0 0 5 48 4 0 0 5 48 4 0
MatAssemblyEnd 2 1.0 9.5987e-04 1.1 0.00e+00 0.0 4.0e+00 6.5e+02 8.0e+00 0 0 6 1 8 0 0 6 1 9 0
MatGetRowIJ 1 1.0 1.9073e-06 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
MatGetOrdering 1 1.0 6.3896e-05 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 2.0e+00 0 0 0 0 2 0 0 0 0 2 0
MatZeroEntries 3 1.0 1.8978e-04 1.2 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
KSPGMRESOrthog 24 1.0 2.4314e-03 2.6 3.14e+06 1.1 0.0e+00 0.0e+00 2.4e+01 0 7 0 0 25 0 7 0 0 26 2425
KSPSetUp 2 1.0 8.0109e-05 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
KSPSolve 1 1.0 8.0268e-02 1.0 4.90e+07 1.3 5.0e+01 2.6e+03 5.7e+01 5100 77 42 60 5100 77 42 61 1087
PCSetUp 2 1.0 6.8043e-02 1.3 2.21e+07 1.4 0.0e+00 0.0e+00 7.0e+00 4 44 0 0 7 4 44 0 0 7 559
PCSetUpOnBlocks 1 1.0 6.7769e-02 1.3 2.21e+07 1.4 0.0e+00 0.0e+00 5.0e+00 4 44 0 0 5 4 44 0 0 5 561
PCApply 26 1.0 8.5933e-03 1.2 2.03e+07 1.2 0.0e+00 0.0e+00 0.0e+00 0 42 0 0 0 0 42 0 0 0 4300
------------------------------------------------------------------------------------------------------------------------
Memory usage is given in bytes:
Object Type Creations Destructions Memory Descendants' Mem.
Reports information only for process 0.
--- Event Stage 0: Main Stage
Vector 43 43 865112 0
Vector Scatter 2 2 2072 0
Index Set 7 7 18248 0
IS L to G Mapping 1 1 564 0
Matrix 4 4 5588632 0
Krylov Solver 2 2 19360 0
Preconditioner 2 2 1784 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 9.53674e-08
Average time for MPI_Barrier(): 1.00136e-06
Average time for zero size MPI_Send(): 1.34706e-05
#PETSc Option Table entries:
-ksp_right_pc
-log_summary
-pc_type bjacobi
-q 0
-sub_pc_factor_levels 4
-sub_pc_factor_zeropivot 0
-sub_pc_type ilu
#End of PETSc Option Table entries
Compiled without FORTRAN kernels
Compiled with full precision matrices (default)
sizeof(short) 2 sizeof(int) 4 sizeof(long) 8 sizeof(void*) 8 sizeof(PetscScalar) 8 sizeof(PetscInt) 4
Configure run at: Thu Nov 8 11:21:02 2012
Configure options: --with-debugging=false --COPTFLAGS=-O3 --CXXOPTFLAGS=-O3 --FOPTFLAGS=-O3 --with-clanguage=C++ --with-shared-libraries=1 --with-mpi-dir=/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1 --with-mumps=true --download-mumps=1 --with-metis=true --download-metis=1 --with-parmetis=true --download-parmetis=1 --with-superlu=true --download-superlu=1 --with-superludir=true --download-superlu_dist=1 --with-blacs=true --download-blacs=1 --with-scalapack=true --download-scalapack=1 --with-hypre=true --download-hypre=1 --with-blas-lib="[/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_intel_lp64.so,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_sequential.so,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_core.so]" --with-lapack-lib="[/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_lapack95_lp64.a]"
-----------------------------------------
Libraries compiled on Thu Nov 8 11:21:02 2012 on daedalus.ices.utexas.edu
Machine characteristics: Linux-2.6.32-279.1.1.el6.x86_64-x86_64-with-redhat-6.3-Carbon
Using PETSc directory: /opt/apps/ossw/libraries/petsc/petsc-3.3-p2
Using PETSc arch: intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt
-----------------------------------------
Using C compiler: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpicxx -wd1572 -O3 -fPIC ${COPTFLAGS} ${CFLAGS}
Using Fortran compiler: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpif90 -fPIC -O3 ${FOPTFLAGS} ${FFLAGS}
-----------------------------------------
Using include paths: -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/include -I/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/include
-----------------------------------------
Using C linker: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpicxx
Using Fortran linker: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpif90
Using libraries: -Wl,-rpath,/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -L/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -lpetsc -lX11 -Wl,-rpath,/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -L/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -lcmumps -ldmumps -lsmumps -lzmumps -lmumps_common -lpord -lHYPRE -lpthread -lsuperlu_dist_3.0 -lparmetis -lmetis -lscalapack -lblacs -lsuperlu_4.3 -Wl,-rpath,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64 -L/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64 -lmkl_lapack95_lp64 -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -Wl,-rpath,/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/lib -L/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/lib -Wl,-rpath,/opt/apps/sysnet/intel/12.1/composer_xe_2011_sp1.7.256/compiler/lib/intel64 -L/opt/apps/sysnet/intel/12.1/composer_xe_2011_sp1.7.256/compiler/lib/intel64 -Wl,-rpath,/usr/lib/gcc/x86_64-redhat-linux/4.4.6 -L/usr/lib/gcc/x86_64-redhat-linux/4.4.6 -lmpichf90 -lifport -lifcore -lm -lm -lmpichcxx -ldl -lmpich -lopa -lmpl -lrt -lpthread -limf -lsvml -lipgo -ldecimal -lcilkrts -lstdc++ -lgcc_s -lirc -lirc_s -ldl
-----------------------------------------
----------------------------------------------------------------------------------------------------------------------
| Processor id: 0 |
| Num Processors: 2 |
| Time: Tue Feb 5 13:38:06 2013 |
| OS: Linux |
| HostName: hbar.ices.utexas.edu |
| OS Release: 2.6.32-279.1.1.el6.x86_64 |
| OS Version: #1 SMP Tue Jul 10 11:24:23 CDT 2012 |
| Machine: x86_64 |
| Username: benkirk |
| Configuration: ./configure '--enable-everything' |
| '--prefix=/workspace/libmesh/install' |
| 'CXX=mpicxx' |
| 'CC=mpicc' |
| 'F77=mpif77' |
| 'FC=mpif90' |
| 'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2' |
| 'PETSC_ARCH=intel-12.1-mkl-intel-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/intel-12.1/mpich2-1.4.1p1/mkl-intel-10.3.12.361'|
| 'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/intel-12.1' |
----------------------------------------------------------------------------------------------------------------------
--------------------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=1.62357, Active time=1.57579 |
--------------------------------------------------------------------------------------------------------------------
| 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.0934 0.093403 0.1240 0.123956 5.93 7.87 |
| build_constraint_matrix_and_vector() 2048 0.0185 0.000009 0.0185 0.000009 1.18 1.18 |
| build_sparsity() 1 0.0858 0.085755 0.1986 0.198610 5.44 12.60 |
| create_dof_constraints() 1 0.0581 0.058149 0.2780 0.278022 3.69 17.64 |
| distribute_dofs() 1 0.0581 0.058134 0.1566 0.156563 3.69 9.94 |
| dof_indices() 10816 0.5785 0.000053 0.5785 0.000053 36.71 36.71 |
| hetero_cnstrn_elem_mat_vec() 2048 0.0191 0.000009 0.0191 0.000009 1.21 1.21 |
| prepare_send_list() 1 0.0005 0.000495 0.0005 0.000495 0.03 0.03 |
| reinit() 1 0.0973 0.097344 0.0973 0.097344 6.18 6.18 |
| |
| EquationSystems |
| build_solution_vector() 1 0.0088 0.008831 0.1172 0.117236 0.56 7.44 |
| |
| ExodusII_IO |
| write_nodal_data() 1 0.0207 0.020730 0.0207 0.020730 1.32 1.32 |
| |
| FE |
| compute_shape_functions() 2048 0.0501 0.000024 0.0501 0.000024 3.18 3.18 |
| init_shape_functions() 1 0.0002 0.000152 0.0002 0.000152 0.01 0.01 |
| |
| FEMap |
| compute_affine_map() 2048 0.0233 0.000011 0.0233 0.000011 1.48 1.48 |
| init_reference_to_physical_map() 1 0.0001 0.000126 0.0001 0.000126 0.01 0.01 |
| |
| Mesh |
| find_neighbors() 1 0.0780 0.078040 0.0781 0.078150 4.95 4.96 |
| renumber_nodes_and_elem() 2 0.0031 0.001546 0.0031 0.001546 0.20 0.20 |
| |
| MeshCommunication |
| compute_hilbert_indices() 2 0.0347 0.017340 0.0347 0.017340 2.20 2.20 |
| find_global_indices() 2 0.0127 0.006369 0.0511 0.025549 0.81 3.24 |
| parallel_sort() 2 0.0032 0.001580 0.0032 0.001619 0.20 0.21 |
| |
| MeshOutput |
| write_equation_systems() 1 0.0001 0.000090 0.1381 0.138107 0.01 8.76 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0123 0.012261 0.0123 0.012261 0.78 0.78 |
| |
| MetisPartitioner |
| partition() 1 0.1333 0.133330 0.1583 0.158297 8.46 10.05 |
| |
| Parallel |
| allgather() 9 0.0007 0.000083 0.0008 0.000085 0.05 0.05 |
| max(bool) 1 0.0000 0.000003 0.0000 0.000003 0.00 0.00 |
| max(scalar) 105 0.0004 0.000003 0.0004 0.000003 0.02 0.02 |
| max(vector) 24 0.0002 0.000006 0.0004 0.000015 0.01 0.02 |
| min(bool) 121 0.0003 0.000003 0.0003 0.000003 0.02 0.02 |
| min(scalar) 99 0.0040 0.000040 0.0040 0.000040 0.25 0.25 |
| min(vector) 24 0.0002 0.000009 0.0004 0.000019 0.01 0.03 |
| probe() 12 0.0007 0.000057 0.0007 0.000057 0.04 0.04 |
| receive() 12 0.0001 0.000012 0.0008 0.000069 0.01 0.05 |
| send() 12 0.0001 0.000005 0.0001 0.000005 0.00 0.00 |
| send_receive() 16 0.0005 0.000028 0.0014 0.000086 0.03 0.09 |
| sum() 20 0.0002 0.000010 0.0003 0.000016 0.01 0.02 |
| |
| Parallel::Request |
| wait() 12 0.0000 0.000003 0.0000 0.000003 0.00 0.00 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0044 0.004352 0.0050 0.004987 0.28 0.32 |
| set_parent_processor_ids() 1 0.0029 0.002885 0.0029 0.002885 0.18 0.18 |
| |
| PetscLinearSolver |
| solve() 1 0.0830 0.083048 0.0830 0.083048 5.27 5.27 |
| |
| System |
| assemble() 1 0.0882 0.088177 0.3123 0.312267 5.60 19.82 |
--------------------------------------------------------------------------------------------------------------------
| Totals: 19501 1.5758 100.00 |
--------------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example introduction_ex5:
* mpirun -np 2 example-devel -q 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with this library; if not, write to the Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
C++ Includes