libMesh
rb_classes.h
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 #ifndef RB_CLASSES_H
21 #define RB_CLASSES_H
22 
23 #if defined(LIBMESH_HAVE_SLEPC) && defined(LIBMESH_HAVE_GLPK)
24 
25 #include "libmesh/rb_construction.h"
26 #include "libmesh/rb_scm_construction.h"
27 #include "libmesh/fe_base.h"
28 #include "libmesh/auto_ptr.h"
29 #include "libmesh/rb_evaluation.h"
30 #include "libmesh/rb_scm_evaluation.h"
31 
32 // Bring in bits from the libMesh namespace.
33 // Just the bits we're using, since this is a header.
44 using libMesh::Real;
45 using libMesh::UniquePtr;
46 
47 // local include
48 #include "assembly.h"
49 
50 // A simple subclass of RBEvaluation. We also store the theta expansion object
51 // for the affine expansion of the PDE as a member variable.
52 class SimpleRBEvaluation : public RBEvaluation
53 {
54 public:
55 
60  : RBEvaluation(comm)
61  {
63  }
64 
70  {
72  return rb_scm_eval->get_SCM_LB();
73  }
74 
79 
85 
86 };
87 
88 // A simple subclass of RBConstruction, which initializes libMesh-related data such
89 // as the number of variables and their finite element type. We also store the objects
90 // that define the affine expansion of the PDE as member variables.
92 {
93 public:
94 
96  const std::string & name,
97  const unsigned int number)
98  : Parent(es, name, number)
99  {}
100 
104  virtual ~SimpleRBConstruction () {}
105 
110 
115 
119  virtual void init_data()
120  {
121  u_var = this->add_variable ("u", libMesh::FIRST);
122 
123  // Generate a DirichletBoundary object
124  dirichlet_bc = build_zero_dirichlet_boundary_object();
125 
126  // Set the Dirichlet boundary IDs
127  // and the Dirichlet boundary variable numbers
128  dirichlet_bc->b.insert(3);
129  dirichlet_bc->variables.push_back(u_var);
130 
131  // Attach dirichlet_bc (must do this _before_ Parent::init_data)
132  get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
133 
134  Parent::init_data();
135 
136  // Set the rb_assembly_expansion for this Construction object.
137  set_rb_assembly_expansion(ex02_rb_assembly_expansion);
138 
139  // We need to define an inner product matrix for this problem
140  set_inner_product_assembly(ex02_rb_assembly_expansion.B_assembly);
141  }
142 
147  virtual void init_context(FEMContext & c)
148  {
149  // For efficiency, we should prerequest all
150  // the data we will need to build the
151  // linear system before doing an element loop.
152  FEBase * elem_fe = libmesh_nullptr;
153  c.get_element_fe(u_var, elem_fe);
154 
155  elem_fe->get_JxW();
156  elem_fe->get_phi();
157  elem_fe->get_dphi();
158  }
159 
163  unsigned int u_var;
164 
171 
175  UniquePtr<DirichletBoundary> dirichlet_bc;
176 
177 };
178 
179 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
180 
181 #endif
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
This is the EquationSystems class.
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
RBConstruction Parent
The type of the parent.
Definition: rb_classes.h:114
const class libmesh_nullptr_t libmesh_nullptr
This class allows one to associate Dirichlet boundary values with a given set of mesh boundary ids an...
virtual ~SimpleRBConstruction()
Destructor.
Definition: rb_classes.h:104
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
Definition: fem_context.h:262
Ex02RBAssemblyExpansion ex02_rb_assembly_expansion
The object that stores the "assembly" expansion of the parameter dependent PDE, i.e.
Definition: rb_classes.h:170
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual Real get_SCM_LB()
Evaluate single SCM lower bound.
virtual void init_context(FEMContext &c)
Pre-request all relevant element data.
Definition: rb_classes.h:147
Ex02RBThetaExpansion ex02_rb_theta_expansion
The object that stores the "theta" expansion of the parameter dependent PDE, i.e. ...
Definition: rb_classes.h:84
This class is part of the rbOOmit framework.
void set_parameters(const RBParameters &params)
Set the current parameters to params.
virtual void init_data()
Initialize data structures.
Definition: rb_classes.h:119
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
RBSCMEvaluation * rb_scm_eval
Pointer to the SCM object that will provide our coercivity constant lower bound.
Definition: rb_classes.h:78
SimpleRBEvaluation(const libMesh::Parallel::Communicator &comm)
Constructor.
Definition: rb_classes.h:59
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:42
FEGenericBase< Real > FEBase
void set_rb_theta_expansion(RBThetaExpansion &rb_theta_expansion_in)
Set the RBThetaExpansion object.
Definition: rb_evaluation.C:84
This class is part of the rbOOmit framework.
Definition: rb_evaluation.h:50
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Parallel::Communicator & comm() const
This class stores the set of ElemAssembly functor objects that define the "parameter-independent expa...
virtual Real get_stability_lower_bound()
We override get_stability_lower_bound so that it calls rb_scm_eval to return a parameter-dependent lo...
Definition: rb_classes.h:69
This class is part of the rbOOmit framework.
const RBParameters & get_parameters() const
Get the current parameters.
SimpleRBConstruction(EquationSystems &es, const std::string &name, const unsigned int number)
Definition: rb_classes.h:95
SimpleRBConstruction sys_type
The type of system.
Definition: rb_classes.h:109
sys get_dof_map()
This class is part of the rbOOmit framework.