The source file assembly.h with comments:
#ifndef __assembly_h__
#define __assembly_h__
#include "libmesh/sparse_matrix.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dense_vector.h"
#include "libmesh/fe.h"
#include "libmesh/fe_interface.h"
#include "libmesh/fe_base.h"
#include "libmesh/elem_assembly.h"
#include "libmesh/quadrature_gauss.h"
rbOOmit includes
#include "libmesh/rb_theta.h"
#include "libmesh/rb_assembly_expansion.h"
Bring in bits from the libMesh namespace.
Just the bits we're using, since this is a header.
using libMesh::ElemAssembly;
using libMesh::FEInterface;
using libMesh::FEMContext;
using libMesh::Number;
using libMesh::Point;
using libMesh::RBAssemblyExpansion;
using libMesh::RBParameters;
using libMesh::RBTheta;
using libMesh::RBThetaExpansion;
using libMesh::Real;
using libMesh::RealGradient;
Functors for the parameter-dependent part of the affine decomposition of the PDE
The RHS and outputs just require a constant value of 1, so use a default RBTheta object there
struct ThetaA0 : RBTheta { virtual Number evaluate(const RBParameters& ) { return 0.05; } };
struct ThetaA1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("x_vel"); } };
struct ThetaA2 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("y_vel"); } };
struct A0 : ElemAssembly
{
Assemble the Laplacian operator
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
The velocity shape function gradients at interior
quadrature points.
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
The number of local degrees of freedom in each variable
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
Now we will build the affine operator
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] * dphi[j][qp]*dphi[i][qp];
}
};
struct A1 : ElemAssembly
{
Convection in the x-direction
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
The number of local degrees of freedom in each variable
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
Now we will build the affine operator
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] *dphi[j][qp](0)*phi[i][qp];
}
};
struct A2 : ElemAssembly
{
Convection in the y-direction
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
The number of local degrees of freedom in each variable
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
Now we will build the affine operator
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] *dphi[j][qp](1)*phi[i][qp];
}
};
struct F0 : ElemAssembly
{
Source term, 1 throughout the domain
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
The number of local degrees of freedom in each variable
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
Now we will build the affine operator
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*phi[i][qp] );
}
};
struct OutputAssembly : ElemAssembly
{
OutputAssembly(Real min_x_in, Real max_x_in,
Real min_y_in, Real max_y_in)
:
min_x(min_x_in),
max_x(max_x_in),
min_y(min_y_in),
max_y(max_y_in)
{}
Output: Average value over the region [min_x,max_x]x[min_y,max_y]
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
The number of local degrees of freedom in each variable
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
Now we will build the affine operator
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
Real output_area = (max_x-min_x) * (max_y-min_y);
Point centroid = c.elem->centroid();
if( (min_x <= centroid(0)) && (centroid(0) <= max_x) &&
(min_y <= centroid(1)) && (centroid(1) <= max_y) )
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*phi[i][qp] ) / output_area;
}
Member variables that define the output region in 2D
Real min_x, max_x, min_y, max_y;
};
Define an RBThetaExpansion class for this PDE
struct CDRBThetaExpansion : RBThetaExpansion
{
/**
* Constructor.
*/
CDRBThetaExpansion()
{
set up the RBThetaExpansion object
attach_A_theta(&theta_a_0); // Attach the lhs theta
attach_A_theta(&theta_a_1);
attach_A_theta(&theta_a_2);
attach_F_theta(&rb_theta); // Attach the rhs theta
attach_output_theta(&rb_theta); // Attach output 0 theta
attach_output_theta(&rb_theta); // Attach output 1 theta
attach_output_theta(&rb_theta); // Attach output 2 theta
attach_output_theta(&rb_theta); // Attach output 3 theta
}
The RBTheta member variables
ThetaA0 theta_a_0;
ThetaA1 theta_a_1;
ThetaA2 theta_a_2;
RBTheta rb_theta; // Default RBTheta object, just returns 1.
};
Define an RBAssemblyExpansion class for this PDE
struct CDRBAssemblyExpansion : RBAssemblyExpansion
{
/**
* Constructor.
*/
CDRBAssemblyExpansion()
:
L0(0.7,0.8,0.7,0.8),
L1(0.2,0.3,0.7,0.8),
L2(0.2,0.3,0.2,0.3),
L3(0.7,0.8,0.2,0.3)
{
And set up the RBAssemblyExpansion object
attach_A_assembly(&A0_assembly); // Attach the lhs assembly
attach_A_assembly(&A1_assembly);
attach_A_assembly(&A2_assembly);
attach_F_assembly(&F0_assembly); // Attach the rhs assembly
attach_output_assembly(&L0); // Attach output 0 assembly
attach_output_assembly(&L1); // Attach output 1 assembly
attach_output_assembly(&L2); // Attach output 2 assembly
attach_output_assembly(&L3); // Attach output 3 assembly
}
The ElemAssembly objects
A0 A0_assembly;
A1 A1_assembly;
A2 A2_assembly;
F0 F0_assembly;
OutputAssembly L0;
OutputAssembly L1;
OutputAssembly L2;
OutputAssembly L3;
};
#endif
The source file rb_classes.h with comments:
rbOOmit 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
#ifndef __rb_classes_h__
#define __rb_classes_h__
#include "libmesh/rb_construction.h"
#include "libmesh/fe_base.h"
local include
#include "assembly.h"
Bring in bits from the libMesh namespace.
Just the bits we're using, since this is a header.
using libMesh::AutoPtr;
using libMesh::DirichletBoundary;
using libMesh::EquationSystems;
using libMesh::FEMContext;
using libMesh::RBConstruction;
using libMesh::RBEvaluation;
using libMesh::Real;
A simple subclass of RBEvaluation, which just needs to specify
(a lower bound for) the coercivity constant for this problem.
For this simple convection-diffusion problem, we can set the
coercivity constant lower bound to 0.05.
class SimpleRBEvaluation : public RBEvaluation
{
public:
/**
* Constructor. Just set the theta expansion.
*/
SimpleRBEvaluation()
{
set_rb_theta_expansion(cd_rb_theta_expansion);
}
/**
* The coercivity constant is bounded below by 0.05.
*/
virtual Real get_stability_lower_bound() { return 0.05; }
/**
* The object that stores the "theta" expansion of the parameter dependent PDE,
* i.e. the set of parameter-dependent functions in the affine expansion of the PDE.
*/
CDRBThetaExpansion cd_rb_theta_expansion;
};
A simple subclass of Construction, which just needs to override build_rb_evaluation
in order to build a SimpleRBEvaluation object, rather than an RBEvaluation object.
class SimpleRBConstruction : public RBConstruction
{
public:
SimpleRBConstruction (EquationSystems& es,
const std::string& name_in,
const unsigned int number_in)
: Parent(es, name_in, number_in),
dirichlet_bc(AutoPtr<DirichletBoundary>(NULL))
{}
/**
* Destructor.
*/
virtual ~SimpleRBConstruction () { }
/**
* The type of system.
*/
typedef SimpleRBConstruction sys_type;
/**
* The type of the parent.
*/
typedef RBConstruction Parent;
/**
* Initialize data structures.
*/
virtual void init_data()
{
u_var = this->add_variable ("u", FIRST);
Generate a DirichletBoundary object
dirichlet_bc = build_zero_dirichlet_boundary_object();
Set the Dirichet boundary IDs
and the Dirichlet boundary variable numbers
dirichlet_bc->b.insert(0);
dirichlet_bc->b.insert(1);
dirichlet_bc->b.insert(2);
dirichlet_bc->b.insert(3);
dirichlet_bc->variables.push_back(u_var);
Attach dirichlet_bc (must do this _before_ Parent::init_data)
get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
Parent::init_data();
Set the rb_assembly_expansion for this Construction object.
The theta expansion comes from the RBEvaluation object.
set_rb_assembly_expansion(cd_rb_assembly_expansion);
We need to define an inner product matrix for this problem
set_inner_product_assembly(cd_rb_assembly_expansion.A0_assembly);
}
/**
* Pre-request all relevant element data.
*/
virtual void init_context(FEMContext &c)
{
For efficiency, we should prerequest all
the data we will need to build the
linear system before doing an element loop.
c.element_fe_var[u_var]->get_JxW();
c.element_fe_var[u_var]->get_phi();
c.element_fe_var[u_var]->get_dphi();
}
/**
* Variable number for u.
*/
unsigned int u_var;
/**
* The object that stores the "assembly" expansion of the parameter dependent PDE,
* i.e. the objects that define how to assemble the set of parameter-independent
* operators in the affine expansion of the PDE.
*/
CDRBAssemblyExpansion cd_rb_assembly_expansion;
/**
* The object that defines which degrees of freedom are on a Dirichlet boundary.
*/
AutoPtr<DirichletBoundary> dirichlet_bc;
};
#endif
The source file reduced_basis_ex1.C with comments:
/* rbOOmit 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 */
Reduced Basis Example 1 - Certified Reduced Basis Method
In this example problem we use the Certified Reduced Basis method to solve a steady convection-diffusion problem on the unit square. The reduced basis method relies on an expansion of the PDE in the form \sum_q=1^Q_a theta_a^q(\mu) * a^q(u,v) = \sum_q=1^Q_f theta_f^q(\mu) f^q(v) where theta_a, theta_f are parameter dependent functions and a^q, f^q are parameter independent operators (\mu denotes a parameter).
We first attach the parameter dependent functions and paramater independent operators to the RBSystem. Then in Offline mode, a reduced basis space is generated and written out to the directory "offline_data". In Online mode, the reduced basis data in "offline_data" is read in and used to solve the reduced problem for the parameters specified in reduced_basis_ex1.in.
We also attach four outputs to the system which are averages over certain subregions of the domain. In Online mode, we print out the values of these outputs as well as rigorous error bounds with respect to the output associated with the "truth" finite element discretization.
C++ include files that we need
#include <iostream>
#include <algorithm>
#include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
#include <cmath>
#include <set>
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/equation_systems.h"
#include "libmesh/dof_map.h"
#include "libmesh/getpot.h"
#include "libmesh/elem.h"
local includes
#include "rb_classes.h"
#include "assembly.h"
Bring in everything from the libMesh namespace
using namespace libMesh;
The main program.
int main (int argc, char** argv)
{
Initialize libMesh.
LibMeshInit init (argc, argv);
#if !defined(LIBMESH_HAVE_XDR)
We need XDR support to write out reduced bases
libmesh_example_assert(false, "--enable-xdr");
#elif defined(LIBMESH_DEFAULT_SINGLE_PRECISION)
XDR binary support requires double precision
libmesh_example_assert(false, "--disable-singleprecision");
#endif
FIXME: This example currently segfaults with Trilinos?
libmesh_example_assert(libMesh::default_solver_package() == PETSC_SOLVERS, "--enable-petsc");
Skip this 2D example if libMesh was compiled as 1D-only.
libmesh_example_assert(2 <= LIBMESH_DIM, "2D support");
Parse the input file (reduced_basis_ex1.in) using GetPot
std::string parameters_filename = "reduced_basis_ex1.in";
GetPot infile(parameters_filename);
unsigned int n_elem = infile("n_elem", 1); // Determines the number of elements in the "truth" mesh
const unsigned int dim = 2; // The number of spatial dimensions
bool store_basis_functions = infile("store_basis_functions", true); // Do we write the RB basis functions to disk?
Read the "online_mode" flag from the command line
GetPot command_line (argc, argv);
int online_mode = 0;
if ( command_line.search(1, "-online_mode") )
online_mode = command_line.next(online_mode);
Build a mesh.
Mesh mesh (dim);
MeshTools::Generation::build_square (mesh,
n_elem, n_elem,
0., 1.,
0., 1.,
QUAD4);
Create an equation systems object.
EquationSystems equation_systems (mesh);
We override RBConstruction with SimpleRBConstruction in order to
specialize a few functions for this particular problem.
SimpleRBConstruction & rb_con =
equation_systems.add_system<SimpleRBConstruction> ("RBConvectionDiffusion");
Initialize the data structures for the equation system.
equation_systems.init ();
Print out some information about the "truth" discretization
equation_systems.print_info();
mesh.print_info();
Build a new RBEvaluation object which will be used to perform
Reduced Basis calculations. This is required in both the
"Offline" and "Online" stages.
SimpleRBEvaluation rb_eval;
We need to give the RBConstruction object a pointer to
our RBEvaluation object
rb_con.set_rb_evaluation(rb_eval);
if(!online_mode) // Perform the Offline stage of the RB method
{
Read in the data that defines this problem from the specified text file
rb_con.process_parameters_file(parameters_filename);
Print out info that describes the current setup of rb_con
rb_con.print_info();
Prepare rb_con for the Construction stage of the RB method.
This sets up the necessary data structures and performs
initial assembly of the "truth" affine expansion of the PDE.
rb_con.initialize_rb_construction();
Compute the reduced basis space by computing "snapshots", i.e.
"truth" solves, at well-chosen parameter values and employing
these snapshots as basis functions.
rb_con.train_reduced_basis();
Write out the data that will subsequently be required for the Evaluation stage
rb_con.get_rb_evaluation().write_offline_data_to_files();
If requested, write out the RB basis functions for visualization purposes
if(store_basis_functions)
{
Write out the basis functions
rb_con.get_rb_evaluation().write_out_basis_functions(rb_con);
}
}
else // Perform the Online stage of the RB method
{
Read in the reduced basis data
rb_eval.read_offline_data_from_files();
Read in online_N and initialize online parameters
unsigned int online_N = infile("online_N",1);
Real online_x_vel = infile("online_x_vel", 0.);
Real online_y_vel = infile("online_y_vel", 0.);
RBParameters online_mu;
online_mu.set_value("x_vel", online_x_vel);
online_mu.set_value("y_vel", online_y_vel);
rb_eval.set_parameters(online_mu);
rb_eval.print_parameters();
Now do the Online solve using the precomputed reduced basis
rb_eval.rb_solve(online_N);
Print out outputs as well as the corresponding output error bounds.
std::cout << "output 1, value = " << rb_eval.RB_outputs[0]
<< ", bound = " << rb_eval.RB_output_error_bounds[0]
<< std::endl;
std::cout << "output 2, value = " << rb_eval.RB_outputs[1]
<< ", bound = " << rb_eval.RB_output_error_bounds[1]
<< std::endl;
std::cout << "output 3, value = " << rb_eval.RB_outputs[2]
<< ", bound = " << rb_eval.RB_output_error_bounds[2]
<< std::endl;
std::cout << "output 4, value = " << rb_eval.RB_outputs[3]
<< ", bound = " << rb_eval.RB_output_error_bounds[3]
<< std::endl << std::endl;
if(store_basis_functions)
{
Read in the basis functions
rb_eval.read_in_basis_functions(rb_con);
Plot the solution
rb_con.load_rb_solution();
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems ("RB_sol.e",equation_systems);
#endif
Plot the first basis function that was generated from the train_reduced_basis
call in the Offline stage
rb_con.load_basis_function(0);
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems ("bf0.e",equation_systems);
#endif
}
}
return 0;
}
The source file assembly.h without comments:
#ifndef __assembly_h__
#define __assembly_h__
#include "libmesh/sparse_matrix.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dense_vector.h"
#include "libmesh/fe.h"
#include "libmesh/fe_interface.h"
#include "libmesh/fe_base.h"
#include "libmesh/elem_assembly.h"
#include "libmesh/quadrature_gauss.h"
#include "libmesh/rb_theta.h"
#include "libmesh/rb_assembly_expansion.h"
using libMesh::ElemAssembly;
using libMesh::FEInterface;
using libMesh::FEMContext;
using libMesh::Number;
using libMesh::Point;
using libMesh::RBAssemblyExpansion;
using libMesh::RBParameters;
using libMesh::RBTheta;
using libMesh::RBThetaExpansion;
using libMesh::Real;
using libMesh::RealGradient;
struct ThetaA0 : RBTheta { virtual Number evaluate(const RBParameters& ) { return 0.05; } };
struct ThetaA1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("x_vel"); } };
struct ThetaA2 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("y_vel"); } };
struct A0 : ElemAssembly
{
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] * dphi[j][qp]*dphi[i][qp];
}
};
struct A1 : ElemAssembly
{
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] *dphi[j][qp](0)*phi[i][qp];
}
};
struct A2 : ElemAssembly
{
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const std::vector<std::vector<RealGradient> >& dphi =
c.element_fe_var[u_var]->get_dphi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
for (unsigned int j=0; j != n_u_dofs; j++)
c.get_elem_jacobian()(i,j) += JxW[qp] *dphi[j][qp](1)*phi[i][qp];
}
};
struct F0 : ElemAssembly
{
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*phi[i][qp] );
}
};
struct OutputAssembly : ElemAssembly
{
OutputAssembly(Real min_x_in, Real max_x_in,
Real min_y_in, Real max_y_in)
:
min_x(min_x_in),
max_x(max_x_in),
min_y(min_y_in),
max_y(max_y_in)
{}
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW =
c.element_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[u_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
Real output_area = (max_x-min_x) * (max_y-min_y);
Point centroid = c.elem->centroid();
if( (min_x <= centroid(0)) && (centroid(0) <= max_x) &&
(min_y <= centroid(1)) && (centroid(1) <= max_y) )
for (unsigned int qp=0; qp != n_qpoints; qp++)
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*phi[i][qp] ) / output_area;
}
Real min_x, max_x, min_y, max_y;
};
struct CDRBThetaExpansion : RBThetaExpansion
{
/**
* Constructor.
*/
CDRBThetaExpansion()
{
attach_A_theta(&theta_a_0); // Attach the lhs theta
attach_A_theta(&theta_a_1);
attach_A_theta(&theta_a_2);
attach_F_theta(&rb_theta); // Attach the rhs theta
attach_output_theta(&rb_theta); // Attach output 0 theta
attach_output_theta(&rb_theta); // Attach output 1 theta
attach_output_theta(&rb_theta); // Attach output 2 theta
attach_output_theta(&rb_theta); // Attach output 3 theta
}
ThetaA0 theta_a_0;
ThetaA1 theta_a_1;
ThetaA2 theta_a_2;
RBTheta rb_theta; // Default RBTheta object, just returns 1.
};
struct CDRBAssemblyExpansion : RBAssemblyExpansion
{
/**
* Constructor.
*/
CDRBAssemblyExpansion()
:
L0(0.7,0.8,0.7,0.8),
L1(0.2,0.3,0.7,0.8),
L2(0.2,0.3,0.2,0.3),
L3(0.7,0.8,0.2,0.3)
{
attach_A_assembly(&A0_assembly); // Attach the lhs assembly
attach_A_assembly(&A1_assembly);
attach_A_assembly(&A2_assembly);
attach_F_assembly(&F0_assembly); // Attach the rhs assembly
attach_output_assembly(&L0); // Attach output 0 assembly
attach_output_assembly(&L1); // Attach output 1 assembly
attach_output_assembly(&L2); // Attach output 2 assembly
attach_output_assembly(&L3); // Attach output 3 assembly
}
A0 A0_assembly;
A1 A1_assembly;
A2 A2_assembly;
F0 F0_assembly;
OutputAssembly L0;
OutputAssembly L1;
OutputAssembly L2;
OutputAssembly L3;
};
#endif
The source file rb_classes.h without comments:
#ifndef __rb_classes_h__
#define __rb_classes_h__
#include "libmesh/rb_construction.h"
#include "libmesh/fe_base.h"
#include "assembly.h"
using libMesh::AutoPtr;
using libMesh::DirichletBoundary;
using libMesh::EquationSystems;
using libMesh::FEMContext;
using libMesh::RBConstruction;
using libMesh::RBEvaluation;
using libMesh::Real;
class SimpleRBEvaluation : public RBEvaluation
{
public:
/**
* Constructor. Just set the theta expansion.
*/
SimpleRBEvaluation()
{
set_rb_theta_expansion(cd_rb_theta_expansion);
}
/**
* The coercivity constant is bounded below by 0.05.
*/
virtual Real get_stability_lower_bound() { return 0.05; }
/**
* The object that stores the "theta" expansion of the parameter dependent PDE,
* i.e. the set of parameter-dependent functions in the affine expansion of the PDE.
*/
CDRBThetaExpansion cd_rb_theta_expansion;
};
class SimpleRBConstruction : public RBConstruction
{
public:
SimpleRBConstruction (EquationSystems& es,
const std::string& name_in,
const unsigned int number_in)
: Parent(es, name_in, number_in),
dirichlet_bc(AutoPtr<DirichletBoundary>(NULL))
{}
/**
* Destructor.
*/
virtual ~SimpleRBConstruction () { }
/**
* The type of system.
*/
typedef SimpleRBConstruction sys_type;
/**
* The type of the parent.
*/
typedef RBConstruction Parent;
/**
* Initialize data structures.
*/
virtual void init_data()
{
u_var = this->add_variable ("u", FIRST);
dirichlet_bc = build_zero_dirichlet_boundary_object();
dirichlet_bc->b.insert(0);
dirichlet_bc->b.insert(1);
dirichlet_bc->b.insert(2);
dirichlet_bc->b.insert(3);
dirichlet_bc->variables.push_back(u_var);
get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
Parent::init_data();
set_rb_assembly_expansion(cd_rb_assembly_expansion);
set_inner_product_assembly(cd_rb_assembly_expansion.A0_assembly);
}
/**
* Pre-request all relevant element data.
*/
virtual void init_context(FEMContext &c)
{
c.element_fe_var[u_var]->get_JxW();
c.element_fe_var[u_var]->get_phi();
c.element_fe_var[u_var]->get_dphi();
}
/**
* Variable number for u.
*/
unsigned int u_var;
/**
* The object that stores the "assembly" expansion of the parameter dependent PDE,
* i.e. the objects that define how to assemble the set of parameter-independent
* operators in the affine expansion of the PDE.
*/
CDRBAssemblyExpansion cd_rb_assembly_expansion;
/**
* The object that defines which degrees of freedom are on a Dirichlet boundary.
*/
AutoPtr<DirichletBoundary> dirichlet_bc;
};
#endif
The source file reduced_basis_ex1.C without comments:
/* rbOOmit 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 */
#include <iostream>
#include <algorithm>
#include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
#include <cmath>
#include <set>
#include "libmesh/libmesh.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/exodusII_io.h"
#include "libmesh/equation_systems.h"
#include "libmesh/dof_map.h"
#include "libmesh/getpot.h"
#include "libmesh/elem.h"
#include "rb_classes.h"
#include "assembly.h"
using namespace libMesh;
int main (int argc, char** argv)
{
LibMeshInit init (argc, argv);
#if !defined(LIBMESH_HAVE_XDR)
libmesh_example_assert(false, "--enable-xdr");
#elif defined(LIBMESH_DEFAULT_SINGLE_PRECISION)
libmesh_example_assert(false, "--disable-singleprecision");
#endif
libmesh_example_assert(libMesh::default_solver_package() == PETSC_SOLVERS, "--enable-petsc");
libmesh_example_assert(2 <= LIBMESH_DIM, "2D support");
std::string parameters_filename = "reduced_basis_ex1.in";
GetPot infile(parameters_filename);
unsigned int n_elem = infile("n_elem", 1); // Determines the number of elements in the "truth" mesh
const unsigned int dim = 2; // The number of spatial dimensions
bool store_basis_functions = infile("store_basis_functions", true); // Do we write the RB basis functions to disk?
GetPot command_line (argc, argv);
int online_mode = 0;
if ( command_line.search(1, "-online_mode") )
online_mode = command_line.next(online_mode);
Mesh mesh (dim);
MeshTools::Generation::build_square (mesh,
n_elem, n_elem,
0., 1.,
0., 1.,
QUAD4);
EquationSystems equation_systems (mesh);
SimpleRBConstruction & rb_con =
equation_systems.add_system<SimpleRBConstruction> ("RBConvectionDiffusion");
equation_systems.init ();
equation_systems.print_info();
mesh.print_info();
SimpleRBEvaluation rb_eval;
rb_con.set_rb_evaluation(rb_eval);
if(!online_mode) // Perform the Offline stage of the RB method
{
rb_con.process_parameters_file(parameters_filename);
rb_con.print_info();
rb_con.initialize_rb_construction();
rb_con.train_reduced_basis();
rb_con.get_rb_evaluation().write_offline_data_to_files();
if(store_basis_functions)
{
rb_con.get_rb_evaluation().write_out_basis_functions(rb_con);
}
}
else // Perform the Online stage of the RB method
{
rb_eval.read_offline_data_from_files();
unsigned int online_N = infile("online_N",1);
Real online_x_vel = infile("online_x_vel", 0.);
Real online_y_vel = infile("online_y_vel", 0.);
RBParameters online_mu;
online_mu.set_value("x_vel", online_x_vel);
online_mu.set_value("y_vel", online_y_vel);
rb_eval.set_parameters(online_mu);
rb_eval.print_parameters();
rb_eval.rb_solve(online_N);
std::cout << "output 1, value = " << rb_eval.RB_outputs[0]
<< ", bound = " << rb_eval.RB_output_error_bounds[0]
<< std::endl;
std::cout << "output 2, value = " << rb_eval.RB_outputs[1]
<< ", bound = " << rb_eval.RB_output_error_bounds[1]
<< std::endl;
std::cout << "output 3, value = " << rb_eval.RB_outputs[2]
<< ", bound = " << rb_eval.RB_output_error_bounds[2]
<< std::endl;
std::cout << "output 4, value = " << rb_eval.RB_outputs[3]
<< ", bound = " << rb_eval.RB_output_error_bounds[3]
<< std::endl << std::endl;
if(store_basis_functions)
{
rb_eval.read_in_basis_functions(rb_con);
rb_con.load_rb_solution();
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems ("RB_sol.e",equation_systems);
#endif
rb_con.load_basis_function(0);
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems ("bf0.e",equation_systems);
#endif
}
}
return 0;
}
The console output of the program:
***************************************************************
* Running Example reduced_basis_ex1:
* mpirun -np 2 example-devel -online_mode 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
*** Warning, This code is untested, experimental, or likely to see future API changes: src/reduced_basis/rb_parametrized.C, line 41, compiled Feb 5 2013 at 13:19:15 ***
EquationSystems
n_systems()=1
System #0, "RBConvectionDiffusion"
Type "RBConstruction"
Variables="u"
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=2601
n_local_dofs()=1330
n_constrained_dofs()=200
n_local_constrained_dofs()=94
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 8.70934
Average Off-Processor Bandwidth <= 0.132257
Maximum On-Processor Bandwidth <= 11
Maximum Off-Processor Bandwidth <= 5
DofMap Constraints
Number of DoF Constraints = 200
Average DoF Constraint Length= 0
Number of Node Constraints = 0
Mesh Information:
mesh_dimension()=2
spatial_dimension()=3
n_nodes()=2601
n_local_nodes()=1330
n_elem()=2500
n_local_elem()=1250
n_active_elem()=2500
n_subdomains()=1
n_partitions()=2
n_processors()=2
n_threads()=1
processor_id()=0
Initializing training parameters with deterministic training set...
Parameter x_vel: log scaling = 0
Parameter y_vel: log scaling = 0
RBConstruction parameters:
system name: RBConvectionDiffusion
constrained_problem: 0
Nmax: 20
Aq operators attached: 3
Fq functions attached: 1
n_outputs: 4
output 0, Q_l = 1
output 1, Q_l = 1
output 2, Q_l = 1
output 3, Q_l = 1
Number of parameters: 2
Parameter x_vel: Min = -2, Max = 2, value = -2
Parameter y_vel: Min = -2, Max = 2, value = 0.79
n_training_samples: 100
single-matrix mode? 0
reuse preconditioner? 1
use a relative error bound in greedy? 0
write out data during basis training? 0
quiet mode? 1
Compute output dual inner products
output_dual_innerprods[0][0] = 0.323675
output_dual_innerprods[1][0] = 0.323675
output_dual_innerprods[2][0] = 0.323675
output_dual_innerprods[3][0] = 0.323675
---- Performing Greedy basis enrichment ----
---- Basis dimension: 0 ----
Performing RB solves on training set
Maximum (absolute) error bound is 3.74824
Performing truth solve at parameter:
x_vel: -2.000000e+00
y_vel: -2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 1 ----
Performing RB solves on training set
Maximum (absolute) error bound is 6.67625
Performing truth solve at parameter:
x_vel: 2.000000e+00
y_vel: 2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 2 ----
Performing RB solves on training set
Maximum (absolute) error bound is 5.18295
Performing truth solve at parameter:
x_vel: 2.000000e+00
y_vel: -2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 3 ----
Performing RB solves on training set
Maximum (absolute) error bound is 3.59028
Performing truth solve at parameter:
x_vel: -2.000000e+00
y_vel: 2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 4 ----
Performing RB solves on training set
Maximum (absolute) error bound is 2.71338
Performing truth solve at parameter:
x_vel: -2.222222e-01
y_vel: -2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 5 ----
Performing RB solves on training set
Maximum (absolute) error bound is 2.50784
Performing truth solve at parameter:
x_vel: 6.666667e-01
y_vel: 2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 6 ----
Performing RB solves on training set
Maximum (absolute) error bound is 2.15653
Performing truth solve at parameter:
x_vel: -2.222222e-01
y_vel: 6.666667e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 7 ----
Performing RB solves on training set
Maximum (absolute) error bound is 1.71352
Performing truth solve at parameter:
x_vel: 2.222222e-01
y_vel: -1.111111e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 8 ----
Performing RB solves on training set
Maximum (absolute) error bound is 1.61442
Performing truth solve at parameter:
x_vel: -1.555556e+00
y_vel: 2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 9 ----
Performing RB solves on training set
Maximum (absolute) error bound is 1.09712
Performing truth solve at parameter:
x_vel: 2.000000e+00
y_vel: -2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 10 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.9794
Performing truth solve at parameter:
x_vel: 2.222222e-01
y_vel: 2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 11 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.763958
Performing truth solve at parameter:
x_vel: 2.222222e-01
y_vel: -2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 12 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.616746
Performing truth solve at parameter:
x_vel: -1.111111e+00
y_vel: -1.111111e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 13 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.530181
Performing truth solve at parameter:
x_vel: -2.222222e-01
y_vel: -2.000000e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 14 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.475932
Performing truth solve at parameter:
x_vel: -6.666667e-01
y_vel: 2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 15 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.440399
Performing truth solve at parameter:
x_vel: 1.111111e+00
y_vel: -6.666667e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 16 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.423537
Performing truth solve at parameter:
x_vel: -1.111111e+00
y_vel: 1.111111e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 17 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.408558
Performing truth solve at parameter:
x_vel: 6.666667e-01
y_vel: 1.111111e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 18 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.333764
Performing truth solve at parameter:
x_vel: 2.222222e-01
y_vel: 2.222222e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 19 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.196867
Performing truth solve at parameter:
x_vel: -2.222222e-01
y_vel: -6.666667e-01
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 20 ----
Performing RB solves on training set
Maximum (absolute) error bound is 0.191475
Maximum number of basis functions reached: Nmax = 20
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/reduced_basis/reduced_basis_ex1/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:46:12 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): 3.726e+00 1.00000 3.726e+00
Objects: 1.148e+03 1.00000 1.148e+03
Flops: 7.805e+08 1.05324 7.607e+08 1.521e+09
Flops/sec: 2.094e+08 1.05324 2.042e+08 4.083e+08
MPI Messages: 8.408e+03 1.00000 8.408e+03 1.682e+04
MPI Message Lengths: 4.408e+06 1.00000 5.243e+02 8.816e+06
MPI Reductions: 2.424e+04 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: 3.7263e+00 100.0% 1.5215e+09 100.0% 1.682e+04 100.0% 5.243e+02 100.0% 2.424e+04 100.0%
------------------------------------------------------------------------------------------------------------------------
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
VecDot 4075 1.0 1.8129e-02 1.4 1.08e+07 1.0 0.0e+00 0.0e+00 4.1e+03 0 1 0 0 17 0 1 0 0 17 1169
VecMDot 4056 1.0 4.9089e-02 1.3 1.46e+08 1.0 0.0e+00 0.0e+00 4.1e+03 1 19 0 0 17 1 19 0 0 17 5826
VecNorm 4308 1.0 1.2588e-02 1.2 1.15e+07 1.0 0.0e+00 0.0e+00 4.3e+03 0 1 0 0 18 0 1 0 0 18 1780
VecScale 4323 1.0 2.5177e-03 1.1 5.72e+06 1.0 0.0e+00 0.0e+00 0.0e+00 0 1 0 0 0 0 1 0 0 0 4445
VecCopy 547 1.0 4.8804e-04 1.1 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 10005 1.0 5.3079e-03 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 549 1.0 6.0558e-04 1.1 1.46e+06 1.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 4716
VecMAXPY 4223 1.0 4.2430e-02 1.1 1.57e+08 1.0 0.0e+00 0.0e+00 0.0e+00 1 20 0 0 0 1 20 0 0 0 7241
VecAssemblyBegin 595 1.0 9.3760e-03 2.6 0.00e+00 0.0 1.0e+01 7.8e+02 1.8e+03 0 0 0 0 7 0 0 0 0 7 0
VecAssemblyEnd 595 1.0 2.3818e-04 1.1 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 8223 1.0 1.0593e-02 1.0 0.00e+00 0.0 1.6e+04 5.3e+02 0.0e+00 0 0 98 99 0 0 0 98 99 0 0
VecScatterEnd 8223 1.0 1.1411e-02 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
VecNormalize 4223 1.0 1.6112e-02 1.2 1.68e+07 1.0 0.0e+00 0.0e+00 4.2e+03 0 2 0 0 17 0 2 0 0 17 2045
MatMult 4223 1.0 7.3325e-02 1.0 9.31e+07 1.1 8.4e+03 5.3e+02 0.0e+00 2 12 50 51 0 2 12 50 51 0 2477
MatMultAdd 3915 1.0 6.6257e-02 1.0 9.15e+07 1.1 7.8e+03 5.3e+02 0.0e+00 2 12 47 47 0 2 12 47 47 0 2695
MatSolve 4308 1.0 1.1205e-01 1.1 2.50e+08 1.1 0.0e+00 0.0e+00 0.0e+00 3 32 0 0 0 3 32 0 0 0 4332
MatLUFactorNum 42 1.0 1.7164e-02 1.1 1.30e+07 1.1 0.0e+00 0.0e+00 0.0e+00 0 2 0 0 0 0 2 0 0 0 1456
MatILUFactorSym 42 1.0 3.2199e-02 1.1 0.00e+00 0.0 0.0e+00 0.0e+00 1.3e+02 1 0 0 0 1 1 0 0 0 1 0
MatAssemblyBegin 4213 1.0 3.2993e-02 1.6 0.00e+00 0.0 1.2e+01 2.8e+03 8.4e+03 1 0 0 0 35 1 0 0 0 35 0
MatAssemblyEnd 4213 1.0 8.0548e-02 1.0 0.00e+00 0.0 3.4e+02 1.3e+02 7.0e+02 2 0 2 1 3 2 0 2 1 3 0
MatGetRow 218120 1.0 1.7584e-02 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
MatGetRowIJ 42 1.0 3.0994e-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 42 1.0 7.6604e-04 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 8.4e+01 0 0 0 0 0 0 0 0 0 0 0
MatZeroEntries 56 1.0 4.2105e-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
MatAXPY 82 1.0 8.3786e-02 1.0 0.00e+00 0.0 3.3e+02 1.3e+02 1.3e+03 2 0 2 0 5 2 0 2 0 5 0
KSPGMRESOrthog 4056 1.0 8.9388e-02 1.1 2.93e+08 1.0 0.0e+00 0.0e+00 4.1e+03 2 38 0 0 17 2 38 0 0 17 6401
KSPSetUp 127 1.0 3.3426e-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 85 1.0 3.6957e-01 1.0 6.77e+08 1.1 8.4e+03 5.3e+02 8.6e+03 10 87 50 51 35 10 87 50 51 35 3573
PCSetUp 84 1.0 5.4420e-02 1.1 1.30e+07 1.1 0.0e+00 0.0e+00 2.1e+02 1 2 0 0 1 1 2 0 0 1 459
PCSetUpOnBlocks 85 1.0 5.3261e-02 1.1 1.30e+07 1.1 0.0e+00 0.0e+00 2.1e+02 1 2 0 0 1 1 2 0 0 1 469
PCApply 4308 1.0 1.3700e-01 1.1 2.50e+08 1.1 0.0e+00 0.0e+00 0.0e+00 4 32 0 0 0 4 32 0 0 0 3543
------------------------------------------------------------------------------------------------------------------------
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 449 449 3585480 0
Vector Scatter 88 88 91168 0
Index Set 302 302 470392 0
IS L to G Mapping 1 1 564 0
Matrix 303 303 31664508 0
Krylov Solver 2 2 19360 0
Preconditioner 2 2 1784 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 0
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
-online_mode 0
-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:46:12 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=3.73708, Active time=3.68534 |
-------------------------------------------------------------------------------------------------------------------
| 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.0279 0.027879 0.0314 0.031380 0.76 0.85 |
| build_constraint_matrix() 11250 0.0501 0.000004 0.0501 0.000004 1.36 1.36 |
| build_sparsity() 1 0.0191 0.019069 0.0535 0.053522 0.52 1.45 |
| cnstrn_elem_mat_vec() 11250 0.0343 0.000003 0.0343 0.000003 0.93 0.93 |
| create_dof_constraints() 1 0.0211 0.021123 0.0889 0.088871 0.57 2.41 |
| distribute_dofs() 1 0.0254 0.025423 0.0687 0.068689 0.69 1.86 |
| dof_indices() 26380 0.8536 0.000032 0.8536 0.000032 23.16 23.16 |
| prepare_send_list() 1 0.0001 0.000058 0.0001 0.000058 0.00 0.00 |
| reinit() 1 0.0430 0.042972 0.0430 0.042972 1.17 1.17 |
| |
| FE |
| compute_shape_functions() 24174 0.3960 0.000016 0.3960 0.000016 10.74 10.74 |
| init_shape_functions() 1692 0.0251 0.000015 0.0251 0.000015 0.68 0.68 |
| inverse_map() 3348 0.0205 0.000006 0.0205 0.000006 0.56 0.56 |
| |
| FEMap |
| compute_affine_map() 24174 0.1502 0.000006 0.1502 0.000006 4.07 4.07 |
| compute_face_map() 1674 0.0199 0.000012 0.0410 0.000025 0.54 1.11 |
| init_face_shape_functions() 18 0.0002 0.000010 0.0002 0.000010 0.00 0.00 |
| init_reference_to_physical_map() 1692 0.0151 0.000009 0.0151 0.000009 0.41 0.41 |
| |
| Mesh |
| find_neighbors() 1 0.0289 0.028879 0.0294 0.029365 0.78 0.80 |
| renumber_nodes_and_elem() 2 0.0010 0.000492 0.0010 0.000492 0.03 0.03 |
| |
| MeshCommunication |
| assign_global_indices() 1 0.0630 0.063016 0.0639 0.063890 1.71 1.73 |
| compute_hilbert_indices() 2 0.0216 0.010815 0.0216 0.010815 0.59 0.59 |
| find_global_indices() 2 0.0086 0.004303 0.0329 0.016445 0.23 0.89 |
| parallel_sort() 2 0.0021 0.001030 0.0023 0.001164 0.06 0.06 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0055 0.005471 0.0055 0.005471 0.15 0.15 |
| |
| MetisPartitioner |
| partition() 1 0.0330 0.033004 0.0489 0.048929 0.90 1.33 |
| |
| Parallel |
| allgather() 14 0.0007 0.000051 0.0007 0.000053 0.02 0.02 |
| barrier() 1 0.0000 0.000014 0.0000 0.000014 0.00 0.00 |
| broadcast() 44 0.0001 0.000003 0.0001 0.000003 0.00 0.00 |
| max(bool) 5 0.0000 0.000002 0.0000 0.000002 0.00 0.00 |
| max(scalar) 119 0.0003 0.000003 0.0003 0.000003 0.01 0.01 |
| max(vector) 28 0.0002 0.000006 0.0004 0.000013 0.00 0.01 |
| maxloc(scalar) 21 0.0663 0.003155 0.0663 0.003155 1.80 1.80 |
| min(bool) 138 0.0003 0.000002 0.0003 0.000002 0.01 0.01 |
| min(scalar) 113 0.0016 0.000014 0.0016 0.000014 0.04 0.04 |
| min(vector) 28 0.0002 0.000008 0.0005 0.000017 0.01 0.01 |
| probe() 32 0.0003 0.000009 0.0003 0.000009 0.01 0.01 |
| receive() 24 0.0006 0.000026 0.0009 0.000037 0.02 0.02 |
| send() 20 0.0001 0.000004 0.0001 0.000004 0.00 0.00 |
| send_receive() 24 0.0004 0.000015 0.0009 0.000037 0.01 0.02 |
| sum() 44 0.0004 0.000010 0.0005 0.000012 0.01 0.01 |
| |
| Parallel::Request |
| wait() 20 0.0000 0.000002 0.0000 0.000002 0.00 0.00 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0023 0.002346 0.0024 0.002424 0.06 0.07 |
| set_parent_processor_ids() 1 0.0018 0.001831 0.0018 0.001831 0.05 0.05 |
| |
| PetscLinearSolver |
| solve() 85 0.3764 0.004429 0.3764 0.004429 10.21 10.21 |
| |
| RBConstruction |
| add_scaled_matrix_and_vector() 9 0.3717 0.041295 1.8486 0.205404 10.08 50.16 |
| clear() 1 0.0007 0.000738 0.0007 0.000738 0.02 0.02 |
| compute_Fq_representor_innerprods() 1 0.0012 0.001226 0.0064 0.006391 0.03 0.17 |
| compute_max_error_bound() 21 0.0126 0.000601 0.4628 0.022037 0.34 12.56 |
| compute_output_dual_innerprods() 1 0.0025 0.002460 0.0203 0.020347 0.07 0.55 |
| enrich_RB_space() 20 0.0142 0.000712 0.0142 0.000712 0.39 0.39 |
| train_reduced_basis() 1 0.0016 0.001640 1.1312 1.131159 0.04 30.69 |
| truth_assembly() 20 0.0644 0.003219 0.0644 0.003219 1.75 1.75 |
| truth_solve() 20 0.0020 0.000102 0.1547 0.007737 0.06 4.20 |
| update_RB_system_matrices() 20 0.0609 0.003044 0.0609 0.003044 1.65 1.65 |
| update_residual_terms() 20 0.1447 0.007235 0.4098 0.020489 3.93 11.12 |
| |
| RBEvaluation |
| clear() 1 0.0002 0.000182 0.0002 0.000182 0.00 0.00 |
| compute_residual_dual_norm() 1050 0.3379 0.000322 0.3379 0.000322 9.17 9.17 |
| rb_solve() 1050 0.0454 0.000043 0.3835 0.000365 1.23 10.41 |
| resize_data_structures() 1 0.0003 0.000314 0.0003 0.000314 0.01 0.01 |
| write_offline_data_to_files() 1 0.0011 0.001085 0.0011 0.001085 0.03 0.03 |
| write_out_basis_functions() 1 0.0001 0.000056 0.3712 0.371209 0.00 10.07 |
| write_out_vectors() 1 0.3065 0.306506 0.3712 0.371153 8.32 10.07 |
-------------------------------------------------------------------------------------------------------------------
| Totals: 108671 3.6853 100.00 |
-------------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example reduced_basis_ex1:
* mpirun -np 2 example-devel -online_mode 0 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
***************************************************************
* Running Example reduced_basis_ex1:
* mpirun -np 2 example-devel -online_mode 1 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
*** Warning, This code is untested, experimental, or likely to see future API changes: src/reduced_basis/rb_parametrized.C, line 41, compiled Feb 5 2013 at 13:19:15 ***
EquationSystems
n_systems()=1
System #0, "RBConvectionDiffusion"
Type "RBConstruction"
Variables="u"
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=2601
n_local_dofs()=1330
n_constrained_dofs()=200
n_local_constrained_dofs()=94
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 8.70934
Average Off-Processor Bandwidth <= 0.132257
Maximum On-Processor Bandwidth <= 11
Maximum Off-Processor Bandwidth <= 5
DofMap Constraints
Number of DoF Constraints = 200
Average DoF Constraint Length= 0
Number of Node Constraints = 0
Mesh Information:
mesh_dimension()=2
spatial_dimension()=3
n_nodes()=2601
n_local_nodes()=1330
n_elem()=2500
n_local_elem()=1250
n_active_elem()=2500
n_subdomains()=1
n_partitions()=2
n_processors()=2
n_threads()=1
processor_id()=0
x_vel: -2.000000e+00
y_vel: 7.900000e-01
output 1, value = 0.12088, bound = 0.0809549
output 2, value = 0.368636, bound = 0.0809549
output 3, value = 0.256529, bound = 0.0809549
output 4, value = 0.116739, bound = 0.0809549
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/reduced_basis/reduced_basis_ex1/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:46:13 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): 8.115e-01 1.00000 8.115e-01
Objects: 3.300e+01 1.00000 3.300e+01
Flops: 5.320e+04 1.04642 5.202e+04 1.040e+05
Flops/sec: 6.556e+04 1.04642 6.410e+04 1.282e+05
MPI Messages: 4.000e+00 1.00000 4.000e+00 8.000e+00
MPI Message Lengths: 1.984e+03 1.00000 4.960e+02 3.968e+03
MPI Reductions: 1.070e+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: 8.1149e-01 100.0% 1.0404e+05 100.0% 8.000e+00 100.0% 4.960e+02 100.0% 1.060e+02 99.1%
------------------------------------------------------------------------------------------------------------------------
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
VecCopy 3 1.0 9.7752e-06 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
VecSet 25 1.0 1.9312e-05 1.1 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 20 1.0 2.2888e-05 1.1 5.32e+04 1.0 0.0e+00 0.0e+00 0.0e+00 0100 0 0 0 0100 0 0 0 4546
VecAssemblyBegin 23 1.0 9.4585e-0367.9 0.00e+00 0.0 0.0e+00 0.0e+00 6.9e+01 1 0 0 0 64 1 0 0 0 65 0
VecAssemblyEnd 23 1.0 2.7180e-05 1.1 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 2 1.0 3.9101e-05 1.0 0.00e+00 0.0 4.0e+00 7.9e+02 0.0e+00 0 0 50 80 0 0 0 50 80 0 0
VecScatterEnd 2 1.0 9.2983e-06 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
MatZeroEntries 2 1.0 3.0994e-05 1.1 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
------------------------------------------------------------------------------------------------------------------------
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 25 25 293224 0
Vector Scatter 1 1 1036 0
Index Set 2 2 1744 0
IS L to G Mapping 1 1 564 0
Matrix 3 3 182024 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 9.53674e-08
Average time for MPI_Barrier(): 6.19888e-07
Average time for zero size MPI_Send(): 1.50204e-05
#PETSc Option Table entries:
-ksp_right_pc
-log_summary
-online_mode 1
-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:46:13 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.823023, Active time=0.790483 |
--------------------------------------------------------------------------------------------------------------
| 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.0276 0.027561 0.0310 0.031017 3.49 3.92 |
| build_sparsity() 1 0.0185 0.018458 0.0532 0.053220 2.34 6.73 |
| create_dof_constraints() 1 0.0207 0.020724 0.0876 0.087614 2.62 11.08 |
| distribute_dofs() 1 0.0253 0.025285 0.0685 0.068541 3.20 8.67 |
| dof_indices() 6380 0.1690 0.000026 0.1690 0.000026 21.39 21.39 |
| prepare_send_list() 1 0.0001 0.000063 0.0001 0.000063 0.01 0.01 |
| reinit() 1 0.0423 0.042271 0.0423 0.042271 5.35 5.35 |
| |
| EquationSystems |
| build_solution_vector() 2 0.0090 0.004490 0.0759 0.037932 1.14 9.60 |
| |
| ExodusII_IO |
| write_nodal_data() 2 0.0188 0.009395 0.0188 0.009395 2.38 2.38 |
| |
| Mesh |
| find_neighbors() 1 0.0287 0.028726 0.0288 0.028778 3.63 3.64 |
| renumber_nodes_and_elem() 2 0.0010 0.000507 0.0010 0.000507 0.13 0.13 |
| |
| MeshCommunication |
| assign_global_indices() 1 0.0626 0.062588 0.0630 0.062961 7.92 7.96 |
| compute_hilbert_indices() 2 0.0217 0.010849 0.0217 0.010849 2.75 2.75 |
| find_global_indices() 2 0.0081 0.004052 0.0323 0.016139 1.03 4.08 |
| parallel_sort() 2 0.0021 0.001028 0.0021 0.001074 0.26 0.27 |
| |
| MeshOutput |
| write_equation_systems() 2 0.0001 0.000050 0.0949 0.047425 0.01 12.00 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0051 0.005138 0.0051 0.005138 0.65 0.65 |
| |
| MetisPartitioner |
| partition() 1 0.0320 0.031957 0.0476 0.047612 4.04 6.02 |
| |
| Parallel |
| allgather() 14 0.0004 0.000027 0.0004 0.000029 0.05 0.05 |
| barrier() 1 0.0000 0.000018 0.0000 0.000018 0.00 0.00 |
| broadcast() 13 0.0002 0.000014 0.0002 0.000012 0.02 0.02 |
| max(bool) 1 0.0000 0.000003 0.0000 0.000003 0.00 0.00 |
| max(scalar) 135 0.0003 0.000002 0.0003 0.000002 0.04 0.04 |
| max(vector) 30 0.0002 0.000006 0.0004 0.000013 0.02 0.05 |
| min(bool) 156 0.0003 0.000002 0.0003 0.000002 0.04 0.04 |
| min(scalar) 129 0.1286 0.000997 0.1286 0.000997 16.27 16.27 |
| min(vector) 30 0.0002 0.000008 0.0005 0.000016 0.03 0.06 |
| probe() 24 0.0003 0.000011 0.0003 0.000011 0.03 0.03 |
| receive() 22 0.0002 0.000010 0.0005 0.000021 0.03 0.06 |
| send() 22 0.0001 0.000007 0.0001 0.000007 0.02 0.02 |
| send_receive() 24 0.0004 0.000015 0.0009 0.000036 0.04 0.11 |
| sum() 31 0.0003 0.000008 0.1254 0.004045 0.03 15.86 |
| |
| Parallel::Request |
| wait() 22 0.0002 0.000009 0.0002 0.000009 0.03 0.03 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0023 0.002341 0.0035 0.003541 0.30 0.45 |
| set_parent_processor_ids() 1 0.0018 0.001784 0.0018 0.001784 0.23 0.23 |
| |
| RBConstruction |
| clear() 1 0.0003 0.000319 0.0003 0.000319 0.04 0.04 |
| load_basis_function() 1 0.0000 0.000049 0.0000 0.000049 0.01 0.01 |
| load_rb_solution() 1 0.0002 0.000176 0.0002 0.000176 0.02 0.02 |
| |
| RBEvaluation |
| clear() 1 0.0001 0.000079 0.0001 0.000079 0.01 0.01 |
| compute_residual_dual_norm() 1 0.0010 0.001017 0.0010 0.001017 0.13 0.13 |
| rb_solve() 1 0.0123 0.012300 0.0133 0.013318 1.56 1.68 |
| read_in_basis_functions() 1 0.0000 0.000027 0.3348 0.334822 0.00 42.36 |
| read_in_vectors() 1 0.1469 0.146875 0.3348 0.334795 18.58 42.35 |
| read_offline_data_from_files() 1 0.0009 0.000862 0.0012 0.001243 0.11 0.16 |
| resize_data_structures() 1 0.0004 0.000381 0.0004 0.000381 0.05 0.05 |
--------------------------------------------------------------------------------------------------------------
| Totals: 7070 0.7905 100.00 |
--------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example reduced_basis_ex1:
* mpirun -np 2 example-devel -online_mode 1 -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************