waves
Basic FE playground
ElementComputation.hpp
Go to the documentation of this file.
1 
5 template <typename scalar, int element_type>
8  const katoptron::ElementsList &elementsList,
9  const katoptron::NodesList &nodesList)
10 {
11  tbox::Cache *cache;
12  size_t element_size;
13 
14  if (element_type == static_cast<int>(tbox::ElType::HEX8))
15  {
16  tbox::CacheHex8 &cachehex8 = trilinosHex8GetCache();
17  cache = &cachehex8;
18  element_size = 8;
19  }
20  else if (element_type == static_cast<int>(tbox::ElType::TETRA4))
21  {
22  tbox::CacheTetra4 &cachetetra4 = trilinosTetra4GetCache();
23  cache = &cachetetra4;
24  element_size = 4;
25  }
26  else if (element_type == static_cast<int>(tbox::ElType::TRI3))
27  {
28  tbox::CacheTri3 &cachetri3 = trilinosTri3GetCache();
29  cache = &cachetri3;
30  element_size = 3;
31 
32  tVector<double, 3> tg1, tg2, n;
33 
34  int i = 0;
35  for (auto n1 = 0; n1 < element_size; ++n1, ++i)
36  {
37  for (auto j = 0; j < 3; ++j)
38  {
39  tg1(j) += (static_cast<tbox::CacheTri3 *>(cache)->getDsf(k))(0, i) *
40  nodesList.nodes(elementsList.getElementNode(e, n1), j);
41  tg2(j) += (static_cast<tbox::CacheTri3 *>(cache)->getDsf(k))(1, i) *
42  nodesList.nodes(elementsList.getElementNode(e, n1), j);
43  }
44  }
45 
46  n(0) = tg1(1) * tg2(2) - tg1(2) * tg2(1);
47  n(1) = tg1(2) * tg2(0) - tg1(0) * tg2(2);
48  n(2) = tg1(0) * tg2(1) - tg1(1) * tg2(0);
49  double detJ = sqrt(pow(n(0), 2) + pow(n(1), 2) + pow(n(2), 2));
50  return detJ;
51  }
52  else if (element_type == static_cast<int>(tbox::ElType::QUAD4))
53  {
54  tbox::CacheQuad4 &cachequad4 = trilinosQuad4GetCache();
55  cache = &cachequad4;
56  element_size = 4;
57 
58  tVector<double, 3> tg1, tg2, n;
59 
60  int i = 0;
61  for (auto n1 = 0; n1 < element_size; ++n1, ++i)
62  {
63  for (auto j = 0; j < 3; ++j)
64  {
65  tg1(j) += (static_cast<tbox::CacheQuad4 *>(cache)->getDsf(k))(0, i) *
66  nodesList.nodes(elementsList.getElementNode(e, n1), j);
67  tg2(j) += (static_cast<tbox::CacheQuad4 *>(cache)->getDsf(k))(1, i) *
68  nodesList.nodes(elementsList.getElementNode(e, n1), j);
69  }
70  }
71 
72  n(0) = tg1(1) * tg2(2) - tg1(2) * tg2(1);
73  n(1) = tg1(2) * tg2(0) - tg1(0) * tg2(2);
74  n(2) = tg1(0) * tg2(1) - tg1(1) * tg2(0);
75  double detJ = sqrt(pow(n(0), 2) + pow(n(1), 2) + pow(n(2), 2));
76  return detJ;
77  }
78 
80  J.clean();
81 
82  tVector<double, 3> tmp, tmp2;
83 
84  int i = 0;
85  for (auto n1 = 0; n1 < elementsList.getElementSize(e); ++n1, ++i)
86  {
87  for (auto j = 0; j < 3; ++j)
88  {
89  tmp(j) = nodesList.nodes(elementsList.getElementNode(e, n1), j);
90  if (element_type == static_cast<int>(tbox::ElType::HEX8))
91  tmp2(j) = (static_cast<tbox::CacheHex8 *>(cache)->getDsf(k))(j, i);
92  else if (element_type == static_cast<int>(tbox::ElType::TETRA4))
93  tmp2(j) = (static_cast<tbox::CacheTetra4 *>(cache)->getDsf(k))(j, i);
94  }
95 
96  Jn = tmp2.ddotproduct(tmp);
97  J += Jn;
98  }
99 
100  double detJ = J.det();
101  return detJ;
102 }
103 
107 template <typename scalar, int element_type>
109 {
110  static tbox::Lazy<tbox::CacheHex8> cache(2);
111  return cache.get();
112 }
113 
117 template <typename scalar, int element_type>
119 {
120  static tbox::Lazy<tbox::CacheTetra4> cache(2);
121  return cache.get();
122 }
123 
127 template <typename scalar, int element_type>
129 {
130  static tbox::Lazy<tbox::CacheTri3> cache(2);
131  return cache.get();
132 }
133 
137 template <typename scalar, int element_type>
139 {
140  static tbox::Lazy<tbox::CacheQuad4> cache(2);
141  return cache.get();
142 }
tVector::ddotproduct
tMatrix< T, rows, rows > ddotproduct(tVector rhs)
Definition: tMatrix.h:207
tMatrix::det
T det()
Definition: tMatrix.h:91
ElementComputation::local_ordinal_type
katoptron::Map::local_ordinal_type local_ordinal_type
Definition: ElementComputation.h:32
ElementComputation::buildJ
double buildJ(int k, tMatrix< double, 3, 3 > &J, local_ordinal_type e, const katoptron::ElementsList &elementsList, const katoptron::NodesList &nodesList)
Compute the Jacobian matrix of an element.
Definition: ElementComputation.hpp:6
tMatrix::clean
void clean()
Definition: tMatrix.h:29
katoptron::NodesList::nodes
Kokkos::View< double **, Kokkos::LayoutRight > nodes
Definition: NodesList.h:22
katoptron::NodesList
Class used to store the node information.
Definition: NodesList.h:19
ElementComputation::trilinosHex8GetCache
tbox::CacheHex8 & trilinosHex8GetCache()
Return the tbox cache of the hexahedron elements.
Definition: ElementComputation.hpp:108
katoptron::ElementsList
Class used to store the element information including:
Definition: ElementsList.h:25
tVector
Definition: tMatrix.h:7
tMatrix
Definition: tMatrix.h:10
ElementComputation::trilinosQuad4GetCache
tbox::CacheQuad4 & trilinosQuad4GetCache()
Return the tbox cache of the quadrangular elements.
Definition: ElementComputation.hpp:138
ElementComputation::trilinosTri3GetCache
tbox::CacheTri3 & trilinosTri3GetCache()
Return the tbox cache of the triangular elements.
Definition: ElementComputation.hpp:128
katoptron::ElementsList::getElementSize
local_ordinal_type getElementSize(const local_ordinal_type e) const
Return the size of the local element of local ID e. Argument:
Definition: ElementsList.h:67
katoptron::ElementsList::getElementNode
local_ordinal_type getElementNode(const local_ordinal_type e, local_ordinal_type i) const
Return the local ID of the ith node of the local element of ID e. Argument:
Definition: ElementsList.h:92
ElementComputation::trilinosTetra4GetCache
tbox::CacheTetra4 & trilinosTetra4GetCache()
Return the tbox cache of the tetrahedron elements.
Definition: ElementComputation.hpp:118