Neko 1.99.7
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
opgrad_kernel.h
Go to the documentation of this file.
1#ifndef __MATH_OPGRAD_KERNEL_H__
2#define __MATH_OPGRAD_KERNEL_H__
3
4#include "elem_block.h"
5/*
6 Copyright (c) 2021-2023, The Neko Authors
7 All rights reserved.
8
9 Redistribution and use in source and binary forms, with or without
10 modification, are permitted provided that the following conditions
11 are met:
12
13 * Redistributions of source code must retain the above copyright
14 notice, this list of conditions and the following disclaimer.
15
16 * Redistributions in binary form must reproduce the above
17 copyright notice, this list of conditions and the following
18 disclaimer in the documentation and/or other materials provided
19 with the distribution.
20
21 * Neither the name of the authors nor the names of its
22 contributors may be used to endorse or promote products derived
23 from this software without specific prior written permission.
24
25 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
26 "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
27 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
28 FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
29 COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
30 INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
31 BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
32 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
33 CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
34 LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
35 ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
36 POSSIBILITY OF SUCH DAMAGE.
37*/
38
43template< typename T, const int LX, const int CHUNKS >
47 const T * __restrict__ u,
48 const T * __restrict__ dx,
49 const T * __restrict__ dy,
50 const T * __restrict__ dz,
51 const T * __restrict__ drdx,
52 const T * __restrict__ dsdx,
53 const T * __restrict__ dtdx,
54 const T * __restrict__ drdy,
55 const T * __restrict__ dsdy,
56 const T * __restrict__ dtdy,
57 const T * __restrict__ drdz,
58 const T * __restrict__ dsdz,
59 const T * __restrict__ dtdz,
60 const T * __restrict__ w3) {
61
62 __shared__ T shu[LX * LX * LX];
63
67
68
69 int i,j,k;
70
71 const int e = blockIdx.x;
72 const int iii = threadIdx.x;
73 const int nchunks = (LX * LX * LX - 1) / CHUNKS + 1;
74
75 if (iii < (LX * LX)) {
76 shdx[iii] = dx[iii];
77 shdy[iii] = dy[iii];
78 shdz[iii] = dz[iii];
79 }
80
81 j = iii;
82 while(j < (LX * LX * LX)) {
83 shu[j] = u[j + e * LX * LX * LX];
84 j = j + CHUNKS;
85 }
86
88
89 for (int n = 0; n < nchunks; n++) {
90 const int ijk = iii + n * CHUNKS;
91 const int jk = ijk / LX;
92 i = ijk - jk * LX;
93 k = jk / LX;
94 j = jk - k * LX;
95 if ( i < LX && j < LX && k < LX ) {
96 T rtmp = 0.0;
97 T stmp = 0.0;
98 T ttmp = 0.0;
99 for (int l = 0; l < LX; l++) {
100 rtmp += shdx[i + l * LX] * shu[l + j * LX + k * LX * LX];
101 stmp += shdy[j + l * LX] * shu[i + l * LX + k * LX * LX];
102 ttmp += shdz[k + l * LX] * shu[i + j * LX + l * LX * LX];
103 }
104
105 ux[ijk + e * LX * LX * LX] = w3[ijk]
106 * (drdx[ijk + e * LX * LX * LX] * rtmp
107 + dsdx[ijk + e * LX * LX * LX] * stmp
108 + dtdx[ijk + e * LX * LX * LX] * ttmp);
109
110 uy[ijk + e * LX * LX * LX] = w3[ijk]
111 * (drdy[ijk + e * LX * LX * LX] * rtmp
112 + dsdy[ijk + e * LX * LX * LX] * stmp
113 + dtdy[ijk + e * LX * LX * LX] * ttmp);
114
115 uz[ijk + e * LX * LX * LX] = w3[ijk]
116 * (drdz[ijk + e * LX * LX * LX] * rtmp
117 + dsdz[ijk + e * LX * LX * LX] * stmp
118 + dtdz[ijk + e * LX * LX * LX] * ttmp);
119
120 }
121 }
122
123}
124
125template< typename T, const int LX, const int EB >
130 const T * __restrict__ u,
144 const int nelv) {
145
146 /* One slice per element in the block */
147 __shared__ T shu[EB * LX * LX];
148
149 /* Element independent, one copy per block */
153
154 static_assert(sizeof(shu) + sizeof(shdx) + sizeof(shdy) + sizeof(shdz)
156 "kstep block exceeds the shared memory budget");
157
158 /* Threads past the last element still have to reach the barriers in the k
159 loop, so clamp their reads and drop their stores rather than returning
160 early. At EB == 1 all of this is constant folded away */
161 const int eb = (EB == 1) ? 0 : threadIdx.z;
162 const int e_blk = blockIdx.x * EB + eb;
163 const bool active = (EB == 1) ? true : (e_blk < nelv);
164 const int e = active ? e_blk : (nelv - 1);
165 const int j = threadIdx.y;
166 const int i = threadIdx.x;
167 const int ij = i + j * LX;
168 const int sh = eb * LX * LX;
169 const int ele = e*LX*LX*LX;
170
171 if (eb == 0) {
172 shdx[ij] = dx[ij];
173 shdy[ij] = dy[ij];
174 shdz[ij] = dz[ij];
175 }
176
178
179#pragma unroll LX
180 for (int k = 0; k < LX; ++k) {
181 ru[k] = u[ij + k*LX*LX + ele];
182 }
183
185
186 #pragma unroll
187 for (int k = 0; k < LX; ++k) {
188 const int ijk = ij + k*LX*LX;
189 const T W3 = w3[ijk];
190 T ttmp = 0.0;
191 shu[sh + ij] = ru[k];
192#pragma unroll
193 for (int l = 0; l < LX; l++) {
194 ttmp += shdz[k+l*LX] * ru[l];
195 }
197
198 T rtmp = 0.0;
199 T stmp = 0.0;
200#pragma unroll
201 for (int l = 0; l < LX; l++) {
202 rtmp += shdx[i+l*LX] * shu[sh + l+j*LX];
203 stmp += shdy[j+l*LX] * shu[sh + i+l*LX];
204 }
205
206 if (active) {
207 ux[ijk + ele] = W3 * (drdx[ijk + ele] * rtmp
208 + dsdx[ijk + ele] * stmp
209 + dtdx[ijk + ele] * ttmp);
210
211 uy[ijk + ele] = W3 * (drdy[ijk + ele] * rtmp
212 + dsdy[ijk + ele] * stmp
213 + dtdy[ijk + ele] * ttmp);
214
215 uz[ijk + ele] = W3 * (drdz[ijk + ele] * rtmp
216 + dsdz[ijk + ele] * stmp
217 + dtdz[ijk + ele] * ttmp);
218 }
220 }
221}
222
223
224#endif // __MATH_OPGRAD_KERNEL_H__
__global__ void ale_add_kinematics_kernel(const int n, T *__restrict__ wx, T *__restrict__ wy, T *__restrict__ wz, const T *__restrict__ x_ref, const T *__restrict__ y_ref, const T *__restrict__ z_ref, const T *__restrict__ phi, const T *__restrict__ x, const T *__restrict__ y, const T *__restrict__ z, const kinematics_params_t kin_params)
__shared__ T shu[LX *LX]
#define NEKO_EB_MAX_SMEM
Definition elem_block.h:60
#define NEKO_EB_BOUNDS(NT)
Definition elem_block.h:95
const bool active
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dsdx
__global__ void T *__restrict__ 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
__global__ void T *__restrict__ uy
__global__ void opgrad_kernel_1d(T *__restrict__ ux, T *__restrict__ uy, T *__restrict__ uz, const T *__restrict__ u, 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__ w3)
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ dx
const int sh
__global__ void T *__restrict__ T *__restrict__ uz
__global__ void T *__restrict__ 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
__global__ void T *__restrict__ 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 int nelv
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dz
T ru[LX]
const int eb
__shared__ T shdy[LX *LX]
__global__ void T *__restrict__ 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
__global__ void T *__restrict__ 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
const int i
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dtdx
const int ij
__shared__ T shdx[LX *LX]
const int e
__global__ void T *__restrict__ 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
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ drdx
__shared__ T shdz[LX *LX]
const int e_blk
__global__ void T *__restrict__ 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
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ const T *__restrict__ const T *__restrict__ dy
__global__ void T *__restrict__ T *__restrict__ const T *__restrict__ u
const int ele
__global__ void T *__restrict__ 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__ w3
const int j
__syncthreads()