Neko  0.9.0
A portable framework for high-order spectral element flow simulations
conv1_kernel.h
Go to the documentation of this file.
1 #ifndef __MATH_CONV1_KERNEL_H__
2 #define __MATH_CONV1_KERNEL_H__
3 /*
4  Copyright (c) 2021-2023, The Neko Authors
5  All rights reserved.
6 
7  Redistribution and use in source and binary forms, with or without
8  modification, are permitted provided that the following conditions
9  are met:
10 
11  * Redistributions of source code must retain the above copyright
12  notice, this list of conditions and the following disclaimer.
13 
14  * Redistributions in binary form must reproduce the above
15  copyright notice, this list of conditions and the following
16  disclaimer in the documentation and/or other materials provided
17  with the distribution.
18 
19  * Neither the name of the authors nor the names of its
20  contributors may be used to endorse or promote products derived
21  from this software without specific prior written permission.
22 
23  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
24  "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
25  LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
26  FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
27  COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
28  INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
29  BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
30  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
31  CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
32  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
33  ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
34  POSSIBILITY OF SUCH DAMAGE.
35 */
36 
41 template< typename T, const int LX, const int CHUNKS >
42 __global__ void conv1_kernel_1d(T * __restrict__ du,
43  const T * __restrict__ u,
44  const T * __restrict__ vx,
45  const T * __restrict__ vy,
46  const T * __restrict__ vz,
47  const T * __restrict__ dx,
48  const T * __restrict__ dy,
49  const T * __restrict__ dz,
50  const T * __restrict__ drdx,
51  const T * __restrict__ dsdx,
52  const T * __restrict__ dtdx,
53  const T * __restrict__ drdy,
54  const T * __restrict__ dsdy,
55  const T * __restrict__ dtdy,
56  const T * __restrict__ drdz,
57  const T * __restrict__ dsdz,
58  const T * __restrict__ dtdz,
59  const T * __restrict__ jacinv) {
60 
61  __shared__ T shu[LX * LX * LX];
62 
63  __shared__ T shvx[LX * LX * LX];
64  __shared__ T shvy[LX * LX * LX];
65  __shared__ T shvz[LX * LX * LX];
66 
67  __shared__ T shdx[LX * LX];
68  __shared__ T shdy[LX * LX];
69  __shared__ T shdz[LX * LX];
70 
71  __shared__ T shjacinv[LX * LX * LX];
72 
73  const int e = blockIdx.x;
74  const int iii = threadIdx.x;
75  const int nchunks = (LX * LX * LX - 1) / CHUNKS + 1;
76  const int ele = e*LX*LX*LX;
77 
78  if (iii < (LX * LX)) {
79  shdx[iii] = dx[iii];
80  shdy[iii] = dy[iii];
81  shdz[iii] = dz[iii];
82  }
83 
84  int l = iii;
85  while(l < (LX * LX * LX)) {
86  shu[l] = u[l + ele];
87 
88  shvx[l] = vx[l + ele];
89  shvy[l] = vy[l + ele];
90  shvz[l] = vz[l + ele];
91 
92  shjacinv[l] = jacinv[l + ele];
93 
94  l = l + CHUNKS;
95  }
96 
97  __syncthreads();
98 
99  for (int n = 0; n < nchunks; n++) {
100  const int ijk = iii + n * CHUNKS;
101  const int jk = ijk / LX;
102  const int i = ijk - jk * LX;
103  const int k = jk / LX;
104  const int j = jk - k * LX;
105  if ( i < LX && j < LX && k < LX) {
106  T rtmp = 0.0;
107  T stmp = 0.0;
108  T ttmp = 0.0;
109  for (int l = 0; l < LX; l++) {
110  rtmp += shdx[i + l * LX] * shu[l + j * LX + k * LX * LX];
111  stmp += shdy[j + l * LX] * shu[i + l * LX + k * LX * LX];
112  ttmp += shdz[k + l * LX] * shu[i + j * LX + l * LX * LX];
113  }
114 
115  du[ijk + e * LX * LX * LX] = shjacinv[ijk] *
116  (shvx[ijk] * (drdx[ijk + ele] * rtmp
117  + dsdx[ijk + ele] * stmp
118  + dtdx[ijk + ele] * ttmp)
119  + shvy[ijk] * (drdy[ijk + ele] * rtmp
120  + dsdy[ijk + ele] * stmp
121  + dtdy[ijk + ele] * ttmp)
122  + shvz[ijk] * (drdz[ijk + ele] * rtmp
123  + dsdz[ijk + ele] * stmp
124  + dtdz[ijk + ele] * ttmp));
125  }
126  }
127 }
128 
129 template< typename T, const int LX >
130 __global__ void __launch_bounds__(LX*LX,3)
131  conv1_kernel_kstep(T * __restrict__ du,
132  const T * __restrict__ u,
133  const T * __restrict__ vx,
134  const T * __restrict__ vy,
135  const T * __restrict__ vz,
136  const T * __restrict__ dx,
137  const T * __restrict__ dy,
138  const T * __restrict__ dz,
139  const T * __restrict__ drdx,
140  const T * __restrict__ dsdx,
141  const T * __restrict__ dtdx,
142  const T * __restrict__ drdy,
143  const T * __restrict__ dsdy,
144  const T * __restrict__ dtdy,
145  const T * __restrict__ drdz,
146  const T * __restrict__ dsdz,
147  const T * __restrict__ dtdz,
148  const T * __restrict__ jacinv) {
149 
150  __shared__ T shu[LX * LX];
151 
152  __shared__ T shdx[LX*LX];
153  __shared__ T shdy[LX*LX];
154  __shared__ T shdz[LX*LX];
155 
156  const int e = blockIdx.x;
157  const int j = threadIdx.y;
158  const int i = threadIdx.x;
159  const int ij = i + j * LX;
160  const int ele = e*LX*LX*LX;
161 
162  shdx[ij] = dx[ij];
163  shdy[ij] = dy[ij];
164  shdz[ij] = dz[ij];
165 
166  T ru[LX];
167  T rvx[LX];
168  T rvy[LX];
169  T rvz[LX];
170  T rjacinv[LX];
171 
172 #pragma unroll LX
173  for (int k = 0; k < LX; ++k) {
174  ru[k] = u[ij + k*LX*LX + ele];
175  rvx[k] = vx[ij + k*LX*LX + ele];
176  rvy[k] = vy[ij + k*LX*LX + ele];
177  rvz[k] = vz[ij + k*LX*LX + ele];
178  rjacinv[k] = jacinv[ij + k*LX*LX + ele];
179  }
180 
182 
183 #pragma unroll
184  for (int k = 0; k < LX; ++k) {
185  const int ijk = ij + k*LX*LX;
186  T ttmp = 0.0;
187  shu[ij] = ru[k];
188  for (int l = 0; l < LX; l++) {
189  ttmp += shdz[k+l*LX] * ru[l];
190  }
191  __syncthreads();
192 
193  T rtmp = 0.0;
194  T stmp = 0.0;
195 #pragma unroll
196  for (int l = 0; l < LX; l++) {
197  rtmp += shdx[i+l*LX] * shu[l+j*LX];
198  stmp += shdy[j+l*LX] * shu[i+l*LX];
199  }
200 
201  du[ijk + ele] = rjacinv[k] *
202  (rvx[k] * (drdx[ijk + ele] * rtmp
203  + dsdx[ijk + ele] * stmp
204  + dtdx[ijk + ele] * ttmp)
205  + rvy[k] * (drdy[ijk + ele] * rtmp
206  + dsdy[ijk + ele] * stmp
207  + dtdy[ijk + ele] * ttmp)
208  + rvz[k] * (drdz[ijk + ele] * rtmp
209  + dsdz[ijk + ele] * stmp
210  + dtdz[ijk + ele] * ttmp));
211  __syncthreads();
212  }
213 }
214 
215 
216 #endif // __MATH_CONV1_KERNEL_H__
__shared__ T shu[LX *LX]
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dz
Definition: conv1_kernel.h:138
T rvy[LX]
Definition: conv1_kernel.h:168
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ jacinv
Definition: conv1_kernel.h:148
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dy
Definition: conv1_kernel.h:137
T ru[LX]
Definition: conv1_kernel.h:166
__shared__ T shdz[LX *LX]
Definition: conv1_kernel.h:154
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ drdx
Definition: conv1_kernel.h:139
const int i
Definition: conv1_kernel.h:158
T rvx[LX]
Definition: conv1_kernel.h:167
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dx
Definition: conv1_kernel.h:136
const int ij
Definition: conv1_kernel.h:159
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dtdx
Definition: conv1_kernel.h:141
T rvz[LX]
Definition: conv1_kernel.h:169
__shared__ T shdy[LX *LX]
Definition: conv1_kernel.h:153
T rjacinv[LX]
Definition: conv1_kernel.h:170
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dtdz
Definition: conv1_kernel.h:147
const int e
Definition: conv1_kernel.h:156
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ vz
Definition: conv1_kernel.h:135
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dsdx
Definition: conv1_kernel.h:140
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dtdy
Definition: conv1_kernel.h:144
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dsdz
Definition: conv1_kernel.h:146
__global__ void const T *__restrict__ const T *__restrict__ vx
Definition: conv1_kernel.h:133
__global__ void conv1_kernel_1d(T *__restrict__ du, const T *__restrict__ u, const T *__restrict__ vx, const T *__restrict__ vy, const T *__restrict__ vz, const T *__restrict__ dx, const T *__restrict__ dy, const T *__restrict__ dz, const T *__restrict__ drdx, const T *__restrict__ dsdx, const T *__restrict__ dtdx, const T *__restrict__ drdy, const T *__restrict__ dsdy, const T *__restrict__ dtdy, const T *__restrict__ drdz, const T *__restrict__ dsdz, const T *__restrict__ dtdz, const T *__restrict__ jacinv)
Definition: conv1_kernel.h:42
__shared__ T shdx[LX *LX]
Definition: conv1_kernel.h:152
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dsdy
Definition: conv1_kernel.h:143
__global__ void __launch_bounds__(LX *LX, 3) conv1_kernel_kstep(T *__restrict__ du
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ vy
Definition: conv1_kernel.h:134
const int ele
Definition: conv1_kernel.h:160
__global__ void const T *__restrict__ u
Definition: conv1_kernel.h:132
const int j
Definition: conv1_kernel.h:157
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ drdy
Definition: conv1_kernel.h:142
__syncthreads()
__global__ void const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ drdz
Definition: conv1_kernel.h:145