MFEM  v4.5.2
Finite element discretization library
grad_by_nodes.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #include "dispatch.hpp"
13 #include "grad.hpp"
14 
15 namespace mfem
16 {
17 
18 namespace internal
19 {
20 
21 namespace quadrature_interpolator
22 {
23 
24 // Tensor-product evaluation of quadrature point derivatives: dispatch function.
25 // Instantiation for the case QVectorLayout::byNODES.
26 template<>
27 void TensorDerivatives<QVectorLayout::byNODES>(const int NE,
28  const int vdim,
29  const DofToQuad &maps,
30  const Vector &e_vec,
31  Vector &q_der)
32 {
33  if (NE == 0) { return; }
34  const int dim = maps.FE->GetDim();
35  const int D1D = maps.ndof;
36  const int Q1D = maps.nqpt;
37  const double *B = maps.B.Read();
38  const double *G = maps.G.Read();
39  const double *J = nullptr; // not used in DERIVATIVES (non-GRAD_PHYS) mode
40  const double *X = e_vec.Read();
41  double *Y = q_der.Write();
42 
44  constexpr bool P = false; // GRAD_PHYS
45 
46  const int id = (vdim<<8) | (D1D<<4) | Q1D;
47 
48  if (dim == 2)
49  {
50  switch (id)
51  {
52  case 0x133: return Derivatives2D<L,P,1,3,3,16>(NE,B,G,J,X,Y);
53  case 0x134: return Derivatives2D<L,P,1,3,4,16>(NE,B,G,J,X,Y);
54  case 0x143: return Derivatives2D<L,P,1,4,3,16>(NE,B,G,J,X,Y);
55  case 0x144: return Derivatives2D<L,P,1,4,4,16>(NE,B,G,J,X,Y);
56 
57  case 0x222: return Derivatives2D<L,P,2,2,2,16>(NE,B,G,J,X,Y);
58  case 0x223: return Derivatives2D<L,P,2,2,3,8>(NE,B,G,J,X,Y);
59  case 0x224: return Derivatives2D<L,P,2,2,4,4>(NE,B,G,J,X,Y);
60  case 0x225: return Derivatives2D<L,P,2,2,5,4>(NE,B,G,J,X,Y);
61  case 0x226: return Derivatives2D<L,P,2,2,6,2>(NE,B,G,J,X,Y);
62 
63  case 0x233: return Derivatives2D<L,P,2,3,3,2>(NE,B,G,J,X,Y);
64  case 0x234: return Derivatives2D<L,P,2,3,4,4>(NE,B,G,J,X,Y);
65  case 0x243: return Derivatives2D<L,P,2,4,3,4>(NE,B,G,J,X,Y);
66  case 0x236: return Derivatives2D<L,P,2,3,6,2>(NE,B,G,J,X,Y);
67 
68  case 0x244: return Derivatives2D<L,P,2,4,4,2>(NE,B,G,J,X,Y);
69  case 0x245: return Derivatives2D<L,P,2,4,5,2>(NE,B,G,J,X,Y);
70  case 0x246: return Derivatives2D<L,P,2,4,6,2>(NE,B,G,J,X,Y);
71  case 0x247: return Derivatives2D<L,P,2,4,7,2>(NE,B,G,J,X,Y);
72 
73  case 0x256: return Derivatives2D<L,P,2,5,6,2>(NE,B,G,J,X,Y);
74  default:
75  {
76  constexpr int MD = MAX_D1D;
77  constexpr int MQ = MAX_Q1D;
78  if (D1D > MD || Q1D > MQ)
79  {
80  MFEM_ABORT("");
81  }
82  Derivatives2D<L,P,0,0,0,0,MD,MQ>(NE,B,G,J,X,Y,vdim,D1D,Q1D);
83  return;
84  }
85  }
86  }
87  if (dim == 3)
88  {
89  switch (id)
90  {
91  case 0x124: return Derivatives3D<L,P,1,2,4>(NE,B,G,J,X,Y);
92  case 0x133: return Derivatives3D<L,P,1,3,3>(NE,B,G,J,X,Y);
93  case 0x134: return Derivatives3D<L,P,1,3,4>(NE,B,G,J,X,Y);
94  case 0x136: return Derivatives3D<L,P,1,3,6>(NE,B,G,J,X,Y);
95  case 0x144: return Derivatives3D<L,P,1,4,4>(NE,B,G,J,X,Y);
96  case 0x148: return Derivatives3D<L,P,1,4,8>(NE,B,G,J,X,Y);
97 
98  case 0x323: return Derivatives3D<L,P,3,2,3>(NE,B,G,J,X,Y);
99  case 0x324: return Derivatives3D<L,P,3,2,4>(NE,B,G,J,X,Y);
100  case 0x325: return Derivatives3D<L,P,3,2,5>(NE,B,G,J,X,Y);
101  case 0x326: return Derivatives3D<L,P,3,2,6>(NE,B,G,J,X,Y);
102 
103  case 0x333: return Derivatives3D<L,P,3,3,3>(NE,B,G,J,X,Y);
104  case 0x334: return Derivatives3D<L,P,3,3,4>(NE,B,G,J,X,Y);
105  case 0x335: return Derivatives3D<L,P,3,3,5>(NE,B,G,J,X,Y);
106  case 0x336: return Derivatives3D<L,P,3,3,6>(NE,B,G,J,X,Y);
107  case 0x344: return Derivatives3D<L,P,3,4,4>(NE,B,G,J,X,Y);
108  case 0x346: return Derivatives3D<L,P,3,4,6>(NE,B,G,J,X,Y);
109  case 0x347: return Derivatives3D<L,P,3,4,7>(NE,B,G,J,X,Y);
110  case 0x348: return Derivatives3D<L,P,3,4,8>(NE,B,G,J,X,Y);
111  default:
112  {
113  constexpr int MD = 8;
114  constexpr int MQ = 8;
115  MFEM_VERIFY(D1D <= MD, "Orders higher than " << MD-1
116  << " are not supported!");
117  MFEM_VERIFY(Q1D <= MQ, "Quadrature rules with more than "
118  << MQ << " 1D points are not supported!");
119  Derivatives3D<L,P,0,0,0,MD,MQ>(NE,B,G,J,X,Y,vdim,D1D,Q1D);
120  return;
121  }
122  }
123  }
124  mfem::out << "Unknown kernel 0x" << std::hex << id << std::endl;
125  MFEM_ABORT("Kernel not supported yet");
126 }
127 
128 } // namespace quadrature_interpolator
129 
130 } // namespace internal
131 
132 } // namespace mfem
const int MAX_Q1D
Definition: forall.hpp:29
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
QVectorLayout
Type describing possible layouts for Q-vectors.
Definition: fespace.hpp:52
int dim
Definition: ex24.cpp:53
const int MAX_D1D
Definition: forall.hpp:28
NQPT x VDIM x NE (values) / NQPT x VDIM x DIM x NE (grads)