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 #include "libmesh/rb_construction.h"
24 #include "libmesh/fe_base.h"
25 #include "libmesh/rb_evaluation.h"
26 
27 // local include
28 #include "assembly.h"
29 
30 // Bring in bits from the libMesh namespace.
31 // Just the bits we're using, since this is a header.
37 using libMesh::Real;
38 using libMesh::UniquePtr;
39 
40 
41 // A simple subclass of RBEvaluation, which just needs to specify
42 // (a lower bound for) the coercivity constant for this problem.
43 // For this simple convection-diffusion problem, we can set the
44 // coercivity constant lower bound to 0.05.
46 {
47 public:
48 
53  RBEvaluation(comm)
54  {
56  }
57 
61  virtual Real get_stability_lower_bound() { return 0.05; }
62 
68 
69 };
70 
71 // A simple subclass of Construction, which just needs to override build_rb_evaluation
72 // in order to build a SimpleRBEvaluation object, rather than an RBEvaluation object.
74 {
75 public:
76 
78  const std::string & name_in,
79  const unsigned int number_in)
80  : Parent(es, name_in, number_in),
81  dirichlet_bc(UniquePtr<DirichletBoundary>())
82  {}
83 
87  virtual ~SimpleRBConstruction () { }
88 
93 
98 
102  virtual void init_data()
103  {
104  u_var = this->add_variable ("u", libMesh::FIRST);
105 
106  // Generate a DirichletBoundary object
107  dirichlet_bc = build_zero_dirichlet_boundary_object();
108 
109  // Set the Dirichlet boundary IDs
110  // and the Dirichlet boundary variable numbers
111  dirichlet_bc->b.insert(0);
112  dirichlet_bc->b.insert(1);
113  dirichlet_bc->b.insert(2);
114  dirichlet_bc->b.insert(3);
115  dirichlet_bc->variables.push_back(u_var);
116 
117  // Attach dirichlet_bc (must do this _before_ Parent::init_data)
118  get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
119 
120  Parent::init_data();
121 
122  // Set the rb_assembly_expansion for this Construction object.
123  // The theta expansion comes from the RBEvaluation object.
124  set_rb_assembly_expansion(cd_rb_assembly_expansion);
125 
126  // We need to define an inner product matrix for this problem
127  set_inner_product_assembly(cd_rb_assembly_expansion.A0_assembly);
128  }
129 
133  virtual void init_context(FEMContext & c)
134  {
135  // For efficiency, we should prerequest all
136  // the data we will need to build the
137  // linear system before doing an element loop.
138  FEBase * elem_fe = libmesh_nullptr;
139  c.get_element_fe(u_var, elem_fe);
140 
141  elem_fe->get_JxW();
142  elem_fe->get_phi();
143  elem_fe->get_dphi();
144  }
145 
149  unsigned int u_var;
150 
157 
161  UniquePtr<DirichletBoundary> dirichlet_bc;
162 };
163 
164 #endif
UniquePtr< DirichletBoundary > dirichlet_bc
The object that defines which degrees of freedom are on a Dirichlet boundary.
Definition: rb_classes.h:161
This is the EquationSystems class.
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
CDRBThetaExpansion cd_rb_theta_expansion
The object that stores the "theta" expansion of the parameter dependent PDE, i.e. ...
Definition: rb_classes.h:67
RBConstruction Parent
The type of the parent.
Definition: rb_classes.h:97
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:87
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
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual void init_context(FEMContext &c)
Pre-request all relevant element data.
Definition: rb_classes.h:133
SimpleRBConstruction(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Definition: rb_classes.h:77
unsigned int u_var
Variable number for u.
Definition: rb_classes.h:149
virtual void init_data()
Initialize data structures.
Definition: rb_classes.h:102
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
SimpleRBEvaluation(const libMesh::Parallel::Communicator &comm)
Constructor.
Definition: rb_classes.h:52
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
virtual Real get_stability_lower_bound()
The coercivity constant is bounded below by 0.05.
Definition: rb_classes.h:61
This class is part of the rbOOmit framework.
CDRBAssemblyExpansion cd_rb_assembly_expansion
The object that stores the "assembly" expansion of the parameter dependent PDE, i.e.
Definition: rb_classes.h:156
SimpleRBConstruction sys_type
The type of system.
Definition: rb_classes.h:92
sys get_dof_map()