libMesh
rb_classes.h
Go to the documentation of this file.
1 #ifndef RB_CLASSES_H
2 #define RB_CLASSES_H
3 
4 // local includes
5 #include "assembly.h"
6 
7 // rbOOmit includes
8 #include "libmesh/rb_construction.h"
9 #include "libmesh/rb_evaluation.h"
10 
11 // libMesh includes
12 #include "libmesh/fe_base.h"
13 #include "libmesh/dof_map.h"
14 
15 using namespace libMesh;
16 
17 
19 {
20 public:
21 
26  : RBEvaluation(comm)
27  {
28  set_rb_theta_expansion(elasticity_theta_expansion);
29  }
30 
35  virtual Real get_stability_lower_bound() { return 1.; }
36 
42 };
43 
44 
46 {
47 public:
48 
50  const std::string & name_in,
51  const unsigned int number_in) :
52  Parent(es, name_in, number_in),
53  elasticity_assembly_expansion(*this),
54  ip_assembly(*this)
55  {}
56 
61 
66 
71 
75  virtual void init_data()
76  {
77  u_var = this->add_variable("u", FIRST);
78  v_var = this->add_variable("v", FIRST);
79  w_var = this->add_variable("w", FIRST);
80 
81  // Generate a DirichletBoundary object
82  dirichlet_bc = build_zero_dirichlet_boundary_object();
83 
84  // Set the Dirichlet boundary condition
85  dirichlet_bc->b.insert(BOUNDARY_ID_MIN_X); // Dirichlet boundary at x=0
86  dirichlet_bc->variables.push_back(u_var);
87  dirichlet_bc->variables.push_back(v_var);
88  dirichlet_bc->variables.push_back(w_var);
89 
90  // Attach dirichlet_bc (must do this _before_ Parent::init_data)
91  get_dof_map().add_dirichlet_boundary(*dirichlet_bc);
92 
93  Parent::init_data();
94 
95  // Set the rb_assembly_expansion for this Construction object
96  set_rb_assembly_expansion(elasticity_assembly_expansion);
97 
98  // We need to define an inner product matrix for this problem
99  set_inner_product_assembly(ip_assembly);
100  }
101 
105  virtual void init_context(FEMContext & c)
106  {
107  // For efficiency, we should prerequest all
108  // the data we will need to build the
109  // linear system before doing an element loop.
110  FEBase * elem_fe = libmesh_nullptr;
111  c.get_element_fe(u_var, elem_fe);
112 
113  elem_fe->get_JxW();
114  elem_fe->get_phi();
115  elem_fe->get_dphi();
116  }
117 
121  unsigned int u_var;
122  unsigned int v_var;
123  unsigned int w_var;
124 
129 
134 
139 };
140 
141 #endif
This is the EquationSystems class.
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
ElasticityRBConstruction(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Definition: rb_classes.h:49
virtual Real get_stability_lower_bound()
Return a "dummy" lower bound for the coercivity constant.
Definition: rb_classes.h:35
ElasticityThetaExpansion elasticity_theta_expansion
The object that stores the "theta" expansion of the parameter dependent PDE, i.e. ...
Definition: rb_classes.h:41
unsigned int u_var
Variable numbers.
Definition: rb_classes.h:121
const class libmesh_nullptr_t libmesh_nullptr
InnerProductAssembly ip_assembly
Object to assemble the inner product matrix.
Definition: rb_classes.h:133
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
The libMesh namespace provides an interface to certain functionality in the library.
virtual void init_data()
Initialize data structures.
Definition: rb_classes.h:75
ElasticityRBConstruction sys_type
The type of system.
Definition: rb_classes.h:65
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual ~ElasticityRBConstruction()
Destructor.
Definition: rb_classes.h:60
const std::vector< Real > & get_JxW() const
Definition: fe_abstract.h:237
UniquePtr< DirichletBoundary > dirichlet_bc
The object that defines which degrees of freedom are on a Dirichlet boundary.
Definition: rb_classes.h:138
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:208
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
RBConstruction Parent
The type of the parent.
Definition: rb_classes.h:70
virtual void init_context(FEMContext &c)
Pre-request all relevant element data.
Definition: rb_classes.h:105
This class is part of the rbOOmit framework.
Definition: rb_evaluation.h:50
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:216
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
This class is part of the rbOOmit framework.
ElasticityRBEvaluation(const Parallel::Communicator &comm)
Constructor.
Definition: rb_classes.h:25
ElasticityAssemblyExpansion elasticity_assembly_expansion
The object that stores the "assembly" expansion of the parameter dependent PDE.
Definition: rb_classes.h:128
sys get_dof_map()
This class forms the foundation from which generic finite elements may be derived.