21 #include "libmesh/dof_map.h" 22 #include "libmesh/elem.h" 23 #include "libmesh/enum_to_string.h" 24 #include "libmesh/fe.h" 25 #include "libmesh/fe_interface.h" 26 #include "libmesh/fe_macro.h" 27 #include "libmesh/tensor_value.h" 42 void nedelec_one_nodal_soln(
const Elem * elem,
44 const std::vector<Number> & elem_soln,
46 std::vector<Number> & nodal_soln,
47 const bool add_p_level)
68 libmesh_assert_equal_to (elem_soln.size(), 3);
70 if (elem_type ==
TRI6)
71 libmesh_assert_equal_to (nodal_soln.size(), 6*2);
73 libmesh_assert_equal_to (nodal_soln.size(), 7*2);
79 libmesh_assert_equal_to (elem_soln.size(), 4);
81 if (elem_type ==
QUAD8)
82 libmesh_assert_equal_to (nodal_soln.size(), 8*2);
84 libmesh_assert_equal_to (nodal_soln.size(), 9*2);
90 libmesh_assert_equal_to (elem_soln.size(), 6);
91 if (elem_type ==
TET10)
92 libmesh_assert_equal_to (nodal_soln.size(), 10*3);
94 libmesh_assert_equal_to (nodal_soln.size(), 14*3);
101 libmesh_assert_equal_to (elem_soln.size(), 12);
103 if (elem_type ==
HEX20)
104 libmesh_assert_equal_to (nodal_soln.size(), 20*3);
106 libmesh_assert_equal_to (nodal_soln.size(), 27*3);
112 libmesh_error_msg(
"ERROR: Invalid ElemType " <<
Utility::enum_to_string(elem_type) <<
" selected for NEDELEC_ONE FE family!");
116 const unsigned int n_sf =
119 std::vector<Point> refspace_nodes;
121 libmesh_assert_equal_to (refspace_nodes.size(),
n_nodes);
128 const std::vector<std::vector<RealGradient>> & vis_phi = vis_fe->get_phi();
130 vis_fe->reinit(elem,&refspace_nodes);
132 for (
unsigned int n = 0; n <
n_nodes; n++)
134 libmesh_assert_equal_to (elem_soln.size(), n_sf);
137 for (
int d = 0; d <
dim; d++)
139 nodal_soln[
dim*n+d] = 0;
143 for (
unsigned int i=0; i<n_sf; i++)
145 for (
int d = 0; d <
dim; d++)
147 nodal_soln[
dim*n+d] += elem_soln[i]*(vis_phi[i][n](d));
156 libmesh_error_msg(
"ERROR: Invalid total order " <<
Utility::enum_to_string(totalorder) <<
" selected for NEDELEC_ONE FE family!");
164 unsigned int nedelec_one_n_dofs(
const ElemType t,
const Order o)
197 libmesh_error_msg(
"ERROR: Invalid Order " <<
Utility::enum_to_string(o) <<
" selected for NEDELEC_ONE FE family!");
204 unsigned int nedelec_one_n_dofs_at_node(
const ElemType t,
206 const unsigned int n)
228 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
246 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
265 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
285 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
306 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
331 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
362 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
400 libmesh_error_msg(
"ERROR: Invalid node ID " << n);
413 libmesh_error_msg(
"ERROR: Invalid Order " <<
Utility::enum_to_string(o) <<
" selected for NEDELEC_ONE FE family!");
419 #ifdef LIBMESH_ENABLE_AMR 423 const Elem * libmesh_dbg_var(elem),
432 libmesh_not_implemented();
434 #endif // #ifdef LIBMESH_ENABLE_AMR 438 #define NEDELEC_LOW_D_ERROR_MESSAGE \ 439 libmesh_error_msg("ERROR: This method makes no sense for low-D elements!"); 447 const std::vector<Number> &,
448 std::vector<Number> &,
450 { NEDELEC_LOW_D_ERROR_MESSAGE }
455 const std::vector<Number> &,
456 std::vector<Number> &,
458 { NEDELEC_LOW_D_ERROR_MESSAGE }
463 const std::vector<Number> & elem_soln,
464 std::vector<Number> & nodal_soln,
465 const bool add_p_level)
466 { nedelec_one_nodal_soln(elem, order, elem_soln, 2 , nodal_soln, add_p_level); }
471 const std::vector<Number> & elem_soln,
472 std::vector<Number> & nodal_soln,
473 const bool add_p_level)
474 { nedelec_one_nodal_soln(elem, order, elem_soln, 3 , nodal_soln, add_p_level); }
521 #ifdef LIBMESH_ENABLE_AMR 527 { NEDELEC_LOW_D_ERROR_MESSAGE }
534 { NEDELEC_LOW_D_ERROR_MESSAGE }
539 const unsigned int variable_number,
541 { nedelec_one_compute_constraints(constraints, dof_map, variable_number, elem, 2); }
546 const unsigned int variable_number,
548 { nedelec_one_compute_constraints(constraints, dof_map, variable_number, elem, 3); }
549 #endif // LIBMESH_ENABLE_AMR 554 { NEDELEC_LOW_D_ERROR_MESSAGE }
557 { NEDELEC_LOW_D_ERROR_MESSAGE }
560 const unsigned int,
const Point &)
561 { NEDELEC_LOW_D_ERROR_MESSAGE }
564 const unsigned int,
const Point &,
const bool)
565 { NEDELEC_LOW_D_ERROR_MESSAGE }
567 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 570 const unsigned int,
const Point &)
571 { NEDELEC_LOW_D_ERROR_MESSAGE }
574 const unsigned int,
const Point &,
const bool)
575 { NEDELEC_LOW_D_ERROR_MESSAGE }
581 { NEDELEC_LOW_D_ERROR_MESSAGE }
584 { NEDELEC_LOW_D_ERROR_MESSAGE }
587 const unsigned int,
const Point &)
588 { NEDELEC_LOW_D_ERROR_MESSAGE }
591 const unsigned int,
const Point &,
const bool)
592 { NEDELEC_LOW_D_ERROR_MESSAGE }
594 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 597 const unsigned int,
const Point &)
598 { NEDELEC_LOW_D_ERROR_MESSAGE }
601 const unsigned int,
const Point &,
const bool)
602 { NEDELEC_LOW_D_ERROR_MESSAGE }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
static unsigned int n_dofs(const ElemType t, const Order o)
ElemType
Defines an enum for geometric element types.
Order
defines an enum for polynomial orders.
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
This is the base class from which all geometric element types are derived.
static unsigned int n_dofs_at_node(const ElemType t, const Order o, const unsigned int n)
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
unsigned int p_level() const
virtual bool shapes_need_reinit() const override
The libMesh namespace provides an interface to certain functionality in the library.
virtual bool is_hierarchic() const override
This class handles the numbering of degrees of freedom on a mesh.
LIBMESH_DEFAULT_VECTORIZED_FE(template<>Real FE< 0, BERNSTEIN)
LIBMESH_FE_SIDE_NODAL_SOLN(HIERARCHIC_VEC)
const dof_id_type n_nodes
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
virtual unsigned int n_nodes() const =0
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
static void nodal_soln(const Elem *elem, const Order o, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, bool add_p_level=true)
Build the nodal soln from the element soln.
static unsigned int n_dofs_per_elem(const ElemType t, const Order o)
virtual FEContinuity get_continuity() const override
std::string enum_to_string(const T e)
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
virtual ElemType type() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
The constraint matrix storage format.
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)