@@ -71,7 +71,9 @@ std::string solver_element_name(consts::ElementType eType)
7171// / choice of basis family and polynomial order for each element type
7272// / (basis_factory::default_basis_request). The switch deliberately has no
7373// / default case so that compilers building with -Wswitch flag any newly added
74- // / solver element type that is missing a mapping here.
74+ // / solver element type that is missing a mapping here. Returns std::nullopt for
75+ // / element types the FE Basis does not implement (NA/PNT/NRB); callers test FE
76+ // / Basis support with has_value().
7577std::optional<fe::ElementType> to_fe_element_type (consts::ElementType eType)
7678{
7779 switch (eType) {
@@ -99,11 +101,12 @@ std::optional<fe::ElementType> to_fe_element_type(consts::ElementType eType)
99101 return std::nullopt ;
100102}
101103
102- bool use_basis_adapter_for (consts::ElementType eType)
103- {
104- return to_fe_element_type (eType).has_value ();
105- }
106-
104+ // / Whether the FE Basis face adapter can evaluate face shape functions for
105+ // / eType. An element face is always a point, line, or surface topology, so the
106+ // / switch restricts support to those types (a volume element never appears as a
107+ // / face); it then defers to to_fe_element_type to confirm the FE Basis library
108+ // / actually provides a mapping for that face type. The face get_gnn uses this
109+ // / to choose between the FE Basis path and the explicit paths.
107110bool supports_face_basis_adapter_for (consts::ElementType eType)
108111{
109112 switch (eType) {
@@ -114,7 +117,7 @@ bool supports_face_basis_adapter_for(consts::ElementType eType)
114117 case consts::ElementType::QUD4 :
115118 case consts::ElementType::QUD8 :
116119 case consts::ElementType::QUD9 :
117- return use_basis_adapter_for (eType);
120+ return to_fe_element_type (eType). has_value ( );
118121 default :
119122 return false ;
120123 }
@@ -145,6 +148,17 @@ const febasis::BasisFunction& basis_for_solver_element(consts::ElementType eType
145148 return *it->second ;
146149}
147150
151+ // / Permutation from a solver element's local node ordering to the FE Basis
152+ // / ReferenceNodeLayout ordering, indexed by the solver-local node number:
153+ // / map[solver_node] is the matching FE Basis node. The solver and the FE Basis
154+ // / library number element nodes with different conventions, so this table
155+ // / reconciles them at the adapter boundary. An empty span means the two
156+ // / orderings already coincide (identity) and no permutation is applied, which
157+ // / holds for every element type not listed below (lines, Quad4/8/9, Hex8/20).
158+ // / Wedge6 (WDG) reuses the Triangle6 table: its two triangular node triples are
159+ // / reordered exactly like a 6-node triangle.
160+ // / \note These tables must stay consistent with the FE Basis lattice ordering;
161+ // / a mismatch would silently assign shape functions to the wrong nodes.
148162std::span<const std::size_t > solver_to_basis_node_map (consts::ElementType eType)
149163{
150164 static constexpr std::array<std::size_t , 3 > tri3{1 , 2 , 0 };
@@ -173,6 +187,10 @@ std::span<const std::size_t> solver_to_basis_node_map(consts::ElementType eType)
173187 }
174188}
175189
190+ // / Map a single solver-local node index to its FE Basis node index for eType by
191+ // / applying solver_to_basis_node_map (identity when no permutation is
192+ // / registered). Throws BasisNodeOrderingException when solver_node is negative
193+ // / or falls outside the element's node map.
176194std::size_t basis_index_for_solver_node (consts::ElementType eType, const int solver_node)
177195{
178196 if (solver_node < 0 ) {
@@ -194,6 +212,10 @@ std::size_t basis_index_for_solver_node(consts::ElementType eType, const int sol
194212 " is outside node map for " + solver_element_name (eType));
195213}
196214
215+ // / Build a 3-component FE Basis reference coordinate from column g of the solver
216+ // / xi array, zero-filling the trailing components that are inactive for
217+ // / lower-dimensional elements. Throws BasisConfigurationException when xi has
218+ // / fewer rows than the basis reference dimension.
197219fe::math::Vector<fe::Real, 3 > make_basis_point (const febasis::BasisFunction& basis,
198220 const int g,
199221 const Array<double >& xi)
@@ -214,6 +236,10 @@ fe::math::Vector<fe::Real, 3> make_basis_point(const febasis::BasisFunction& bas
214236 return point;
215237}
216238
239+ // / Scatter FE Basis values and gradients (in ReferenceNodeLayout order) into the
240+ // / solver N and Nx arrays at Gauss point g, permuting into solver node order via
241+ // / basis_index_for_solver_node. Validates the value and gradient counts against
242+ // / eNoN and zeroes unused gradient rows.
217243void copy_basis_values_to_solver_arrays (consts::ElementType eType,
218244 const int eNoN,
219245 const int g,
@@ -254,6 +280,9 @@ void copy_basis_values_to_solver_arrays(consts::ElementType eType,
254280 }
255281}
256282
283+ // / Evaluate the cached FE Basis for eType at Gauss point g and write the solver
284+ // / N and Nx arrays. Nx holds reference-space gradients only; physical-coordinate
285+ // / derivatives are formed later by the solver from the mapping Jacobian.
257286void evaluate_basis_values_and_gradients (const int insd,
258287 consts::ElementType eType,
259288 const int eNoN,
@@ -279,6 +308,8 @@ void evaluate_basis_values_and_gradients(const int insd,
279308 copy_basis_values_to_solver_arrays (eType, eNoN, g, values, gradients, N, Nx);
280309}
281310
311+ // / evaluate_basis_values_and_gradients specialized to a faceType, using the
312+ // / face's own reference dimension (xi rows) and N/Nx storage.
282313void evaluate_face_basis_values_and_gradients (const int gaus_pt, faceType& face)
283314{
284315 evaluate_basis_values_and_gradients (
@@ -291,6 +322,9 @@ void evaluate_face_basis_values_and_gradients(const int gaus_pt, faceType& face)
291322 face.Nx );
292323}
293324
325+ // / Number of packed second-derivative components the solver Nxx stores for a
326+ // / given reference dimension: 1 in 1D, 3 in 2D, 6 in 3D. Throws
327+ // / BasisConfigurationException for any other dimension.
294328int required_nxx_components_for_dimension (const int dimension)
295329{
296330 switch (dimension) {
@@ -306,6 +340,10 @@ int required_nxx_components_for_dimension(const int dimension)
306340 }
307341}
308342
343+ // / Scatter FE Basis Hessians (in ReferenceNodeLayout order) into the packed
344+ // / solver Nxx array at Gauss point g, permuting into solver node order. Packing
345+ // / is [dxx, dyy, dxy] in 2D and [dxx, dyy, dzz, dxy, dyz, dxz] in 3D. Validates
346+ // / the Hessian count against eNoN and the Nxx row count against the dimension.
309347void copy_basis_hessians_to_solver_nxx (consts::ElementType eType,
310348 const int eNoN,
311349 const int g,
@@ -354,6 +392,9 @@ void copy_basis_hessians_to_solver_nxx(consts::ElementType eType,
354392 }
355393}
356394
395+ // / Evaluate the cached FE Basis Hessians for eType at Gauss point gaus_pt and
396+ // / write the packed solver Nxx array. Validates insd and ind2 against the basis
397+ // / reference dimension and the required packed-component count.
357398void evaluate_basis_hessians (const int insd,
358399 const int ind2,
359400 consts::ElementType eType,
@@ -384,6 +425,8 @@ void evaluate_basis_hessians(const int insd,
384425 copy_basis_hessians_to_solver_nxx (eType, eNoN, gaus_pt, basis.dimension (), hessians, Nxx);
385426}
386427
428+ // / Shape data for a point (0-D) face: a single unit basis value with zero
429+ // / derivatives. Used for the PNT face case, which has no FE Basis evaluator.
387430void set_point_face_shape_data (const int gaus_pt, faceType& face)
388431{
389432 face.N (0 , gaus_pt) = 1.0 ;
@@ -438,7 +481,7 @@ void get_gip(Simulation* simulation, faceType& face)
438481void get_gnn (const int insd, consts::ElementType eType, const int eNoN, const int g, Array<double >& xi,
439482 Array<double >& N, Array3<double >& Nx)
440483{
441- if (!use_basis_adapter_for (eType)) {
484+ if (!to_fe_element_type (eType). has_value ( )) {
442485 fe::raise<febasis::BasisElementCompatibilityException>(SVMP_HERE ,
443486 " [get_gnn] FE Basis does not support solver element " + solver_element_name (eType));
444487 }
@@ -505,7 +548,7 @@ void get_gn_nxx(const int insd, const int ind2, consts::ElementType eType, const
505548 return ;
506549 }
507550
508- if (!use_basis_adapter_for (eType)) {
551+ if (!to_fe_element_type (eType). has_value ( )) {
509552 fe::raise<febasis::BasisElementCompatibilityException>(SVMP_HERE ,
510553 " [get_gn_nxx] FE Basis Hessian evaluation does not support solver element " +
511554 solver_element_name (eType));
0 commit comments