Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 :
19 :
20 : #ifndef LIBMESH_FE_BASE_H
21 : #define LIBMESH_FE_BASE_H
22 :
23 : // Local includes
24 : #include "libmesh/libmesh_common.h"
25 : #include "libmesh/compare_types.h"
26 : #include "libmesh/fe_abstract.h"
27 : #include "libmesh/fe_transformation_base.h"
28 : #include "libmesh/point.h"
29 : #include "libmesh/reference_counted_object.h"
30 : #include "libmesh/tensor_tools.h"
31 : #include "libmesh/type_n_tensor.h"
32 : #include "libmesh/vector_value.h"
33 : #include "libmesh/dense_matrix.h"
34 :
35 : // C++ includes
36 : #include <cstddef>
37 : #include <vector>
38 : #include <memory>
39 :
40 : namespace libMesh
41 : {
42 :
43 :
44 : // forward declarations
45 : template <typename T> class DenseMatrix;
46 : template <typename T> class DenseVector;
47 : class BoundaryInfo;
48 : class DofConstraints;
49 : class DofMap;
50 : class Elem;
51 : class MeshBase;
52 : template <typename T> class NumericVector;
53 : class QBase;
54 : template <typename T> class FETransformationBase;
55 : class FEType;
56 :
57 : #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
58 : class NodeConstraints;
59 : #endif
60 :
61 : #ifdef LIBMESH_ENABLE_PERIODIC
62 : class PeriodicBoundaries;
63 : class PointLocatorBase;
64 : #endif
65 :
66 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
67 : template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
68 : class InfFE;
69 : #endif
70 :
71 : /**
72 : * This class forms the foundation from which generic finite
73 : * elements may be derived. In the current implementation the
74 : * templated derived class \p FE offers a wide variety of commonly
75 : * used finite element concepts. Check there for details.
76 : *
77 : * Use the \p FEGenericBase<OutputType>::build() method to create an
78 : * object of any of the derived classes which is compatible with
79 : * OutputType.
80 : *
81 : * \author Benjamin S. Kirk
82 : * \date 2002
83 : */
84 : template <typename OutputType>
85 : class FEGenericBase : public FEAbstract
86 : {
87 : protected:
88 :
89 : /**
90 : * Constructor. Optionally initializes required data
91 : * structures. Protected so that this base class
92 : * cannot be explicitly instantiated.
93 : */
94 : FEGenericBase (const unsigned int dim,
95 : const FEType & fet);
96 :
97 : public:
98 :
99 : /**
100 : * Destructor.
101 : */
102 : virtual ~FEGenericBase();
103 :
104 : /**
105 : * Builds a specific finite element type. A \p
106 : * std::unique_ptr<FEGenericBase> is returned to prevent a memory leak. This
107 : * way the user need not remember to delete the object.
108 : *
109 : * The build call will fail if the OutputType of this class is not
110 : * compatible with the output required for the requested \p type
111 : */
112 : static std::unique_ptr<FEGenericBase> build (const unsigned int dim,
113 : const FEType & type);
114 :
115 : /**
116 : * Convenient typedefs for gradients of output, hessians of output,
117 : * and potentially-complex-valued versions of same.
118 : */
119 : typedef OutputType OutputShape;
120 : typedef typename TensorTools::IncrementRank<OutputShape>::type OutputGradient;
121 : typedef typename TensorTools::IncrementRank<OutputGradient>::type OutputTensor;
122 : typedef typename TensorTools::DecrementRank<OutputShape>::type OutputDivergence;
123 : typedef typename TensorTools::MakeNumber<OutputShape>::type OutputNumber;
124 : typedef typename TensorTools::IncrementRank<OutputNumber>::type OutputNumberGradient;
125 : typedef typename TensorTools::IncrementRank<OutputNumberGradient>::type OutputNumberTensor;
126 : typedef typename TensorTools::DecrementRank<OutputNumber>::type OutputNumberDivergence;
127 :
128 :
129 :
130 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
131 :
132 : /**
133 : * Builds a specific infinite element type. A \p
134 : * std::unique_ptr<FEGenericBase> is returned to prevent a memory leak. This
135 : * way the user need not remember to delete the object.
136 : *
137 : * The build call will fail if the OutputShape of this class is not
138 : * compatible with the output required for the requested \p type
139 : */
140 : static std::unique_ptr<FEGenericBase> build_InfFE (const unsigned int dim,
141 : const FEType & type);
142 :
143 : #endif
144 :
145 : #ifdef LIBMESH_ENABLE_AMR
146 :
147 : /**
148 : * Computes the constraint matrix contributions (for
149 : * non-conforming adapted meshes) corresponding to
150 : * variable number \p var_number, using generic
151 : * projections.
152 : */
153 : static void compute_proj_constraints (DofConstraints & constraints,
154 : DofMap & dof_map,
155 : const unsigned int variable_number,
156 : const Elem * elem);
157 :
158 : /**
159 : * Creates a local projection on \p coarse_elem, based on the
160 : * DoF values in \p global_vector for it's children. Computes a
161 : * vector of coefficients corresponding to dof_indices for only the
162 : * single given \p var
163 : */
164 :
165 : static void coarsened_dof_values(const NumericVector<Number> & global_vector,
166 : const DofMap & dof_map,
167 : const Elem * coarse_elem,
168 : DenseVector<Number> & coarse_dofs,
169 : const unsigned int var,
170 : const bool use_old_dof_indices = false);
171 :
172 : /**
173 : * Creates a local projection on \p coarse_elem, based on the
174 : * DoF values in \p global_vector for it's children. Computes a
175 : * vector of coefficients corresponding to all dof_indices.
176 : */
177 :
178 : static void coarsened_dof_values(const NumericVector<Number> & global_vector,
179 : const DofMap & dof_map,
180 : const Elem * coarse_elem,
181 : DenseVector<Number> & coarse_dofs,
182 : const bool use_old_dof_indices = false);
183 :
184 : #endif // #ifdef LIBMESH_ENABLE_AMR
185 :
186 : #ifdef LIBMESH_ENABLE_PERIODIC
187 :
188 : /**
189 : * Computes the constraint matrix contributions (for
190 : * meshes with periodic boundary conditions) corresponding to
191 : * variable number \p var_number, using generic projections.
192 : */
193 : static void compute_periodic_constraints (DofConstraints & constraints,
194 : DofMap & dof_map,
195 : const PeriodicBoundaries & boundaries,
196 : const MeshBase & mesh,
197 : const PointLocatorBase * point_locator,
198 : const unsigned int variable_number,
199 : const Elem * elem);
200 :
201 : #endif // LIBMESH_ENABLE_PERIODIC
202 :
203 : /**
204 : * \returns The shape function values at the quadrature points
205 : * on the element.
206 : */
207 21096 : const std::vector<std::vector<OutputShape>> & get_phi() const
208 21096 : { libmesh_assert(!calculations_started || calculate_phi);
209 315998006 : calculate_phi = true; return phi; }
210 :
211 0 : const std::vector<std::vector<OutputShape>> & get_dual_phi() const
212 : {
213 0 : libmesh_assert(!calculations_started || calculate_dual);
214 0 : calculate_dual = true;
215 : // Dual phi computation relies on primal phi computation
216 0 : this->request_phi();
217 0 : return dual_phi;
218 : }
219 :
220 835107 : virtual void request_phi() const override
221 835107 : { get_phi(); }
222 :
223 0 : virtual void request_dual_phi() const override
224 0 : { get_dual_phi(); }
225 :
226 : /**
227 : * \returns The shape function derivatives at the quadrature
228 : * points.
229 : */
230 120021 : const std::vector<std::vector<OutputGradient>> & get_dphi() const
231 120021 : { libmesh_assert(!calculations_started || calculate_dphi);
232 274662976 : calculate_dphi = calculate_dphiref = true; return dphi; }
233 :
234 0 : const std::vector<std::vector<OutputGradient>> & get_dual_dphi() const
235 0 : { libmesh_assert(!calculations_started || calculate_dphi);
236 0 : calculate_dphi = calculate_dual = calculate_dphiref = true; return dual_dphi; }
237 :
238 98759 : virtual void request_dphi() const override
239 98759 : { get_dphi(); }
240 :
241 0 : virtual void request_dual_dphi() const override
242 0 : { get_dual_dphi(); }
243 :
244 0 : const DenseMatrix<Real> & get_dual_coeff() const
245 0 : { return dual_coeff; }
246 :
247 : /**
248 : * \returns The curl of the shape function at the quadrature
249 : * points.
250 : */
251 : virtual_for_inffe
252 642716 : const std::vector<std::vector<OutputShape>> & get_curl_phi() const
253 245740 : { libmesh_assert(!calculations_started || calculate_curl_phi);
254 2607164 : calculate_curl_phi = calculate_dphiref = true; return curl_phi; }
255 :
256 : /**
257 : * \returns The divergence of the shape function at the quadrature
258 : * points.
259 : */
260 : virtual_for_inffe
261 177876 : const std::vector<std::vector<OutputDivergence>> & get_div_phi() const
262 59412 : { libmesh_assert(!calculations_started || calculate_div_phi);
263 710964 : calculate_div_phi = calculate_dphiref = true; return div_phi; }
264 :
265 : /**
266 : * \returns The shape function x-derivative at the quadrature
267 : * points.
268 : */
269 0 : const std::vector<std::vector<OutputShape>> & get_dphidx() const
270 0 : { libmesh_assert(!calculations_started || calculate_dphi);
271 0 : calculate_dphi = calculate_dphiref = true; return dphidx; }
272 :
273 : /**
274 : * \returns The shape function y-derivative at the quadrature
275 : * points.
276 : */
277 0 : const std::vector<std::vector<OutputShape>> & get_dphidy() const
278 0 : { libmesh_assert(!calculations_started || calculate_dphi);
279 0 : calculate_dphi = calculate_dphiref = true; return dphidy; }
280 :
281 : /**
282 : * \returns The shape function z-derivative at the quadrature
283 : * points.
284 : */
285 0 : const std::vector<std::vector<OutputShape>> & get_dphidz() const
286 0 : { libmesh_assert(!calculations_started || calculate_dphi);
287 0 : calculate_dphi = calculate_dphiref = true; return dphidz; }
288 :
289 : /**
290 : * \returns The shape function xi-derivative at the quadrature
291 : * points.
292 : */
293 8026275 : const std::vector<std::vector<OutputShape>> & get_dphidxi() const
294 8026275 : { libmesh_assert(!calculations_started || calculate_dphiref);
295 87428895 : calculate_dphiref = true; return dphidxi; }
296 :
297 : /**
298 : * \returns The shape function eta-derivative at the quadrature
299 : * points.
300 : */
301 8005115 : const std::vector<std::vector<OutputShape>> & get_dphideta() const
302 8005115 : { libmesh_assert(!calculations_started || calculate_dphiref);
303 8005115 : calculate_dphiref = true; return dphideta; }
304 :
305 : /**
306 : * \returns The shape function zeta-derivative at the quadrature
307 : * points.
308 : */
309 1246460 : const std::vector<std::vector<OutputShape>> & get_dphidzeta() const
310 1246460 : { libmesh_assert(!calculations_started || calculate_dphiref);
311 1246460 : calculate_dphiref = true; return dphidzeta; }
312 :
313 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
314 :
315 : /**
316 : * \returns The shape function second derivatives at the quadrature
317 : * points.
318 : */
319 176228 : const std::vector<std::vector<OutputTensor>> & get_d2phi() const
320 176228 : { libmesh_assert(!calculations_started || calculate_d2phi);
321 2226351 : calculate_d2phi = calculate_dphiref = true; return d2phi; }
322 :
323 0 : const std::vector<std::vector<OutputTensor>> & get_dual_d2phi() const
324 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
325 0 : calculate_d2phi = calculate_dual = calculate_dphiref = true; return dual_d2phi; }
326 :
327 : /**
328 : * \returns The shape function second derivatives at the quadrature
329 : * points.
330 : */
331 0 : const std::vector<std::vector<OutputShape>> & get_d2phidx2() const
332 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
333 0 : calculate_d2phi = calculate_dphiref = true; return d2phidx2; }
334 :
335 : /**
336 : * \returns The shape function second derivatives at the quadrature
337 : * points.
338 : */
339 0 : const std::vector<std::vector<OutputShape>> & get_d2phidxdy() const
340 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
341 0 : calculate_d2phi = calculate_dphiref = true; return d2phidxdy; }
342 :
343 : /**
344 : * \returns The shape function second derivatives at the quadrature
345 : * points.
346 : */
347 0 : const std::vector<std::vector<OutputShape>> & get_d2phidxdz() const
348 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
349 0 : calculate_d2phi = calculate_dphiref = true; return d2phidxdz; }
350 :
351 : /**
352 : * \returns The shape function second derivatives at the quadrature
353 : * points.
354 : */
355 0 : const std::vector<std::vector<OutputShape>> & get_d2phidy2() const
356 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
357 0 : calculate_d2phi = calculate_dphiref = true; return d2phidy2; }
358 :
359 : /**
360 : * \returns The shape function second derivatives at the quadrature
361 : * points.
362 : */
363 0 : const std::vector<std::vector<OutputShape>> & get_d2phidydz() const
364 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
365 0 : calculate_d2phi = calculate_dphiref = true; return d2phidydz; }
366 :
367 : /**
368 : * \returns The shape function second derivatives at the quadrature
369 : * points.
370 : */
371 0 : const std::vector<std::vector<OutputShape>> & get_d2phidz2() const
372 0 : { libmesh_assert(!calculations_started || calculate_d2phi);
373 0 : calculate_d2phi = calculate_dphiref = true; return d2phidz2; }
374 :
375 : /**
376 : * \returns The shape function second derivatives at the quadrature
377 : * points, in reference coordinates
378 : */
379 521543 : const std::vector<std::vector<OutputShape>> & get_d2phidxi2() const
380 521543 : { libmesh_assert(!calculations_started || calculate_d2phi);
381 6568893 : calculate_d2phi = calculate_dphiref = true; return d2phidxi2; }
382 :
383 : /**
384 : * \returns The shape function second derivatives at the quadrature
385 : * points, in reference coordinates
386 : */
387 513883 : const std::vector<std::vector<OutputShape>> & get_d2phidxideta() const
388 513883 : { libmesh_assert(!calculations_started || calculate_d2phi);
389 513883 : calculate_d2phi = calculate_dphiref = true; return d2phidxideta; }
390 :
391 : /**
392 : * \returns The shape function second derivatives at the quadrature
393 : * points, in reference coordinates
394 : */
395 432120 : const std::vector<std::vector<OutputShape>> & get_d2phidxidzeta() const
396 432120 : { libmesh_assert(!calculations_started || calculate_d2phi);
397 432120 : calculate_d2phi = calculate_dphiref = true; return d2phidxidzeta; }
398 :
399 : /**
400 : * \returns The shape function second derivatives at the quadrature
401 : * points, in reference coordinates
402 : */
403 513883 : const std::vector<std::vector<OutputShape>> & get_d2phideta2() const
404 513883 : { libmesh_assert(!calculations_started || calculate_d2phi);
405 513883 : calculate_d2phi = calculate_dphiref = true; return d2phideta2; }
406 :
407 : /**
408 : * \returns The shape function second derivatives at the quadrature
409 : * points, in reference coordinates
410 : */
411 432120 : const std::vector<std::vector<OutputShape>> & get_d2phidetadzeta() const
412 432120 : { libmesh_assert(!calculations_started || calculate_d2phi);
413 432120 : calculate_d2phi = calculate_dphiref = true; return d2phidetadzeta; }
414 :
415 : /**
416 : * \returns The shape function second derivatives at the quadrature
417 : * points, in reference coordinates
418 : */
419 432120 : const std::vector<std::vector<OutputShape>> & get_d2phidzeta2() const
420 432120 : { libmesh_assert(!calculations_started || calculate_d2phi);
421 432120 : calculate_d2phi = calculate_dphiref = true; return d2phidzeta2; }
422 :
423 : #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
424 :
425 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
426 :
427 : /**
428 : * \returns The global first derivative of the phase term
429 : * which is used in infinite elements, evaluated at the
430 : * quadrature points.
431 : *
432 : * In case of the general finite element class \p FE this
433 : * field is initialized to all zero, so that the variational
434 : * formulation for an infinite element produces correct element
435 : * matrices for a mesh using both finite and infinite elements.
436 : */
437 0 : const std::vector<OutputGradient> & get_dphase() const
438 0 : { return dphase; }
439 :
440 :
441 : /**
442 : * \returns The multiplicative weight at each quadrature point.
443 : * This weight is used for certain infinite element weak
444 : * formulations, so that weighted Sobolev spaces are
445 : * used for the trial function space. This renders the
446 : * variational form easily computable.
447 : *
448 : * In case of the general finite element class \p FE this
449 : * field is initialized to all ones, so that the variational
450 : * formulation for an infinite element produces correct element
451 : * matrices for a mesh using both finite and infinite elements.
452 : */
453 2 : virtual const std::vector<Real> & get_Sobolev_weight() const
454 2 : { return weight; }
455 :
456 : /**
457 : * \returns The first global derivative of the multiplicative
458 : * weight at each quadrature point. See \p get_Sobolev_weight()
459 : * for details. In case of \p FE initialized to all zero.
460 : */
461 0 : virtual const std::vector<RealGradient> & get_Sobolev_dweight() const
462 0 : { return dweight; }
463 :
464 : /**
465 : * \returns The multiplicative weight (see \p get_Sobolev_weight)
466 : * but weighted with the radial coordinate square.
467 : *
468 : * In finite elements, this gives just 1, similar to \p get_Sobolev_Weight()
469 : */
470 1340 : virtual const std::vector<Real> & get_Sobolev_weightxR_sq() const
471 1340 : { return weight; }
472 :
473 : /**
474 : * \returns The first global derivative of the multiplicative weight
475 : * (see \p dget_Sobolev_weight) but weighted with the square of the
476 : * radial coordinate.
477 : *
478 : * In finite elements, this is 0.
479 : */
480 1340 : virtual const std::vector<RealGradient> & get_Sobolev_dweightxR_sq() const
481 1340 : { return dweight; }
482 :
483 : /**
484 : * \returns The shape function \p phi (for FE) and \p phi weighted by r/decay
485 : * for InfFE.
486 : *
487 : * To compensate for the decay function applied to the Jacobian (see \p get_JxWxdecay_sq),
488 : * the wave function \p phi should be divided by this function.
489 : *
490 : * The factor r must be compensated for by the Sobolev \p weight.
491 : * (i.e. by using \p get_Sobolev_weightxR_sq())
492 : **/
493 1340 : virtual const std::vector<std::vector<OutputShape>> & get_phi_over_decayxR () const
494 1340 : { return get_phi();}
495 :
496 : /**
497 : * \returns the gradient of the shape function (see \p get_dphi()),
498 : * but in case of \p InfFE, weighted with r/decay.
499 : * See \p get_phi_over_decayxR() for details.
500 : */
501 1340 : virtual const std::vector<std::vector<OutputGradient>> & get_dphi_over_decayxR () const
502 1340 : { return get_dphi();}
503 :
504 : /**
505 : * \returns the gradient of the shape function (see \p get_dphi()),
506 : * but in case of \p InfFE, weighted with 1/decay.
507 : *
508 : * In contrast to the shape function, its gradient stays finite
509 : * when divided by the decay function.
510 : */
511 0 : virtual const std::vector<std::vector<OutputGradient>> & get_dphi_over_decay () const
512 0 : { return get_dphi();}
513 :
514 : #endif
515 :
516 : /**
517 : * Prints the value of each shape function at each quadrature point.
518 : */
519 : virtual void print_phi(std::ostream & os) const override;
520 : virtual void print_dual_phi(std::ostream & os) const override;
521 :
522 : /**
523 : * Prints the value of each shape function's derivative
524 : * at each quadrature point.
525 : */
526 : virtual void print_dphi(std::ostream & os) const override;
527 : virtual void print_dual_dphi(std::ostream & os) const override;
528 :
529 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
530 :
531 : /**
532 : * Prints the value of each shape function's second derivatives
533 : * at each quadrature point.
534 : */
535 : virtual void print_d2phi(std::ostream & os) const override;
536 : virtual void print_dual_d2phi(std::ostream & os) const override;
537 :
538 : #endif
539 :
540 :
541 : protected:
542 :
543 :
544 :
545 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
546 :
547 : /**
548 : * Initialize the data fields for the base of an
549 : * an infinite element. Implement this in the derived
550 : * class \p FE<Dim,T>.
551 : */
552 : virtual void init_base_shape_functions(const std::vector<Point> & qp,
553 : const Elem * e) = 0;
554 :
555 : #endif
556 :
557 : /**
558 : * Determine which values are to be calculated, for both the FE
559 : * itself and for the FEMap.
560 : */
561 : virtual_for_inffe
562 : void determine_calculations();
563 :
564 : /**
565 : * \returns true iff no calculations have been requested of this
566 : * FE object or of its associated FEMap
567 : */
568 26448917 : bool calculating_nothing() const
569 : {
570 28251018 : return calculate_nothing &&
571 1802101 : !this->calculate_phi && !this->calculate_dphi &&
572 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
573 1802101 : !this->calculate_d2phi &&
574 : #endif
575 28409303 : !this->calculate_curl_phi && !this->calculate_div_phi &&
576 4190748 : !this->calculate_map;
577 : }
578 :
579 : /**
580 : * After having updated the jacobian and the transformation
581 : * from local to global coordinates in \p FEMap::compute_map(),
582 : * the first derivatives of the shape functions are
583 : * transformed to global coordinates, giving \p dphi,
584 : * \p dphidx, \p dphidy, and \p dphidz. This method
585 : * should rarely be re-defined in derived classes, but
586 : * still should be usable for children. Therefore, keep
587 : * it protected.
588 : */
589 : virtual void compute_shape_functions(const Elem * elem, const std::vector<Point> & qp) override;
590 :
591 : /**
592 : * Compute the dual basis coefficients \p dual_coeff
593 : * we rely on the \p JxW (or weights) and the \p phi values,
594 : * which can come from default or customized qrule
595 : */
596 : void compute_dual_shape_coeffs(const std::vector<Real> & JxW, const std::vector<std::vector<OutputShape>> & phi);
597 :
598 : /**
599 : * Compute \p dual_phi, \p dual_dphi, \p dual_d2phi
600 : * It is only valid for this to be called after reinit has occurred with a
601 : * quadrature rule
602 : */
603 : void compute_dual_shape_functions();
604 :
605 : /**
606 : * Object that handles computing shape function values, gradients, etc
607 : * in the physical domain.
608 : */
609 : std::unique_ptr<FETransformationBase<OutputType>> _fe_trans;
610 :
611 : /**
612 : * Shape function values.
613 : */
614 : std::vector<std::vector<OutputShape>> phi;
615 : std::vector<std::vector<OutputShape>> dual_phi;
616 :
617 : /**
618 : * Shape function derivative values.
619 : */
620 : std::vector<std::vector<OutputGradient>> dphi;
621 : std::vector<std::vector<OutputGradient>> dual_dphi;
622 :
623 : /**
624 : * Coefficient matrix for the dual basis.
625 : */
626 : mutable DenseMatrix<Real> dual_coeff;
627 :
628 : /**
629 : * Shape function curl values. Only defined for vector types.
630 : */
631 : std::vector<std::vector<OutputShape>> curl_phi;
632 :
633 : /**
634 : * Shape function divergence values. Only defined for vector types.
635 : */
636 : std::vector<std::vector<OutputDivergence>> div_phi;
637 :
638 : /**
639 : * Shape function derivatives in the xi direction.
640 : */
641 : std::vector<std::vector<OutputShape>> dphidxi;
642 :
643 : /**
644 : * Shape function derivatives in the eta direction.
645 : */
646 : std::vector<std::vector<OutputShape>> dphideta;
647 :
648 : /**
649 : * Shape function derivatives in the zeta direction.
650 : */
651 : std::vector<std::vector<OutputShape>> dphidzeta;
652 :
653 : /**
654 : * Shape function derivatives in the x direction.
655 : */
656 : std::vector<std::vector<OutputShape>> dphidx;
657 :
658 : /**
659 : * Shape function derivatives in the y direction.
660 : */
661 : std::vector<std::vector<OutputShape>> dphidy;
662 :
663 : /**
664 : * Shape function derivatives in the z direction.
665 : */
666 : std::vector<std::vector<OutputShape>> dphidz;
667 :
668 :
669 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
670 :
671 : /**
672 : * Shape function second derivative values.
673 : */
674 : std::vector<std::vector<OutputTensor>> d2phi;
675 : std::vector<std::vector<OutputTensor>> dual_d2phi;
676 :
677 : /**
678 : * Shape function second derivatives in the xi direction.
679 : */
680 : std::vector<std::vector<OutputShape>> d2phidxi2;
681 :
682 : /**
683 : * Shape function second derivatives in the xi-eta direction.
684 : */
685 : std::vector<std::vector<OutputShape>> d2phidxideta;
686 :
687 : /**
688 : * Shape function second derivatives in the xi-zeta direction.
689 : */
690 : std::vector<std::vector<OutputShape>> d2phidxidzeta;
691 :
692 : /**
693 : * Shape function second derivatives in the eta direction.
694 : */
695 : std::vector<std::vector<OutputShape>> d2phideta2;
696 :
697 : /**
698 : * Shape function second derivatives in the eta-zeta direction.
699 : */
700 : std::vector<std::vector<OutputShape>> d2phidetadzeta;
701 :
702 : /**
703 : * Shape function second derivatives in the zeta direction.
704 : */
705 : std::vector<std::vector<OutputShape>> d2phidzeta2;
706 :
707 : /**
708 : * Shape function second derivatives in the x direction.
709 : */
710 : std::vector<std::vector<OutputShape>> d2phidx2;
711 :
712 : /**
713 : * Shape function second derivatives in the x-y direction.
714 : */
715 : std::vector<std::vector<OutputShape>> d2phidxdy;
716 :
717 : /**
718 : * Shape function second derivatives in the x-z direction.
719 : */
720 : std::vector<std::vector<OutputShape>> d2phidxdz;
721 :
722 : /**
723 : * Shape function second derivatives in the y direction.
724 : */
725 : std::vector<std::vector<OutputShape>> d2phidy2;
726 :
727 : /**
728 : * Shape function second derivatives in the y-z direction.
729 : */
730 : std::vector<std::vector<OutputShape>> d2phidydz;
731 :
732 : /**
733 : * Shape function second derivatives in the z direction.
734 : */
735 : std::vector<std::vector<OutputShape>> d2phidz2;
736 :
737 : #endif
738 :
739 :
740 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
741 :
742 : //--------------------------------------------------------------
743 : /* protected members for infinite elements, which are accessed
744 : * from the outside through some inline functions
745 : */
746 :
747 :
748 : /**
749 : * Used for certain infinite element families:
750 : * the first derivatives of the phase term in global coordinates,
751 : * over all quadrature points.
752 : */
753 : std::vector<OutputGradient> dphase;
754 :
755 : /**
756 : * Used for certain infinite element families:
757 : * the global derivative of the additional radial weight \f$ 1/{r^2} \f$,
758 : * over all quadrature points.
759 : */
760 : std::vector<RealGradient> dweight;
761 :
762 : /**
763 : * Used for certain infinite element families:
764 : * the additional radial weight \f$ 1/{r^2} \f$ in local coordinates,
765 : * over all quadrature points.
766 : */
767 : std::vector<Real> weight;
768 :
769 : #endif
770 :
771 : private:
772 :
773 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
774 :
775 : /**
776 : * Make all \p InfFE<Dim,T_radial,T_map> classes friends
777 : * so that they can safely used \p FE<Dim-1,T_base> through
778 : * a \p FEGenericBase * as base approximation.
779 : */
780 : template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
781 : friend class InfFE;
782 :
783 : #endif
784 :
785 :
786 : };
787 :
788 : // --------------------------------------------------------------------
789 : // Generic templates. We specialize for OutputType = Real, so these are
790 : // only used for OutputType = RealVectorValue
791 : template <typename OutputType>
792 0 : void FEGenericBase<OutputType>::compute_dual_shape_functions ()
793 : {
794 0 : libmesh_error_msg(
795 : "Computation of dual shape functions for vector finite element "
796 : "families is not currently implemented");
797 : }
798 :
799 : template <typename OutputType>
800 0 : void FEGenericBase<OutputType>::compute_dual_shape_coeffs(const std::vector<Real> & /*JxW*/, const std::vector<std::vector<OutputShape>> & /*phi_vals*/)
801 : {
802 0 : libmesh_error_msg(
803 : "Computation of dual shape functions for vector finite element "
804 : "families is not currently implemented");
805 : }
806 :
807 : // -----------------------------------------------------------
808 : // Forward declaration of specialization
809 : template <>
810 : void FEGenericBase<Real>::compute_dual_shape_functions();
811 :
812 : template <>
813 : void FEGenericBase<Real>::compute_dual_shape_coeffs(const std::vector<Real> & /*JxW*/, const std::vector<std::vector<OutputShape>> & /*phi_vals*/);
814 :
815 :
816 : // Typedefs for convenience and backwards compatibility
817 : typedef FEGenericBase<Real> FEBase;
818 : typedef FEGenericBase<RealGradient> FEVectorBase;
819 :
820 :
821 :
822 :
823 : // ------------------------------------------------------------
824 : // FEGenericBase class inline members
825 : template <typename OutputType>
826 : inline
827 15326404 : FEGenericBase<OutputType>::FEGenericBase(const unsigned int d,
828 : const FEType & fet) :
829 : FEAbstract(d,fet),
830 13587717 : _fe_trans( FETransformationBase<OutputType>::build(fet) ),
831 13587717 : phi(),
832 13587717 : dual_phi(),
833 13587717 : dphi(),
834 13587717 : dual_dphi(),
835 : curl_phi(),
836 : div_phi(),
837 : dphidxi(),
838 : dphideta(),
839 : dphidzeta(),
840 : dphidx(),
841 : dphidy(),
842 : dphidz()
843 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
844 : ,d2phi(),
845 : dual_d2phi(),
846 : d2phidxi2(),
847 : d2phidxideta(),
848 : d2phidxidzeta(),
849 : d2phideta2(),
850 : d2phidetadzeta(),
851 : d2phidzeta2(),
852 : d2phidx2(),
853 : d2phidxdy(),
854 : d2phidxdz(),
855 : d2phidy2(),
856 : d2phidydz(),
857 13018191 : d2phidz2()
858 : #endif
859 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
860 : ,dphase(),
861 : dweight(),
862 2308213 : weight()
863 : #endif
864 : {
865 15326404 : }
866 :
867 :
868 :
869 : template <typename OutputType>
870 : inline
871 15326404 : FEGenericBase<OutputType>::~FEGenericBase()
872 : {
873 17634024 : }
874 :
875 : } // namespace libMesh
876 :
877 : #endif // LIBMESH_FE_BASE_H
|
