xref: /libCEED/include/ceed.h (revision 6f5dc8ba493eb609f13fe336f2a5e06988b93c70)
1 /// Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at
2 /// the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights
3 /// reserved. See files LICENSE and NOTICE for details.
4 ///
5 /// This file is part of CEED, a collection of benchmarks, miniapps, software
6 /// libraries and APIs for efficient high-order finite element and spectral
7 /// element discretizations for exascale applications. For more information and
8 /// source code availability see http://github.com/ceed.
9 ///
10 /// The CEED research is supported by the Exascale Computing Project 17-SC-20-SC,
11 /// a collaborative effort of two U.S. Department of Energy organizations (Office
12 /// of Science and the National Nuclear Security Administration) responsible for
13 /// the planning and preparation of a capable exascale ecosystem, including
14 /// software, applications, hardware, advanced system engineering and early
15 /// testbed platforms, in support of the nation's exascale computing imperative.
16 
17 /// @file
18 /// Public header for user and utility components of libCEED
19 #ifndef _ceed_h
20 #define _ceed_h
21 
22 /// @defgroup Ceed Ceed: core components
23 /// @defgroup CeedVector CeedVector: storing and manipulating vectors
24 /// @defgroup CeedElemRestriction CeedElemRestriction: restriction from local vectors to elements
25 /// @defgroup CeedBasis CeedBasis: fully discrete finite element-like objects
26 /// @defgroup CeedQFunction CeedQFunction: independent operations at quadrature points
27 /// @defgroup CeedOperator CeedOperator: composed FE-type operations on vectors
28 ///
29 /// @page FunctionCategories libCEED: Types of Functions
30 ///    libCEED provides three different header files depending upon the type of
31 ///    functions a user requires.
32 /// @section Utility Utility Functions
33 ///    These functions are intended general utilities that may be useful to
34 ///    libCEED developers and users. These functions can generally be found in
35 ///    "ceed.h".
36 /// @section User User Functions
37 ///    These functions are intended to be used by general users of libCEED
38 ///    and can generally be found in "ceed.h".
39 /// @section Backend Backend Developer Functions
40 ///    These functions are intended to be used by backend developers of
41 ///    libCEED and can generally be found in "ceed-backend.h".
42 /// @section Developer Library Developer Functions
43 ///    These functions are intended to be used by library developers of
44 ///    libCEED and can generally be found in "ceed-impl.h".
45 
46 /**
47   CEED_EXTERN is used in this header to denote all publicly visible symbols.
48 
49   No other file should declare publicly visible symbols, thus it should never be
50   used outside ceed.h.
51  */
52 #ifdef __cplusplus
53 #  define CEED_EXTERN extern "C"
54 #else
55 #  define CEED_EXTERN extern
56 #endif
57 
58 /**
59   @ingroup CeedQFunction
60   This macro populates the correct function annotations for User QFunction
61     source for code generation backends or populates default values for CPU
62     backends.
63 **/
64 #ifndef CEED_QFUNCTION
65 #define CEED_QFUNCTION(name) \
66   static const char name ## _loc[] = __FILE__ ":" #name;        \
67   static int name
68 #endif
69 
70 /**
71   @ingroup CeedQFunction
72   Using VLA syntax to reshape User QFunction inputs and outputs can make
73     user code more readable. VLA is a C99 feature that is not supported by
74     the C++ dialect used by CUDA. This macro allows users to use the VLA
75     syntax with the CUDA backends.
76 **/
77 #ifndef CEED_Q_VLA
78 #  define CEED_Q_VLA Q
79 #endif
80 
81 /**
82   @ingroup Ceed
83   This macro provides the appropriate SIMD Pragma for the compilation
84     environment. Code generation backends may redefine this macro, as needed.
85 **/
86 #ifndef CeedPragmaSIMD
87 #  if defined(__INTEL_COMPILER)
88 #    define CeedPragmaSIMD _Pragma("vector")
89 // Cannot use Intel pragma ivdep because it miscompiles unpacking symmetric tensors, as in
90 // Poisson2DApply, where the SIMD loop body contains temporaries such as the following.
91 //
92 //     const CeedScalar dXdxdXdxT[2][2] = {{qd[i+0*Q], qd[i+2*Q]},
93 //                                         {qd[i+2*Q], qd[i+1*Q]}};
94 //     for (int j=0; j<2; j++)
95 //        vg[i+j*Q] = (du[0] * dXdxdXdxT[0][j] + du[1] * dXdxdXdxT[1][j]);
96 //
97 // Miscompilation with pragma ivdep observed with icc (ICC) 19.0.5.281 20190815
98 // at -O2 and above.
99 #  elif defined(__GNUC__) && __GNUC__ >= 5
100 #    define CeedPragmaSIMD _Pragma("GCC ivdep")
101 #  elif defined(_OPENMP) && _OPENMP >= 201307 // OpenMP-4.0 (July, 2013)
102 #    define CeedPragmaSIMD _Pragma("omp simd")
103 #  else
104 #    define CeedPragmaSIMD
105 #  endif
106 #endif
107 
108 #include <assert.h>
109 #include <stdint.h>
110 #include <stddef.h>
111 #include <stdarg.h>
112 #include <stdio.h>
113 #include <stdbool.h>
114 
115 /// Integer type, used for indexing
116 /// @ingroup Ceed
117 typedef int32_t CeedInt;
118 /// Scalar (floating point) type
119 /// @ingroup Ceed
120 typedef double CeedScalar;
121 
122 /// Library context created by CeedInit()
123 /// @ingroup CeedUser
124 typedef struct Ceed_private *Ceed;
125 /// Non-blocking Ceed interfaces return a CeedRequest.
126 /// To perform an operation immediately, pass \ref CEED_REQUEST_IMMEDIATE instead.
127 /// @ingroup CeedUser
128 typedef struct CeedRequest_private *CeedRequest;
129 /// Handle for vectors over the field \ref CeedScalar
130 /// @ingroup CeedVectorUser
131 typedef struct CeedVector_private *CeedVector;
132 /// Handle for object describing restriction to elements
133 /// @ingroup CeedElemRestrictionUser
134 typedef struct CeedElemRestriction_private *CeedElemRestriction;
135 /// Handle for object describing discrete finite element evaluations
136 /// @ingroup CeedBasisUser
137 typedef struct CeedBasis_private *CeedBasis;
138 /// Handle for object describing functions evaluated independently at quadrature points
139 /// @ingroup CeedQFunctionUser
140 typedef struct CeedQFunction_private *CeedQFunction;
141 /// Handle for object describing context data for CeedQFunctions
142 /// @ingroup CeedQFunctionUser
143 typedef struct CeedQFunctionContext_private *CeedQFunctionContext;
144 /// Handle for object describing FE-type operators acting on vectors
145 ///
146 /// Given an element restriction \f$E\f$, basis evaluator \f$B\f$, and
147 ///   quadrature function\f$f\f$, a CeedOperator expresses operations of the form
148 ///   $$ E^T B^T f(B E u) $$
149 ///   acting on the vector \f$u\f$.
150 /// @ingroup CeedOperatorUser
151 typedef struct CeedOperator_private *CeedOperator;
152 
153 CEED_EXTERN int CeedInit(const char *resource, Ceed *ceed);
154 CEED_EXTERN int CeedGetResource(Ceed ceed, const char **resource);
155 CEED_EXTERN int CeedIsDeterministic(Ceed ceed, bool *isDeterministic);
156 CEED_EXTERN int CeedView(Ceed ceed, FILE *stream);
157 CEED_EXTERN int CeedDestroy(Ceed *ceed);
158 
159 CEED_EXTERN int CeedErrorImpl(Ceed, const char *, int, const char *, int,
160                               const char *, ...);
161 /// Raise an error on ceed object
162 ///
163 /// @param ceed Ceed library context or NULL
164 /// @param ecode Error code (int)
165 /// @param ... printf-style format string followed by arguments as needed
166 ///
167 /// @ingroup Ceed
168 /// @sa CeedSetErrorHandler()
169 #if defined(__clang__)
170 /// Use nonstandard ternary to convince the compiler/clang-tidy that this
171 /// function never returns zero.
172 #  define CeedError(ceed, ecode, ...)                                     \
173   (CeedErrorImpl((ceed), __FILE__, __LINE__, __func__, (ecode), __VA_ARGS__) ?: (ecode))
174 #else
175 #  define CeedError(ceed, ecode, ...)                                     \
176   CeedErrorImpl((ceed), __FILE__, __LINE__, __func__, (ecode), __VA_ARGS__) ?: (ecode)
177 #endif
178 
179 /// Ceed error handlers
180 CEED_EXTERN int CeedErrorReturn(Ceed, const char *, int, const char *, int,
181                                 const char *, va_list *);
182 CEED_EXTERN int CeedErrorStore(Ceed, const char *, int, const char *, int,
183                                const char *, va_list *);
184 CEED_EXTERN int CeedErrorAbort(Ceed, const char *, int, const char *, int,
185                                const char *, va_list *);
186 CEED_EXTERN int CeedErrorExit(Ceed, const char *, int, const char *, int,
187                               const char *, va_list *);
188 CEED_EXTERN int CeedSetErrorHandler(Ceed ceed,
189                                     int (*eh)(Ceed, const char *, int,
190                                         const char *, int, const char *,
191                                         va_list *));
192 CEED_EXTERN int CeedGetErrorMessage(Ceed, const char **errmsg);
193 CEED_EXTERN int CeedResetErrorMessage(Ceed, const char **errmsg);
194 
195 /// Specify memory type
196 ///
197 /// Many Ceed interfaces take or return pointers to memory.  This enum is used to
198 /// specify where the memory being provided or requested must reside.
199 /// @ingroup Ceed
200 typedef enum {
201   /// Memory resides on the host
202   CEED_MEM_HOST,
203   /// Memory resides on a device (corresponding to \ref Ceed resource)
204   CEED_MEM_DEVICE,
205 } CeedMemType;
206 
207 CEED_EXTERN const char *const CeedMemTypes[];
208 
209 CEED_EXTERN int CeedGetPreferredMemType(Ceed ceed, CeedMemType *type);
210 
211 /// Conveys ownership status of arrays passed to Ceed interfaces.
212 /// @ingroup Ceed
213 typedef enum {
214   /// Implementation will copy the values and not store the passed pointer.
215   CEED_COPY_VALUES,
216   /// Implementation can use and modify the data provided by the user, but does
217   /// not take ownership.
218   CEED_USE_POINTER,
219   /// Implementation takes ownership of the pointer and will free using
220   /// CeedFree() when done using it.  The user should not assume that the
221   /// pointer remains valid after ownership has been transferred.  Note that
222   /// arrays allocated using C++ operator new or other allocators cannot
223   /// generally be freed using CeedFree().  CeedFree() is capable of freeing any
224   /// memory that can be freed using free(3).
225   CEED_OWN_POINTER,
226 } CeedCopyMode;
227 
228 /// Denotes type of vector norm to be computed
229 /// @ingroup CeedVector
230 typedef enum {
231   /// L_1 norm: sum_i |x_i|
232   CEED_NORM_1,
233   /// L_2 norm: sqrt(sum_i |x_i|^2)
234   CEED_NORM_2,
235   /// L_Infinity norm: max_i |x_i|
236   CEED_NORM_MAX,
237 } CeedNormType;
238 
239 CEED_EXTERN const char *const CeedCopyModes[];
240 
241 CEED_EXTERN int CeedVectorCreate(Ceed ceed, CeedInt len, CeedVector *vec);
242 CEED_EXTERN int CeedVectorSetArray(CeedVector vec, CeedMemType mtype,
243                                    CeedCopyMode cmode, CeedScalar *array);
244 CEED_EXTERN int CeedVectorSetValue(CeedVector vec, CeedScalar value);
245 CEED_EXTERN int CeedVectorSyncArray(CeedVector vec, CeedMemType mtype);
246 CEED_EXTERN int CeedVectorTakeArray(CeedVector vec, CeedMemType mtype,
247                                     CeedScalar **array);
248 CEED_EXTERN int CeedVectorGetArray(CeedVector vec, CeedMemType mtype,
249                                    CeedScalar **array);
250 CEED_EXTERN int CeedVectorGetArrayRead(CeedVector vec, CeedMemType mtype,
251                                        const CeedScalar **array);
252 CEED_EXTERN int CeedVectorRestoreArray(CeedVector vec, CeedScalar **array);
253 CEED_EXTERN int CeedVectorRestoreArrayRead(CeedVector vec,
254     const CeedScalar **array);
255 CEED_EXTERN int CeedVectorNorm(CeedVector vec, CeedNormType type,
256                                CeedScalar *norm);
257 CEED_EXTERN int CeedVectorReciprocal(CeedVector vec);
258 CEED_EXTERN int CeedVectorView(CeedVector vec, const char *fpfmt, FILE *stream);
259 CEED_EXTERN int CeedVectorGetLength(CeedVector vec, CeedInt *length);
260 CEED_EXTERN int CeedVectorDestroy(CeedVector *vec);
261 
262 CEED_EXTERN CeedRequest *const CEED_REQUEST_IMMEDIATE;
263 CEED_EXTERN CeedRequest *const CEED_REQUEST_ORDERED;
264 CEED_EXTERN int CeedRequestWait(CeedRequest *req);
265 
266 /// Argument for CeedOperatorSetField that vector is collocated with
267 /// quadrature points, used with QFunction eval mode CEED_EVAL_NONE
268 /// or CEED_EVAL_INTERP only, not with CEED_EVAL_GRAD, CEED_EVAL_DIV,
269 /// or CEED_EVAL_CURL
270 /// @ingroup CeedBasis
271 CEED_EXTERN const CeedBasis CEED_BASIS_COLLOCATED;
272 
273 /// Argument for CeedOperatorSetField to use active input or output
274 /// @ingroup CeedVector
275 CEED_EXTERN const CeedVector CEED_VECTOR_ACTIVE;
276 
277 /// Argument for CeedOperatorSetField to use no vector, used with
278 /// qfunction input with eval mode CEED_EVAL_WEIGHT
279 /// @ingroup CeedVector
280 CEED_EXTERN const CeedVector CEED_VECTOR_NONE;
281 
282 /// Argument for CeedOperatorSetField to use no ElemRestriction, only used with
283 /// eval mode CEED_EVAL_WEIGHT.
284 /// @ingroup CeedElemRestriction
285 CEED_EXTERN const CeedElemRestriction CEED_ELEMRESTRICTION_NONE;
286 
287 /// Argument for CeedOperatorCreate that QFunction is not created by user.
288 /// Only used for QFunctions dqf and dqfT. If implemented, a backend may
289 /// attempt to provide the action of these QFunctions.
290 /// @ingroup CeedQFunction
291 CEED_EXTERN const CeedQFunction CEED_QFUNCTION_NONE;
292 
293 /// Denotes whether a linear transformation or its transpose should be applied
294 /// @ingroup CeedBasis
295 typedef enum {
296   /// Apply the linear transformation
297   CEED_NOTRANSPOSE,
298   /// Apply the transpose
299   CEED_TRANSPOSE
300 } CeedTransposeMode;
301 
302 CEED_EXTERN const char *const CeedTransposeModes[];
303 
304 /// Argument for CeedElemRestrictionCreateStrided that L-vector is in
305 /// the Ceed backend's preferred layout. This argument should only be used
306 /// with vectors created by a Ceed backend.
307 /// @ingroup CeedElemRestriction
308 CEED_EXTERN const CeedInt CEED_STRIDES_BACKEND[3];
309 
310 CEED_EXTERN int CeedElemRestrictionCreate(Ceed ceed, CeedInt nelem,
311     CeedInt elemsize, CeedInt ncomp, CeedInt compstride, CeedInt lsize,
312     CeedMemType mtype, CeedCopyMode cmode, const CeedInt *offsets,
313     CeedElemRestriction *rstr);
314 CEED_EXTERN int CeedElemRestrictionCreateStrided(Ceed ceed,
315     CeedInt nelem, CeedInt elemsize, CeedInt ncomp, CeedInt lsize,
316     const CeedInt strides[3], CeedElemRestriction *rstr);
317 CEED_EXTERN int CeedElemRestrictionCreateBlocked(Ceed ceed, CeedInt nelem,
318     CeedInt elemsize, CeedInt blksize, CeedInt ncomp, CeedInt compstride,
319     CeedInt lsize, CeedMemType mtype, CeedCopyMode cmode,
320     const CeedInt *offsets, CeedElemRestriction *rstr);
321 CEED_EXTERN int CeedElemRestrictionCreateBlockedStrided(Ceed ceed,
322     CeedInt nelem, CeedInt elemsize, CeedInt blksize, CeedInt ncomp,
323     CeedInt lsize, const CeedInt strides[3], CeedElemRestriction *rstr);
324 CEED_EXTERN int CeedElemRestrictionCreateVector(CeedElemRestriction rstr,
325     CeedVector *lvec, CeedVector *evec);
326 CEED_EXTERN int CeedElemRestrictionApply(CeedElemRestriction rstr,
327     CeedTransposeMode tmode, CeedVector u, CeedVector ru, CeedRequest *request);
328 CEED_EXTERN int CeedElemRestrictionApplyBlock(CeedElemRestriction rstr,
329     CeedInt block, CeedTransposeMode tmode, CeedVector u, CeedVector ru,
330     CeedRequest *request);
331 CEED_EXTERN int CeedElemRestrictionGetCompStride(CeedElemRestriction rstr,
332     CeedInt *compstride);
333 CEED_EXTERN int CeedElemRestrictionGetNumElements(CeedElemRestriction rstr,
334     CeedInt *numelem);
335 CEED_EXTERN int CeedElemRestrictionGetElementSize(CeedElemRestriction rstr,
336     CeedInt *elemsize);
337 CEED_EXTERN int CeedElemRestrictionGetLVectorSize(CeedElemRestriction rstr,
338     CeedInt *lsize);
339 CEED_EXTERN int CeedElemRestrictionGetNumComponents(CeedElemRestriction rstr,
340     CeedInt *numcomp);
341 CEED_EXTERN int CeedElemRestrictionGetNumBlocks(CeedElemRestriction rstr,
342     CeedInt *numblk);
343 CEED_EXTERN int CeedElemRestrictionGetBlockSize(CeedElemRestriction rstr,
344     CeedInt *blksize);
345 CEED_EXTERN int CeedElemRestrictionGetMultiplicity(CeedElemRestriction rstr,
346     CeedVector mult);
347 CEED_EXTERN int CeedElemRestrictionView(CeedElemRestriction rstr, FILE *stream);
348 CEED_EXTERN int CeedElemRestrictionDestroy(CeedElemRestriction *rstr);
349 
350 // The formalism here is that we have the structure
351 //  \int_\Omega v^T f_0(u, \nabla u, qdata) + (\nabla v)^T f_1(u, \nabla u, qdata)
352 // where gradients are with respect to the reference element.
353 
354 /// Basis evaluation mode
355 ///
356 /// Modes can be bitwise ORed when passing to most functions.
357 /// @ingroup CeedBasis
358 typedef enum {
359   /// Perform no evaluation (either because there is no data or it is already at
360   /// quadrature points)
361   CEED_EVAL_NONE   = 0,
362   /// Interpolate from nodes to quadrature points
363   CEED_EVAL_INTERP = 1,
364   /// Evaluate gradients at quadrature points from input in a nodal basis
365   CEED_EVAL_GRAD   = 2,
366   /// Evaluate divergence at quadrature points from input in a nodal basis
367   CEED_EVAL_DIV    = 4,
368   /// Evaluate curl at quadrature points from input in a nodal basis
369   CEED_EVAL_CURL   = 8,
370   /// Using no input, evaluate quadrature weights on the reference element
371   CEED_EVAL_WEIGHT = 16,
372 } CeedEvalMode;
373 
374 CEED_EXTERN const char *const CeedEvalModes[];
375 
376 /// Type of quadrature; also used for location of nodes
377 /// @ingroup CeedBasis
378 typedef enum {
379   /// Gauss-Legendre quadrature
380   CEED_GAUSS = 0,
381   /// Gauss-Legendre-Lobatto quadrature
382   CEED_GAUSS_LOBATTO = 1,
383 } CeedQuadMode;
384 
385 CEED_EXTERN const char *const CeedQuadModes[];
386 
387 /// Type of basis shape to create non-tensor H1 element basis
388 ///
389 /// Dimension can be extracted with bitwise AND
390 /// (CeedElemTopology & 2**(dim + 2)) == TRUE
391 /// @ingroup CeedBasis
392 typedef enum {
393   /// Line
394   CEED_LINE = 1 << 16 | 0,
395   /// Triangle - 2D shape
396   CEED_TRIANGLE = 2 << 16 | 1,
397   /// Quadralateral - 2D shape
398   CEED_QUAD = 2 << 16 | 2,
399   /// Tetrahedron - 3D shape
400   CEED_TET = 3 << 16 | 3,
401   /// Pyramid - 3D shape
402   CEED_PYRAMID = 3 << 16 | 4,
403   /// Prism - 3D shape
404   CEED_PRISM = 3 << 16 | 5,
405   /// Hexehedron - 3D shape
406   CEED_HEX = 3 << 16 | 6,
407 } CeedElemTopology;
408 
409 CEED_EXTERN const char *const CeedElemTopologies[];
410 
411 CEED_EXTERN int CeedBasisCreateTensorH1Lagrange(Ceed ceed, CeedInt dim,
412     CeedInt ncomp, CeedInt P, CeedInt Q, CeedQuadMode qmode, CeedBasis *basis);
413 CEED_EXTERN int CeedBasisCreateTensorH1(Ceed ceed, CeedInt dim, CeedInt ncomp,
414                                         CeedInt P1d, CeedInt Q1d,
415                                         const CeedScalar *interp1d,
416                                         const CeedScalar *grad1d,
417                                         const CeedScalar *qref1d,
418                                         const CeedScalar *qweight1d,
419                                         CeedBasis *basis);
420 CEED_EXTERN int CeedBasisCreateH1(Ceed ceed, CeedElemTopology topo,
421                                   CeedInt ncomp,
422                                   CeedInt nnodes, CeedInt nqpts,
423                                   const CeedScalar *interp,
424                                   const CeedScalar *grad,
425                                   const CeedScalar *qref,
426                                   const CeedScalar *qweight, CeedBasis *basis);
427 CEED_EXTERN int CeedBasisView(CeedBasis basis, FILE *stream);
428 CEED_EXTERN int CeedBasisApply(CeedBasis basis, CeedInt nelem,
429                                CeedTransposeMode tmode,
430                                CeedEvalMode emode, CeedVector u, CeedVector v);
431 CEED_EXTERN int CeedBasisGetDimension(CeedBasis basis, CeedInt *dim);
432 CEED_EXTERN int CeedBasisGetTopology(CeedBasis basis, CeedElemTopology *topo);
433 CEED_EXTERN int CeedBasisGetNumComponents(CeedBasis basis, CeedInt *numcomp);
434 CEED_EXTERN int CeedBasisGetNumNodes(CeedBasis basis, CeedInt *P);
435 CEED_EXTERN int CeedBasisGetNumNodes1D(CeedBasis basis, CeedInt *P1d);
436 CEED_EXTERN int CeedBasisGetNumQuadraturePoints(CeedBasis basis, CeedInt *Q);
437 CEED_EXTERN int CeedBasisGetNumQuadraturePoints1D(CeedBasis basis,
438     CeedInt *Q1d);
439 CEED_EXTERN int CeedBasisGetQRef(CeedBasis basis, const CeedScalar **qref);
440 CEED_EXTERN int CeedBasisGetQWeights(CeedBasis basis,
441                                      const CeedScalar **qweight);
442 CEED_EXTERN int CeedBasisGetInterp(CeedBasis basis, const CeedScalar **interp);
443 CEED_EXTERN int CeedBasisGetInterp1D(CeedBasis basis,
444                                      const CeedScalar **interp1d);
445 CEED_EXTERN int CeedBasisGetGrad(CeedBasis basis, const CeedScalar **grad);
446 CEED_EXTERN int CeedBasisGetGrad1D(CeedBasis basis, const CeedScalar **grad1d);
447 CEED_EXTERN int CeedBasisDestroy(CeedBasis *basis);
448 
449 CEED_EXTERN int CeedGaussQuadrature(CeedInt Q, CeedScalar *qref1d,
450                                     CeedScalar *qweight1d);
451 CEED_EXTERN int CeedLobattoQuadrature(CeedInt Q, CeedScalar *qref1d,
452                                       CeedScalar *qweight1d);
453 CEED_EXTERN int CeedQRFactorization(Ceed ceed, CeedScalar *mat, CeedScalar *tau,
454                                     CeedInt m, CeedInt n);
455 CEED_EXTERN int CeedSymmetricSchurDecomposition(Ceed ceed, CeedScalar *mat,
456     CeedScalar *lambda, CeedInt n);
457 CEED_EXTERN int CeedSimultaneousDiagonalization(Ceed ceed, CeedScalar *matA,
458     CeedScalar *matB, CeedScalar *x, CeedScalar *lambda, CeedInt n);
459 
460 /** Handle for the object describing the user CeedQFunction
461 
462  @param ctx user-defined context set using CeedQFunctionSetContext() or NULL
463 
464  @param Q   number of quadrature points at which to evaluate
465 
466  @param in  array of pointers to each input argument in the order provided
467               by the user in CeedQFunctionAddInput().  Each array has shape
468               `[dim, ncomp, Q]` where `dim` is the geometric dimension for
469               \ref CEED_EVAL_GRAD (`dim=1` for \ref CEED_EVAL_INTERP) and
470               `ncomp` is the number of field components (`ncomp=1` for
471               scalar fields).  This results in indexing the `i`th input at
472               quadrature point `j` as `in[i][(d*ncomp + c)*Q + j]`.
473 
474  @param out array of pointers to each output array in the order provided
475               using CeedQFunctionAddOutput().  The shapes are as above for
476               \a in.
477 
478  @return An error code: 0 - success, otherwise - failure
479 
480  @ingroup CeedQFunction
481 **/
482 typedef int (*CeedQFunctionUser)(void *ctx, const CeedInt Q,
483                                  const CeedScalar *const *in,
484                                  CeedScalar *const *out);
485 
486 CEED_EXTERN int CeedQFunctionCreateInterior(Ceed ceed, CeedInt vlength,
487     CeedQFunctionUser f, const char *source, CeedQFunction *qf);
488 CEED_EXTERN int CeedQFunctionCreateInteriorByName(Ceed ceed, const char *name,
489     CeedQFunction *qf);
490 CEED_EXTERN int CeedQFunctionCreateIdentity(Ceed ceed, CeedInt size,
491     CeedEvalMode inmode, CeedEvalMode outmode, CeedQFunction *qf);
492 CEED_EXTERN int CeedQFunctionAddInput(CeedQFunction qf, const char *fieldname,
493                                       CeedInt size, CeedEvalMode emode);
494 CEED_EXTERN int CeedQFunctionAddOutput(CeedQFunction qf, const char *fieldname,
495                                        CeedInt size, CeedEvalMode emode);
496 CEED_EXTERN int CeedQFunctionSetContext(CeedQFunction qf,
497                                         CeedQFunctionContext ctx);
498 CEED_EXTERN int CeedQFunctionView(CeedQFunction qf, FILE *stream);
499 CEED_EXTERN int CeedQFunctionApply(CeedQFunction qf, CeedInt Q,
500                                    CeedVector *u, CeedVector *v);
501 CEED_EXTERN int CeedQFunctionDestroy(CeedQFunction *qf);
502 
503 CEED_EXTERN int CeedQFunctionContextCreate(Ceed ceed,
504     CeedQFunctionContext *ctx);
505 CEED_EXTERN int CeedQFunctionContextSetData(CeedQFunctionContext ctx,
506     CeedMemType mtype, CeedCopyMode cmode, size_t size, void *data);
507 CEED_EXTERN int CeedQFunctionContextGetData(CeedQFunctionContext ctx,
508     CeedMemType mtype,
509     void *data);
510 CEED_EXTERN int CeedQFunctionContextRestoreData(CeedQFunctionContext ctx,
511     void *data);
512 CEED_EXTERN int CeedQFunctionContextView(CeedQFunctionContext ctx,
513     FILE *stream);
514 CEED_EXTERN int CeedQFunctionContextDestroy(CeedQFunctionContext *ctx);
515 
516 CEED_EXTERN int CeedOperatorCreate(Ceed ceed, CeedQFunction qf,
517                                    CeedQFunction dqf, CeedQFunction dqfT,
518                                    CeedOperator *op);
519 CEED_EXTERN int CeedCompositeOperatorCreate(Ceed ceed, CeedOperator *op);
520 CEED_EXTERN int CeedOperatorSetField(CeedOperator op, const char *fieldname,
521                                      CeedElemRestriction r, CeedBasis b,
522                                      CeedVector v);
523 CEED_EXTERN int CeedCompositeOperatorAddSub(CeedOperator compositeop,
524     CeedOperator subop);
525 CEED_EXTERN int CeedOperatorLinearAssembleQFunction(CeedOperator op,
526     CeedVector *assembled, CeedElemRestriction *rstr, CeedRequest *request);
527 CEED_EXTERN int CeedOperatorLinearAssembleDiagonal(CeedOperator op,
528     CeedVector assembled, CeedRequest *request);
529 CEED_EXTERN int CeedOperatorLinearAssembleAddDiagonal(CeedOperator op,
530     CeedVector assembled, CeedRequest *request);
531 CEED_EXTERN int CeedOperatorLinearAssemblePointBlockDiagonal(CeedOperator op,
532     CeedVector assembled, CeedRequest *request);
533 CEED_EXTERN int CeedOperatorLinearAssembleAddPointBlockDiagonal(CeedOperator op,
534     CeedVector assembled, CeedRequest *request);
535 CEED_EXTERN int CeedOperatorMultigridLevelCreate(CeedOperator opFine,
536     CeedVector PMultFine, CeedElemRestriction rstrCoarse, CeedBasis basisCoarse,
537     CeedOperator *opCoarse, CeedOperator *opProlong, CeedOperator *opRestrict);
538 CEED_EXTERN int CeedOperatorMultigridLevelCreateTensorH1(
539   CeedOperator opFine, CeedVector PMultFine, CeedElemRestriction rstrCoarse,
540   CeedBasis basisCoarse, const CeedScalar *interpCtoF, CeedOperator *opCoarse,
541   CeedOperator *opProlong, CeedOperator *opRestrict);
542 CEED_EXTERN int CeedOperatorMultigridLevelCreateH1(CeedOperator opFine,
543     CeedVector PMultFine, CeedElemRestriction rstrCoarse, CeedBasis basisCoarse,
544     const CeedScalar *interpCtoF, CeedOperator *opCoarse,
545     CeedOperator *opProlong, CeedOperator *opRestrict);
546 CEED_EXTERN int CeedOperatorCreateFDMElementInverse(CeedOperator op,
547     CeedOperator *fdminv, CeedRequest *request);
548 CEED_EXTERN int CeedOperatorView(CeedOperator op, FILE *stream);
549 CEED_EXTERN int CeedOperatorApply(CeedOperator op, CeedVector in,
550                                   CeedVector out, CeedRequest *request);
551 CEED_EXTERN int CeedOperatorApplyAdd(CeedOperator op, CeedVector in,
552                                      CeedVector out, CeedRequest *request);
553 CEED_EXTERN int CeedOperatorDestroy(CeedOperator *op);
554 
555 /**
556   @brief Return integer power
557 
558   @param[in] base   The base to exponentiate
559   @param[in] power  The power to raise the base to
560 
561   @return base^power
562 
563   @ref Utility
564 **/
565 static inline CeedInt CeedIntPow(CeedInt base, CeedInt power) {
566   CeedInt result = 1;
567   while (power) {
568     if (power & 1) result *= base;
569     power >>= 1;
570     base *= base;
571   }
572   return result;
573 }
574 
575 /**
576   @brief Return minimum of two integers
577 
578   @param[in] a  The first integer to compare
579   @param[in] b  The second integer to compare
580 
581   @return The minimum of the two integers
582 
583   @ref Utility
584 **/
585 static inline CeedInt CeedIntMin(CeedInt a, CeedInt b) { return a < b ? a : b; }
586 
587 /**
588   @brief Return maximum of two integers
589 
590   @param[in] a  The first integer to compare
591   @param[in] b  The second integer to compare
592 
593   @return The maximum of the two integers
594 
595   @ref Utility
596 **/
597 static inline CeedInt CeedIntMax(CeedInt a, CeedInt b) { return a > b ? a : b; }
598 
599 #endif
600