xref: /libCEED/backends/cuda-gen/ceed-cuda-gen-operator.c (revision 621cd4619de62aa879af10f6c433e7e1cbbb1a1b)
1 // Copyright (c) 2017-2018, Lawrence Livermore National Security, LLC.
2 // Produced at the Lawrence Livermore National Laboratory. LLNL-CODE-734707.
3 // All Rights 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 #include <ceed-backend.h>
18 #include "ceed-cuda-gen.h"
19 #include "ceed-cuda-gen-operator-build.h"
20 #include "../cuda/ceed-cuda.h"
21 
22 static int CeedOperatorDestroy_Cuda_gen(CeedOperator op) {
23   int ierr;
24   CeedOperator_Cuda_gen *impl;
25   ierr = CeedOperatorGetData(op, (void *)&impl); CeedChk(ierr);
26   ierr = CeedFree(&impl); CeedChk(ierr);
27   return 0;
28 }
29 
30 static int CeedOperatorApplyAdd_Cuda_gen(CeedOperator op, CeedVector invec,
31     CeedVector outvec, CeedRequest *request) {
32   int ierr;
33   Ceed ceed;
34   ierr = CeedOperatorGetCeed(op, &ceed); CeedChk(ierr);
35   CeedOperator_Cuda_gen *data;
36   ierr = CeedOperatorGetData(op, (void *)&data); CeedChk(ierr);
37   CeedQFunction qf;
38   CeedQFunction_Cuda_gen *qf_data;
39   ierr = CeedOperatorGetQFunction(op, &qf); CeedChk(ierr);
40   ierr = CeedQFunctionGetData(qf, (void **)&qf_data); CeedChk(ierr);
41   CeedInt nelem, numinputfields, numoutputfields;
42   ierr = CeedOperatorGetNumElements(op, &nelem); CeedChk(ierr);
43   ierr = CeedQFunctionGetNumArgs(qf, &numinputfields, &numoutputfields);
44   CeedChk(ierr);
45   CeedOperatorField *opinputfields, *opoutputfields;
46   ierr = CeedOperatorGetFields(op, &opinputfields, &opoutputfields);
47   CeedChk(ierr);
48   CeedQFunctionField *qfinputfields, *qfoutputfields;
49   ierr = CeedQFunctionGetFields(qf, &qfinputfields, &qfoutputfields);
50   CeedChk(ierr);
51   CeedEvalMode emode;
52   CeedVector vec;
53 
54   //Creation of the operator
55   ierr = CeedCudaGenOperatorBuild(op); CeedChk(ierr);
56 
57   // Input vectors
58   for (CeedInt i = 0; i < numinputfields; i++) {
59     ierr = CeedQFunctionFieldGetEvalMode(qfinputfields[i], &emode);
60     CeedChk(ierr);
61     if (emode == CEED_EVAL_WEIGHT) { // Skip
62       data->fields.in[i] = NULL;
63     } else {
64       // Get input vector
65       ierr = CeedOperatorFieldGetVector(opinputfields[i], &vec); CeedChk(ierr);
66       if (vec == CEED_VECTOR_ACTIVE) vec = invec;
67       ierr = CeedVectorGetArrayRead(vec, CEED_MEM_DEVICE, &data->fields.in[i]);
68       CeedChk(ierr);
69     }
70   }
71 
72   // Output vectors
73   for (CeedInt i = 0; i < numoutputfields; i++) {
74     ierr = CeedQFunctionFieldGetEvalMode(qfoutputfields[i], &emode);
75     CeedChk(ierr);
76     if (emode == CEED_EVAL_WEIGHT) { // Skip
77       data->fields.out[i] = NULL;
78     } else {
79       // Get output vector
80       ierr = CeedOperatorFieldGetVector(opoutputfields[i], &vec); CeedChk(ierr);
81       if (vec == CEED_VECTOR_ACTIVE) vec = outvec;
82       ierr = CeedVectorGetArray(vec, CEED_MEM_DEVICE, &data->fields.out[i]);
83       CeedChk(ierr);
84     }
85   }
86 
87   // Copy the context
88   size_t ctxsize;
89   ierr = CeedQFunctionGetContextSize(qf, &ctxsize); CeedChk(ierr);
90   if (ctxsize > 0) {
91     if (!qf_data->d_c) {
92       ierr = cudaMalloc(&qf_data->d_c, ctxsize); CeedChk_Cu(ceed, ierr);
93     }
94     void *ctx;
95     ierr = CeedQFunctionGetInnerContext(qf, &ctx); CeedChk(ierr);
96     ierr = cudaMemcpy(qf_data->d_c, ctx, ctxsize, cudaMemcpyHostToDevice);
97     CeedChk_Cu(ceed, ierr);
98   }
99 
100   // Apply operator
101   void *opargs[] = {(void *) &nelem, &qf_data->d_c, &data->indices,
102                     &data->fields, &data->B, &data->G, &data->W
103                    };
104   const CeedInt dim = data->dim;
105   const CeedInt Q1d = data->Q1d;
106   if (dim==1) {
107     const CeedInt elemsPerBlock = 32;
108     CeedInt grid = nelem/elemsPerBlock + ( (nelem/elemsPerBlock*elemsPerBlock<nelem)
109                                            ? 1 : 0 );
110     CeedInt sharedMem = elemsPerBlock*Q1d*sizeof(CeedScalar);
111     ierr = CeedRunKernelDimSharedCuda(ceed, data->op, grid, Q1d, 1, elemsPerBlock,
112                                       sharedMem, opargs);
113   } else if (dim==2) {
114     const CeedInt elemsPerBlock = Q1d<4? 16 : 2;
115     CeedInt grid = nelem/elemsPerBlock + ( (nelem/elemsPerBlock*elemsPerBlock<nelem)
116                                            ? 1 : 0 );
117     CeedInt sharedMem = elemsPerBlock*Q1d*Q1d*sizeof(CeedScalar);
118     ierr = CeedRunKernelDimSharedCuda(ceed, data->op, grid, Q1d, Q1d,
119                                       elemsPerBlock, sharedMem, opargs);
120   } else if (dim==3) {
121     const CeedInt elemsPerBlock = Q1d<6? 4 : (Q1d<8? 2 : 1);
122     CeedInt grid = nelem/elemsPerBlock + ( (nelem/elemsPerBlock*elemsPerBlock<nelem)
123                                            ? 1 : 0 );
124     CeedInt sharedMem = elemsPerBlock*Q1d*Q1d*sizeof(CeedScalar);
125     ierr = CeedRunKernelDimSharedCuda(ceed, data->op, grid, Q1d, Q1d,
126                                       elemsPerBlock, sharedMem, opargs);
127   }
128   CeedChk(ierr);
129 
130   // Restore input arrays
131   for (CeedInt i = 0; i < numinputfields; i++) {
132     ierr = CeedQFunctionFieldGetEvalMode(qfinputfields[i], &emode);
133     CeedChk(ierr);
134     if (emode == CEED_EVAL_WEIGHT) { // Skip
135     } else {
136       ierr = CeedOperatorFieldGetVector(opinputfields[i], &vec); CeedChk(ierr);
137       if (vec == CEED_VECTOR_ACTIVE) vec = invec;
138       ierr = CeedVectorRestoreArrayRead(vec, &data->fields.in[i]);
139       CeedChk(ierr);
140     }
141   }
142 
143   // Restore output arrays
144   for (CeedInt i = 0; i < numoutputfields; i++) {
145     ierr = CeedQFunctionFieldGetEvalMode(qfoutputfields[i], &emode);
146     CeedChk(ierr);
147     if (emode == CEED_EVAL_WEIGHT) { // Skip
148     } else {
149       ierr = CeedOperatorFieldGetVector(opoutputfields[i], &vec); CeedChk(ierr);
150       if (vec == CEED_VECTOR_ACTIVE) vec = outvec;
151       ierr = CeedVectorRestoreArray(vec, &data->fields.out[i]);
152       CeedChk(ierr);
153     }
154   }
155 
156   return 0;
157 }
158 
159 static int CeedOperatorAssembleLinearQFunction_Cuda(CeedOperator op) {
160   int ierr;
161   Ceed ceed;
162   ierr = CeedOperatorGetCeed(op, &ceed); CeedChk(ierr);
163   return CeedError(ceed, 1, "Backend does not implement QFunction assembly");
164 }
165 
166 int CeedOperatorCreate_Cuda_gen(CeedOperator op) {
167   int ierr;
168   Ceed ceed;
169   ierr = CeedOperatorGetCeed(op, &ceed); CeedChk(ierr);
170   CeedOperator_Cuda_gen *impl;
171 
172   ierr = CeedCalloc(1, &impl); CeedChk(ierr);
173   ierr = CeedOperatorSetData(op, (void *)&impl);
174 
175   ierr = CeedSetBackendFunction(ceed, "Operator", op, "AssembleLinearQFunction",
176                                 CeedOperatorAssembleLinearQFunction_Cuda);
177   CeedChk(ierr);
178   ierr = CeedSetBackendFunction(ceed, "Operator", op, "ApplyAdd",
179                                 CeedOperatorApplyAdd_Cuda_gen); CeedChk(ierr);
180   ierr = CeedSetBackendFunction(ceed, "Operator", op, "Destroy",
181                                 CeedOperatorDestroy_Cuda_gen); CeedChk(ierr);
182   return 0;
183 }
184