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"
#include "libmesh/boundary_info.h"
rbOOmit includes
#include "libmesh/rb_theta.h"
#include "libmesh/rb_assembly_expansion.h"
#include "libmesh/rb_parametrized_function.h"
#include "libmesh/rb_eim_construction.h"
#include "libmesh/rb_eim_theta.h"
Bring in bits from the libMesh namespace.
Just the bits we're using, since this is a header.
using libMesh::boundary_id_type;
using libMesh::DenseSubMatrix;
using libMesh::ElemAssembly;
using libMesh::FEInterface;
using libMesh::FEMContext;
using libMesh::Number;
using libMesh::Point;
using libMesh::RBAssemblyExpansion;
using libMesh::RBConstruction;
using libMesh::RBParameters;
using libMesh::RBParametrizedFunction;
using libMesh::RBTheta;
using libMesh::RBThetaExpansion;
using libMesh::RBEIMAssembly;
using libMesh::RBEIMConstruction;
using libMesh::RBEIMEvaluation;
using libMesh::RBEIMTheta;
using libMesh::Real;
using libMesh::RealGradient;
struct ElemAssemblyWithConstruction : ElemAssembly
{
RBConstruction* rb_con;
};
The "x component" of the function we're approximating with EIM
struct Gx : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1. + curvature*p(0);
}
};
The "y component" of the function we're approximating with EIM
struct Gy : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1. + curvature*p(0);
}
};
The "z component" of the function we're approximating with EIM
struct Gz : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1./(1. + curvature*p(0));
}
};
struct ThetaA0 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * mu.get_value("Bi");
}
};
struct AssemblyA0 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 1 || *b == 2 || *b == 3 || *b == 4 )
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_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_sidepoints = c.side_qrule->n_points();
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
break;
}
}
};
struct ThetaA1 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * mu.get_value("Bi") * mu.get_value("curvature");
}
};
struct AssemblyA1 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 1 || *b == 3 ) // y == -0.2, y == 0.2
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_fe_var[u_var]->get_phi();
const std::vector<Point>& xyz =
c.side_fe_var[u_var]->get_xyz();
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_sidepoints = c.side_qrule->n_points();
for (unsigned int qp=0; qp != n_sidepoints; qp++)
{
Real x_hat = xyz[qp](0);
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_side[qp] * x_hat * phi_side[j][qp]*phi_side[i][qp];
}
break;
}
}
};
struct ThetaA2 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return 0.2*mu.get_value("kappa") * mu.get_value("Bi") * mu.get_value("curvature");
}
};
struct AssemblyA2 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 2 || *b == 4) // x == 0.2, x == -0.2
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_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_sidepoints = c.side_qrule->n_points();
if(*b==2)
{
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
}
}
if(*b==4)
{
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
}
}
}
}
};
struct ThetaEIM : RBEIMTheta {
ThetaEIM(RBEIMEvaluation& rb_eim_eval_in, unsigned int index_in)
:
RBEIMTheta(rb_eim_eval_in, index_in)
{}
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * RBEIMTheta::evaluate(mu);
}
};
struct AssemblyEIM : RBEIMAssembly
{
AssemblyEIM(RBEIMConstruction& rb_eim_con_in,
unsigned int basis_function_index_in)
: RBEIMAssembly(rb_eim_con_in,
basis_function_index_in)
{}
virtual void interior_assembly(FEMContext &c)
{
PDE variable numbers
const unsigned int u_var = 0;
EIM variable numbers
const unsigned int Gx_var = 0;
const unsigned int Gy_var = 1;
const unsigned int Gz_var = 2;
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 std::vector<Point>& qpoints =
c.element_fe_var[u_var]->get_xyz();
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();
std::vector<Number> eim_values_Gx;
evaluate_basis_function(Gx_var,
*c.elem,
qpoints,
eim_values_Gx);
std::vector<Number> eim_values_Gy;
evaluate_basis_function(Gy_var,
*c.elem,
qpoints,
eim_values_Gy);
std::vector<Number> eim_values_Gz;
evaluate_basis_function(Gz_var,
*c.elem,
qpoints,
eim_values_Gz);
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] * ( eim_values_Gx[qp]*dphi[i][qp](0)*dphi[j][qp](0) +
eim_values_Gy[qp]*dphi[i][qp](1)*dphi[j][qp](1) +
eim_values_Gz[qp]*dphi[i][qp](2)*dphi[j][qp](2) );
}
}
}
};
struct ThetaF0 : RBTheta { virtual Number evaluate(const RBParameters& ) { return 1.; } };
struct AssemblyF0 : 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();
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 ThetaF1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("curvature"); } };
struct AssemblyF1 : 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<Point>& xyz =
c.element_fe_var[u_var]->get_xyz();
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++)
{
Real x_hat = xyz[qp](0);
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*x_hat*phi[i][qp] );
}
}
};
struct Ex6InnerProduct : 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 Ex6EIMInnerProduct : ElemAssembly
{
Use the L2 norm to find the best fit
virtual void interior_assembly(FEMContext &c)
{
const unsigned int Gx_var = 0;
const unsigned int Gy_var = 1;
const unsigned int Gz_var = 2;
const std::vector<Real> &JxW =
c.element_fe_var[Gx_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[Gx_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[Gx_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
DenseSubMatrix<Number>& Kxx = c.get_elem_jacobian(Gx_var,Gx_var);
DenseSubMatrix<Number>& Kyy = c.get_elem_jacobian(Gy_var,Gy_var);
DenseSubMatrix<Number>& Kzz = c.get_elem_jacobian(Gz_var,Gz_var);
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++)
{
Kxx(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
Kyy(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
Kzz(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
}
}
};
Define an RBThetaExpansion class for this PDE
The A terms depend on EIM, so we deal with them later
struct Ex6ThetaExpansion : RBThetaExpansion
{
/**
* Constructor.
*/
Ex6ThetaExpansion()
{
attach_A_theta(&theta_a0);
attach_A_theta(&theta_a1);
attach_A_theta(&theta_a2);
attach_F_theta(&theta_f0); // Attach the rhs theta
attach_F_theta(&theta_f1);
}
The RBTheta member variables
ThetaA0 theta_a0;
ThetaA1 theta_a1;
ThetaA2 theta_a2;
ThetaF0 theta_f0;
ThetaF1 theta_f1;
};
Define an RBAssemblyExpansion class for this PDE
The A terms depend on EIM, so we deal with them later
struct Ex6AssemblyExpansion : RBAssemblyExpansion
{
/**
* Constructor.
*/
Ex6AssemblyExpansion(RBConstruction& rb_con)
{
Point to the RBConstruction object
assembly_a0.rb_con = &rb_con;
assembly_a1.rb_con = &rb_con;
assembly_a2.rb_con = &rb_con;
attach_A_assembly(&assembly_a0);
attach_A_assembly(&assembly_a1);
attach_A_assembly(&assembly_a2);
attach_F_assembly(&assembly_f0); // Attach the rhs assembly
attach_F_assembly(&assembly_f1);
}
The ElemAssembly objects
AssemblyA0 assembly_a0;
AssemblyA1 assembly_a1;
AssemblyA2 assembly_a2;
AssemblyF0 assembly_f0;
AssemblyF1 assembly_f1;
};
#endif
The source file eim_classes.h with comments:
#ifndef __eim_classes_h__
#define __eim_classes_h__
local includes
#include "libmesh/rb_eim_construction.h"
#include "libmesh/rb_eim_evaluation.h"
#include "assembly.h"
A simple subclass of RBEIMEvaluation. Overload
evaluate_parametrized_function to define the
function that we "empirically" interpolate.
class SimpleEIMEvaluation : public RBEIMEvaluation
{
public:
SimpleEIMEvaluation()
{
attach_parametrized_function(&g_x);
attach_parametrized_function(&g_y);
attach_parametrized_function(&g_z);
}
/**
* Build a ThetaEIM rather than an RBEIMTheta.
*/
virtual AutoPtr<RBTheta> build_eim_theta(unsigned int index)
{
return AutoPtr<RBTheta>(new ThetaEIM(*this, index));
}
/**
* Parametrized functions that we approximate with EIM
*/
Gx g_x;
Gy g_y;
Gz g_z;
};
A simple subclass of RBEIMConstruction.
class SimpleEIMConstruction : public RBEIMConstruction
{
public:
/**
* Constructor.
*/
SimpleEIMConstruction (EquationSystems& es,
const std::string& name_in,
const unsigned int number_in)
: Parent(es, name_in, number_in)
{
}
/**
* The type of the parent.
*/
typedef RBEIMConstruction Parent;
/**
* Provide an implementation of build_eim_assembly
*/
virtual AutoPtr<ElemAssembly> build_eim_assembly(unsigned int index)
{
return AutoPtr<ElemAssembly>(new AssemblyEIM(*this, index));
}
/**
* Initialize data structures.
*/
virtual void init_data()
{
Gx_var = this->add_variable ("x_comp_of_G", FIRST);
Gy_var = this->add_variable ("y_comp_of_G", FIRST);
Gz_var = this->add_variable ("z_comp_of_G", FIRST);
Parent::init_data();
set_inner_product_assembly(eim_ip);
}
/**
* Variable numbers.
*/
unsigned int Gx_var;
unsigned int Gy_var;
unsigned int Gz_var;
/**
* Inner product assembly object
*/
Ex6EIMInnerProduct eim_ip;
};
#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(ex6_theta_expansion);
}
/**
* Return a "dummy" lower bound for the coercivity constant.
* To do this rigorously we should use the SCM classes.
*/
virtual Real get_stability_lower_bound() { return 1.; }
/**
* 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.
*/
Ex6ThetaExpansion ex6_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),
ex6_assembly_expansion(*this),
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(5);
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(ex6_assembly_expansion);
We need to define an inner product matrix for this problem
set_inner_product_assembly(ex6_ip);
}
/**
* 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.
*/
Ex6AssemblyExpansion ex6_assembly_expansion;
/**
* The inner product assembly object
*/
Ex6InnerProduct ex6_ip;
/**
* The object that defines which degrees of freedom are on a Dirichlet boundary.
*/
AutoPtr<DirichletBoundary> dirichlet_bc;
};
#endif
The source file reduced_basis_ex6.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 6 - Heat transfer on a curved domain in 3D
In this example we consider heat transfer modeled by a Poisson equation with Robin boundary condition: -kappa \Laplacian u = 1, on \Omega -kappa du\dn = kappa Bi u, on \partial\Omega_Biot, u = 0 on \partial\Omega_Dirichlet,
We consider a reference domain \Omega_hat = [-0.2,0.2]x[-0.2,0.2]x[0,3], and the physical domain is then obtain via the parametrized mapping: x = -1/mu + (1/mu+x_hat)*cos(mu*z_hat) y = y_hat z = (1/mu+x_hat)*sin(mu*z_hat) for (x_hat,y_hat,z_hat) \in \Omega_hat. (Here "hats" denotes reference domain.) Also, the "reference Dirichlet boundaries" are [-0.2,0.2]x[-0.2,0.2]x{0} and [-0.2,0.2]x[-0.2,0.2]x{3}, and the remaining boundaries are the "Biot" boundaries.
Then, after putting the PDE into weak form and mapping it to the reference domain, we obtain: \kappa \int_\Omega_hat [ (1+mu*x_hat) v_x w_x + (1+mu*x_hat) v_y w_y + 1/(1+mu*x_hat) v_z w_z ] + \kappa Bi \int_\partial\Omega_hat_Biot1 (1-0.2mu) u v + \kappa Bi \int_\partial\Omega_hat_Biot2 (1+mu x_hat) u v + \kappa Bi \int_\partial\Omega_hat_Biot3 (1+0.2mu) u v = \int_\Omega_hat (1+mu x_hat) v where \partial\Omega_hat_Biot1 = [-0.2] x [-0.2,0.2] x [0,3] \partial\Omega_hat_Biot2 = [-0.2,0.2] x {-0.2} x [0,3] \UNION [-0.2,0.2] x {0.2} x [0,3] \partial\Omega_hat_Biot3 = [0.2] x [-0.2,0.2] x [0,3]
The term \kappa \int_\Omega_hat 1/(1+mu*x_hat) v_z w_z is "non-affine" (in the Reduced Basis sense), since we can't express it in the form \sum theta_q(kappa,mu) a(v,w). As a result, (as in reduced_basis_ex4) we must employ the Empirical Interpolation Method (EIM) in order to apply the Reduced Basis method here.
The approach we use is to construct an EIM approximation, G_EIM, to the vector-valued function G(x_hat,y_hat;mu) = (1 + mu*x_hat, 1 + mu*x_hat, 1/(1+mu*x_hat)) and then we express the "volumetric integral part" of the left-hand side operator as a(v,w;mu) = \int_\hat\Omega G_EIM(x_hat,y_hat;mu) \dot (v_x w_x, v_y w_y, v_z w_z). (We actually only need EIM for the third component of G_EIM, but it's helpful to demonstrate "vector-valued" EIM here.)
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 "eim_classes.h"
#include "assembly.h"
Bring in everything from the libMesh namespace
using namespace libMesh;
Define a function to scale the mesh according to the parameter.
void transform_mesh_and_plot(EquationSystems& es, Real curvature, const std::string& filename);
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
This is a 3D example
libmesh_example_assert(3 == LIBMESH_DIM, "3D support");
Parse the input file using GetPot
std::string eim_parameters = "eim.in";
std::string rb_parameters = "rb.in";
std::string main_parameters = "reduced_basis_ex6.in";
GetPot infile(main_parameters);
unsigned int n_elem_xy = infile("n_elem_xy", 1);
unsigned int n_elem_z = infile("n_elem_z", 1);
Do we write the RB basis functions to disk?
bool store_basis_functions = infile("store_basis_functions", true);
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;
MeshTools::Generation::build_cube (mesh,
n_elem_xy, n_elem_xy, n_elem_z,
-0.2, 0.2,
-0.2, 0.2,
0., 3.,
HEX8);
Create an equation systems object.
EquationSystems equation_systems (mesh);
SimpleEIMConstruction & eim_construction =
equation_systems.add_system<SimpleEIMConstruction> ("EIM");
SimpleRBConstruction & rb_construction =
equation_systems.add_system<SimpleRBConstruction> ("RB");
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();
Initialize the standard RBEvaluation object
SimpleRBEvaluation rb_eval;
Initialize the EIM RBEvaluation object
SimpleEIMEvaluation eim_rb_eval;
Set the rb_eval objects for the RBConstructions
eim_construction.set_rb_evaluation(eim_rb_eval);
rb_construction.set_rb_evaluation(rb_eval);
if(!online_mode) // Perform the Offline stage of the RB method
{
Read data from input file and print state
eim_construction.process_parameters_file(eim_parameters);
eim_construction.print_info();
Perform the EIM Greedy and write out the data
eim_construction.initialize_rb_construction();
eim_construction.train_reduced_basis();
eim_construction.get_rb_evaluation().write_offline_data_to_files("eim_data");
Read data from input file and print state
rb_construction.process_parameters_file(rb_parameters);
attach the EIM theta objects to the RBEvaluation
eim_rb_eval.initialize_eim_theta_objects();
rb_eval.get_rb_theta_expansion().attach_multiple_A_theta(eim_rb_eval.get_eim_theta_objects());
attach the EIM assembly objects to the RBConstruction
eim_construction.initialize_eim_assembly_objects();
rb_construction.get_rb_assembly_expansion().attach_multiple_A_assembly(eim_construction.get_eim_assembly_objects());
Print out the state of rb_construction now that the EIM objects have been attached
rb_construction.print_info();
Need to initialize _after_ EIM greedy so that
the system knows how many affine terms there are
rb_construction.initialize_rb_construction();
rb_construction.train_reduced_basis();
rb_construction.get_rb_evaluation().write_offline_data_to_files("rb_data");
Write out the basis functions, if requested
if(store_basis_functions)
{
Write out the basis functions
eim_construction.get_rb_evaluation().write_out_basis_functions(eim_construction,"eim_data");
rb_construction.get_rb_evaluation().write_out_basis_functions(rb_construction,"rb_data");
}
}
else // Perform the Online stage of the RB method
{
eim_rb_eval.read_offline_data_from_files("eim_data");
attach the EIM theta objects to rb_eval objects
eim_rb_eval.initialize_eim_theta_objects();
rb_eval.get_rb_theta_expansion().attach_multiple_A_theta(eim_rb_eval.get_eim_theta_objects());
Read in the offline data for rb_eval
rb_eval.read_offline_data_from_files("rb_data");
Get the parameters at which we will do a reduced basis solve
Real online_curvature = infile("online_curvature", 0.);
Real online_Bi = infile("online_Bi", 0.);
Real online_kappa = infile("online_kappa", 0.);
RBParameters online_mu;
online_mu.set_value("curvature", online_curvature);
online_mu.set_value("Bi", online_Bi);
online_mu.set_value("kappa", online_kappa);
rb_eval.set_parameters(online_mu);
rb_eval.print_parameters();
rb_eval.rb_solve( rb_eval.get_n_basis_functions() );
plot the solution, if requested
if(store_basis_functions)
{
read in the data from files
eim_rb_eval.read_in_basis_functions(eim_construction,"eim_data");
rb_eval.read_in_basis_functions(rb_construction,"rb_data");
eim_construction.load_rb_solution();
rb_construction.load_rb_solution();
transform_mesh_and_plot(equation_systems,online_curvature,"RB_sol.e");
}
}
return 0;
}
void transform_mesh_and_plot(EquationSystems& es, Real curvature, const std::string& filename)
{
Loop over the mesh nodes and move them!
MeshBase& mesh = es.get_mesh();
MeshBase::node_iterator node_it = mesh.nodes_begin();
const MeshBase::node_iterator node_end = mesh.nodes_end();
for( ; node_it != node_end; node_it++)
{
Node* node = *node_it;
Real x = (*node)(0);
Real z = (*node)(2);
(*node)(0) = -1./curvature + (1./curvature + x)*cos(curvature*z);
(*node)(2) = (1./curvature + x)*sin(curvature*z);
}
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems(filename, es);
#endif
}
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/boundary_info.h"
#include "libmesh/rb_theta.h"
#include "libmesh/rb_assembly_expansion.h"
#include "libmesh/rb_parametrized_function.h"
#include "libmesh/rb_eim_construction.h"
#include "libmesh/rb_eim_theta.h"
using libMesh::boundary_id_type;
using libMesh::DenseSubMatrix;
using libMesh::ElemAssembly;
using libMesh::FEInterface;
using libMesh::FEMContext;
using libMesh::Number;
using libMesh::Point;
using libMesh::RBAssemblyExpansion;
using libMesh::RBConstruction;
using libMesh::RBParameters;
using libMesh::RBParametrizedFunction;
using libMesh::RBTheta;
using libMesh::RBThetaExpansion;
using libMesh::RBEIMAssembly;
using libMesh::RBEIMConstruction;
using libMesh::RBEIMEvaluation;
using libMesh::RBEIMTheta;
using libMesh::Real;
using libMesh::RealGradient;
struct ElemAssemblyWithConstruction : ElemAssembly
{
RBConstruction* rb_con;
};
struct Gx : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1. + curvature*p(0);
}
};
struct Gy : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1. + curvature*p(0);
}
};
struct Gz : public RBParametrizedFunction
{
virtual Number evaluate(const RBParameters& mu,
const Point& p)
{
Real curvature = mu.get_value("curvature");
return 1./(1. + curvature*p(0));
}
};
struct ThetaA0 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * mu.get_value("Bi");
}
};
struct AssemblyA0 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 1 || *b == 2 || *b == 3 || *b == 4 )
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_fe_var[u_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_sidepoints = c.side_qrule->n_points();
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
break;
}
}
};
struct ThetaA1 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * mu.get_value("Bi") * mu.get_value("curvature");
}
};
struct AssemblyA1 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 1 || *b == 3 ) // y == -0.2, y == 0.2
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_fe_var[u_var]->get_phi();
const std::vector<Point>& xyz =
c.side_fe_var[u_var]->get_xyz();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_sidepoints = c.side_qrule->n_points();
for (unsigned int qp=0; qp != n_sidepoints; qp++)
{
Real x_hat = xyz[qp](0);
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_side[qp] * x_hat * phi_side[j][qp]*phi_side[i][qp];
}
break;
}
}
};
struct ThetaA2 : RBTheta {
virtual Number evaluate(const RBParameters& mu)
{
return 0.2*mu.get_value("kappa") * mu.get_value("Bi") * mu.get_value("curvature");
}
};
struct AssemblyA2 : ElemAssemblyWithConstruction
{
virtual void boundary_assembly(FEMContext &c)
{
const std::vector<boundary_id_type> bc_ids =
rb_con->get_mesh().boundary_info->boundary_ids (c.elem,c.side);
for (std::vector<boundary_id_type>::const_iterator b =
bc_ids.begin(); b != bc_ids.end(); ++b)
if( *b == 2 || *b == 4) // x == 0.2, x == -0.2
{
const unsigned int u_var = 0;
const std::vector<Real> &JxW_side =
c.side_fe_var[u_var]->get_JxW();
const std::vector<std::vector<Real> >& phi_side =
c.side_fe_var[u_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_sidepoints = c.side_qrule->n_points();
if(*b==2)
{
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
}
}
if(*b==4)
{
for (unsigned int qp=0; qp != n_sidepoints; 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_side[qp] * phi_side[j][qp]*phi_side[i][qp];
}
}
}
}
};
struct ThetaEIM : RBEIMTheta {
ThetaEIM(RBEIMEvaluation& rb_eim_eval_in, unsigned int index_in)
:
RBEIMTheta(rb_eim_eval_in, index_in)
{}
virtual Number evaluate(const RBParameters& mu)
{
return mu.get_value("kappa") * RBEIMTheta::evaluate(mu);
}
};
struct AssemblyEIM : RBEIMAssembly
{
AssemblyEIM(RBEIMConstruction& rb_eim_con_in,
unsigned int basis_function_index_in)
: RBEIMAssembly(rb_eim_con_in,
basis_function_index_in)
{}
virtual void interior_assembly(FEMContext &c)
{
const unsigned int u_var = 0;
const unsigned int Gx_var = 0;
const unsigned int Gy_var = 1;
const unsigned int Gz_var = 2;
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 std::vector<Point>& qpoints =
c.element_fe_var[u_var]->get_xyz();
const unsigned int n_u_dofs = c.dof_indices_var[u_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
std::vector<Number> eim_values_Gx;
evaluate_basis_function(Gx_var,
*c.elem,
qpoints,
eim_values_Gx);
std::vector<Number> eim_values_Gy;
evaluate_basis_function(Gy_var,
*c.elem,
qpoints,
eim_values_Gy);
std::vector<Number> eim_values_Gz;
evaluate_basis_function(Gz_var,
*c.elem,
qpoints,
eim_values_Gz);
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] * ( eim_values_Gx[qp]*dphi[i][qp](0)*dphi[j][qp](0) +
eim_values_Gy[qp]*dphi[i][qp](1)*dphi[j][qp](1) +
eim_values_Gz[qp]*dphi[i][qp](2)*dphi[j][qp](2) );
}
}
}
};
struct ThetaF0 : RBTheta { virtual Number evaluate(const RBParameters& ) { return 1.; } };
struct AssemblyF0 : 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 ThetaF1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("curvature"); } };
struct AssemblyF1 : 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<Point>& xyz =
c.element_fe_var[u_var]->get_xyz();
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++)
{
Real x_hat = xyz[qp](0);
for (unsigned int i=0; i != n_u_dofs; i++)
c.get_elem_residual()(i) += JxW[qp] * ( 1.*x_hat*phi[i][qp] );
}
}
};
struct Ex6InnerProduct : 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 Ex6EIMInnerProduct : ElemAssembly
{
virtual void interior_assembly(FEMContext &c)
{
const unsigned int Gx_var = 0;
const unsigned int Gy_var = 1;
const unsigned int Gz_var = 2;
const std::vector<Real> &JxW =
c.element_fe_var[Gx_var]->get_JxW();
const std::vector<std::vector<Real> >& phi =
c.element_fe_var[Gx_var]->get_phi();
const unsigned int n_u_dofs = c.dof_indices_var[Gx_var].size();
unsigned int n_qpoints = (c.get_element_qrule())->n_points();
DenseSubMatrix<Number>& Kxx = c.get_elem_jacobian(Gx_var,Gx_var);
DenseSubMatrix<Number>& Kyy = c.get_elem_jacobian(Gy_var,Gy_var);
DenseSubMatrix<Number>& Kzz = c.get_elem_jacobian(Gz_var,Gz_var);
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++)
{
Kxx(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
Kyy(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
Kzz(i,j) += JxW[qp] * phi[j][qp]*phi[i][qp];
}
}
};
struct Ex6ThetaExpansion : RBThetaExpansion
{
/**
* Constructor.
*/
Ex6ThetaExpansion()
{
attach_A_theta(&theta_a0);
attach_A_theta(&theta_a1);
attach_A_theta(&theta_a2);
attach_F_theta(&theta_f0); // Attach the rhs theta
attach_F_theta(&theta_f1);
}
ThetaA0 theta_a0;
ThetaA1 theta_a1;
ThetaA2 theta_a2;
ThetaF0 theta_f0;
ThetaF1 theta_f1;
};
struct Ex6AssemblyExpansion : RBAssemblyExpansion
{
/**
* Constructor.
*/
Ex6AssemblyExpansion(RBConstruction& rb_con)
{
assembly_a0.rb_con = &rb_con;
assembly_a1.rb_con = &rb_con;
assembly_a2.rb_con = &rb_con;
attach_A_assembly(&assembly_a0);
attach_A_assembly(&assembly_a1);
attach_A_assembly(&assembly_a2);
attach_F_assembly(&assembly_f0); // Attach the rhs assembly
attach_F_assembly(&assembly_f1);
}
AssemblyA0 assembly_a0;
AssemblyA1 assembly_a1;
AssemblyA2 assembly_a2;
AssemblyF0 assembly_f0;
AssemblyF1 assembly_f1;
};
#endif
The source file eim_classes.h without comments:
#ifndef __eim_classes_h__
#define __eim_classes_h__
#include "libmesh/rb_eim_construction.h"
#include "libmesh/rb_eim_evaluation.h"
#include "assembly.h"
class SimpleEIMEvaluation : public RBEIMEvaluation
{
public:
SimpleEIMEvaluation()
{
attach_parametrized_function(&g_x);
attach_parametrized_function(&g_y);
attach_parametrized_function(&g_z);
}
/**
* Build a ThetaEIM rather than an RBEIMTheta.
*/
virtual AutoPtr<RBTheta> build_eim_theta(unsigned int index)
{
return AutoPtr<RBTheta>(new ThetaEIM(*this, index));
}
/**
* Parametrized functions that we approximate with EIM
*/
Gx g_x;
Gy g_y;
Gz g_z;
};
class SimpleEIMConstruction : public RBEIMConstruction
{
public:
/**
* Constructor.
*/
SimpleEIMConstruction (EquationSystems& es,
const std::string& name_in,
const unsigned int number_in)
: Parent(es, name_in, number_in)
{
}
/**
* The type of the parent.
*/
typedef RBEIMConstruction Parent;
/**
* Provide an implementation of build_eim_assembly
*/
virtual AutoPtr<ElemAssembly> build_eim_assembly(unsigned int index)
{
return AutoPtr<ElemAssembly>(new AssemblyEIM(*this, index));
}
/**
* Initialize data structures.
*/
virtual void init_data()
{
Gx_var = this->add_variable ("x_comp_of_G", FIRST);
Gy_var = this->add_variable ("y_comp_of_G", FIRST);
Gz_var = this->add_variable ("z_comp_of_G", FIRST);
Parent::init_data();
set_inner_product_assembly(eim_ip);
}
/**
* Variable numbers.
*/
unsigned int Gx_var;
unsigned int Gy_var;
unsigned int Gz_var;
/**
* Inner product assembly object
*/
Ex6EIMInnerProduct eim_ip;
};
#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(ex6_theta_expansion);
}
/**
* Return a "dummy" lower bound for the coercivity constant.
* To do this rigorously we should use the SCM classes.
*/
virtual Real get_stability_lower_bound() { return 1.; }
/**
* 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.
*/
Ex6ThetaExpansion ex6_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),
ex6_assembly_expansion(*this),
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(5);
dirichlet_bc->variables.push_back(u_var);
get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
Parent::init_data();
set_rb_assembly_expansion(ex6_assembly_expansion);
set_inner_product_assembly(ex6_ip);
}
/**
* 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.
*/
Ex6AssemblyExpansion ex6_assembly_expansion;
/**
* The inner product assembly object
*/
Ex6InnerProduct ex6_ip;
/**
* The object that defines which degrees of freedom are on a Dirichlet boundary.
*/
AutoPtr<DirichletBoundary> dirichlet_bc;
};
#endif
The source file reduced_basis_ex6.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 "eim_classes.h"
#include "assembly.h"
using namespace libMesh;
void transform_mesh_and_plot(EquationSystems& es, Real curvature, const std::string& filename);
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(3 == LIBMESH_DIM, "3D support");
std::string eim_parameters = "eim.in";
std::string rb_parameters = "rb.in";
std::string main_parameters = "reduced_basis_ex6.in";
GetPot infile(main_parameters);
unsigned int n_elem_xy = infile("n_elem_xy", 1);
unsigned int n_elem_z = infile("n_elem_z", 1);
bool store_basis_functions = infile("store_basis_functions", true);
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;
MeshTools::Generation::build_cube (mesh,
n_elem_xy, n_elem_xy, n_elem_z,
-0.2, 0.2,
-0.2, 0.2,
0., 3.,
HEX8);
EquationSystems equation_systems (mesh);
SimpleEIMConstruction & eim_construction =
equation_systems.add_system<SimpleEIMConstruction> ("EIM");
SimpleRBConstruction & rb_construction =
equation_systems.add_system<SimpleRBConstruction> ("RB");
equation_systems.init ();
equation_systems.print_info();
mesh.print_info();
SimpleRBEvaluation rb_eval;
SimpleEIMEvaluation eim_rb_eval;
eim_construction.set_rb_evaluation(eim_rb_eval);
rb_construction.set_rb_evaluation(rb_eval);
if(!online_mode) // Perform the Offline stage of the RB method
{
eim_construction.process_parameters_file(eim_parameters);
eim_construction.print_info();
eim_construction.initialize_rb_construction();
eim_construction.train_reduced_basis();
eim_construction.get_rb_evaluation().write_offline_data_to_files("eim_data");
rb_construction.process_parameters_file(rb_parameters);
eim_rb_eval.initialize_eim_theta_objects();
rb_eval.get_rb_theta_expansion().attach_multiple_A_theta(eim_rb_eval.get_eim_theta_objects());
eim_construction.initialize_eim_assembly_objects();
rb_construction.get_rb_assembly_expansion().attach_multiple_A_assembly(eim_construction.get_eim_assembly_objects());
rb_construction.print_info();
rb_construction.initialize_rb_construction();
rb_construction.train_reduced_basis();
rb_construction.get_rb_evaluation().write_offline_data_to_files("rb_data");
if(store_basis_functions)
{
eim_construction.get_rb_evaluation().write_out_basis_functions(eim_construction,"eim_data");
rb_construction.get_rb_evaluation().write_out_basis_functions(rb_construction,"rb_data");
}
}
else // Perform the Online stage of the RB method
{
eim_rb_eval.read_offline_data_from_files("eim_data");
eim_rb_eval.initialize_eim_theta_objects();
rb_eval.get_rb_theta_expansion().attach_multiple_A_theta(eim_rb_eval.get_eim_theta_objects());
rb_eval.read_offline_data_from_files("rb_data");
Real online_curvature = infile("online_curvature", 0.);
Real online_Bi = infile("online_Bi", 0.);
Real online_kappa = infile("online_kappa", 0.);
RBParameters online_mu;
online_mu.set_value("curvature", online_curvature);
online_mu.set_value("Bi", online_Bi);
online_mu.set_value("kappa", online_kappa);
rb_eval.set_parameters(online_mu);
rb_eval.print_parameters();
rb_eval.rb_solve( rb_eval.get_n_basis_functions() );
if(store_basis_functions)
{
eim_rb_eval.read_in_basis_functions(eim_construction,"eim_data");
rb_eval.read_in_basis_functions(rb_construction,"rb_data");
eim_construction.load_rb_solution();
rb_construction.load_rb_solution();
transform_mesh_and_plot(equation_systems,online_curvature,"RB_sol.e");
}
}
return 0;
}
void transform_mesh_and_plot(EquationSystems& es, Real curvature, const std::string& filename)
{
MeshBase& mesh = es.get_mesh();
MeshBase::node_iterator node_it = mesh.nodes_begin();
const MeshBase::node_iterator node_end = mesh.nodes_end();
for( ; node_it != node_end; node_it++)
{
Node* node = *node_it;
Real x = (*node)(0);
Real z = (*node)(2);
(*node)(0) = -1./curvature + (1./curvature + x)*cos(curvature*z);
(*node)(2) = (1./curvature + x)*sin(curvature*z);
}
#ifdef LIBMESH_HAVE_EXODUS_API
ExodusII_IO(mesh).write_equation_systems(filename, es);
#endif
}
The console output of the program:
***************************************************************
* Running Example reduced_basis_ex6:
* 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()=2
System #0, "EIM"
Type "RBConstruction"
Variables={ "x_comp_of_G" "y_comp_of_G" "z_comp_of_G" }
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=18513
n_local_dofs()=9438
n_constrained_dofs()=0
n_local_constrained_dofs()=0
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 23.3593
Average Off-Processor Bandwidth <= 0.311457
Maximum On-Processor Bandwidth <= 27
Maximum Off-Processor Bandwidth <= 9
DofMap Constraints
Number of DoF Constraints = 0
Number of Node Constraints = 0
System #1, "RB"
Type "RBConstruction"
Variables="u"
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=6171
n_local_dofs()=3146
n_constrained_dofs()=242
n_local_constrained_dofs()=121
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 23.3593
Average Off-Processor Bandwidth <= 0.311457
Maximum On-Processor Bandwidth <= 27
Maximum Off-Processor Bandwidth <= 9
DofMap Constraints
Number of DoF Constraints = 242
Average DoF Constraint Length= 0
Number of Node Constraints = 0
Mesh Information:
mesh_dimension()=3
spatial_dimension()=3
n_nodes()=6171
n_local_nodes()=3146
n_elem()=5000
n_local_elem()=2500
n_active_elem()=5000
n_subdomains()=1
n_partitions()=2
n_processors()=2
n_threads()=1
processor_id()=0
Initializing training parameters with deterministic training set...
Parameter curvature: log scaling = 0
RBConstruction parameters:
system name: EIM
constrained_problem: 0
Nmax: 20
Basis training error tolerance: 0.001
Aq operators attached: 0
Fq functions attached: 0
n_outputs: 0
Number of parameters: 1
Parameter curvature: Min = 0.1, Max = 1.0472, value = 0.5
n_training_samples: 25
single-matrix mode? 0
reuse preconditioner? 1
use a relative error bound in greedy? 1
write out data during basis training? 0
quiet mode? 1
RBEIMConstruction parameters:
best fit type: projection
Initializing parametrized functions in training set...
Completed solve for training sample 1 of 25
Completed solve for training sample 2 of 25
Completed solve for training sample 3 of 25
Completed solve for training sample 4 of 25
Completed solve for training sample 5 of 25
Completed solve for training sample 6 of 25
Completed solve for training sample 7 of 25
Completed solve for training sample 8 of 25
Completed solve for training sample 9 of 25
Completed solve for training sample 10 of 25
Completed solve for training sample 11 of 25
Completed solve for training sample 12 of 25
Completed solve for training sample 13 of 25
Completed solve for training sample 14 of 25
Completed solve for training sample 15 of 25
Completed solve for training sample 16 of 25
Completed solve for training sample 17 of 25
Completed solve for training sample 18 of 25
Completed solve for training sample 19 of 25
Completed solve for training sample 20 of 25
Completed solve for training sample 21 of 25
Completed solve for training sample 22 of 25
Completed solve for training sample 23 of 25
Completed solve for training sample 24 of 25
Completed solve for training sample 25 of 25
Parametrized functions in training set initialized
---- Performing Greedy basis enrichment ----
---- Basis dimension: 0 ----
Performing truth solve at parameter:
curvature: 1.047200e+00
Enriching the RB space
Updating RB matrices
---- Basis dimension: 1 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.221031
Performing truth solve at parameter:
curvature: 1.000000e-01
Enriching the RB space
Updating RB matrices
---- Basis dimension: 2 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.0106172
Performing truth solve at parameter:
curvature: 6.130667e-01
Enriching the RB space
Updating RB matrices
---- Basis dimension: 3 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.000307645
Specified error tolerance reached.
Perform one more Greedy iteration for error bounds.
Performing truth solve at parameter:
curvature: 2.973333e-01
Enriching the RB space
Updating RB matrices
---- Basis dimension: 3 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.000307645
Extra Greedy iteration finished.
Initializing training parameters with random training set...
Parameter Bi: log scaling = 0
Parameter curvature: log scaling = 1
Parameter kappa: log scaling = 1
RBConstruction parameters:
system name: RB
constrained_problem: 0
Nmax: 15
Basis training error tolerance: 0.001
Aq operators attached: 6
Fq functions attached: 2
n_outputs: 0
Number of parameters: 3
Parameter Bi: Min = 0.001, Max = 0.01, value = 0.005
Parameter curvature: Min = 0.1, Max = 1.0472, value = 0.5
Parameter kappa: Min = 0.5, Max = 2, value = 1
n_training_samples: 1000
single-matrix mode? 0
reuse preconditioner? 1
use a relative error bound in greedy? 1
write out data during basis training? 0
quiet mode? 1
---- Performing Greedy basis enrichment ----
---- Basis dimension: 0 ----
Performing RB solves on training set
Maximum (relative) error bound is inf
Performing truth solve at parameter:
Bi: 4.064144e-03
curvature: 4.308919e-01
kappa: 1.193010e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 1 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.035047
Performing truth solve at parameter:
Bi: 2.992444e-03
curvature: 1.033701e+00
kappa: 1.746540e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 2 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.0148518
Performing truth solve at parameter:
Bi: 9.780694e-03
curvature: 9.289260e-01
kappa: 1.826589e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 3 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.00191765
Performing truth solve at parameter:
Bi: 8.246993e-03
curvature: 1.003271e-01
kappa: 1.813348e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 4 ----
Performing RB solves on training set
Maximum (relative) error bound is 0.00106327
Performing truth solve at parameter:
Bi: 1.134944e-03
curvature: 1.039075e-01
kappa: 1.915826e+00
Enriching the RB space
Updating RB matrices
Updating RB residual terms
---- Basis dimension: 5 ----
Performing RB solves on training set
Maximum (relative) error bound is 4.45341e-05
Specified error tolerance reached.
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/reduced_basis/reduced_basis_ex6/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:50:59 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): 6.504e+01 1.00000 6.504e+01
Objects: 6.780e+02 1.00000 6.780e+02
Flops: 2.374e+09 1.03767 2.331e+09 4.661e+09
Flops/sec: 3.650e+07 1.03767 3.584e+07 7.167e+07
MPI Messages: 2.388e+03 1.00000 2.388e+03 4.776e+03
MPI Message Lengths: 3.772e+06 1.00000 1.580e+03 7.544e+06
MPI Reductions: 7.975e+03 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: 6.5040e+01 100.0% 4.6614e+09 100.0% 4.776e+03 100.0% 1.580e+03 100.0% 7.974e+03 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 1002 1.0 8.4479e-03 1.8 9.25e+06 1.0 0.0e+00 0.0e+00 1.0e+03 0 0 0 0 13 0 0 0 0 13 2147
VecMDot 1026 1.0 4.7181e-02 1.5 8.55e+07 1.0 0.0e+00 0.0e+00 1.0e+03 0 4 0 0 13 0 4 0 0 13 3556
VecNorm 1250 1.0 1.5167e-02 1.9 8.49e+06 1.0 0.0e+00 0.0e+00 1.2e+03 0 0 0 0 16 0 0 0 0 16 1099
VecScale 1137 1.0 1.9779e-03 1.2 4.06e+06 1.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 4025
VecCopy 437 1.0 2.2285e-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
VecSet 2569 1.0 5.1396e-03 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 377 1.0 8.7936e-03 1.0 5.91e+06 1.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 1318
VecMAXPY 1088 1.0 2.4819e-02 1.0 9.26e+07 1.0 0.0e+00 0.0e+00 0.0e+00 0 4 0 0 0 0 4 0 0 0 7321
VecAssemblyBegin 511 1.0 3.6898e-02 7.8 0.00e+00 0.0 5.4e+01 6.8e+03 1.5e+03 0 0 1 5 19 0 0 1 5 19 0
VecAssemblyEnd 511 1.0 3.5214e-04 1.2 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00 0 0 0 0 0 0 0 0 0 0 0
VecScatterBegin 2245 1.0 6.2368e-03 1.0 0.00e+00 0.0 4.5e+03 1.5e+03 0.0e+00 0 0 94 87 0 0 0 94 87 0 0
VecScatterEnd 2245 1.0 4.6679e-0212.5 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 1088 1.0 1.1799e-02 1.4 1.17e+07 1.0 0.0e+00 0.0e+00 1.1e+03 0 0 0 0 14 0 0 0 0 14 1942
MatMult 1088 1.0 1.4656e-01 1.4 1.79e+08 1.0 2.2e+03 1.1e+03 0.0e+00 0 8 46 32 0 0 8 46 32 0 2399
MatMultAdd 971 1.0 9.2335e-02 1.0 2.13e+08 1.0 1.9e+03 1.4e+03 0.0e+00 0 9 41 37 0 0 9 41 37 0 4523
MatSolve 1150 1.0 5.9445e-01 1.0 1.26e+09 1.0 0.0e+00 0.0e+00 0.0e+00 1 53 0 0 0 1 53 0 0 0 4163
MatLUFactorNum 12 1.0 3.0877e-01 1.0 5.16e+08 1.0 0.0e+00 0.0e+00 0.0e+00 0 22 0 0 0 0 22 0 0 0 3284
MatILUFactorSym 12 1.0 8.9547e-01 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 3.6e+01 1 0 0 0 0 1 0 0 0 0 0
MatAssemblyBegin 1154 1.0 5.4400e-01 5.9 0.00e+00 0.0 2.4e+01 2.1e+04 2.3e+03 0 0 1 7 29 0 0 1 7 29 0
MatAssemblyEnd 1154 1.0 1.0272e-01 1.0 0.00e+00 0.0 1.8e+02 2.8e+02 3.8e+02 0 0 4 1 5 0 0 4 1 5 0
MatGetRow 264264 1.0 2.7475e-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 12 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 12 1.0 3.3140e-04 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 2.4e+01 0 0 0 0 0 0 0 0 0 0 0
MatZeroEntries 41 1.0 3.8140e-03 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
MatAXPY 38 1.0 1.5989e-01 1.0 0.00e+00 0.0 1.5e+02 2.7e+02 6.1e+02 0 0 3 1 8 0 0 3 1 8 0
KSPGMRESOrthog 1026 1.0 7.0882e-02 1.3 1.71e+08 1.0 0.0e+00 0.0e+00 1.0e+03 0 7 0 0 13 0 7 0 0 13 4735
KSPSetUp 74 1.0 4.0960e-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 62 1.0 2.0017e+00 1.0 2.15e+09 1.0 2.2e+03 1.1e+03 2.2e+03 3 90 46 32 28 3 90 46 32 28 2106
PCSetUp 24 1.0 1.2062e+00 1.0 5.16e+08 1.0 0.0e+00 0.0e+00 6.4e+01 2 22 0 0 1 2 22 0 0 1 841
PCSetUpOnBlocks 62 1.0 1.2056e+00 1.0 5.16e+08 1.0 0.0e+00 0.0e+00 6.0e+01 2 22 0 0 1 2 22 0 0 1 841
PCApply 1150 1.0 6.0656e-01 1.0 1.26e+09 1.0 0.0e+00 0.0e+00 0.0e+00 1 53 0 0 0 1 53 0 0 0 4079
------------------------------------------------------------------------------------------------------------------------
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 309 309 9068464 0
Vector Scatter 54 54 55944 0
Index Set 144 144 315792 0
IS L to G Mapping 6 6 3384 0
Matrix 156 156 134583592 0
Krylov Solver 4 4 38720 0
Preconditioner 4 4 3568 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 0
Average time for MPI_Barrier(): 1.23978e-06
Average time for zero size MPI_Send(): 1.2517e-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:50:59 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=65.0525, Active time=64.7939 |
-------------------------------------------------------------------------------------------------------------------
| 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() 2 0.3989 0.199453 0.4391 0.219570 0.62 0.68 |
| build_constraint_matrix() 22500 0.1135 0.000005 0.1135 0.000005 0.18 0.18 |
| build_sparsity() 2 0.4126 0.206300 0.9299 0.464946 0.64 1.44 |
| cnstrn_elem_mat_vec() 22500 0.0437 0.000002 0.0437 0.000002 0.07 0.07 |
| create_dof_constraints() 2 0.0655 0.032731 0.3236 0.161778 0.10 0.50 |
| distribute_dofs() 2 0.1382 0.069094 0.3780 0.188979 0.21 0.58 |
| dof_indices() 382924 27.4234 0.000072 27.4234 0.000072 42.32 42.32 |
| prepare_send_list() 2 0.0011 0.000575 0.0011 0.000575 0.00 0.00 |
| reinit() 2 0.2364 0.118191 0.2364 0.118191 0.36 0.36 |
| |
| FE |
| compute_shape_functions() 217000 14.8113 0.000068 14.8113 0.000068 22.86 22.86 |
| init_shape_functions() 22078 0.6926 0.000031 0.6926 0.000031 1.07 1.07 |
| inverse_map() 180030 2.1291 0.000012 2.1291 0.000012 3.29 3.29 |
| |
| FEMap |
| compute_affine_map() 217000 2.3921 0.000011 2.3921 0.000011 3.69 3.69 |
| compute_face_map() 22000 0.2490 0.000011 0.2490 0.000011 0.38 0.38 |
| init_face_shape_functions() 20 0.0005 0.000027 0.0005 0.000027 0.00 0.00 |
| init_reference_to_physical_map() 22078 0.6356 0.000029 0.6356 0.000029 0.98 0.98 |
| |
| Mesh |
| find_neighbors() 1 0.0938 0.093783 0.0938 0.093835 0.14 0.14 |
| renumber_nodes_and_elem() 2 0.0038 0.001895 0.0038 0.001895 0.01 0.01 |
| |
| MeshCommunication |
| assign_global_indices() 2 0.2847 0.142371 0.2885 0.144236 0.44 0.45 |
| compute_hilbert_indices() 2 0.0464 0.023209 0.0464 0.023209 0.07 0.07 |
| find_global_indices() 2 0.0162 0.008117 0.0674 0.033684 0.03 0.10 |
| parallel_sort() 2 0.0040 0.002002 0.0042 0.002115 0.01 0.01 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0165 0.016540 0.0165 0.016540 0.03 0.03 |
| |
| MetisPartitioner |
| partition() 1 0.1198 0.119848 0.1528 0.152842 0.18 0.24 |
| |
| Parallel |
| allgather() 22 0.0043 0.000195 0.0043 0.000196 0.01 0.01 |
| barrier() 2 0.0000 0.000013 0.0000 0.000013 0.00 0.00 |
| broadcast() 47 0.0002 0.000005 0.0002 0.000005 0.00 0.00 |
| max(bool) 10 0.0000 0.000002 0.0000 0.000002 0.00 0.00 |
| max(scalar) 206 0.0017 0.000008 0.0017 0.000008 0.00 0.00 |
| max(vector) 45 0.0003 0.000006 0.0006 0.000014 0.00 0.00 |
| maxloc(scalar) 14 0.0410 0.002927 0.0410 0.002927 0.06 0.06 |
| min(bool) 220 0.0005 0.000002 0.0005 0.000002 0.00 0.00 |
| min(scalar) 180 0.0073 0.000041 0.0073 0.000041 0.01 0.01 |
| min(vector) 45 0.0004 0.000008 0.0007 0.000017 0.00 0.00 |
| probe() 58 0.0016 0.000027 0.0016 0.000027 0.00 0.00 |
| receive() 42 0.0008 0.000020 0.0024 0.000057 0.00 0.00 |
| send() 34 0.0002 0.000005 0.0002 0.000005 0.00 0.00 |
| send_receive() 38 0.0010 0.000026 0.0031 0.000082 0.00 0.00 |
| sum() 46 0.0008 0.000018 0.0009 0.000020 0.00 0.00 |
| |
| Parallel::Request |
| wait() 34 0.0001 0.000003 0.0001 0.000003 0.00 0.00 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0054 0.005379 0.0075 0.007497 0.01 0.01 |
| set_parent_processor_ids() 1 0.0035 0.003454 0.0035 0.003454 0.01 0.01 |
| |
| PetscLinearSolver |
| solve() 62 2.5177 0.040609 2.5177 0.040609 3.89 3.89 |
| |
| PointLocatorTree |
| init(no master) 1 0.0390 0.039044 0.0391 0.039131 0.06 0.06 |
| operator() 15 0.0018 0.000121 0.0021 0.000139 0.00 0.00 |
| |
| RBConstruction |
| add_scaled_matrix_and_vector() 10 3.0056 0.300556 17.1904 1.719036 4.64 26.53 |
| clear() 3 0.0014 0.000459 0.0014 0.000459 0.00 0.00 |
| compute_Fq_representor_innerprods() 2 0.0160 0.008013 0.1358 0.067897 0.02 0.21 |
| compute_max_error_bound() 10 0.0435 0.004351 2.5751 0.257510 0.07 3.97 |
| enrich_RB_space() 5 0.0029 0.000571 0.0029 0.000571 0.00 0.00 |
| train_reduced_basis() 2 0.0020 0.000980 9.4049 4.702470 0.00 14.52 |
| truth_assembly() 5 0.1196 0.023924 0.1199 0.023974 0.18 0.19 |
| truth_solve() 5 0.0013 0.000262 0.6304 0.126075 0.00 0.97 |
| update_RB_system_matrices() 9 0.0272 0.003024 0.0272 0.003024 0.04 0.04 |
| update_residual_terms() 5 0.1013 0.020251 1.0404 0.208084 0.16 1.61 |
| |
| RBEIMConstruction |
| compute_best_fit_error() 100 0.0866 0.000866 0.0926 0.000926 0.13 0.14 |
| enrich_RB_space() 4 0.3429 0.085737 4.9870 1.246762 0.53 7.70 |
| truth_solve() 129 5.1139 0.039643 35.1052 0.272133 7.89 54.18 |
| update_RB_system_matrices() 4 0.0007 0.000182 0.0074 0.001843 0.00 0.01 |
| |
| RBEIMEvaluation |
| rb_solve() 3007 0.0342 0.000011 0.0342 0.000011 0.05 0.05 |
| write_offline_data_to_files() 1 0.0001 0.000107 0.0006 0.000646 0.00 0.00 |
| |
| RBEvaluation |
| clear() 3 0.0001 0.000034 0.0001 0.000034 0.00 0.00 |
| compute_residual_dual_norm() 3000 2.2869 0.000762 2.2869 0.000762 3.53 3.53 |
| rb_solve() 3000 0.0978 0.000033 2.4197 0.000807 0.15 3.73 |
| resize_data_structures() 2 0.0005 0.000257 0.0005 0.000257 0.00 0.00 |
| write_offline_data_to_files() 2 0.0018 0.000909 0.0018 0.000909 0.00 0.00 |
| write_out_basis_functions() 2 0.0001 0.000060 0.8403 0.420133 0.00 1.30 |
| write_out_vectors() 2 0.5509 0.275473 0.8401 0.420073 0.85 1.30 |
-------------------------------------------------------------------------------------------------------------------
| Totals: 1118590 64.7939 100.00 |
-------------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example reduced_basis_ex6:
* 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_ex6:
* 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()=2
System #0, "EIM"
Type "RBConstruction"
Variables={ "x_comp_of_G" "y_comp_of_G" "z_comp_of_G" }
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=18513
n_local_dofs()=9438
n_constrained_dofs()=0
n_local_constrained_dofs()=0
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 23.3593
Average Off-Processor Bandwidth <= 0.311457
Maximum On-Processor Bandwidth <= 27
Maximum Off-Processor Bandwidth <= 9
DofMap Constraints
Number of DoF Constraints = 0
Number of Node Constraints = 0
System #1, "RB"
Type "RBConstruction"
Variables="u"
Finite Element Types="LAGRANGE", "JACOBI_20_00"
Infinite Element Mapping="CARTESIAN"
Approximation Orders="FIRST", "THIRD"
n_dofs()=6171
n_local_dofs()=3146
n_constrained_dofs()=242
n_local_constrained_dofs()=121
n_vectors()=1
n_matrices()=1
DofMap Sparsity
Average On-Processor Bandwidth <= 23.3593
Average Off-Processor Bandwidth <= 0.311457
Maximum On-Processor Bandwidth <= 27
Maximum Off-Processor Bandwidth <= 9
DofMap Constraints
Number of DoF Constraints = 242
Average DoF Constraint Length= 0
Number of Node Constraints = 0
Mesh Information:
mesh_dimension()=3
spatial_dimension()=3
n_nodes()=6171
n_local_nodes()=3146
n_elem()=5000
n_local_elem()=2500
n_active_elem()=5000
n_subdomains()=1
n_partitions()=2
n_processors()=2
n_threads()=1
processor_id()=0
Bi: 5.000000e-03
curvature: 1.047200e+00
kappa: 1.300000e+00
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/reduced_basis/reduced_basis_ex6/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:51:04 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): 4.409e+00 1.00000 4.409e+00
Objects: 3.300e+01 1.00000 3.300e+01
Flops: 8.809e+04 1.04000 8.639e+04 1.728e+05
Flops/sec: 1.998e+04 1.04000 1.959e+04 3.919e+04
MPI Messages: 6.000e+00 1.00000 6.000e+00 1.200e+01
MPI Message Lengths: 8.720e+03 1.00000 1.453e+03 1.744e+04
MPI Reductions: 7.300e+01 1.00000
Flop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)
e.g., VecAXPY() for real vectors of length N --> 2N flops
and VecAXPY() for complex vectors of length N --> 8N flops
Summary of Stages: ----- Time ------ ----- Flops ----- --- Messages --- -- Message Lengths -- -- Reductions --
Avg %Total Avg %Total counts %Total Avg %Total counts %Total
0: Main Stage: 4.4091e+00 100.0% 1.7279e+05 100.0% 1.200e+01 100.0% 1.453e+03 100.0% 7.200e+01 98.6%
------------------------------------------------------------------------------------------------------------------------
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 2 1.0 2.0027e-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
VecSet 18 1.0 5.0545e-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 8 1.0 7.0095e-05 1.2 8.81e+04 1.0 0.0e+00 0.0e+00 0.0e+00 0100 0 0 0 0100 0 0 0 2465
VecAssemblyBegin 10 1.0 1.4758e-04 1.1 0.00e+00 0.0 0.0e+00 0.0e+00 3.0e+01 0 0 0 0 41 0 0 0 0 42 0
VecAssemblyEnd 10 1.0 3.8147e-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 2 1.0 7.3195e-05 1.0 0.00e+00 0.0 4.0e+00 2.9e+03 0.0e+00 0 0 33 67 0 0 0 33 67 0 0
VecScatterEnd 2 1.0 7.8678e-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
MatZeroEntries 4 1.0 7.7105e-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
------------------------------------------------------------------------------------------------------------------------
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 18 18 785712 0
Vector Scatter 2 2 2072 0
Index Set 4 4 4896 0
IS L to G Mapping 2 2 1128 0
Matrix 6 6 3915088 0
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 9.53674e-08
Average time for MPI_Barrier(): 8.10623e-07
Average time for zero size MPI_Send(): 1.5974e-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:51:04 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=4.41784, Active time=4.33897 |
--------------------------------------------------------------------------------------------------------------
| 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() 2 0.3941 0.197066 0.4336 0.216822 9.08 9.99 |
| build_sparsity() 2 0.4049 0.202453 0.9255 0.462728 9.33 21.33 |
| create_dof_constraints() 2 0.0642 0.032095 0.3218 0.160876 1.48 7.42 |
| distribute_dofs() 2 0.1397 0.069852 0.3763 0.188134 3.22 8.67 |
| dof_indices() 25400 1.3241 0.000052 1.3241 0.000052 30.52 30.52 |
| prepare_send_list() 2 0.0012 0.000576 0.0012 0.000576 0.03 0.03 |
| reinit() 2 0.2357 0.117866 0.2357 0.117866 5.43 5.43 |
| |
| EquationSystems |
| build_solution_vector() 1 0.0376 0.037630 0.5532 0.553196 0.87 12.75 |
| |
| ExodusII_IO |
| write_nodal_data() 1 0.0267 0.026654 0.0267 0.026654 0.61 0.61 |
| |
| Mesh |
| find_neighbors() 1 0.0912 0.091225 0.0933 0.093333 2.10 2.15 |
| renumber_nodes_and_elem() 2 0.0037 0.001868 0.0037 0.001868 0.09 0.09 |
| |
| MeshCommunication |
| assign_global_indices() 2 0.4816 0.240818 0.4839 0.241952 11.10 11.15 |
| compute_hilbert_indices() 2 0.0472 0.023620 0.0472 0.023620 1.09 1.09 |
| find_global_indices() 2 0.0159 0.007926 0.0676 0.033792 0.37 1.56 |
| parallel_sort() 2 0.0039 0.001962 0.0040 0.002004 0.09 0.09 |
| |
| MeshOutput |
| write_equation_systems() 1 0.0001 0.000093 0.5800 0.580018 0.00 13.37 |
| |
| MeshTools::Generation |
| build_cube() 1 0.0156 0.015646 0.0156 0.015646 0.36 0.36 |
| |
| MetisPartitioner |
| partition() 1 0.1197 0.119664 0.1529 0.152856 2.76 3.52 |
| |
| Parallel |
| allgather() 22 0.0007 0.000031 0.0007 0.000033 0.02 0.02 |
| barrier() 2 0.0000 0.000010 0.0000 0.000010 0.00 0.00 |
| broadcast() 48 0.0004 0.000008 0.0003 0.000007 0.01 0.01 |
| max(bool) 2 0.0000 0.000003 0.0000 0.000003 0.00 0.00 |
| max(scalar) 177 0.0005 0.000003 0.0005 0.000003 0.01 0.01 |
| max(vector) 39 0.0003 0.000007 0.0006 0.000015 0.01 0.01 |
| min(bool) 201 0.0005 0.000003 0.0005 0.000003 0.01 0.01 |
| min(scalar) 166 0.3853 0.002321 0.3853 0.002321 8.88 8.88 |
| min(vector) 39 0.0003 0.000009 0.0007 0.000018 0.01 0.02 |
| probe() 42 0.0011 0.000026 0.0011 0.000026 0.03 0.03 |
| receive() 38 0.0007 0.000017 0.0017 0.000045 0.02 0.04 |
| send() 38 0.0002 0.000006 0.0002 0.000006 0.00 0.00 |
| send_receive() 38 0.0011 0.000029 0.0028 0.000074 0.03 0.07 |
| sum() 48 0.0015 0.000031 0.3811 0.007939 0.03 8.78 |
| |
| Parallel::Request |
| wait() 38 0.0003 0.000009 0.0003 0.000009 0.01 0.01 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0053 0.005304 0.0054 0.005385 0.12 0.12 |
| set_parent_processor_ids() 1 0.0038 0.003775 0.0038 0.003775 0.09 0.09 |
| |
| RBConstruction |
| clear() 3 0.0007 0.000220 0.0007 0.000220 0.02 0.02 |
| load_rb_solution() 2 0.0004 0.000186 0.0004 0.000186 0.01 0.01 |
| |
| RBEIMEvaluation |
| rb_solve() 1 0.0101 0.010114 0.0101 0.010114 0.23 0.23 |
| read_offline_data_from_files() 1 0.0001 0.000087 0.0007 0.000662 0.00 0.02 |
| |
| RBEvaluation |
| clear() 3 0.0001 0.000032 0.0001 0.000032 0.00 0.00 |
| compute_residual_dual_norm() 1 0.0038 0.003769 0.0038 0.003769 0.09 0.09 |
| rb_solve() 1 0.0001 0.000139 0.0140 0.014022 0.00 0.32 |
| read_in_basis_functions() 2 0.0000 0.000021 1.3774 0.688724 0.00 31.75 |
| read_in_vectors() 2 0.5134 0.256689 1.3774 0.688703 11.83 31.75 |
| read_offline_data_from_files() 2 0.0008 0.000396 0.0012 0.000584 0.02 0.03 |
| resize_data_structures() 2 0.0004 0.000187 0.0004 0.000187 0.01 0.01 |
--------------------------------------------------------------------------------------------------------------
| Totals: 26388 4.3390 100.00 |
--------------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example reduced_basis_ex6:
* 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
***************************************************************