The source file adaptivity_ex1.C with comments:
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/edge_edge3.h"
#include "libmesh/gnuplot_io.h"
#include "libmesh/equation_systems.h"
#include "libmesh/linear_implicit_system.h"
#include "libmesh/fe.h"
#include "libmesh/getpot.h"
#include "libmesh/quadrature_gauss.h"
#include "libmesh/sparse_matrix.h"
#include "libmesh/dof_map.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dense_vector.h"
#include "libmesh/error_vector.h"
#include "libmesh/kelly_error_estimator.h"
#include "libmesh/mesh_refinement.h"
Bring in everything from the libMesh namespace
using namespace libMesh;
void assemble_1D(EquationSystems& es, const std::string& system_name);
int main(int argc, char** argv)
{
Initialize the library. This is necessary because the library
may depend on a number of other libraries (i.e. MPI and PETSc)
that require initialization before use. When the LibMeshInit
object goes out of scope, other libraries and resources are
finalized.
LibMeshInit init (argc, argv);
Skip adaptive examples on a non-adaptive libMesh build
#ifndef LIBMESH_ENABLE_AMR
libmesh_example_assert(false, "--enable-amr");
#else
Create a new mesh
Mesh mesh;
GetPot command_line (argc, argv);
int n = 4;
if ( command_line.search(1, "-n") )
n = command_line.next(n);
Build a 1D mesh with 4 elements from x=0 to x=1, using
EDGE3 (i.e. quadratic) 1D elements. They are called EDGE3 elements
because a quadratic element contains 3 nodes.
MeshTools::Generation::build_line(mesh,n,0.,1.,EDGE3);
Define the equation systems object and the system we are going
to solve. See Introduction Example 2 for more details.
EquationSystems equation_systems(mesh);
LinearImplicitSystem& system = equation_systems.add_system
<LinearImplicitSystem>("1D");
Add a variable "u" to the system, using second-order approximation
system.add_variable("u",SECOND);
Give the system a pointer to the matrix assembly function. This
will be called when needed by the library.
system.attach_assemble_function(assemble_1D);
Define the mesh refinement object that takes care of adaptively
refining the mesh.
MeshRefinement mesh_refinement(mesh);
These parameters determine the proportion of elements that will
be refined and coarsened. Any element within 30% of the maximum
error on any element will be refined, and any element within 30%
of the minimum error on any element might be coarsened
mesh_refinement.refine_fraction() = 0.7;
mesh_refinement.coarsen_fraction() = 0.3;
We won't refine any element more than 5 times in total
mesh_refinement.max_h_level() = 5;
Initialize the data structures for the equation system.
equation_systems.init();
Refinement parameters
const unsigned int max_r_steps = 5; // Refine the mesh 5 times
Define the refinement loop
for(unsigned int r_step=0; r_step<=max_r_steps; r_step++)
{
Solve the equation system
equation_systems.get_system("1D").solve();
We need to ensure that the mesh is not refined on the last iteration
of this loop, since we do not want to refine the mesh unless we are
going to solve the equation system for that refined mesh.
if(r_step != max_r_steps)
{
Error estimation objects, see Adaptivity Example 2 for details
ErrorVector error;
KellyErrorEstimator error_estimator;
Compute the error for each active element
error_estimator.estimate_error(system, error);
Flag elements to be refined and coarsened
mesh_refinement.flag_elements_by_error_fraction (error);
Perform refinement and coarsening
mesh_refinement.refine_and_coarsen_elements();
Reinitialize the equation_systems object for the newly refined
mesh. One of the steps in this is project the solution onto the
new mesh
equation_systems.reinit();
}
}
Construct gnuplot plotting object, pass in mesh, title of plot
and boolean to indicate use of grid in plot. The grid is used to
show the edges of each element in the mesh.
GnuPlotIO plot(mesh,"Adaptivity Example 1", GnuPlotIO::GRID_ON);
Write out script to be called from within gnuplot:
Load gnuplot, then type "call 'gnuplot_script'" from gnuplot prompt
plot.write_equation_systems("gnuplot_script",equation_systems);
#endif // #ifndef LIBMESH_ENABLE_AMR
All done. libMesh objects are destroyed here. Because the
LibMeshInit object was created first, its destruction occurs
last, and it's destructor finalizes any external libraries and
checks for leaked memory.
return 0;
}
Define the matrix assembly function for the 1D PDE we are solving
void assemble_1D(EquationSystems& es, const std::string& system_name)
{
#ifdef LIBMESH_ENABLE_AMR
It is a good idea to check we are solving the correct system
libmesh_assert_equal_to (system_name, "1D");
Get a reference to the mesh object
const MeshBase& mesh = es.get_mesh();
The dimension we are using, i.e. dim==1
const unsigned int dim = mesh.mesh_dimension();
Get a reference to the system we are solving
LinearImplicitSystem& system = es.get_system<LinearImplicitSystem>("1D");
Get a reference to the DofMap object for this system. The DofMap object
handles the index translation from node and element numbers to degree of
freedom numbers. DofMap's are discussed in more detail in future examples.
const DofMap& dof_map = system.get_dof_map();
Get a constant reference to the Finite Element type for the first
(and only) variable in the system.
FEType fe_type = dof_map.variable_type(0);
Build a finite element object of the specified type. The build
function dynamically allocates memory so we use an AutoPtr in this case.
An AutoPtr is a pointer that cleans up after itself. See examples 3 and 4
for more details on AutoPtr.
AutoPtr<FEBase> fe(FEBase::build(dim, fe_type));
Tell the finite element object to use fifth order Gaussian quadrature
QGauss qrule(dim,FIFTH);
fe->attach_quadrature_rule(&qrule);
Here we define some references to cell-specific data that will be used to
assemble the linear system.
The element Jacobian * quadrature weight at each integration point.
The element Jacobian * quadrature weight at each integration point.
const std::vector<Real>& JxW = fe->get_JxW();
The element shape functions evaluated at the quadrature points.
const std::vector<std::vector<Real> >& phi = fe->get_phi();
The element shape function gradients evaluated at the quadrature points.
const std::vector<std::vector<RealGradient> >& dphi = fe->get_dphi();
Declare a dense matrix and dense vector to hold the element matrix
and right-hand-side contribution
DenseMatrix<Number> Ke;
DenseVector<Number> Fe;
This vector will hold the degree of freedom indices for the element.
These define where in the global system the element degrees of freedom
get mapped.
std::vector<dof_id_type> dof_indices;
We now loop over all the active elements in the mesh in order to calculate
the matrix and right-hand-side contribution from each element. Use a
const_element_iterator to loop over the elements. We make
el_end const as it is used only for the stopping condition of the loop.
MeshBase::const_element_iterator el = mesh.active_local_elements_begin();
const MeshBase::const_element_iterator el_end = mesh.active_local_elements_end();
Note that ++el is preferred to el++ when using loops with iterators
for( ; el != el_end; ++el)
{
It is convenient to store a pointer to the current element
const Elem* elem = *el;
Get the degree of freedom indices for the current element.
These define where in the global matrix and right-hand-side this
element will contribute to.
dof_map.dof_indices(elem, dof_indices);
Compute the element-specific data for the current element. This
involves computing the location of the quadrature points (q_point)
and the shape functions (phi, dphi) for the current element.
fe->reinit(elem);
Store the number of local degrees of freedom contained in this element
const int n_dofs = dof_indices.size();
We resize and zero out Ke and Fe (resize() also clears the matrix and
vector). In this example, all elements in the mesh are EDGE3's, so
Ke will always be 3x3, and Fe will always be 3x1. If the mesh contained
different element types, then the size of Ke and Fe would change.
Ke.resize(n_dofs, n_dofs);
Fe.resize(n_dofs);
Now loop over quadrature points to handle numerical integration
for(unsigned int qp=0; qp<qrule.n_points(); qp++)
{
Now build the element matrix and right-hand-side using loops to
integrate the test functions (i) against the trial functions (j).
for(unsigned int i=0; i<phi.size(); i++)
{
Fe(i) += JxW[qp]*phi[i][qp];
for(unsigned int j=0; j<phi.size(); j++)
{
Ke(i,j) += JxW[qp]*(1.e-3*dphi[i][qp]*dphi[j][qp] +
phi[i][qp]*phi[j][qp]);
}
}
}
At this point we have completed the matrix and RHS summation. The
final step is to apply boundary conditions, which in this case are
simple Dirichlet conditions with u(0) = u(1) = 0.
Define the penalty parameter used to enforce the BC's
Define the penalty parameter used to enforce the BC's
double penalty = 1.e10;
Loop over the sides of this element. For a 1D element, the "sides"
are defined as the nodes on each edge of the element, i.e. 1D elements
have 2 sides.
for(unsigned int s=0; s<elem->n_sides(); s++)
{
If this element has a NULL neighbor, then it is on the edge of the
mesh and we need to enforce a boundary condition using the penalty
method.
if(elem->neighbor(s) == NULL)
{
Ke(s,s) += penalty;
Fe(s) += 0*penalty;
}
}
This is a function call that is necessary when using adaptive
mesh refinement. See Adaptivity Example 2 for more details.
dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
Add Ke and Fe to the global matrix and right-hand-side.
system.matrix->add_matrix(Ke, dof_indices);
system.rhs->add_vector(Fe, dof_indices);
}
#endif // #ifdef LIBMESH_ENABLE_AMR
}
The source file adaptivity_ex1.C without comments:
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/edge_edge3.h"
#include "libmesh/gnuplot_io.h"
#include "libmesh/equation_systems.h"
#include "libmesh/linear_implicit_system.h"
#include "libmesh/fe.h"
#include "libmesh/getpot.h"
#include "libmesh/quadrature_gauss.h"
#include "libmesh/sparse_matrix.h"
#include "libmesh/dof_map.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dense_vector.h"
#include "libmesh/error_vector.h"
#include "libmesh/kelly_error_estimator.h"
#include "libmesh/mesh_refinement.h"
using namespace libMesh;
void assemble_1D(EquationSystems& es, const std::string& system_name);
int main(int argc, char** argv)
{
LibMeshInit init (argc, argv);
#ifndef LIBMESH_ENABLE_AMR
libmesh_example_assert(false, "--enable-amr");
#else
Mesh mesh;
GetPot command_line (argc, argv);
int n = 4;
if ( command_line.search(1, "-n") )
n = command_line.next(n);
MeshTools::Generation::build_line(mesh,n,0.,1.,EDGE3);
EquationSystems equation_systems(mesh);
LinearImplicitSystem& system = equation_systems.add_system
<LinearImplicitSystem>("1D");
system.add_variable("u",SECOND);
system.attach_assemble_function(assemble_1D);
MeshRefinement mesh_refinement(mesh);
mesh_refinement.refine_fraction() = 0.7;
mesh_refinement.coarsen_fraction() = 0.3;
mesh_refinement.max_h_level() = 5;
equation_systems.init();
const unsigned int max_r_steps = 5; // Refine the mesh 5 times
for(unsigned int r_step=0; r_step<=max_r_steps; r_step++)
{
equation_systems.get_system("1D").solve();
if(r_step != max_r_steps)
{
ErrorVector error;
KellyErrorEstimator error_estimator;
error_estimator.estimate_error(system, error);
mesh_refinement.flag_elements_by_error_fraction (error);
mesh_refinement.refine_and_coarsen_elements();
equation_systems.reinit();
}
}
GnuPlotIO plot(mesh,"Adaptivity Example 1", GnuPlotIO::GRID_ON);
plot.write_equation_systems("gnuplot_script",equation_systems);
#endif // #ifndef LIBMESH_ENABLE_AMR
return 0;
}
void assemble_1D(EquationSystems& es, const std::string& system_name)
{
#ifdef LIBMESH_ENABLE_AMR
libmesh_assert_equal_to (system_name, "1D");
const MeshBase& mesh = es.get_mesh();
const unsigned int dim = mesh.mesh_dimension();
LinearImplicitSystem& system = es.get_system<LinearImplicitSystem>("1D");
const DofMap& dof_map = system.get_dof_map();
FEType fe_type = dof_map.variable_type(0);
AutoPtr<FEBase> fe(FEBase::build(dim, fe_type));
QGauss qrule(dim,FIFTH);
fe->attach_quadrature_rule(&qrule);
const std::vector<Real>& JxW = fe->get_JxW();
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 el_end = mesh.active_local_elements_end();
for( ; el != el_end; ++el)
{
const Elem* elem = *el;
dof_map.dof_indices(elem, dof_indices);
fe->reinit(elem);
const int n_dofs = dof_indices.size();
Ke.resize(n_dofs, n_dofs);
Fe.resize(n_dofs);
for(unsigned int qp=0; qp<qrule.n_points(); qp++)
{
for(unsigned int i=0; i<phi.size(); i++)
{
Fe(i) += JxW[qp]*phi[i][qp];
for(unsigned int j=0; j<phi.size(); j++)
{
Ke(i,j) += JxW[qp]*(1.e-3*dphi[i][qp]*dphi[j][qp] +
phi[i][qp]*phi[j][qp]);
}
}
}
double penalty = 1.e10;
for(unsigned int s=0; s<elem->n_sides(); s++)
{
if(elem->neighbor(s) == NULL)
{
Ke(s,s) += penalty;
Fe(s) += 0*penalty;
}
}
dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
system.matrix->add_matrix(Ke, dof_indices);
system.rhs->add_vector(Fe, dof_indices);
}
#endif // #ifdef LIBMESH_ENABLE_AMR
}
The console output of the program:
***************************************************************
* Running Example adaptivity_ex1:
* mpirun -np 2 example-devel -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/adaptivity/adaptivity_ex1/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:38:07 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.116e-01 1.00001 1.116e-01
Objects: 2.850e+02 1.00000 2.850e+02
Flops: 1.613e+04 1.06717 1.562e+04 3.124e+04
Flops/sec: 1.445e+05 1.06716 1.400e+05 2.799e+05
MPI Messages: 9.300e+01 1.00000 9.300e+01 1.860e+02
MPI Message Lengths: 1.280e+03 1.00000 1.376e+01 2.560e+03
MPI Reductions: 4.300e+02 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.1157e-01 100.0% 3.1239e+04 100.0% 1.860e+02 100.0% 1.376e+01 100.0% 4.290e+02 99.8%
------------------------------------------------------------------------------------------------------------------------
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 18 1.0 3.1853e-04 6.5 1.28e+03 1.1 0.0e+00 0.0e+00 1.8e+01 0 8 0 0 4 0 8 0 0 4 8
VecNorm 30 1.0 5.0516e-0312.9 1.10e+03 1.1 0.0e+00 0.0e+00 3.0e+01 2 7 0 0 7 2 7 0 0 7 0
VecScale 24 1.0 1.5736e-05 1.0 4.40e+02 1.1 0.0e+00 0.0e+00 0.0e+00 0 3 0 0 0 0 3 0 0 0 54
VecCopy 27 1.0 1.3113e-05 1.3 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 86 1.0 2.3603e-05 1.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
VecAXPY 12 1.0 1.0234e-02 1.8 4.40e+02 1.1 0.0e+00 0.0e+00 0.0e+00 7 3 0 0 0 7 3 0 0 0 0
VecMAXPY 24 1.0 9.2983e-06 1.2 1.98e+03 1.1 0.0e+00 0.0e+00 0.0e+00 0 12 0 0 0 0 12 0 0 0 414
VecAssemblyBegin 53 1.0 1.8680e-03 1.2 0.00e+00 0.0 2.2e+01 6.0e+00 1.3e+02 2 0 12 5 30 2 0 12 5 30 0
VecAssemblyEnd 53 1.0 4.0770e-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
VecScatterBegin 45 1.0 9.8944e-05 1.1 0.00e+00 0.0 8.6e+01 1.9e+01 0.0e+00 0 0 46 63 0 0 0 46 63 0 0
VecScatterEnd 45 1.0 1.6713e-03 3.3 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 1 0 0 0 0 1 0 0 0 0 0
VecNormalize 24 1.0 5.0077e-0316.3 1.32e+03 1.1 0.0e+00 0.0e+00 2.4e+01 2 8 0 0 6 2 8 0 0 6 1
MatMult 24 1.0 6.6614e-04 5.1 3.03e+03 1.1 4.8e+01 1.2e+01 0.0e+00 0 19 26 22 0 0 19 26 22 0 9
MatSolve 30 1.0 2.2888e-05 1.1 5.93e+03 1.1 0.0e+00 0.0e+00 0.0e+00 0 37 0 0 0 0 37 0 0 0 502
MatLUFactorNum 6 1.0 4.6015e-05 1.0 1.92e+03 1.1 0.0e+00 0.0e+00 0.0e+00 0 12 0 0 0 0 12 0 0 0 80
MatILUFactorSym 6 1.0 1.1897e-04 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 1.8e+01 0 0 0 0 4 0 0 0 0 4 0
MatAssemblyBegin 12 1.0 8.0061e-04 2.5 0.00e+00 0.0 1.8e+01 1.7e+01 2.4e+01 1 0 10 12 6 1 0 10 12 6 0
MatAssemblyEnd 12 1.0 7.4410e-04 1.1 0.00e+00 0.0 2.4e+01 5.0e+00 4.8e+01 1 0 13 5 11 1 0 13 5 11 0
MatGetRowIJ 6 1.0 4.0531e-06 1.9 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 6 1.0 1.3995e-04 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 1.4e+01 0 0 0 0 3 0 0 0 0 3 0
MatZeroEntries 18 1.0 1.5974e-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
KSPGMRESOrthog 18 1.0 3.5739e-04 4.2 2.60e+03 1.1 0.0e+00 0.0e+00 1.8e+01 0 16 0 0 4 0 16 0 0 4 14
KSPSetUp 12 1.0 1.1420e-04 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 6 1.0 1.3312e-02 1.0 1.61e+04 1.1 4.8e+01 1.2e+01 9.2e+01 12100 26 22 21 12100 26 22 21 2
PCSetUp 12 1.0 1.1818e-03 1.0 1.92e+03 1.1 0.0e+00 0.0e+00 4.4e+01 1 12 0 0 10 1 12 0 0 10 3
PCSetUpOnBlocks 6 1.0 6.2490e-04 1.0 1.92e+03 1.1 0.0e+00 0.0e+00 3.2e+01 1 12 0 0 7 1 12 0 0 7 6
PCApply 30 1.0 2.2101e-04 1.0 5.93e+03 1.1 0.0e+00 0.0e+00 0.0e+00 0 37 0 0 0 0 37 0 0 0 52
------------------------------------------------------------------------------------------------------------------------
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 164 164 264128 0
Vector Scatter 17 17 17612 0
Index Set 49 49 36932 0
IS L to G Mapping 6 6 3384 0
Matrix 24 24 77416 0
Krylov Solver 12 12 116160 0
Preconditioner 12 12 10704 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 0
Average time for MPI_Barrier(): 2.00272e-06
Average time for zero size MPI_Send(): 1.50204e-05
#PETSc Option Table entries:
-ksp_right_pc
-log_summary
-pc_type bjacobi
-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:07 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=0.122248, Active time=0.096812 |
----------------------------------------------------------------------------------------------------------------
| 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() 6 0.0013 0.000219 0.0016 0.000271 1.36 1.68 |
| build_sparsity() 6 0.0013 0.000225 0.0041 0.000676 1.39 4.19 |
| create_dof_constraints() 6 0.0000 0.000004 0.0000 0.000004 0.02 0.02 |
| distribute_dofs() 6 0.0028 0.000474 0.0091 0.001509 2.94 9.36 |
| dof_indices() 249 0.0058 0.000023 0.0058 0.000023 6.00 6.00 |
| old_dof_indices() 100 0.0025 0.000025 0.0025 0.000025 2.60 2.60 |
| prepare_send_list() 6 0.0000 0.000006 0.0000 0.000006 0.04 0.04 |
| reinit() 6 0.0033 0.000547 0.0033 0.000547 3.39 3.39 |
| |
| EquationSystems |
| build_solution_vector() 1 0.0002 0.000248 0.0009 0.000858 0.26 0.89 |
| |
| FE |
| compute_shape_functions() 114 0.0004 0.000003 0.0004 0.000003 0.37 0.37 |
| init_shape_functions() 68 0.0005 0.000008 0.0005 0.000008 0.54 0.54 |
| inverse_map() 108 0.0004 0.000003 0.0004 0.000003 0.36 0.36 |
| |
| FEMap |
| compute_affine_map() 114 0.0005 0.000004 0.0005 0.000004 0.48 0.48 |
| compute_face_map() 31 0.0001 0.000003 0.0001 0.000003 0.09 0.09 |
| init_face_shape_functions() 5 0.0000 0.000005 0.0000 0.000005 0.02 0.02 |
| init_reference_to_physical_map() 68 0.0003 0.000005 0.0003 0.000005 0.33 0.33 |
| |
| GnuPlotIO |
| write_nodal_data() 1 0.0011 0.001065 0.0011 0.001065 1.10 1.10 |
| |
| JumpErrorEstimator |
| estimate_error() 5 0.0017 0.000336 0.0057 0.001140 1.74 5.89 |
| |
| LocationMap |
| find() 76 0.0004 0.000005 0.0004 0.000005 0.42 0.42 |
| init() 10 0.0005 0.000050 0.0005 0.000050 0.52 0.52 |
| |
| Mesh |
| contract() 5 0.0001 0.000015 0.0001 0.000026 0.08 0.13 |
| find_neighbors() 6 0.0024 0.000398 0.0033 0.000556 2.47 3.44 |
| renumber_nodes_and_elem() 17 0.0002 0.000014 0.0002 0.000014 0.25 0.25 |
| |
| MeshCommunication |
| compute_hilbert_indices() 7 0.0007 0.000095 0.0007 0.000095 0.69 0.69 |
| find_global_indices() 7 0.0008 0.000114 0.0041 0.000589 0.83 4.26 |
| parallel_sort() 7 0.0011 0.000151 0.0016 0.000233 1.09 1.69 |
| |
| MeshOutput |
| write_equation_systems() 1 0.0001 0.000109 0.0021 0.002084 0.11 2.15 |
| |
| MeshRefinement |
| _coarsen_elements() 10 0.0001 0.000012 0.0003 0.000027 0.12 0.28 |
| _refine_elements() 10 0.0012 0.000120 0.0029 0.000285 1.24 2.95 |
| add_point() 76 0.0004 0.000005 0.0008 0.000011 0.41 0.85 |
| make_coarsening_compatible() 21 0.0019 0.000091 0.0019 0.000091 1.98 1.98 |
| make_refinement_compatible() 21 0.0002 0.000011 0.0005 0.000023 0.24 0.50 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0003 0.000294 0.0003 0.000294 0.30 0.30 |
| |
| MetisPartitioner |
| partition() 6 0.0021 0.000351 0.0065 0.001076 2.17 6.67 |
| |
| Parallel |
| allgather() 32 0.0009 0.000029 0.0011 0.000033 0.97 1.09 |
| max(bool) 52 0.0016 0.000031 0.0016 0.000031 1.65 1.65 |
| max(scalar) 921 0.0086 0.000009 0.0086 0.000009 8.91 8.91 |
| max(vector) 224 0.0030 0.000013 0.0093 0.000041 3.08 9.58 |
| min(bool) 1163 0.0108 0.000009 0.0108 0.000009 11.18 11.18 |
| min(scalar) 909 0.0099 0.000011 0.0099 0.000011 10.22 10.22 |
| min(vector) 224 0.0032 0.000014 0.0098 0.000044 3.26 10.11 |
| probe() 62 0.0010 0.000016 0.0010 0.000016 1.03 1.03 |
| receive() 62 0.0003 0.000005 0.0013 0.000021 0.33 1.37 |
| send() 62 0.0002 0.000003 0.0002 0.000003 0.18 0.18 |
| send_receive() 76 0.0006 0.000007 0.0022 0.000029 0.58 2.26 |
| sum() 34 0.0005 0.000015 0.0010 0.000031 0.53 1.08 |
| |
| Parallel::Request |
| wait() 62 0.0001 0.000002 0.0001 0.000002 0.11 0.11 |
| |
| Partitioner |
| set_node_processor_ids() 6 0.0009 0.000146 0.0019 0.000325 0.90 2.01 |
| set_parent_processor_ids() 6 0.0002 0.000037 0.0002 0.000037 0.23 0.23 |
| |
| PetscLinearSolver |
| solve() 6 0.0161 0.002681 0.0161 0.002681 16.62 16.62 |
| |
| ProjectVector |
| operator() 5 0.0009 0.000176 0.0038 0.000769 0.91 3.97 |
| |
| System |
| assemble() 6 0.0009 0.000151 0.0029 0.000479 0.93 2.97 |
| project_vector() 5 0.0024 0.000473 0.0091 0.001816 2.44 9.38 |
----------------------------------------------------------------------------------------------------------------
| Totals: 5098 0.0968 100.00 |
----------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example adaptivity_ex1:
* mpirun -np 2 example-devel -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
Adaptivity Example 1 - Solving 1D PDE Using Adaptive Mesh Refinement
This example demonstrates how to solve a simple 1D problem using adaptive mesh refinement. The PDE that is solved is: -epsilon*u''(x) + u(x) = 1, on the domain [0,1] with boundary conditions u(0) = u(1) = 0 and where epsilon << 1.
The approach used to solve 1D problems in libMesh is virtually identical to solving 2D or 3D problems, so in this sense this example represents a good starting point for new users. Note that many concepts are used in this example which are explained more fully in subsequent examples.
Libmesh includes