4 #include "libmesh/sparse_matrix.h" 5 #include "libmesh/numeric_vector.h" 6 #include "libmesh/dense_matrix.h" 7 #include "libmesh/dense_vector.h" 8 #include "libmesh/fe.h" 9 #include "libmesh/fe_interface.h" 10 #include "libmesh/fe_base.h" 11 #include "libmesh/elem_assembly.h" 12 #include "libmesh/quadrature_gauss.h" 15 #include "libmesh/rb_theta.h" 16 #include "libmesh/rb_assembly_expansion.h" 18 #define damping_epsilon 0.001 21 #ifdef LIBMESH_USE_COMPLEX_NUMBERS 92 const unsigned int p_var = 0;
94 FEBase * elem_fe =
nullptr;
97 const std::vector<Real> & JxW = elem_fe->get_JxW();
99 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
101 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
110 for (
unsigned int qp=0; qp != n_qpoints; qp++)
111 for (
unsigned int i=0; i != n_p_dofs; i++)
112 for (
unsigned int j=0; j != n_p_dofs; j++)
114 (phi[j][qp]*phi[i][qp]));
122 const unsigned int p_var = 0;
124 FEBase * elem_fe =
nullptr;
127 const std::vector<Real> & JxW = elem_fe->get_JxW();
131 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
139 for (
unsigned int qp=0; qp != n_qpoints; qp++)
140 for (
unsigned int i=0; i != n_p_dofs; i++)
141 for (
unsigned int j=0; j != n_p_dofs; j++)
151 const unsigned int p_var = 0;
153 FEBase * elem_fe =
nullptr;
156 const std::vector<Real> & JxW = elem_fe->get_JxW();
158 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
166 for (
unsigned int qp=0; qp != n_qpoints; qp++)
167 for (
unsigned int i=0; i != n_p_dofs; i++)
168 for (
unsigned int j=0; j != n_p_dofs; j++)
179 const unsigned int p_var = 0;
181 FEBase * side_fe =
nullptr;
184 const std::vector<Real> & JxW_face = side_fe->get_JxW();
186 const std::vector<std::vector<Real>> & phi_face = side_fe->get_phi();
194 for (
unsigned int qp=0; qp != n_sidepoints; qp++)
195 for (
unsigned int i=0; i != n_p_dofs; i++)
196 for (
unsigned int j=0; j != n_p_dofs; j++)
208 const unsigned int p_var = 0;
210 FEBase * side_fe =
nullptr;
213 const std::vector<Real> & JxW_face = side_fe->get_JxW();
215 const std::vector<std::vector<Real>> & phi_face = side_fe->get_phi();
223 for (
unsigned int qp=0; qp != n_sidepoints; qp++)
224 for (
unsigned int i=0; i != n_p_dofs; i++)
225 for (
unsigned int j=0; j != n_p_dofs; j++)
237 const unsigned int p_var = 0;
239 FEBase * side_fe =
nullptr;
242 const std::vector<Real> & JxW_face = side_fe->get_JxW();
244 const std::vector<std::vector<Real>> & phi_face = side_fe->get_phi();
252 for (
unsigned int qp=0; qp != n_sidepoints; qp++)
253 for (
unsigned int i=0; i != n_p_dofs; i++)
265 const unsigned int p_var = 0;
267 FEBase * side_fe =
nullptr;
270 const std::vector<Real> & JxW_face = side_fe->get_JxW();
272 const std::vector<std::vector<Real>> & phi_face = side_fe->get_phi();
280 for (
unsigned int qp=0; qp != n_sidepoints; qp++)
281 for (
unsigned int i=0; i != n_p_dofs; i++)
Real get_value(const std::string ¶m_name) const
Get the value of the specified parameter, throw an error if it does not exist.
RealVectorValue RealGradient
const DenseMatrix< Number > & get_elem_jacobian() const
Const accessor for element Jacobian.
void get_side_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for edge/face (2D/3D) finite element object for variable var for the largest dimension in th...
ThetaOutput0 theta_output_0
virtual void boundary_assembly(FEMContext &c)
Perform the element boundary assembly.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
AcousticsRBThetaExpansion()
Constructor.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
bool has_side_boundary_id(boundary_id_type id) const
Reports if the boundary id is found on the current side.
void attach_F_assembly(ElemAssembly *Fq_assembly)
Attach ElemAssembly object for the right-hand side (both interior and boundary assembly).
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
This is the MeshBase class.
virtual void attach_F_theta(RBTheta *theta_q_f)
Attach a pointer to a functor object that defines one of the theta_q_a terms.
virtual Number evaluate(const RBParameters &)
Evaluate the functor object for the given parameter.
virtual void interior_assembly(FEMContext &c)
Perform the element interior assembly.
virtual void attach_output_theta(std::vector< std::unique_ptr< RBTheta >> &theta_q_l)
Attach a vector of pointers to functor objects that define one of the outputs.
AcousticsRBAssemblyExpansion()
Constructor.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
This class provides all data required for a physics package (e.g.
const std::vector< dof_id_type > & get_dof_indices() const
Accessor for element dof indices.
unsigned int n_points() const
This class is part of the rbOOmit framework.
FEGenericBase< Real > FEBase
This class provides an encapsulated access to all static public member functions of finite element cl...
virtual void interior_assembly(FEMContext &c)
Perform the element interior assembly.
void attach_A_assembly(ElemAssembly *Aq_assembly)
Attach ElemAssembly object for the left-hand side (both interior and boundary assembly).
const DenseVector< Number > & get_elem_residual() const
Const accessor for element residual.
virtual void boundary_assembly(FEMContext &c)
Perform the element boundary assembly.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
This class stores the set of ElemAssembly functor objects that define the "parameter-independent expa...
ElemAssembly provides a per-element (interior and boundary) assembly functionality.
AcousticsInnerProduct acoustics_inner_product
virtual void boundary_assembly(FEMContext &c)
Perform the element boundary assembly.
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...
virtual void attach_output_assembly(std::vector< std::unique_ptr< ElemAssembly >> &output_assembly)
Attach ElemAssembly object for an output (both interior and boundary assembly).
virtual void interior_assembly(FEMContext &c)
Perform the element interior assembly.
This class is part of the rbOOmit framework.
virtual void boundary_assembly(FEMContext &c)
Perform the element boundary assembly.
A Point defines a location in LIBMESH_DIM dimensional Real space.
const QBase & get_element_qrule() const
Accessor for element interior quadrature rule for the dimension of the current _elem.
const QBase & get_side_qrule() const
Accessor for element side quadrature rule for the dimension of the current _elem. ...
virtual void attach_A_theta(RBTheta *theta_q_a)
Attach a pointer to a functor object that defines one of the theta_q_a terms.