libMesh
rb_eim_assembly.C
Go to the documentation of this file.
1 // rbOOmit: An implementation of the Certified Reduced Basis method.
2 // Copyright (C) 2009, 2010 David J. Knezevic
3 
4 // This file is part of rbOOmit.
5 
6 // rbOOmit is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
10 
11 // rbOOmit is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
15 
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20 // rbOOmit includes
21 #include "libmesh/rb_eim_assembly.h"
22 #include "libmesh/rb_eim_construction.h"
23 #include "libmesh/rb_evaluation.h"
24 
25 // libMesh includes
26 #include "libmesh/fem_context.h"
27 #include "libmesh/dof_map.h"
28 #include "libmesh/quadrature.h"
29 #include "libmesh/libmesh_logging.h"
30 
31 namespace libMesh
32 {
33 
35  unsigned int basis_function_index_in)
36  : _rb_eim_con(rb_eim_con_in),
37  _basis_function_index(basis_function_index_in),
38  _ghosted_basis_function(NumericVector<Number>::build(rb_eim_con_in.get_explicit_system().comm()))
39 {
40  // localize the vector that stores the basis function for this assembly object,
41  // i.e. the vector that is used in evaluate_basis_function_at_quad_pts
42 #ifdef LIBMESH_ENABLE_GHOSTED
46  false,
47  GHOSTED);
49  localize(*_ghosted_basis_function,
51 #else
54  localize(*_ghosted_basis_function);
55 #endif
56 
57  initialize_fe();
58 }
59 
61 {
62 }
63 
65  const Elem & element,
66  const QBase & element_qrule,
67  std::vector<Number> & values)
68 {
69  LOG_SCOPE("evaluate_basis_function", "RBEIMAssembly");
70 
71  bool repeated_qrule = false;
72  if (_qrule.get() != libmesh_nullptr)
73  {
74  repeated_qrule =
75  ( (element_qrule.type() == _qrule->type()) &&
76  (element_qrule.get_dim() == _qrule->get_dim()) &&
77  (element_qrule.get_order() == _qrule->get_order()) );
78  }
79 
80  // If the qrule is not repeated, then we need to make a new copy of element_qrule.
81  if (!repeated_qrule)
82  {
83  _qrule.reset(QBase::build(element_qrule.type(),
84  element_qrule.get_dim(),
85  element_qrule.get_order()).release());
86 
88  }
89 
90  const std::vector<std::vector<Real>> & phi = get_fe().get_phi();
91 
92  // The FE object caches data, hence we recompute as little as
93  // possible on the call to reinit.
94  get_fe().reinit (&element);
95 
96  std::vector<dof_id_type> dof_indices_var;
97 
99  dof_map.dof_indices (&element, dof_indices_var, var);
100 
101  libmesh_assert(dof_indices_var.size() == phi.size());
102 
103  unsigned int n_qpoints = _qrule->n_points();
104  values.resize(n_qpoints);
105 
106  for (unsigned int qp=0; qp<n_qpoints; qp++)
107  {
108  values[qp] = 0.;
109  for (std::size_t i=0; i<dof_indices_var.size(); i++)
110  values[qp] += (*_ghosted_basis_function)(dof_indices_var[i]) * phi[i][qp];
111  }
112 }
113 
115 {
116  return _rb_eim_con;
117 }
118 
120 {
121  return *_ghosted_basis_function;
122 }
123 
125 {
126  return *_fe;
127 }
128 
130 {
132 
133  const unsigned int dim =
135 
136  FEType fe_type = dof_map.variable_type(0);
137  _fe.reset(FEBase::build(dim, fe_type).release());
138 
139  // Pre-request the shape function for efficieny's sake
140  _fe->get_phi();
141 }
142 
143 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
Order get_order() const
Definition: quadrature.h:189
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
unsigned int dim
This class is part of the rbOOmit framework.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
const class libmesh_nullptr_t libmesh_nullptr
virtual ~RBEIMAssembly()
Destructor.
Numeric vector.
Definition: dof_map.h:66
UniquePtr< QBase > _qrule
The libMesh namespace provides an interface to certain functionality in the library.
libmesh_assert(j)
RBEIMConstruction & get_rb_eim_construction()
Get a reference to the RBEIMConstruction object.
void initialize_fe()
Initialize the FE object.
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
RBEIMAssembly(RBEIMConstruction &rb_eim_con_in, unsigned int basis_function_index_in)
Constructor.
const MeshBase & get_mesh() const
Definition: system.h:2014
const DofMap & get_dof_map() const
Definition: system.h:2030
NumericVector< Number > & get_ghosted_basis_function()
Get a reference to the ghosted_basis_function.
ExplicitSystem & get_explicit_system()
Get the ExplicitSystem associated with this system.
NumericVector< Number > & get_basis_function(unsigned int i)
Get a reference to the i^th basis function.
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:208
unsigned int _basis_function_index
The EIM basis function index (from rb_eim_eval) for this assembly object.
virtual void evaluate_basis_function(unsigned int var, const Elem &element, const QBase &element_qrule, std::vector< Number > &values)
Evaluate variable var_number of this object&#39;s EIM basis function at the points qpoints.
FEBase & get_fe()
Retrieve the FE object.
virtual QuadratureType type() const =0
dof_id_type n_local_dofs() const
Definition: system.C:185
unsigned int get_dim() const
Definition: quadrature.h:122
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
dof_id_type n_dofs() const
Definition: system.C:148
RBEvaluation & get_rb_evaluation()
Get a reference to the RBEvaluation object.
RBEIMConstruction & _rb_eim_con
The RBEIMConstruction object that this RBEIMAssembly is based on.
static UniquePtr< QBase > build(const std::string &name, const unsigned int dim, const Order order=INVALID_ORDER)
Builds a specific quadrature rule based on the name string.
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:394
UniquePtr< FEBase > _fe
We store an FE object and an associated quadrature rule.
virtual void attach_quadrature_rule(QBase *q)=0
Provides the class with the quadrature rule.
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53
This class forms the foundation from which generic finite elements may be derived.
UniquePtr< NumericVector< Number > > _ghosted_basis_function
The basis function that we sample to evaluate the empirical interpolation approximation.
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr)=0
This is at the core of this class.
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.