Reference documentation for deal.II version 8.1.0
fe_base.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_base.h 30271 2013-08-09 22:37:31Z bangerth @f$
3 //
4 // Copyright (C) 2000 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__fe_base_h
18 #define __deal2__fe_base_h
19 
20 #include <deal.II/base/config.h>
22 #include <deal.II/base/subscriptor.h>
23 #include <deal.II/base/point.h>
24 #include <deal.II/base/tensor.h>
25 #include <deal.II/base/table.h>
26 #include <deal.II/base/vector_slice.h>
27 #include <deal.II/base/geometry_info.h>
28 #include <deal.II/lac/full_matrix.h>
29 #include <deal.II/fe/fe_update_flags.h>
30 #include <deal.II/fe/mapping.h>
31 
32 #include <string>
33 #include <vector>
34 
36 
37 template<int dim, int spacedim> class FESystem;
38 
39 
45 {
96  {
97  this_element_dominates,
98  other_element_dominates,
99  neither_element_dominates,
100  either_element_can_dominate,
101  no_requirements
102  };
103 
104 
123  inline Domination operator & (const Domination d1,
124  const Domination d2);
125 }
126 
127 
136 template <int dim>
137 class FiniteElementData
138 {
139 public:
170  {
174  unknown = 0x00,
175 
179  L2 = 0x01,
180 
185  Hcurl = 0x02,
186 
191  Hdiv = 0x04,
192 
196  H1 = Hcurl | Hdiv,
197 
202  H2 = 0x0e
203  };
204 
209  static const unsigned int dimension = dim;
210 
214  const unsigned int dofs_per_vertex;
215 
220  const unsigned int dofs_per_line;
221 
226  const unsigned int dofs_per_quad;
227 
232  const unsigned int dofs_per_hex;
233 
237  const unsigned int first_line_index;
238 
242  const unsigned int first_quad_index;
243 
247  const unsigned int first_hex_index;
248 
252  const unsigned int first_face_line_index;
253 
257  const unsigned int first_face_quad_index;
258 
264  const unsigned int dofs_per_face;
265 
271  const unsigned int dofs_per_cell;
272 
281  const unsigned int components;
282 
287  const unsigned int degree;
288 
293 
300 
307 
320  FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
321  const unsigned int n_components,
322  const unsigned int degree,
323  const Conformity conformity = unknown,
324  const unsigned int n_blocks = numbers::invalid_unsigned_int);
325 
329  unsigned int n_dofs_per_vertex () const;
330 
334  unsigned int n_dofs_per_line () const;
335 
339  unsigned int n_dofs_per_quad () const;
340 
344  unsigned int n_dofs_per_hex () const;
345 
350  unsigned int n_dofs_per_face () const;
351 
356  unsigned int n_dofs_per_cell () const;
357 
366  template <int structdim>
367  unsigned int n_dofs_per_object () const;
368 
373  unsigned int n_components () const;
374 
379  unsigned int n_blocks () const;
380 
384  const BlockIndices &block_indices() const;
385 
394  bool is_primitive () const;
395 
402  unsigned int tensor_degree () const;
403 
410  bool conforms (const Conformity) const;
411 
415  bool operator == (const FiniteElementData &) const;
416 
417 protected:
418 
424  void set_primitivity(const bool value);
425 
426 private:
433 };
434 
435 
436 
437 // --------- inline and template functions ---------------
438 
439 
440 #ifndef DOXYGEN
441 
442 namespace FiniteElementDomination
443 {
444  inline
446  const Domination d2)
447  {
448  // go through the entire list of possibilities. note that if we were into
449  // speed, obfuscation and cared enough, we could implement this operator
450  // by doing a bitwise & (and) if we gave these values to the enum values:
451  // neither_element_dominates=0, this_element_dominates=1,
452  // other_element_dominates=2, either_element_can_dominate=3
453  // =this_element_dominates|other_element_dominates
454  switch (d1)
455  {
456  case this_element_dominates:
457  if ((d2 == this_element_dominates) ||
458  (d2 == either_element_can_dominate) ||
459  (d2 == no_requirements))
460  return this_element_dominates;
461  else
462  return neither_element_dominates;
463 
464  case other_element_dominates:
465  if ((d2 == other_element_dominates) ||
466  (d2 == either_element_can_dominate) ||
467  (d2 == no_requirements))
468  return other_element_dominates;
469  else
470  return neither_element_dominates;
471 
472  case neither_element_dominates:
473  return neither_element_dominates;
474 
475  case either_element_can_dominate:
476  if (d2 == no_requirements)
477  return either_element_can_dominate;
478  else
479  return d2;
480 
481  case no_requirements:
482  return d2;
483 
484  default:
485  // shouldn't get here
486  Assert (false, ExcInternalError());
487  }
488 
489  return neither_element_dominates;
490  }
491 }
492 
493 
494 template <int dim>
495 inline
496 unsigned int
498 {
499  return dofs_per_vertex;
500 }
501 
502 
503 
504 template <int dim>
505 inline
506 unsigned int
508 {
509  return dofs_per_line;
510 }
511 
512 
513 
514 template <int dim>
515 inline
516 unsigned int
518 {
519  return dofs_per_quad;
520 }
521 
522 
523 
524 template <int dim>
525 inline
526 unsigned int
528 {
529  return dofs_per_hex;
530 }
531 
532 
533 
534 template <int dim>
535 inline
536 unsigned int
538 {
539  return dofs_per_face;
540 }
541 
542 
543 
544 template <int dim>
545 inline
546 unsigned int
548 {
549  return dofs_per_cell;
550 }
551 
552 
553 
554 template <int dim>
555 template <int structdim>
556 inline
557 unsigned int
559 {
560  switch (structdim)
561  {
562  case 0:
563  return dofs_per_vertex;
564  case 1:
565  return dofs_per_line;
566  case 2:
567  return dofs_per_quad;
568  case 3:
569  return dofs_per_hex;
570  default:
571  Assert (false, ExcInternalError());
572  }
574 }
575 
576 
577 
578 template <int dim>
579 inline
580 unsigned int
582 {
583  return components;
584 }
585 
586 
587 
588 template <int dim>
589 inline
590 bool
592 {
593  return cached_primitivity;
594 }
595 
596 
597 template <int dim>
598 inline
599 void
601 {
602  cached_primitivity = value;
603 }
604 
605 
606 template <int dim>
607 inline
608 const BlockIndices &
610 {
611  return block_indices_data;
612 }
613 
614 
615 
616 template <int dim>
617 inline
618 unsigned int
620 {
621  return block_indices_data.size();
622 }
623 
624 
625 
626 template <int dim>
627 inline
628 unsigned int
630 {
631  return degree;
632 }
633 
634 
635 template <int dim>
636 inline
637 bool
639 {
640  return ((space & conforming_space) == space);
641 }
642 
643 
644 #endif // DOXYGEN
645 
646 
647 DEAL_II_NAMESPACE_CLOSE
648 
649 #endif
const unsigned int first_hex_index
Definition: fe_base.h:247
static const unsigned int invalid_unsigned_int
Definition: types.h:191
const unsigned int components
Definition: fe_base.h:281
Auxiliary class aiding in the handling of block structures like in BlockVector or FESystem...
Definition: block_indices.h:55
const unsigned int dofs_per_quad
Definition: fe_base.h:226
const unsigned int degree
Definition: fe_base.h:287
bool is_primitive() const
const unsigned int dofs_per_line
Definition: fe_base.h:220
const unsigned int first_face_line_index
Definition: fe_base.h:252
const BlockIndices & block_indices() const
unsigned int tensor_degree() const
unsigned int n_dofs_per_face() const
const unsigned int first_quad_index
Definition: fe_base.h:242
const unsigned int dofs_per_hex
Definition: fe_base.h:232
#define Assert(cond, exc)
Definition: exceptions.h:299
bool cached_primitivity
Definition: fe_base.h:432
unsigned int n_components() const
unsigned int n_dofs_per_vertex() const
const unsigned int dofs_per_cell
Definition: fe_base.h:271
Domination operator&(const Domination d1, const Domination d2)
bool conforms(const Conformity) const
unsigned int n_dofs_per_object() const
const unsigned int first_face_quad_index
Definition: fe_base.h:257
unsigned int n_dofs_per_line() const
const Conformity conforming_space
Definition: fe_base.h:292
unsigned int n_blocks() const
const unsigned int dofs_per_face
Definition: fe_base.h:264
const unsigned int first_line_index
Definition: fe_base.h:237
unsigned int n_dofs_per_cell() const
const unsigned int dofs_per_vertex
Definition: fe_base.h:214
void set_primitivity(const bool value)
::ExceptionBase & ExcInternalError()
unsigned int n_dofs_per_hex() const
unsigned int n_dofs_per_quad() const
BlockIndices block_indices_data
Definition: fe_base.h:299