5 #include "libmesh/libmesh.h" 6 #include "libmesh/mesh.h" 7 #include "libmesh/equation_systems.h" 8 #include "libmesh/fe.h" 9 #include "libmesh/quadrature.h" 10 #include "libmesh/dof_map.h" 11 #include "libmesh/dense_matrix.h" 12 #include "libmesh/dense_vector.h" 13 #include "libmesh/fe_interface.h" 14 #include "libmesh/elem.h" 15 #include "libmesh/utility.h" 16 #include "libmesh/print_trace.h" 17 #include "libmesh/int_range.h" 18 #include "libmesh/bounding_box.h" 21 #include "libmesh/rb_assembly_expansion.h" 22 #include "libmesh/rb_eim_theta.h" 23 #include "libmesh/rb_parametrized_function.h" 58 virtual std::vector<Number>
64 const std::vector<Point> & ,
65 const std::vector<Real> & )
override 81 ret[i] = std::exp(-2. * (pow<2>(center_x - p(0)) + pow<2>(center_y - p(1))));
107 "You should only call the evaluate_vec() API when using multi-sample RBParameters objects.");
119 virtual std::vector<Number>
evaluate_vec(
const std::vector<RBParameters> & mus)
override 127 unsigned int count = 0;
128 for (
const auto & mu : mus)
129 count += mu.n_samples();
131 return std::vector<Number>(count,
_val);
143 const unsigned int u_var = 0;
145 FEBase * elem_fe =
nullptr;
148 const std::vector<Real> & JxW = elem_fe->get_JxW();
152 const std::vector<std::vector<RealGradient>> & dphi = elem_fe->get_dphi();
160 for (
unsigned int qp=0; qp != n_qpoints; qp++)
161 for (
unsigned int i=0; i != n_u_dofs; i++)
162 for (
unsigned int j=0; j != n_u_dofs; j++)
173 const unsigned int u_var = 0;
175 FEBase * elem_fe =
nullptr;
178 const std::vector<Real> & JxW = elem_fe->get_JxW();
180 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
186 for (
unsigned int qp=0; qp != n_qpoints; qp++)
187 for (
unsigned int i=0; i != n_u_dofs; i++)
188 for (
unsigned int j=0; j != n_u_dofs; j++)
196 unsigned int basis_function_index_in) :
198 basis_function_index_in)
204 const unsigned int u_var = 0;
206 FEBase * elem_fe =
nullptr;
210 const unsigned int eim_var = 0;
212 const std::vector<Real> & JxW = elem_fe->get_JxW();
214 const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
219 std::vector<Number> eim_values;
225 for (
unsigned int i=0; i != n_u_dofs; i++)
244 const unsigned int u_var = 0;
249 const std::vector<Real> & JxW = fe->get_JxW();
250 const std::vector<std::vector<Real>> & phi = fe->get_phi();
262 for (
unsigned int qp=0; qp != n_qpoints; qp++)
263 for (
unsigned int i=0; i != n_u_dofs; i++)
RealVectorValue RealGradient
const DenseMatrix< Number > & get_elem_jacobian() const
Const accessor for element Jacobian.
bool contains_point(const Point &) const
const Elem & get_elem() const
Accessor for current Elem object.
A simple functor class that provides a RBParameter-dependent function.
TestClass subdomain_id_type
Based on the 4-byte comment warning above, this probably doesn't work with exodusII at all...
This class is part of the rbOOmit framework.
This is the base class from which all geometric element types are derived.
void evaluate_basis_function(dof_id_type elem_id, unsigned int var, std::vector< Number > &values)
Return the basis function values for all quadrature points for variable var on element elem_id...
EimTestRBAssemblyExpansion()
Constructor.
virtual Number evaluate(const RBParameters &mu) override
Evaluate theta for a single scalar-valued RBParameters object.
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
EimTestRBThetaExpansion()
Constructor.
This class provides functionality required to define an assembly object that arises from an "Empirica...
virtual void interior_assembly(FEMContext &c)
Perform the element interior assembly.
EIM_F(RBEIMConstruction &rb_eim_con_in, unsigned int basis_function_index_in)
ThetaConstant output_theta
const Point & min() const
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.
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
OutputAssembly(const BoundingBox &bbox_in)
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).
Defines a Cartesian bounding box by the two corner extremum.
unsigned int n_samples() const
Returns the number of samples stored for all parameters.
const DenseVector< Number > & get_elem_residual() const
Const accessor for element residual.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real get_sample_value(const std::string ¶m_name, std::size_t sample_idx) const
Get the value of the specified parameter at the specified sample, throwing an error if it does not ex...
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.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
virtual std::vector< Number > evaluate_vec(const std::vector< RBParameters > &mus) override
Evaluate theta for multiple mu values, each of which may have multiple "samples". ...
unsigned int get_n_components() const override
Specify the number of components in this parametrized function.
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.
ThetaConstant(Number val)
Constructor.
const Point & max() const
virtual std::vector< Number > evaluate(const RBParameters &mu, const Point &p, dof_id_type, unsigned int, subdomain_id_type, const std::vector< Point > &, const std::vector< Real > &) override
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.
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 interior_assembly(FEMContext &c)
Perform the element interior assembly.
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.