Line data    Source code

       1              : 
       2              : namespace ippl {
       3              :     template <typename T>
       4              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::vertex_points_t
       5        13232 :     QuadrilateralElement<T>::getLocalVertices() const {
       6              :         // For the ordering of local vertices, see section 3.3.1:
       7              :         // https://amas.web.psi.ch/people/aadelmann/ETH-Accel-Lecture-1/projectscompleted/phys/bachelor_thesis_buehlluk.pdf
       8        66160 :         QuadrilateralElement::vertex_points_t vertices = {
       9              :             {0.0, 0.0}, {1.0, 0.0}, {1.0, 1.0}, {0.0, 1.0}};
      10              : 
      11        26432 :         return vertices;
      12        26464 :     }
      13              : 
      14              :     template <typename T>
      15              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::point_t
      16           56 :     QuadrilateralElement<T>::getTransformationJacobian(
      17              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices) const {
      18           56 :         QuadrilateralElement::point_t jacobian;
      19              : 
      20           56 :         jacobian[0] = (global_vertices[1][0] - global_vertices[0][0]);
      21           56 :         jacobian[1] = (global_vertices[2][1] - global_vertices[0][1]);
      22              : 
      23           56 :         return jacobian;
      24              :     }
      25              : 
      26              :     template <typename T>
      27              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::point_t
      28           52 :     QuadrilateralElement<T>::getInverseTransformationJacobian(
      29              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices) const {
      30           52 :         QuadrilateralElement::point_t inv_jacobian;
      31              : 
      32           52 :         inv_jacobian[0] = 1.0 / (global_vertices[1][0] - global_vertices[0][0]);
      33           52 :         inv_jacobian[1] = 1.0 / (global_vertices[2][1] - global_vertices[0][1]);
      34              : 
      35           52 :         return inv_jacobian;
      36              :     }
      37              : 
      38              :     template <typename T>
      39              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::point_t
      40           48 :     QuadrilateralElement<T>::globalToLocal(
      41              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices,
      42              :         const QuadrilateralElement<T>::point_t& global_point) const {
      43              :         // This is actually not a matrix, but an IPPL vector that represents a diagonal matrix
      44           48 :         const QuadrilateralElement<T>::point_t glob2loc_matrix =
      45              :             getInverseTransformationJacobian(global_vertices);
      46              : 
      47           48 :         QuadrilateralElement<T>::point_t local_point =
      48           96 :             glob2loc_matrix * (global_point - global_vertices[0]);
      49              : 
      50           48 :         return local_point;
      51           48 :     }
      52              : 
      53              :     template <typename T>
      54              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::point_t
      55           48 :     QuadrilateralElement<T>::localToGlobal(
      56              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices,
      57              :         const QuadrilateralElement<T>::point_t& local_point) const {
      58              :         // This is actually not a matrix but an IPPL vector that represents a diagonal matrix
      59           48 :         const QuadrilateralElement<T>::point_t loc2glob_matrix =
      60              :             getTransformationJacobian(global_vertices);
      61              : 
      62           48 :         QuadrilateralElement<T>::point_t global_point =
      63           96 :             (loc2glob_matrix * local_point) + global_vertices[0];
      64              : 
      65           48 :         return global_point;
      66           48 :     }
      67              : 
      68              :     template <typename T>
      69            8 :     KOKKOS_FUNCTION T QuadrilateralElement<T>::getDeterminantOfTransformationJacobian(
      70              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices) const {
      71            8 :         T determinant = 1.0;
      72              : 
      73              :         // Since the jacobian is a diagonal matrix in our case the determinant is the product of the
      74              :         // diagonal elements
      75           24 :         for (const T& jacobian_val : getTransformationJacobian(global_vertices)) {
      76           16 :             determinant *= jacobian_val;
      77              :         }
      78              : 
      79            8 :         return determinant;
      80              :     }
      81              : 
      82              :     template <typename T>
      83              :     KOKKOS_FUNCTION typename QuadrilateralElement<T>::point_t
      84            4 :     QuadrilateralElement<T>::getInverseTransposeTransformationJacobian(
      85              :         const QuadrilateralElement<T>::vertex_points_t& global_vertices) const {
      86              :         // Simply return the inverse transformation jacobian since it is a diagonal matrix
      87            4 :         return getInverseTransformationJacobian(global_vertices);
      88              :     }
      89              : 
      90              :     template <typename T>
      91        13200 :     KOKKOS_FUNCTION bool QuadrilateralElement<T>::isPointInRefElement(
      92              :         const Vector<T, 2>& point) const {
      93              :         // check if the local coordinates are inside the reference element
      94              : 
      95        39600 :         for (size_t d = 0; d < 2; d++) {
      96        26400 :             if (point[d] > 1.0 || point[d] < 0.0) {
      97              :                 // The global coordinates are outside of the support.
      98            0 :                 return false;
      99              :             }
     100              :         }
     101              : 
     102        13200 :         return true;
     103              :     }
     104              : 
     105              : }  // namespace ippl