The source file assembly.h with comments:

        #ifndef __assembly_h__
        #define __assembly_h__
        
        #if defined(LIBMESH_HAVE_SLEPC) && defined(LIBMESH_HAVE_GLPK)
        
        #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::RBTheta;
        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& mu) { return mu.get_value("mu_0"); } };
        struct ThetaA1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("mu_1"); } };
        struct ThetaA2 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("mu_2"); } };
        
        struct B : ElemAssembly
        {
Assemble the H1 inner product. This will be used as the inner product for this problem.
          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 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] * (phi[j][qp]*phi[i][qp] + dphi[j][qp]*dphi[i][qp]);
          }
        };
        
        
        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();
        
            Real min_x=0.,max_x=0.5;
        
            Point centroid = c.elem->centroid();
            if( (min_x <= centroid(0)) && (centroid(0) < max_x) )
            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
        {
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();
        
            Real min_x=0.5, max_x=1.;
        
            Point centroid = c.elem->centroid();
            if( (min_x <= centroid(0)) && (centroid(0) <= max_x) )
            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 A2 : 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[i][qp](0)*phi[j][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 Ex02RBThetaExpansion : RBThetaExpansion
        {
        
          /**
           * Constructor.
           */
          Ex02RBThetaExpansion()
          {
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 Ex02RBAssemblyExpansion : RBAssemblyExpansion
        {
        
          /**
           * Constructor.
           */
          Ex02RBAssemblyExpansion()
            :
            L0(0.72,0.88,0.72,0.88), // We make sure these output regions conform to the mesh
            L1(0.12,0.28,0.72,0.88),
            L2(0.12,0.28,0.12,0.28),
            L3(0.72,0.88,0.12,0.28)
          {
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
          B B_assembly;
          A0 A0_assembly;
          A1 A1_assembly;
          A2 A2_assembly;
          F0 F0_assembly;
          OutputAssembly L0;
          OutputAssembly L1;
          OutputAssembly L2;
          OutputAssembly L3;
        };
        
        #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
        
        #endif



The source file rb_classes.h with comments:

        #ifndef __rb_classes_h__
        #define __rb_classes_h__
        
        #if defined(LIBMESH_HAVE_SLEPC) && defined(LIBMESH_HAVE_GLPK)
        
        #include "libmesh/rb_construction.h"
        #include "libmesh/rb_scm_construction.h"
        #include "libmesh/fe_base.h"
        
Bring in bits from the libMesh namespace. Just the bits we're using, since this is a header.
        using libMesh::EquationSystems;
        using libMesh::FEMContext;
        using libMesh::RBConstruction;
        using libMesh::RBSCMConstruction;
        using libMesh::RBEvaluation;
        using libMesh::RBSCMEvaluation;
        using libMesh::RBParameters;
        using libMesh::RBThetaExpansion;
        using libMesh::RBAssemblyExpansion;
        using libMesh::AutoPtr;
        using libMesh::DirichletBoundary;
        using libMesh::Real;
        
local include
        #include "assembly.h"
        
A simple subclass of RBEvaluation. We also store the theta expansion object for the affine expansion of the PDE as a member variable.
        class SimpleRBEvaluation : public RBEvaluation
        {
        public:
        
          /**
           * Constructor. Just set the theta expansion.
           */
          SimpleRBEvaluation(const Parallel::Communicator& comm)
            : RBEvaluation(comm)
          {
            set_rb_theta_expansion(ex02_rb_theta_expansion);
          }
        
          /**
           * We override get_stability_lower_bound so that it calls rb_scm_eval to return
           * a parameter-dependent lower bound for the coercivity constant.
           */
          virtual Real get_stability_lower_bound()
          {
            rb_scm_eval->set_parameters( get_parameters() );
            return rb_scm_eval->get_SCM_LB() ;
          }
        
          /**
           * Pointer to the SCM object that will provide our coercivity constant lower bound.
           */
          RBSCMEvaluation* rb_scm_eval;
        
          /**
           * 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.
           */
          Ex02RBThetaExpansion ex02_rb_theta_expansion;
        
        };
        
A simple subclass of RBConstruction, which initializes libMesh-related data such as the number of variables and their finite element type. We also store the objects that define the affine expansion of the PDE as member variables.
        class SimpleRBConstruction : public RBConstruction
        {
        public:
        
          SimpleRBConstruction (EquationSystems& es,
                                const std::string& name,
                                const unsigned int number)
          : Parent(es, name, number)
          {}
        
          /**
           * 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(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.
            set_rb_assembly_expansion(ex02_rb_assembly_expansion);
        
We need to define an inner product matrix for this problem
            set_inner_product_assembly(ex02_rb_assembly_expansion.B_assembly);
          }
        
          /**
           * Pre-request all relevant element data. (This is not essential, but it
           * allows libMesh to cache data and hence run faster.)
           */
          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.
           */
          Ex02RBAssemblyExpansion ex02_rb_assembly_expansion;
        
          /**
           * The object that defines which degrees of freedom are on a Dirichlet boundary.
           */
          AutoPtr<DirichletBoundary> dirichlet_bc;
        
        };
        
        #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
        
        #endif



The source file reduced_basis_ex2.C with comments:

        #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);
        
This example requires SLEPc and GLPK
        #if !defined(LIBMESH_HAVE_SLEPC) || !defined(LIBMESH_HAVE_GLPK)
          libmesh_example_assert(false, "--enable-slepc --enable-glpk");
        #else
        
        #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_ex2.in) using GetPot
          std::string parameters_filename = "reduced_basis_ex2.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 on the default MPI communicator.
          Mesh mesh (init.comm(), 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 SCM Construction object
          RBSCMConstruction & rb_scm_con =
            equation_systems.add_system<RBSCMConstruction> ("RBSCMConvectionDiffusion");
          rb_scm_con.set_RB_system_name("RBConvectionDiffusion");
          rb_scm_con.add_variable("p", FIRST);
        
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();
        
Set parameters for the eigenvalue problems that will be solved by rb_scm_con
          equation_systems.parameters.set<unsigned int>("eigenpairs")    = 1;
          equation_systems.parameters.set<unsigned int>("basis vectors") = 3;
          equation_systems.parameters.set<unsigned int>
            ("linear solver maximum iterations") = 1000;
        
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(mesh.comm());
        
We need to give the RBConstruction object a pointer to our RBEvaluation object
          rb_con.set_rb_evaluation(rb_eval);
        
We also need a SCM evaluation object to perform SCM calculations
          RBSCMEvaluation rb_scm_eval(mesh.comm());
          rb_scm_eval.set_rb_theta_expansion( rb_eval.get_rb_theta_expansion() );
        
Tell rb_eval about rb_scm_eval
          rb_eval.rb_scm_eval = &rb_scm_eval;
        
Finally, need to give rb_scm_con and rb_eval a pointer to the SCM evaluation object, rb_scm_eval
          rb_scm_con.set_rb_scm_evaluation(rb_scm_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);
            rb_scm_con.process_parameters_file(parameters_filename);
        
Print out info that describes the current setup of rb_con
            rb_con.print_info();
            rb_scm_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();
        
Perform the SCM Greedy algorithm to derive the data required for rb_scm_eval to provide a coercivity lower bound.
            rb_scm_con.perform_SCM_greedy();
        
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("rb_data");
            rb_scm_con.get_rb_scm_evaluation().write_offline_data_to_files("scm_data");
        
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,"rb_data");
            }
          }
          else // Perform the Online stage of the RB method
          {
        
Read in the reduced basis data
            rb_eval.read_offline_data_from_files("rb_data");
            rb_scm_eval.read_offline_data_from_files("scm_data");
        
Read in online_N and initialize online parameters
            unsigned int online_N = infile("online_N",1);
            Real online_mu_0 = infile("online_mu_0", 0.);
            Real online_mu_1 = infile("online_mu_1", 0.);
            Real online_mu_2 = infile("online_mu_2", 0.);
            RBParameters online_mu;
            online_mu.set_value("mu_0", online_mu_0);
            online_mu.set_value("mu_1", online_mu_1);
            online_mu.set_value("mu_2", online_mu_2);
            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,"rb_data");
        
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;
        
        #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
        }



The source file assembly.h without comments:

 
  #ifndef __assembly_h__
  #define __assembly_h__
  
  #if defined(LIBMESH_HAVE_SLEPC) && defined(LIBMESH_HAVE_GLPK)
  
  #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::RBTheta;
  using libMesh::Real;
  using libMesh::RealGradient;
  
  struct ThetaA0 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("mu_0"); } };
  struct ThetaA1 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("mu_1"); } };
  struct ThetaA2 : RBTheta { virtual Number evaluate(const RBParameters& mu) { return mu.get_value("mu_2"); } };
  
  struct B : 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] * (phi[j][qp]*phi[i][qp] + dphi[j][qp]*dphi[i][qp]);
    }
  };
  
  
  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();
  
      Real min_x=0.,max_x=0.5;
  
      Point centroid = c.elem->centroid();
      if( (min_x <= centroid(0)) && (centroid(0) < max_x) )
      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<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();
  
      Real min_x=0.5, max_x=1.;
  
      Point centroid = c.elem->centroid();
      if( (min_x <= centroid(0)) && (centroid(0) <= max_x) )
      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 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[i][qp](0)*phi[j][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 Ex02RBThetaExpansion : RBThetaExpansion
  {
  
    /**
     * Constructor.
     */
    Ex02RBThetaExpansion()
    {
      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 Ex02RBAssemblyExpansion : RBAssemblyExpansion
  {
  
    /**
     * Constructor.
     */
    Ex02RBAssemblyExpansion()
      :
      L0(0.72,0.88,0.72,0.88), // We make sure these output regions conform to the mesh
      L1(0.12,0.28,0.72,0.88),
      L2(0.12,0.28,0.12,0.28),
      L3(0.72,0.88,0.12,0.28)
    {
      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
    }
  
    B B_assembly;
    A0 A0_assembly;
    A1 A1_assembly;
    A2 A2_assembly;
    F0 F0_assembly;
    OutputAssembly L0;
    OutputAssembly L1;
    OutputAssembly L2;
    OutputAssembly L3;
  };
  
  #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
  
  #endif



The source file rb_classes.h without comments:

 
  #ifndef __rb_classes_h__
  #define __rb_classes_h__
  
  #if defined(LIBMESH_HAVE_SLEPC) && defined(LIBMESH_HAVE_GLPK)
  
  #include "libmesh/rb_construction.h"
  #include "libmesh/rb_scm_construction.h"
  #include "libmesh/fe_base.h"
  
  using libMesh::EquationSystems;
  using libMesh::FEMContext;
  using libMesh::RBConstruction;
  using libMesh::RBSCMConstruction;
  using libMesh::RBEvaluation;
  using libMesh::RBSCMEvaluation;
  using libMesh::RBParameters;
  using libMesh::RBThetaExpansion;
  using libMesh::RBAssemblyExpansion;
  using libMesh::AutoPtr;
  using libMesh::DirichletBoundary;
  using libMesh::Real;
  
  #include "assembly.h"
  
  class SimpleRBEvaluation : public RBEvaluation
  {
  public:
  
    /**
     * Constructor. Just set the theta expansion.
     */
    SimpleRBEvaluation(const Parallel::Communicator& comm)
      : RBEvaluation(comm)
    {
      set_rb_theta_expansion(ex02_rb_theta_expansion);
    }
  
    /**
     * We override get_stability_lower_bound so that it calls rb_scm_eval to return
     * a parameter-dependent lower bound for the coercivity constant.
     */
    virtual Real get_stability_lower_bound()
    {
      rb_scm_eval->set_parameters( get_parameters() );
      return rb_scm_eval->get_SCM_LB() ;
    }
  
    /**
     * Pointer to the SCM object that will provide our coercivity constant lower bound.
     */
    RBSCMEvaluation* rb_scm_eval;
  
    /**
     * 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.
     */
    Ex02RBThetaExpansion ex02_rb_theta_expansion;
  
  };
  
  class SimpleRBConstruction : public RBConstruction
  {
  public:
  
    SimpleRBConstruction (EquationSystems& es,
                          const std::string& name,
                          const unsigned int number)
    : Parent(es, name, number)
    {}
  
    /**
     * 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(3);
      dirichlet_bc->variables.push_back(u_var);
  
      get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
  
      Parent::init_data();
  
      set_rb_assembly_expansion(ex02_rb_assembly_expansion);
  
      set_inner_product_assembly(ex02_rb_assembly_expansion.B_assembly);
    }
  
    /**
     * Pre-request all relevant element data. (This is not essential, but it
     * allows libMesh to cache data and hence run faster.)
     */
    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.
     */
    Ex02RBAssemblyExpansion ex02_rb_assembly_expansion;
  
    /**
     * The object that defines which degrees of freedom are on a Dirichlet boundary.
     */
    AutoPtr<DirichletBoundary> dirichlet_bc;
  
  };
  
  #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
  
  #endif



The source file reduced_basis_ex2.C without comments:

 
  #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_SLEPC) || !defined(LIBMESH_HAVE_GLPK)
    libmesh_example_assert(false, "--enable-slepc --enable-glpk");
  #else
  
  #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_ex2.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 (init.comm(), 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");
  
    RBSCMConstruction & rb_scm_con =
      equation_systems.add_system<RBSCMConstruction> ("RBSCMConvectionDiffusion");
    rb_scm_con.set_RB_system_name("RBConvectionDiffusion");
    rb_scm_con.add_variable("p", FIRST);
  
    equation_systems.init ();
  
    equation_systems.print_info();
    mesh.print_info();
  
    equation_systems.parameters.set<unsigned int>("eigenpairs")    = 1;
    equation_systems.parameters.set<unsigned int>("basis vectors") = 3;
    equation_systems.parameters.set<unsigned int>
      ("linear solver maximum iterations") = 1000;
  
    SimpleRBEvaluation rb_eval(mesh.comm());
  
    rb_con.set_rb_evaluation(rb_eval);
  
    RBSCMEvaluation rb_scm_eval(mesh.comm());
    rb_scm_eval.set_rb_theta_expansion( rb_eval.get_rb_theta_expansion() );
  
    rb_eval.rb_scm_eval = &rb_scm_eval;
  
    rb_scm_con.set_rb_scm_evaluation(rb_scm_eval);
  
    if(!online_mode) // Perform the Offline stage of the RB method
    {
      rb_con.process_parameters_file(parameters_filename);
      rb_scm_con.process_parameters_file(parameters_filename);
  
      rb_con.print_info();
      rb_scm_con.print_info();
  
      rb_con.initialize_rb_construction();
  
      rb_scm_con.perform_SCM_greedy();
  
      rb_con.train_reduced_basis();
  
      rb_con.get_rb_evaluation().write_offline_data_to_files("rb_data");
      rb_scm_con.get_rb_scm_evaluation().write_offline_data_to_files("scm_data");
  
      if(store_basis_functions)
      {
        rb_con.get_rb_evaluation().write_out_basis_functions(rb_con,"rb_data");
      }
    }
    else // Perform the Online stage of the RB method
    {
  
      rb_eval.read_offline_data_from_files("rb_data");
      rb_scm_eval.read_offline_data_from_files("scm_data");
  
      unsigned int online_N = infile("online_N",1);
      Real online_mu_0 = infile("online_mu_0", 0.);
      Real online_mu_1 = infile("online_mu_1", 0.);
      Real online_mu_2 = infile("online_mu_2", 0.);
      RBParameters online_mu;
      online_mu.set_value("mu_0", online_mu_0);
      online_mu.set_value("mu_1", online_mu_1);
      online_mu.set_value("mu_2", online_mu_2);
      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_data");
  
        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;
  
  #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
  }



The console output of the program:

make[4]: Entering directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/reduced_basis/reduced_basis_ex2'
***************************************************************
* Running Example reduced_basis_ex2:
*  mpirun -np 4 example-devel -online_mode 0 -eps_type lapack -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
 
*** Warning, This code is untested, experimental, or likely to see future API changes: ../src/reduced_basis/rb_parametrized.C, line 41, compiled Apr 19 2013 at 11:42:51 ***
 EquationSystems
  n_systems()=2
   System #0, "RBConvectionDiffusion"
    Type "RBConstruction"
    Variables="u" 
    Finite Element Types="LAGRANGE", "JACOBI_20_00" 
    Infinite Element Mapping="CARTESIAN" 
    Approximation Orders="FIRST", "THIRD" 
    n_dofs()=676
    n_local_dofs()=186
    n_constrained_dofs()=26
    n_local_constrained_dofs()=11
    n_vectors()=1
    n_matrices()=1
    DofMap Sparsity
      Average  On-Processor Bandwidth <= 8.33136
      Average Off-Processor Bandwidth <= 0.488166
      Maximum  On-Processor Bandwidth <= 11
      Maximum Off-Processor Bandwidth <= 8
    DofMap Constraints
      Number of DoF Constraints = 26
      Average DoF Constraint Length= 0
      Number of Node Constraints = 0
   System #1, "RBSCMConvectionDiffusion"
    Type "Eigen"
    Variables="p" 
    Finite Element Types="LAGRANGE", "JACOBI_20_00" 
    Infinite Element Mapping="CARTESIAN" 
    Approximation Orders="FIRST", "THIRD" 
    n_dofs()=676
    n_local_dofs()=186
    n_constrained_dofs()=0
    n_local_constrained_dofs()=0
    n_vectors()=0
    n_matrices()=2
    DofMap Sparsity
      Average  On-Processor Bandwidth <= 8.33136
      Average Off-Processor Bandwidth <= 0.488166
      Maximum  On-Processor Bandwidth <= 11
      Maximum Off-Processor Bandwidth <= 8
    DofMap Constraints
      Number of DoF Constraints = 0
      Number of Node Constraints = 0

 Mesh Information:
  mesh_dimension()=2
  spatial_dimension()=3
  n_nodes()=676
    n_local_nodes()=186
  n_elem()=625
    n_local_elem()=156
    n_active_elem()=625
  n_subdomains()=1
  n_partitions()=4
  n_processors()=4
  n_threads()=1
  processor_id()=0

Initializing training parameters with random training set...
Parameter mu_0: log scaling = 1
Parameter mu_1: log scaling = 1
Parameter mu_2: log scaling = 1

Initializing training parameters with random training set...
Parameter mu_0: log scaling = 1
Parameter mu_1: log scaling = 1
Parameter mu_2: log scaling = 1


RBConstruction parameters:
system name: RBConvectionDiffusion
constrained_problem: 0
Nmax: 10
Basis training error tolerance: 1e-05
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: 3
Parameter mu_0: Min = 0.1, Max = 1, value = 0.2
Parameter mu_1: Min = 0.1, Max = 1, value = 0.7
Parameter mu_2: Min = 0.01, Max = 0.1, value = 0.1
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


RBSCMConstruction parameters:
system name: RBSCMConvectionDiffusion
SCM Greedy tolerance: 0.1
A_q operators attached: 3
Number of parameters: 3
Parameter mu_0: Min = 0.1, Max = 1, value = 0.2
Parameter mu_1: Min = 0.1, Max = 1, value = 0.7
Parameter mu_2: Min = 0.01, Max = 0.1, value = 0.1
n_training_samples: 100


B_min(0) = -1.0134e-15
B_max(0) = 0.999932

B_min(1) = -1.43169e-15
B_max(1) = 0.999933

B_min(2) = -0.380786
B_max(2) = 5.58608e-17

SCM: Added mu = (0.69213,0.251735,0.0335734)

Stability constant for C_J(0) = 0.216559

SCM iteration 0, max_SCM_error = 0.980015

SCM: Added mu = (0.143443,0.889871,0.0852727)

-----------------------------------


Stability constant for C_J(1) = 0.0914409

SCM iteration 1, max_SCM_error = 0.737033

SCM: Added mu = (0.156446,0.113811,0.0968659)

-----------------------------------


Stability constant for C_J(2) = 0.0620143

SCM iteration 2, max_SCM_error = 0.442637

SCM: Added mu = (0.467928,0.460038,0.0337675)

-----------------------------------


Stability constant for C_J(3) = 0.321689

SCM iteration 3, max_SCM_error = 0.16307

SCM: Added mu = (0.130217,0.260824,0.0709816)

-----------------------------------


Stability constant for C_J(4) = 0.0811065

SCM iteration 4, max_SCM_error = 0.147242

SCM: Added mu = (0.895839,0.590687,0.0184785)

-----------------------------------


Stability constant for C_J(5) = 0.503506

SCM iteration 5, max_SCM_error = 0.141612

SCM: Added mu = (0.637055,0.117424,0.0854742)

-----------------------------------


Stability constant for C_J(6) = 0.0849774

SCM iteration 6, max_SCM_error = 0.101869

SCM: Added mu = (0.136139,0.181766,0.0787585)

-----------------------------------


Stability constant for C_J(7) = 0.0789329

SCM iteration 7, max_SCM_error = 0.0851546

SCM tolerance of 0.1 reached.

Compute output dual inner products
output_dual_innerprods[0][0] = 0.839698
output_dual_innerprods[1][0] = 0.318298
output_dual_innerprods[2][0] = 0.318298
output_dual_innerprods[3][0] = 0.839698

---- Performing Greedy basis enrichment ----

---- Basis dimension: 0 ----
Performing RB solves on training set
Maximum (absolute) error bound is 7.87217

Performing truth solve at parameter:
mu_0: 1.564464e-01
mu_1: 1.138114e-01
mu_2: 9.686588e-02

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 8.23467

Performing truth solve at parameter:
mu_0: 1.023809e-01
mu_1: 5.313481e-01
mu_2: 2.107070e-02

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 1.93695

Performing truth solve at parameter:
mu_0: 9.735663e-01
mu_1: 1.093829e-01
mu_2: 3.718569e-02

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 0.477013

Performing truth solve at parameter:
mu_0: 6.370547e-01
mu_1: 1.174236e-01
mu_2: 8.547420e-02

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 0.387463

Performing truth solve at parameter:
mu_0: 1.434429e-01
mu_1: 8.898710e-01
mu_2: 8.527266e-02

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 0.00459194

Performing truth solve at parameter:
mu_0: 1.145036e-01
mu_1: 7.649858e-01
mu_2: 1.009016e-02

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 0.000485597

Performing truth solve at parameter:
mu_0: 1.009280e-01
mu_1: 2.392213e-01
mu_2: 1.097415e-02

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 4.75365e-05

Performing truth solve at parameter:
mu_0: 1.322754e-01
mu_1: 1.267559e-01
mu_2: 5.962927e-02

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.07344e-05

Performing truth solve at parameter:
mu_0: 1.165157e-01
mu_1: 7.154460e-01
mu_2: 4.061487e-02

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 6.24584e-07

Specified error tolerance reached.
In RBSCMEvaluation::write_offline_data_to_files, directory scm_data already exists, overwriting contents.

 -------------------------------------------------------------------------------------------------------------------
| Processor id:   0                                                                                                 |
| Num Processors: 4                                                                                                 |
| Time:           Fri Apr 19 11:55:00 2013                                                                          |
| OS:             Linux                                                                                             |
| HostName:       spark.ices.utexas.edu                                                                             |
| OS Release:     2.6.32-279.22.1.el6.x86_64                                                                        |
| OS Version:     #1 SMP Tue Feb 5 14:33:39 CST 2013                                                                |
| Machine:        x86_64                                                                                            |
| Username:       roystgnr                                                                                          |
| Configuration:  ../configure  '--enable-everything'                                                               |
|  'METHODS=devel'                                                                                                  |
|  '--prefix=/h2/roystgnr/libmesh-test'                                                                             |
|  'CXX=distcc /usr/bin/g++'                                                                                        |
|  'CC=distcc /usr/bin/gcc'                                                                                         |
|  'FC=distcc /usr/bin/gfortran'                                                                                    |
|  'F77=distcc /usr/bin/gfortran'                                                                                   |
|  'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2'                                                          |
|  'PETSC_ARCH=gcc-system-mkl-gf-10.3.12.361-mpich2-1.4.1p1-cxx-opt'                                                |
|  'SLEPC_DIR=/opt/apps/ossw/libraries/slepc/slepc-3.3-p2-petsc-3.3-p2-cxx-opt'                                     |
|  'TRILINOS_DIR=/opt/apps/ossw/libraries/trilinos/trilinos-10.12.2/sl6/gcc-system/mpich2-1.4.1p1/mkl-gf-10.3.12.361'|
|  'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/gcc-system'                                                 |
|  'HDF5_DIR=/opt/apps/ossw/libraries/hdf5/hdf5-1.8.9/sl6/gcc-system'                                               |
 -------------------------------------------------------------------------------------------------------------------
 --------------------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=7.47322, Active time=7.38946                                                       |
 --------------------------------------------------------------------------------------------------------------------
| 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   |
|--------------------------------------------------------------------------------------------------------------------|
|                                                                                                                    |
|                                                                                                                    |
| CondensedEigenSystem                                                                                               |
|   get_eigenpair()                      14        0.0046      0.000325    0.4832      0.034514    0.06     6.54     |
|   solve()                              14        0.0190      0.001356    5.9937      0.428124    0.26     81.11    |
|                                                                                                                    |
| DofMap                                                                                                             |
|   add_neighbors_to_send_list()         2         0.0016      0.000815    0.0023      0.001163    0.02     0.03     |
|   build_constraint_matrix()            9360      0.0121      0.000001    0.0121      0.000001    0.16     0.16     |
|   build_sparsity()                     2         0.0013      0.000631    0.0034      0.001725    0.02     0.05     |
|   cnstrn_elem_mat_vec()                9360      0.0205      0.000002    0.0205      0.000002    0.28     0.28     |
|   create_dof_constraints()             2         0.0019      0.000965    0.0049      0.002434    0.03     0.07     |
|   distribute_dofs()                    2         0.0033      0.001649    0.0093      0.004650    0.04     0.13     |
|   dof_indices()                        19795     0.0952      0.000005    0.0952      0.000005    1.29     1.29     |
|   prepare_send_list()                  2         0.0000      0.000007    0.0000      0.000007    0.00     0.00     |
|   reinit()                             2         0.0053      0.002657    0.0053      0.002657    0.07     0.07     |
|                                                                                                                    |
| FE                                                                                                                 |
|   compute_shape_functions()            21600     0.0995      0.000005    0.0995      0.000005    1.35     1.35     |
|   init_shape_functions()               3000      0.0174      0.000006    0.0174      0.000006    0.24     0.24     |
|   inverse_map()                        5760      0.0216      0.000004    0.0216      0.000004    0.29     0.29     |
|                                                                                                                    |
| FEMap                                                                                                              |
|   compute_affine_map()                 21600     0.0573      0.000003    0.0573      0.000003    0.78     0.78     |
|   compute_face_map()                   2880      0.0170      0.000006    0.0393      0.000014    0.23     0.53     |
|   init_face_shape_functions()          120       0.0004      0.000003    0.0004      0.000003    0.01     0.01     |
|   init_reference_to_physical_map()     3000      0.0133      0.000004    0.0133      0.000004    0.18     0.18     |
|                                                                                                                    |
| Mesh                                                                                                               |
|   find_neighbors()                     1         0.0017      0.001736    0.0020      0.002021    0.02     0.03     |
|   renumber_nodes_and_elem()            2         0.0002      0.000087    0.0002      0.000087    0.00     0.00     |
|                                                                                                                    |
| MeshCommunication                                                                                                  |
|   assign_global_indices()              1         0.0126      0.012564    0.0129      0.012944    0.17     0.18     |
|   compute_hilbert_indices()            2         0.0055      0.002750    0.0055      0.002750    0.07     0.07     |
|   find_global_indices()                2         0.0008      0.000395    0.0071      0.003572    0.01     0.10     |
|   parallel_sort()                      2         0.0004      0.000209    0.0005      0.000262    0.01     0.01     |
|                                                                                                                    |
| MeshTools::Generation                                                                                              |
|   build_cube()                         1         0.0005      0.000547    0.0005      0.000547    0.01     0.01     |
|                                                                                                                    |
| MetisPartitioner                                                                                                   |
|   partition()                          1         0.0055      0.005515    0.0090      0.008996    0.07     0.12     |
|                                                                                                                    |
| Parallel                                                                                                           |
|   allgather()                          17        0.0002      0.000011    0.0002      0.000013    0.00     0.00     |
|   barrier()                            1         0.0000      0.000014    0.0000      0.000014    0.00     0.00     |
|   broadcast()                          38        0.0002      0.000005    0.0002      0.000005    0.00     0.00     |
|   max(bool)                            7         0.0000      0.000005    0.0000      0.000005    0.00     0.00     |
|   max(scalar)                          168       0.0009      0.000005    0.0009      0.000005    0.01     0.01     |
|   max(vector)                          37        0.0003      0.000008    0.0008      0.000022    0.00     0.01     |
|   maxloc(scalar)                       18        0.0008      0.000047    0.0008      0.000047    0.01     0.01     |
|   min(bool)                            183       0.0008      0.000004    0.0008      0.000004    0.01     0.01     |
|   min(scalar)                          149       0.0014      0.000009    0.0014      0.000009    0.02     0.02     |
|   min(vector)                          37        0.0003      0.000009    0.0009      0.000024    0.00     0.01     |
|   probe()                              98        0.0003      0.000003    0.0003      0.000003    0.00     0.00     |
|   receive()                            82        0.0003      0.000004    0.0005      0.000007    0.00     0.01     |
|   send()                               70        0.0002      0.000002    0.0002      0.000002    0.00     0.00     |
|   send_receive()                       74        0.0004      0.000005    0.0011      0.000014    0.01     0.01     |
|   sum()                                58        0.4789      0.008258    0.4790      0.008259    6.48     6.48     |
|                                                                                                                    |
| Parallel::Request                                                                                                  |
|   wait()                               70        0.0001      0.000001    0.0001      0.000001    0.00     0.00     |
|                                                                                                                    |
| Partitioner                                                                                                        |
|   set_node_processor_ids()             1         0.0006      0.000570    0.0007      0.000735    0.01     0.01     |
|   set_parent_processor_ids()           1         0.0002      0.000183    0.0002      0.000183    0.00     0.00     |
|                                                                                                                    |
| PetscLinearSolver                                                                                                  |
|   solve()                              41        0.1385      0.003379    0.1385      0.003379    1.87     1.87     |
|                                                                                                                    |
| RBConstruction                                                                                                     |
|   add_scaled_Aq()                      51        0.0014      0.000028    0.5024      0.009850    0.02     6.80     |
|   add_scaled_matrix_and_vector()       60        0.2320      0.003866    0.5918      0.009864    3.14     8.01     |
|   clear()                              1         0.0005      0.000524    0.0005      0.000524    0.01     0.01     |
|   compute_Fq_representor_innerprods()  1         0.0007      0.000666    0.0042      0.004213    0.01     0.06     |
|   compute_max_error_bound()            10        0.0029      0.000290    0.0506      0.005055    0.04     0.68     |
|   compute_output_dual_innerprods()     1         0.0013      0.001330    0.0147      0.014684    0.02     0.20     |
|   enrich_RB_space()                    9         0.0037      0.000408    0.0037      0.000408    0.05     0.05     |
|   train_reduced_basis()                1         0.0012      0.001153    0.2516      0.251554    0.02     3.40     |
|   truth_assembly()                     9         0.0147      0.001628    0.0147      0.001628    0.20     0.20     |
|   truth_solve()                        9         0.0011      0.000125    0.0487      0.005414    0.02     0.66     |
|   update_RB_system_matrices()          9         0.0111      0.001233    0.0111      0.001233    0.15     0.15     |
|   update_residual_terms()              9         0.0287      0.003191    0.1174      0.013046    0.39     1.59     |
|                                                                                                                    |
| RBEvaluation                                                                                                       |
|   clear()                              1         0.0001      0.000142    0.0001      0.000142    0.00     0.00     |
|   compute_residual_dual_norm()         250       0.0212      0.000085    0.0212      0.000085    0.29     0.29     |
|   rb_solve()                           250       0.0070      0.000028    0.0473      0.000189    0.09     0.64     |
|   resize_data_structures()             1         0.0000      0.000040    0.0000      0.000040    0.00     0.00     |
|   write_offline_data_to_files()        1         0.0016      0.001617    0.0016      0.001617    0.02     0.02     |
|   write_out_basis_functions()          1         0.0001      0.000139    0.0173      0.017258    0.00     0.23     |
|   write_out_vectors()                  1         0.0039      0.003911    0.0171      0.017119    0.05     0.23     |
|                                                                                                                    |
| RBSCMConstruction                                                                                                  |
|   add_scaled_symm_Aq()                 51        0.0002      0.000004    0.5026      0.009854    0.00     6.80     |
|   compute_SCM_bounding_box()           1         0.0003      0.000320    2.7776      2.777602    0.00     37.59    |
|   compute_SCM_bounds_on_training_set() 8         0.0021      0.000262    0.0150      0.001872    0.03     0.20     |
|   enrich_C_J()                         8         0.0003      0.000038    0.0004      0.000052    0.00     0.01     |
|   evaluate_stability_constant()        8         0.0028      0.000352    4.2050      0.525631    0.04     56.91    |
|   perform_SCM_greedy()                 1         0.0011      0.001058    6.9991      6.999101    0.01     94.72    |
|                                                                                                                    |
| RBSCMEvaluation                                                                                                    |
|   get_SCM_LB()                         450       0.0299      0.000067    0.0299      0.000067    0.41     0.41     |
|   get_SCM_UB()                         200       0.0012      0.000006    0.0012      0.000006    0.02     0.02     |
|   write_offline_data_to_files()        1         0.0006      0.000576    0.0006      0.000576    0.01     0.01     |
|                                                                                                                    |
| SlepcEigenSolver                                                                                                   |
|   solve_generalized()                  14        5.9748      0.426768    5.9748      0.426768    80.86    80.86    |
 --------------------------------------------------------------------------------------------------------------------
| Totals:                                99094     7.3895                                          100.00            |
 --------------------------------------------------------------------------------------------------------------------

 
***************************************************************
* Done Running Example reduced_basis_ex2:
*  mpirun -np 4 example-devel -online_mode 0 -eps_type lapack -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
***************************************************************
* Running Example reduced_basis_ex2:
*  mpirun -np 4 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
***************************************************************
 
*** Warning, This code is untested, experimental, or likely to see future API changes: ../src/reduced_basis/rb_parametrized.C, line 41, compiled Apr 19 2013 at 11:42:51 ***
 EquationSystems
  n_systems()=2
   System #0, "RBConvectionDiffusion"
    Type "RBConstruction"
    Variables="u" 
    Finite Element Types="LAGRANGE", "JACOBI_20_00" 
    Infinite Element Mapping="CARTESIAN" 
    Approximation Orders="FIRST", "THIRD" 
    n_dofs()=676
    n_local_dofs()=186
    n_constrained_dofs()=26
    n_local_constrained_dofs()=11
    n_vectors()=1
    n_matrices()=1
    DofMap Sparsity
      Average  On-Processor Bandwidth <= 8.33136
      Average Off-Processor Bandwidth <= 0.488166
      Maximum  On-Processor Bandwidth <= 11
      Maximum Off-Processor Bandwidth <= 8
    DofMap Constraints
      Number of DoF Constraints = 26
      Average DoF Constraint Length= 0
      Number of Node Constraints = 0
   System #1, "RBSCMConvectionDiffusion"
    Type "Eigen"
    Variables="p" 
    Finite Element Types="LAGRANGE", "JACOBI_20_00" 
    Infinite Element Mapping="CARTESIAN" 
    Approximation Orders="FIRST", "THIRD" 
    n_dofs()=676
    n_local_dofs()=186
    n_constrained_dofs()=0
    n_local_constrained_dofs()=0
    n_vectors()=0
    n_matrices()=2
    DofMap Sparsity
      Average  On-Processor Bandwidth <= 8.33136
      Average Off-Processor Bandwidth <= 0.488166
      Maximum  On-Processor Bandwidth <= 11
      Maximum Off-Processor Bandwidth <= 8
    DofMap Constraints
      Number of DoF Constraints = 0
      Number of Node Constraints = 0

 Mesh Information:
  mesh_dimension()=2
  spatial_dimension()=3
  n_nodes()=676
    n_local_nodes()=186
  n_elem()=625
    n_local_elem()=156
    n_active_elem()=625
  n_subdomains()=1
  n_partitions()=4
  n_processors()=4
  n_threads()=1
  processor_id()=0

mu_0: 2.000000e-01
mu_1: 7.000000e-01
mu_2: 1.000000e-01

output 1, value = 2.35239, bound = 2.79023e-05
output 2, value = 0.944939, bound = 1.71789e-05
output 3, value = 0.944939, bound = 1.71789e-05
output 4, value = 2.35239, bound = 2.79023e-05


 -------------------------------------------------------------------------------------------------------------------
| Processor id:   0                                                                                                 |
| Num Processors: 4                                                                                                 |
| Time:           Fri Apr 19 11:55:00 2013                                                                          |
| OS:             Linux                                                                                             |
| HostName:       spark.ices.utexas.edu                                                                             |
| OS Release:     2.6.32-279.22.1.el6.x86_64                                                                        |
| OS Version:     #1 SMP Tue Feb 5 14:33:39 CST 2013                                                                |
| Machine:        x86_64                                                                                            |
| Username:       roystgnr                                                                                          |
| Configuration:  ../configure  '--enable-everything'                                                               |
|  'METHODS=devel'                                                                                                  |
|  '--prefix=/h2/roystgnr/libmesh-test'                                                                             |
|  'CXX=distcc /usr/bin/g++'                                                                                        |
|  'CC=distcc /usr/bin/gcc'                                                                                         |
|  'FC=distcc /usr/bin/gfortran'                                                                                    |
|  'F77=distcc /usr/bin/gfortran'                                                                                   |
|  'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2'                                                          |
|  'PETSC_ARCH=gcc-system-mkl-gf-10.3.12.361-mpich2-1.4.1p1-cxx-opt'                                                |
|  'SLEPC_DIR=/opt/apps/ossw/libraries/slepc/slepc-3.3-p2-petsc-3.3-p2-cxx-opt'                                     |
|  'TRILINOS_DIR=/opt/apps/ossw/libraries/trilinos/trilinos-10.12.2/sl6/gcc-system/mpich2-1.4.1p1/mkl-gf-10.3.12.361'|
|  'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/gcc-system'                                                 |
|  'HDF5_DIR=/opt/apps/ossw/libraries/hdf5/hdf5-1.8.9/sl6/gcc-system'                                               |
 -------------------------------------------------------------------------------------------------------------------
 --------------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=0.312184, Active time=0.232118                                               |
 --------------------------------------------------------------------------------------------------------------
| 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.0009      0.000447    0.0013      0.000635    0.38     0.55     |
|   build_sparsity()               2         0.0007      0.000365    0.0031      0.001560    0.31     1.34     |
|   create_dof_constraints()       2         0.0012      0.000585    0.0028      0.001393    0.50     1.20     |
|   distribute_dofs()              2         0.0018      0.000890    0.0085      0.004226    0.77     3.64     |
|   dof_indices()                  1699      0.0051      0.000003    0.0051      0.000003    2.18     2.18     |
|   prepare_send_list()            2         0.0000      0.000005    0.0000      0.000005    0.00     0.00     |
|   reinit()                       2         0.0028      0.001391    0.0028      0.001391    1.20     1.20     |
|                                                                                                              |
| EquationSystems                                                                                              |
|   build_solution_vector()        2         0.0012      0.000619    0.0051      0.002542    0.53     2.19     |
|                                                                                                              |
| ExodusII_IO                                                                                                  |
|   write_nodal_data()             2         0.1689      0.084441    0.1689      0.084441    72.76    72.76    |
|                                                                                                              |
| Mesh                                                                                                         |
|   find_neighbors()               1         0.0011      0.001107    0.0020      0.001951    0.48     0.84     |
|   renumber_nodes_and_elem()      2         0.0001      0.000058    0.0001      0.000058    0.05     0.05     |
|                                                                                                              |
| MeshCommunication                                                                                            |
|   assign_global_indices()        1         0.0079      0.007865    0.0114      0.011392    3.39     4.91     |
|   compute_hilbert_indices()      2         0.0033      0.001660    0.0033      0.001660    1.43     1.43     |
|   find_global_indices()          2         0.0005      0.000231    0.0065      0.003268    0.20     2.82     |
|   parallel_sort()                2         0.0003      0.000161    0.0024      0.001177    0.14     1.01     |
|                                                                                                              |
| MeshOutput                                                                                                   |
|   write_equation_systems()       2         0.0001      0.000042    0.1742      0.087090    0.04     75.04    |
|                                                                                                              |
| MeshTools::Generation                                                                                        |
|   build_cube()                   1         0.0004      0.000425    0.0004      0.000425    0.18     0.18     |
|                                                                                                              |
| MetisPartitioner                                                                                             |
|   partition()                    1         0.0034      0.003441    0.0068      0.006835    1.48     2.94     |
|                                                                                                              |
| Parallel                                                                                                     |
|   allgather()                    17        0.0053      0.000309    0.0053      0.000314    2.27     2.30     |
|   barrier()                      1         0.0007      0.000715    0.0007      0.000715    0.31     0.31     |
|   broadcast()                    13        0.0000      0.000003    0.0000      0.000003    0.02     0.01     |
|   max(bool)                      3         0.0000      0.000004    0.0000      0.000004    0.00     0.00     |
|   max(scalar)                    186       0.0008      0.000004    0.0008      0.000004    0.36     0.36     |
|   max(vector)                    41        0.0003      0.000006    0.0008      0.000018    0.11     0.33     |
|   min(bool)                      213       0.0008      0.000004    0.0008      0.000004    0.35     0.35     |
|   min(scalar)                    175       0.0106      0.000061    0.0106      0.000061    4.57     4.57     |
|   min(vector)                    41        0.0003      0.000008    0.0011      0.000026    0.15     0.46     |
|   probe()                        82        0.0012      0.000014    0.0012      0.000014    0.50     0.50     |
|   receive()                      76        0.0002      0.000003    0.0013      0.000018    0.08     0.58     |
|   send()                         76        0.0001      0.000002    0.0001      0.000002    0.06     0.06     |
|   send_receive()                 74        0.0002      0.000003    0.0017      0.000023    0.11     0.73     |
|   sum()                          44        0.0027      0.000062    0.0037      0.000084    1.17     1.59     |
|                                                                                                              |
| Parallel::Request                                                                                            |
|   wait()                         76        0.0001      0.000001    0.0001      0.000001    0.03     0.03     |
|                                                                                                              |
| Partitioner                                                                                                  |
|   set_node_processor_ids()       1         0.0003      0.000307    0.0024      0.002428    0.13     1.05     |
|   set_parent_processor_ids()     1         0.0001      0.000103    0.0001      0.000103    0.04     0.04     |
|                                                                                                              |
| RBConstruction                                                                                               |
|   clear()                        1         0.0002      0.000224    0.0002      0.000224    0.10     0.10     |
|   load_basis_function()          1         0.0001      0.000119    0.0001      0.000119    0.05     0.05     |
|   load_rb_solution()             1         0.0001      0.000135    0.0001      0.000135    0.06     0.06     |
|                                                                                                              |
| RBEvaluation                                                                                                 |
|   clear()                        1         0.0001      0.000070    0.0001      0.000070    0.03     0.03     |
|   compute_residual_dual_norm()   1         0.0001      0.000080    0.0001      0.000080    0.03     0.03     |
|   rb_solve()                     1         0.0039      0.003922    0.0042      0.004236    1.69     1.82     |
|   read_in_basis_functions()      1         0.0000      0.000012    0.0159      0.015892    0.01     6.85     |
|   read_in_vectors()              1         0.0032      0.003162    0.0159      0.015880    1.36     6.84     |
|   read_offline_data_from_files() 1         0.0005      0.000469    0.0005      0.000525    0.20     0.23     |
|   resize_data_structures()       1         0.0001      0.000056    0.0001      0.000056    0.02     0.02     |
|                                                                                                              |
| RBSCMEvaluation                                                                                              |
|   get_SCM_LB()                   1         0.0002      0.000234    0.0002      0.000234    0.10     0.10     |
|   read_offline_data_from_files() 1         0.0001      0.000131    0.0001      0.000131    0.06     0.06     |
 --------------------------------------------------------------------------------------------------------------
| Totals:                          2861      0.2321                                          100.00            |
 --------------------------------------------------------------------------------------------------------------

 
***************************************************************
* Done Running Example reduced_basis_ex2:
*  mpirun -np 4 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
***************************************************************
make[4]: Leaving directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/reduced_basis/reduced_basis_ex2'

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

Hosted By:
SourceForge.net Logo