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" 44 const unsigned int u_var = 0;
46 FEBase * elem_fe =
nullptr;
49 const std::vector<Real> & JxW = elem_fe->get_JxW();
53 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
61 for (
unsigned int qp=0; qp != n_qpoints; qp++)
62 for (
unsigned int i=0; i != n_u_dofs; i++)
63 for (
unsigned int j=0; j != n_u_dofs; j++)
74 const unsigned int u_var = 0;
76 FEBase * elem_fe =
nullptr;
79 const std::vector<Real> & JxW = elem_fe->get_JxW();
81 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
83 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
91 for (
unsigned int qp=0; qp != n_qpoints; qp++)
92 for (
unsigned int i=0; i != n_u_dofs; i++)
93 for (
unsigned int j=0; j != n_u_dofs; j++)
103 const unsigned int u_var = 0;
105 FEBase * elem_fe =
nullptr;
108 const std::vector<Real> & JxW = elem_fe->get_JxW();
110 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
112 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
120 for (
unsigned int qp=0; qp != n_qpoints; qp++)
121 for (
unsigned int i=0; i != n_u_dofs; i++)
122 for (
unsigned int j=0; j != n_u_dofs; j++)
132 const unsigned int u_var = 0;
134 FEBase * elem_fe =
nullptr;
137 const std::vector<Real> & JxW = elem_fe->get_JxW();
139 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
147 for (
unsigned int qp=0; qp != n_qpoints; qp++)
148 for (
unsigned int i=0; i != n_u_dofs; i++)
167 const unsigned int u_var = 0;
169 FEBase * elem_fe =
nullptr;
172 const std::vector<Real> & JxW = elem_fe->get_JxW();
174 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
187 for (
unsigned int qp=0; qp != n_qpoints; qp++)
188 for (
unsigned int i=0; i != n_u_dofs; i++)
233 L0(0.7, 0.8, 0.7, 0.8),
234 L1(0.2, 0.3, 0.7, 0.8),
235 L2(0.2, 0.3, 0.2, 0.3),
236 L3(0.7, 0.8, 0.2, 0.3)
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.
const Elem & get_elem() const
Accessor for current Elem object.
virtual Number evaluate(const RBParameters &mu)
Evaluate the functor object for the given parameter.
CDRBAssemblyExpansion()
Constructor.
OutputAssembly(Real min_x_in, Real max_x_in, Real min_y_in, Real max_y_in)
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" ...
CDRBThetaExpansion()
Constructor.
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 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.
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.
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.
virtual void interior_assembly(FEMContext &c)
Perform the element interior 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 interior_assembly(FEMContext &c)
Perform the element interior assembly.
Output assembly object which computes the average value of the solution variable inside a user-provid...
A Point defines a location in LIBMESH_DIM dimensional Real space.
virtual Number evaluate(const RBParameters &)
Evaluate the functor object for the given parameter.
const QBase & get_element_qrule() const
Accessor for element interior quadrature rule for the dimension of the current _elem.
Point vertex_average() const
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.