MFEM  v4.6.0
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 
43  constexpr QVectorLayout L = QVectorLayout::byNODES;
44  constexpr bool P = false; // GRAD_PHYS
45 
46  const int id = (vdim<<8) | (D1D<<4) | Q1D;
47 
48  if (dim == 1)
49  {
50  return Derivatives1D<L,P>(NE,G,J,X,Y,dim,vdim,D1D,Q1D);
51  }
52  if (dim == 2)
53  {
54  switch (id)
55  {
56  case 0x133: return Derivatives2D<L,P,1,3,3,16>(NE,B,G,J,X,Y);
57  case 0x134: return Derivatives2D<L,P,1,3,4,16>(NE,B,G,J,X,Y);
58  case 0x143: return Derivatives2D<L,P,1,4,3,16>(NE,B,G,J,X,Y);
59  case 0x144: return Derivatives2D<L,P,1,4,4,16>(NE,B,G,J,X,Y);
60 
61  case 0x222: return Derivatives2D<L,P,2,2,2,16>(NE,B,G,J,X,Y);
62  case 0x223: return Derivatives2D<L,P,2,2,3,8>(NE,B,G,J,X,Y);
63  case 0x224: return Derivatives2D<L,P,2,2,4,4>(NE,B,G,J,X,Y);
64  case 0x225: return Derivatives2D<L,P,2,2,5,4>(NE,B,G,J,X,Y);
65  case 0x226: return Derivatives2D<L,P,2,2,6,2>(NE,B,G,J,X,Y);
66 
67  case 0x233: return Derivatives2D<L,P,2,3,3,2>(NE,B,G,J,X,Y);
68  case 0x234: return Derivatives2D<L,P,2,3,4,4>(NE,B,G,J,X,Y);
69  case 0x243: return Derivatives2D<L,P,2,4,3,4>(NE,B,G,J,X,Y);
70  case 0x236: return Derivatives2D<L,P,2,3,6,2>(NE,B,G,J,X,Y);
71 
72  case 0x244: return Derivatives2D<L,P,2,4,4,2>(NE,B,G,J,X,Y);
73  case 0x245: return Derivatives2D<L,P,2,4,5,2>(NE,B,G,J,X,Y);
74  case 0x246: return Derivatives2D<L,P,2,4,6,2>(NE,B,G,J,X,Y);
75  case 0x247: return Derivatives2D<L,P,2,4,7,2>(NE,B,G,J,X,Y);
76 
77  case 0x256: return Derivatives2D<L,P,2,5,6,2>(NE,B,G,J,X,Y);
78  default:
79  {
80  const int MD = DeviceDofQuadLimits::Get().MAX_D1D;
81  const int MQ = DeviceDofQuadLimits::Get().MAX_Q1D;
82  if (D1D > MD || Q1D > MQ)
83  {
84  MFEM_ABORT("");
85  }
86  Derivatives2D<L,P>(NE,B,G,J,X,Y,dim,vdim,D1D,Q1D);
87  return;
88  }
89  }
90  }
91  if (dim == 3)
92  {
93  switch (id)
94  {
95  case 0x124: return Derivatives3D<L,P,1,2,4>(NE,B,G,J,X,Y);
96  case 0x133: return Derivatives3D<L,P,1,3,3>(NE,B,G,J,X,Y);
97  case 0x134: return Derivatives3D<L,P,1,3,4>(NE,B,G,J,X,Y);
98  case 0x136: return Derivatives3D<L,P,1,3,6>(NE,B,G,J,X,Y);
99  case 0x144: return Derivatives3D<L,P,1,4,4>(NE,B,G,J,X,Y);
100  case 0x148: return Derivatives3D<L,P,1,4,8>(NE,B,G,J,X,Y);
101 
102  case 0x323: return Derivatives3D<L,P,3,2,3>(NE,B,G,J,X,Y);
103  case 0x324: return Derivatives3D<L,P,3,2,4>(NE,B,G,J,X,Y);
104  case 0x325: return Derivatives3D<L,P,3,2,5>(NE,B,G,J,X,Y);
105  case 0x326: return Derivatives3D<L,P,3,2,6>(NE,B,G,J,X,Y);
106 
107  case 0x333: return Derivatives3D<L,P,3,3,3>(NE,B,G,J,X,Y);
108  case 0x334: return Derivatives3D<L,P,3,3,4>(NE,B,G,J,X,Y);
109  case 0x335: return Derivatives3D<L,P,3,3,5>(NE,B,G,J,X,Y);
110  case 0x336: return Derivatives3D<L,P,3,3,6>(NE,B,G,J,X,Y);
111  case 0x344: return Derivatives3D<L,P,3,4,4>(NE,B,G,J,X,Y);
112  case 0x346: return Derivatives3D<L,P,3,4,6>(NE,B,G,J,X,Y);
113  case 0x347: return Derivatives3D<L,P,3,4,7>(NE,B,G,J,X,Y);
114  case 0x348: return Derivatives3D<L,P,3,4,8>(NE,B,G,J,X,Y);
115  default:
116  {
117  const int MD = DeviceDofQuadLimits::Get().MAX_INTERP_1D;
118  const int MQ = DeviceDofQuadLimits::Get().MAX_INTERP_1D;
119  MFEM_VERIFY(D1D <= MD, "Orders higher than " << MD-1
120  << " are not supported!");
121  MFEM_VERIFY(Q1D <= MQ, "Quadrature rules with more than "
122  << MQ << " 1D points are not supported!");
123  Derivatives3D<L,P>(NE,B,G,J,X,Y,vdim,D1D,Q1D);
124  return;
125  }
126  }
127  }
128  mfem::out << "Unknown kernel 0x" << std::hex << id << std::endl;
129  MFEM_ABORT("Kernel not supported yet");
130 }
131 
132 } // namespace quadrature_interpolator
133 
134 } // namespace internal
135 
136 } // namespace mfem
mfem::DeviceDofQuadLimits::MAX_Q1D
int MAX_Q1D
Maximum number of 1D quadrature points.
Definition: forall.hpp:113
mfem::DeviceDofQuadLimits::MAX_D1D
int MAX_D1D
Maximum number of 1D nodal points.
Definition: forall.hpp:112
mfem::out
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
mfem::DeviceDofQuadLimits::Get
static const DeviceDofQuadLimits & Get()
Return a const reference to the DeviceDofQuadLimits singleton.
Definition: forall.hpp:122
mfem::QVectorLayout
QVectorLayout
Type describing possible layouts for Q-vectors.
Definition: fespace.hpp:52
dim
int dim
Definition: ex24.cpp:53
mfem::DeviceDofQuadLimits::MAX_INTERP_1D
int MAX_INTERP_1D
Maximum number of points for use in QuadratureInterpolator.
Definition: forall.hpp:118
mfem::QVectorLayout::byNODES
NQPT x VDIM x NE (values) / NQPT x VDIM x DIM x NE (grads)
Generated on Wed Sep 27 2023 14:52:30 for MFEM by   doxygen 1.8.14