Neko 0.9.99
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
cpu_dudxyz.f90
Go to the documentation of this file.
1! Copyright (c) 2021, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
34submodule(opr_cpu) cpu_dudxyz
35 implicit none
36
37contains
38
39 module subroutine opr_cpu_dudxyz(du, u, dr, ds, dt, coef)
40 type(coef_t), intent(in), target :: coef
41 real(kind=rp), dimension(coef%Xh%lx, coef%Xh%ly, & coef%Xh%lz, coef%msh%nelv), intent(inout) :: du
42 real(kind=rp), dimension(coef%Xh%lx, coef%Xh%ly, & coef%Xh%lz, coef%msh%nelv), intent(in) :: u, dr, ds, dt
43
44 associate(xh => coef%Xh, msh => coef%msh, dof => coef%dof)
45 select case (coef%Xh%lx)
46 case (14)
47 call cpu_dudxyz_lx14(du, u, dr, ds, dt, &
48 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
49 case (13)
50 call cpu_dudxyz_lx13(du, u, dr, ds, dt, &
51 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
52 case (12)
53 call cpu_dudxyz_lx12(du, u, dr, ds, dt, &
54 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
55 case (11)
56 call cpu_dudxyz_lx11(du, u, dr, ds, dt, &
57 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
58 case (10)
59 call cpu_dudxyz_lx10(du, u, dr, ds, dt, &
60 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
61 case (9)
62 call cpu_dudxyz_lx9(du, u, dr, ds, dt, &
63 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
64 case (8)
65 call cpu_dudxyz_lx8(du, u, dr, ds, dt, &
66 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
67 case (7)
68 call cpu_dudxyz_lx7(du, u, dr, ds, dt, &
69 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
70 case (6)
71 call cpu_dudxyz_lx6(du, u, dr, ds, dt, &
72 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
73 case (5)
74 call cpu_dudxyz_lx5(du, u, dr, ds, dt, &
75 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
76 case (4)
77 call cpu_dudxyz_lx4(du, u, dr, ds, dt, &
78 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
79 case (3)
80 call cpu_dudxyz_lx3(du, u, dr, ds, dt, &
81 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
82 case (2)
83 call cpu_dudxyz_lx2(du, u, dr, ds, dt, &
84 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv)
85 case default
86 call cpu_dudxyz_lx(du, u, dr, ds, dt, &
87 xh%dx, xh%dy, xh%dz, coef%jacinv, msh%nelv, xh%lx)
88 end select
89
90 end associate
91
92 end subroutine opr_cpu_dudxyz
93
94 subroutine cpu_dudxyz_lx(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel, lx)
95 integer, intent(in) :: nel, lx
96 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
97 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
98 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
99 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
100 real(kind=rp), dimension(lx, lx, lx) :: drst
101 real(kind=rp) :: tmp
102 integer :: e, i, j, k, l
103
104 do e = 1, nel
105 do j = 1, lx * lx
106 do i = 1, lx
107 tmp = 0.0_rp
108 do k = 1, lx
109 tmp = tmp + dx(i,k) * u(k,j,1,e)
110 end do
111 du(i,j,1,e) = tmp
112 end do
113 end do
114
115 do i = 1, lx * lx * lx
116 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
117 end do
118
119 do k = 1, lx
120 do j = 1, lx
121 do i = 1, lx
122 tmp = 0.0_rp
123 do l = 1, lx
124 tmp = tmp + dy(j,l) * u(i,l,k,e)
125 end do
126 drst(i,j,k) = tmp
127 end do
128 end do
129 end do
130
131 do i = 1, lx * lx * lx
132 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
133 end do
134
135 do k = 1, lx
136 do i = 1, lx*lx
137 tmp = 0.0_rp
138 do l = 1, lx
139 tmp = tmp + dz(k,l) * u(i,1,l,e)
140 end do
141 drst(i,1,k) = tmp
142 end do
143 end do
144
145 do i = 1, lx * lx * lx
146 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
147 end do
148
149 do i = 1, lx * lx * lx
150 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
151 end do
152
153 end do
154
155 end subroutine cpu_dudxyz_lx
156
157 subroutine cpu_dudxyz_lx14(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
158 integer, parameter :: lx = 14
159 integer, intent(in) :: nel
160 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
161 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
162 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
163 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
164 real(kind=rp), dimension(lx, lx, lx) :: drst
165 integer :: e, i, j, k
166
167 do e = 1, nel
168 do j = 1, lx * lx
169 do i = 1, lx
170 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
171 + dx(i,2) * u(2,j,1,e) &
172 + dx(i,3) * u(3,j,1,e) &
173 + dx(i,4) * u(4,j,1,e) &
174 + dx(i,5) * u(5,j,1,e) &
175 + dx(i,6) * u(6,j,1,e) &
176 + dx(i,7) * u(7,j,1,e) &
177 + dx(i,8) * u(8,j,1,e) &
178 + dx(i,9) * u(9,j,1,e) &
179 + dx(i,10) * u(10,j,1,e) &
180 + dx(i,11) * u(11,j,1,e) &
181 + dx(i,12) * u(12,j,1,e) &
182 + dx(i,13) * u(13,j,1,e) &
183 + dx(i,14) * u(14,j,1,e)
184 end do
185 end do
186
187 do i = 1, lx * lx * lx
188 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
189 end do
190
191 do k = 1, lx
192 do j = 1, lx
193 do i = 1, lx
194 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
195 + dy(j,2) * u(i,2,k,e) &
196 + dy(j,3) * u(i,3,k,e) &
197 + dy(j,4) * u(i,4,k,e) &
198 + dy(j,5) * u(i,5,k,e) &
199 + dy(j,6) * u(i,6,k,e) &
200 + dy(j,7) * u(i,7,k,e) &
201 + dy(j,8) * u(i,8,k,e) &
202 + dy(j,9) * u(i,9,k,e) &
203 + dy(j,10) * u(i,10,k,e) &
204 + dy(j,11) * u(i,11,k,e) &
205 + dy(j,12) * u(i,12,k,e) &
206 + dy(j,13) * u(i,13,k,e) &
207 + dy(j,14) * u(i,14,k,e)
208 end do
209 end do
210 end do
211
212 do i = 1, lx * lx * lx
213 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
214 end do
215
216 do k = 1, lx
217 do i = 1, lx*lx
218 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
219 + dz(k,2) * u(i,1,2,e) &
220 + dz(k,3) * u(i,1,3,e) &
221 + dz(k,4) * u(i,1,4,e) &
222 + dz(k,5) * u(i,1,5,e) &
223 + dz(k,6) * u(i,1,6,e) &
224 + dz(k,7) * u(i,1,7,e) &
225 + dz(k,8) * u(i,1,8,e) &
226 + dz(k,9) * u(i,1,9,e) &
227 + dz(k,10) * u(i,1,10,e) &
228 + dz(k,11) * u(i,1,11,e) &
229 + dz(k,12) * u(i,1,12,e) &
230 + dz(k,13) * u(i,1,13,e) &
231 + dz(k,14) * u(i,1,14,e)
232 end do
233 end do
234
235 do i = 1, lx * lx * lx
236 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
237 end do
238
239 do i = 1, lx * lx * lx
240 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
241 end do
242
243 end do
244
245 end subroutine cpu_dudxyz_lx14
246
247 subroutine cpu_dudxyz_lx13(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
248 integer, parameter :: lx = 13
249 integer, intent(in) :: nel
250 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
251 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
252 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
253 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
254 real(kind=rp), dimension(lx, lx, lx) :: drst
255 integer :: e, i, j, k
256
257 do e = 1, nel
258 do j = 1, lx * lx
259 do i = 1, lx
260 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
261 + dx(i,2) * u(2,j,1,e) &
262 + dx(i,3) * u(3,j,1,e) &
263 + dx(i,4) * u(4,j,1,e) &
264 + dx(i,5) * u(5,j,1,e) &
265 + dx(i,6) * u(6,j,1,e) &
266 + dx(i,7) * u(7,j,1,e) &
267 + dx(i,8) * u(8,j,1,e) &
268 + dx(i,9) * u(9,j,1,e) &
269 + dx(i,10) * u(10,j,1,e) &
270 + dx(i,11) * u(11,j,1,e) &
271 + dx(i,12) * u(12,j,1,e) &
272 + dx(i,13) * u(13,j,1,e)
273 end do
274 end do
275
276 do i = 1, lx * lx * lx
277 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
278 end do
279
280 do k = 1, lx
281 do j = 1, lx
282 do i = 1, lx
283 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
284 + dy(j,2) * u(i,2,k,e) &
285 + dy(j,3) * u(i,3,k,e) &
286 + dy(j,4) * u(i,4,k,e) &
287 + dy(j,5) * u(i,5,k,e) &
288 + dy(j,6) * u(i,6,k,e) &
289 + dy(j,7) * u(i,7,k,e) &
290 + dy(j,8) * u(i,8,k,e) &
291 + dy(j,9) * u(i,9,k,e) &
292 + dy(j,10) * u(i,10,k,e) &
293 + dy(j,11) * u(i,11,k,e) &
294 + dy(j,12) * u(i,12,k,e) &
295 + dy(j,13) * u(i,13,k,e)
296 end do
297 end do
298 end do
299
300 do i = 1, lx * lx * lx
301 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
302 end do
303
304 do k = 1, lx
305 do i = 1, lx*lx
306 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
307 + dz(k,2) * u(i,1,2,e) &
308 + dz(k,3) * u(i,1,3,e) &
309 + dz(k,4) * u(i,1,4,e) &
310 + dz(k,5) * u(i,1,5,e) &
311 + dz(k,6) * u(i,1,6,e) &
312 + dz(k,7) * u(i,1,7,e) &
313 + dz(k,8) * u(i,1,8,e) &
314 + dz(k,9) * u(i,1,9,e) &
315 + dz(k,10) * u(i,1,10,e) &
316 + dz(k,11) * u(i,1,11,e) &
317 + dz(k,12) * u(i,1,12,e) &
318 + dz(k,13) * u(i,1,13,e)
319 end do
320 end do
321
322 do i = 1, lx * lx * lx
323 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
324 end do
325
326 do i = 1, lx * lx * lx
327 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
328 end do
329
330 end do
331
332 end subroutine cpu_dudxyz_lx13
333
334 subroutine cpu_dudxyz_lx12(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
335 integer, parameter :: lx = 12
336 integer, intent(in) :: nel
337 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
338 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
339 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
340 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
341 real(kind=rp), dimension(lx, lx, lx) :: drst
342 integer :: e, i, j, k
343
344 do e = 1, nel
345 do j = 1, lx * lx
346 do i = 1, lx
347 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
348 + dx(i,2) * u(2,j,1,e) &
349 + dx(i,3) * u(3,j,1,e) &
350 + dx(i,4) * u(4,j,1,e) &
351 + dx(i,5) * u(5,j,1,e) &
352 + dx(i,6) * u(6,j,1,e) &
353 + dx(i,7) * u(7,j,1,e) &
354 + dx(i,8) * u(8,j,1,e) &
355 + dx(i,9) * u(9,j,1,e) &
356 + dx(i,10) * u(10,j,1,e) &
357 + dx(i,11) * u(11,j,1,e) &
358 + dx(i,12) * u(12,j,1,e)
359 end do
360 end do
361
362 do i = 1, lx * lx * lx
363 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
364 end do
365
366 do k = 1, lx
367 do j = 1, lx
368 do i = 1, lx
369 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
370 + dy(j,2) * u(i,2,k,e) &
371 + dy(j,3) * u(i,3,k,e) &
372 + dy(j,4) * u(i,4,k,e) &
373 + dy(j,5) * u(i,5,k,e) &
374 + dy(j,6) * u(i,6,k,e) &
375 + dy(j,7) * u(i,7,k,e) &
376 + dy(j,8) * u(i,8,k,e) &
377 + dy(j,9) * u(i,9,k,e) &
378 + dy(j,10) * u(i,10,k,e) &
379 + dy(j,11) * u(i,11,k,e) &
380 + dy(j,12) * u(i,12,k,e)
381 end do
382 end do
383 end do
384
385 do i = 1, lx * lx * lx
386 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
387 end do
388
389 do k = 1, lx
390 do i = 1, lx*lx
391 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
392 + dz(k,2) * u(i,1,2,e) &
393 + dz(k,3) * u(i,1,3,e) &
394 + dz(k,4) * u(i,1,4,e) &
395 + dz(k,5) * u(i,1,5,e) &
396 + dz(k,6) * u(i,1,6,e) &
397 + dz(k,7) * u(i,1,7,e) &
398 + dz(k,8) * u(i,1,8,e) &
399 + dz(k,9) * u(i,1,9,e) &
400 + dz(k,10) * u(i,1,10,e) &
401 + dz(k,11) * u(i,1,11,e) &
402 + dz(k,12) * u(i,1,12,e)
403 end do
404 end do
405
406 do i = 1, lx * lx * lx
407 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
408 end do
409
410 do i = 1, lx * lx * lx
411 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
412 end do
413
414 end do
415
416 end subroutine cpu_dudxyz_lx12
417
418 subroutine cpu_dudxyz_lx11(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
419 integer, parameter :: lx = 11
420 integer, intent(in) :: nel
421 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
422 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
423 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
424 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
425 real(kind=rp), dimension(lx, lx, lx) :: drst
426 integer :: e, i, j, k
427
428 do e = 1, nel
429 do j = 1, lx * lx
430 do i = 1, lx
431 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
432 + dx(i,2) * u(2,j,1,e) &
433 + dx(i,3) * u(3,j,1,e) &
434 + dx(i,4) * u(4,j,1,e) &
435 + dx(i,5) * u(5,j,1,e) &
436 + dx(i,6) * u(6,j,1,e) &
437 + dx(i,7) * u(7,j,1,e) &
438 + dx(i,8) * u(8,j,1,e) &
439 + dx(i,9) * u(9,j,1,e) &
440 + dx(i,10) * u(10,j,1,e) &
441 + dx(i,11) * u(11,j,1,e)
442 end do
443 end do
444
445 do i = 1, lx * lx * lx
446 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
447 end do
448
449 do k = 1, lx
450 do j = 1, lx
451 do i = 1, lx
452 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
453 + dy(j,2) * u(i,2,k,e) &
454 + dy(j,3) * u(i,3,k,e) &
455 + dy(j,4) * u(i,4,k,e) &
456 + dy(j,5) * u(i,5,k,e) &
457 + dy(j,6) * u(i,6,k,e) &
458 + dy(j,7) * u(i,7,k,e) &
459 + dy(j,8) * u(i,8,k,e) &
460 + dy(j,9) * u(i,9,k,e) &
461 + dy(j,10) * u(i,10,k,e) &
462 + dy(j,11) * u(i,11,k,e)
463 end do
464 end do
465 end do
466
467 do i = 1, lx * lx * lx
468 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
469 end do
470
471 do k = 1, lx
472 do i = 1, lx*lx
473 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
474 + dz(k,2) * u(i,1,2,e) &
475 + dz(k,3) * u(i,1,3,e) &
476 + dz(k,4) * u(i,1,4,e) &
477 + dz(k,5) * u(i,1,5,e) &
478 + dz(k,6) * u(i,1,6,e) &
479 + dz(k,7) * u(i,1,7,e) &
480 + dz(k,8) * u(i,1,8,e) &
481 + dz(k,9) * u(i,1,9,e) &
482 + dz(k,10) * u(i,1,10,e) &
483 + dz(k,11) * u(i,1,11,e)
484 end do
485 end do
486
487 do i = 1, lx * lx * lx
488 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
489 end do
490
491 do i = 1, lx * lx * lx
492 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
493 end do
494
495 end do
496
497 end subroutine cpu_dudxyz_lx11
498
499 subroutine cpu_dudxyz_lx10(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
500 integer, parameter :: lx = 10
501 integer, intent(in) :: nel
502 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
503 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
504 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
505 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
506 real(kind=rp), dimension(lx, lx, lx) :: drst
507 integer :: e, i, j, k
508
509 do e = 1, nel
510 do j = 1, lx * lx
511 do i = 1, lx
512 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
513 + dx(i,2) * u(2,j,1,e) &
514 + dx(i,3) * u(3,j,1,e) &
515 + dx(i,4) * u(4,j,1,e) &
516 + dx(i,5) * u(5,j,1,e) &
517 + dx(i,6) * u(6,j,1,e) &
518 + dx(i,7) * u(7,j,1,e) &
519 + dx(i,8) * u(8,j,1,e) &
520 + dx(i,9) * u(9,j,1,e) &
521 + dx(i,10) * u(10,j,1,e)
522 end do
523 end do
524
525 do i = 1, lx * lx * lx
526 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
527 end do
528
529 do k = 1, lx
530 do j = 1, lx
531 do i = 1, lx
532 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
533 + dy(j,2) * u(i,2,k,e) &
534 + dy(j,3) * u(i,3,k,e) &
535 + dy(j,4) * u(i,4,k,e) &
536 + dy(j,5) * u(i,5,k,e) &
537 + dy(j,6) * u(i,6,k,e) &
538 + dy(j,7) * u(i,7,k,e) &
539 + dy(j,8) * u(i,8,k,e) &
540 + dy(j,9) * u(i,9,k,e) &
541 + dy(j,10) * u(i,10,k,e)
542 end do
543 end do
544 end do
545
546 do i = 1, lx * lx * lx
547 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
548 end do
549
550 do k = 1, lx
551 do i = 1, lx*lx
552 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
553 + dz(k,2) * u(i,1,2,e) &
554 + dz(k,3) * u(i,1,3,e) &
555 + dz(k,4) * u(i,1,4,e) &
556 + dz(k,5) * u(i,1,5,e) &
557 + dz(k,6) * u(i,1,6,e) &
558 + dz(k,7) * u(i,1,7,e) &
559 + dz(k,8) * u(i,1,8,e) &
560 + dz(k,9) * u(i,1,9,e) &
561 + dz(k,10) * u(i,1,10,e)
562 end do
563 end do
564
565 do i = 1, lx * lx * lx
566 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
567 end do
568
569 do i = 1, lx * lx * lx
570 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
571 end do
572
573 end do
574
575 end subroutine cpu_dudxyz_lx10
576
577 subroutine cpu_dudxyz_lx9(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
578 integer, parameter :: lx = 9
579 integer, intent(in) :: nel
580 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
581 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
582 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
583 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
584 real(kind=rp), dimension(lx, lx, lx) :: drst
585 integer :: e, i, j, k
586
587 do e = 1, nel
588 do j = 1, lx * lx
589 do i = 1, lx
590 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
591 + dx(i,2) * u(2,j,1,e) &
592 + dx(i,3) * u(3,j,1,e) &
593 + dx(i,4) * u(4,j,1,e) &
594 + dx(i,5) * u(5,j,1,e) &
595 + dx(i,6) * u(6,j,1,e) &
596 + dx(i,7) * u(7,j,1,e) &
597 + dx(i,8) * u(8,j,1,e) &
598 + dx(i,9) * u(9,j,1,e)
599 end do
600 end do
601
602 do i = 1, lx * lx * lx
603 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
604 end do
605
606 do k = 1, lx
607 do j = 1, lx
608 do i = 1, lx
609 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
610 + dy(j,2) * u(i,2,k,e) &
611 + dy(j,3) * u(i,3,k,e) &
612 + dy(j,4) * u(i,4,k,e) &
613 + dy(j,5) * u(i,5,k,e) &
614 + dy(j,6) * u(i,6,k,e) &
615 + dy(j,7) * u(i,7,k,e) &
616 + dy(j,8) * u(i,8,k,e) &
617 + dy(j,9) * u(i,9,k,e)
618 end do
619 end do
620 end do
621
622 do i = 1, lx * lx * lx
623 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
624 end do
625
626 do k = 1, lx
627 do i = 1, lx*lx
628 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
629 + dz(k,2) * u(i,1,2,e) &
630 + dz(k,3) * u(i,1,3,e) &
631 + dz(k,4) * u(i,1,4,e) &
632 + dz(k,5) * u(i,1,5,e) &
633 + dz(k,6) * u(i,1,6,e) &
634 + dz(k,7) * u(i,1,7,e) &
635 + dz(k,8) * u(i,1,8,e) &
636 + dz(k,9) * u(i,1,9,e)
637 end do
638 end do
639
640 do i = 1, lx * lx * lx
641 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
642 end do
643
644 do i = 1, lx * lx * lx
645 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
646 end do
647
648 end do
649
650 end subroutine cpu_dudxyz_lx9
651
652 subroutine cpu_dudxyz_lx8(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
653 integer, parameter :: lx = 8
654 integer, intent(in) :: nel
655 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
656 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
657 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
658 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
659 real(kind=rp), dimension(lx, lx, lx) :: drst
660 integer :: e, i, j, k
661
662 do e = 1, nel
663 do j = 1, lx * lx
664 do i = 1, lx
665 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
666 + dx(i,2) * u(2,j,1,e) &
667 + dx(i,3) * u(3,j,1,e) &
668 + dx(i,4) * u(4,j,1,e) &
669 + dx(i,5) * u(5,j,1,e) &
670 + dx(i,6) * u(6,j,1,e) &
671 + dx(i,7) * u(7,j,1,e) &
672 + dx(i,8) * u(8,j,1,e)
673 end do
674 end do
675
676 do i = 1, lx * lx * lx
677 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
678 end do
679
680 do k = 1, lx
681 do j = 1, lx
682 do i = 1, lx
683 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
684 + dy(j,2) * u(i,2,k,e) &
685 + dy(j,3) * u(i,3,k,e) &
686 + dy(j,4) * u(i,4,k,e) &
687 + dy(j,5) * u(i,5,k,e) &
688 + dy(j,6) * u(i,6,k,e) &
689 + dy(j,7) * u(i,7,k,e) &
690 + dy(j,8) * u(i,8,k,e)
691 end do
692 end do
693 end do
694
695 do i = 1, lx * lx * lx
696 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
697 end do
698
699 do k = 1, lx
700 do i = 1, lx*lx
701 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
702 + dz(k,2) * u(i,1,2,e) &
703 + dz(k,3) * u(i,1,3,e) &
704 + dz(k,4) * u(i,1,4,e) &
705 + dz(k,5) * u(i,1,5,e) &
706 + dz(k,6) * u(i,1,6,e) &
707 + dz(k,7) * u(i,1,7,e) &
708 + dz(k,8) * u(i,1,8,e)
709 end do
710 end do
711
712 do i = 1, lx * lx * lx
713 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
714 end do
715
716 do i = 1, lx * lx * lx
717 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
718 end do
719
720 end do
721
722 end subroutine cpu_dudxyz_lx8
723
724 subroutine cpu_dudxyz_lx7(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
725 integer, parameter :: lx = 7
726 integer, intent(in) :: nel
727 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
728 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
729 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
730 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
731 real(kind=rp), dimension(lx, lx, lx) :: drst
732 integer :: e, i, j, k
733
734 do e = 1, nel
735 do j = 1, lx * lx
736 do i = 1, lx
737 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
738 + dx(i,2) * u(2,j,1,e) &
739 + dx(i,3) * u(3,j,1,e) &
740 + dx(i,4) * u(4,j,1,e) &
741 + dx(i,5) * u(5,j,1,e) &
742 + dx(i,6) * u(6,j,1,e) &
743 + dx(i,7) * u(7,j,1,e)
744 end do
745 end do
746
747 do i = 1, lx * lx * lx
748 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
749 end do
750
751 do k = 1, lx
752 do j = 1, lx
753 do i = 1, lx
754 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
755 + dy(j,2) * u(i,2,k,e) &
756 + dy(j,3) * u(i,3,k,e) &
757 + dy(j,4) * u(i,4,k,e) &
758 + dy(j,5) * u(i,5,k,e) &
759 + dy(j,6) * u(i,6,k,e) &
760 + dy(j,7) * u(i,7,k,e)
761 end do
762 end do
763 end do
764
765 do i = 1, lx * lx * lx
766 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
767 end do
768
769 do k = 1, lx
770 do i = 1, lx*lx
771 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
772 + dz(k,2) * u(i,1,2,e) &
773 + dz(k,3) * u(i,1,3,e) &
774 + dz(k,4) * u(i,1,4,e) &
775 + dz(k,5) * u(i,1,5,e) &
776 + dz(k,6) * u(i,1,6,e) &
777 + dz(k,7) * u(i,1,7,e)
778 end do
779 end do
780
781 do i = 1, lx * lx * lx
782 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
783 end do
784
785 do i = 1, lx * lx * lx
786 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
787 end do
788
789 end do
790
791 end subroutine cpu_dudxyz_lx7
792
793 subroutine cpu_dudxyz_lx6(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
794 integer, parameter :: lx = 6
795 integer, intent(in) :: nel
796 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
797 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
798 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
799 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
800 real(kind=rp), dimension(lx, lx, lx) :: drst
801 integer :: e, i, j, k
802
803 do e = 1, nel
804 do j = 1, lx * lx
805 do i = 1, lx
806 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
807 + dx(i,2) * u(2,j,1,e) &
808 + dx(i,3) * u(3,j,1,e) &
809 + dx(i,4) * u(4,j,1,e) &
810 + dx(i,5) * u(5,j,1,e) &
811 + dx(i,6) * u(6,j,1,e)
812 end do
813 end do
814
815 do i = 1, lx * lx * lx
816 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
817 end do
818
819 do k = 1, lx
820 do j = 1, lx
821 do i = 1, lx
822 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
823 + dy(j,2) * u(i,2,k,e) &
824 + dy(j,3) * u(i,3,k,e) &
825 + dy(j,4) * u(i,4,k,e) &
826 + dy(j,5) * u(i,5,k,e) &
827 + dy(j,6) * u(i,6,k,e)
828 end do
829 end do
830 end do
831
832 do i = 1, lx * lx * lx
833 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
834 end do
835
836 do k = 1, lx
837 do i = 1, lx*lx
838 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
839 + dz(k,2) * u(i,1,2,e) &
840 + dz(k,3) * u(i,1,3,e) &
841 + dz(k,4) * u(i,1,4,e) &
842 + dz(k,5) * u(i,1,5,e) &
843 + dz(k,6) * u(i,1,6,e)
844 end do
845 end do
846
847 do i = 1, lx * lx * lx
848 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
849 end do
850
851 do i = 1, lx * lx * lx
852 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
853 end do
854
855 end do
856
857 end subroutine cpu_dudxyz_lx6
858
859 subroutine cpu_dudxyz_lx5(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
860 integer, parameter :: lx = 5
861 integer, intent(in) :: nel
862 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
863 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
864 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
865 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
866 real(kind=rp), dimension(lx, lx, lx) :: drst
867 integer :: e, i, j, k
868
869 do e = 1, nel
870 do j = 1, lx * lx
871 do i = 1, lx
872 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
873 + dx(i,2) * u(2,j,1,e) &
874 + dx(i,3) * u(3,j,1,e) &
875 + dx(i,4) * u(4,j,1,e) &
876 + dx(i,5) * u(5,j,1,e)
877 end do
878 end do
879
880 do i = 1, lx * lx * lx
881 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
882 end do
883
884 do k = 1, lx
885 do j = 1, lx
886 do i = 1, lx
887 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
888 + dy(j,2) * u(i,2,k,e) &
889 + dy(j,3) * u(i,3,k,e) &
890 + dy(j,4) * u(i,4,k,e) &
891 + dy(j,5) * u(i,5,k,e)
892 end do
893 end do
894 end do
895
896 do i = 1, lx * lx * lx
897 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
898 end do
899
900 do k = 1, lx
901 do i = 1, lx*lx
902 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
903 + dz(k,2) * u(i,1,2,e) &
904 + dz(k,3) * u(i,1,3,e) &
905 + dz(k,4) * u(i,1,4,e) &
906 + dz(k,5) * u(i,1,5,e)
907 end do
908 end do
909
910 do i = 1, lx * lx * lx
911 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
912 end do
913
914 do i = 1, lx * lx * lx
915 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
916 end do
917
918 end do
919
920 end subroutine cpu_dudxyz_lx5
921
922 subroutine cpu_dudxyz_lx4(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
923 integer, parameter :: lx = 4
924 integer, intent(in) :: nel
925 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
926 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
927 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
928 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
929 real(kind=rp), dimension(lx, lx, lx) :: drst
930 integer :: e, i, j, k
931
932 do e = 1, nel
933 do j = 1, lx * lx
934 do i = 1, lx
935 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
936 + dx(i,2) * u(2,j,1,e) &
937 + dx(i,3) * u(3,j,1,e) &
938 + dx(i,4) * u(4,j,1,e)
939 end do
940 end do
941
942 do i = 1, lx * lx * lx
943 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
944 end do
945
946 do k = 1, lx
947 do j = 1, lx
948 do i = 1, lx
949 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
950 + dy(j,2) * u(i,2,k,e) &
951 + dy(j,3) * u(i,3,k,e) &
952 + dy(j,4) * u(i,4,k,e)
953 end do
954 end do
955 end do
956
957 do i = 1, lx * lx * lx
958 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
959 end do
960
961 do k = 1, lx
962 do i = 1, lx*lx
963 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
964 + dz(k,2) * u(i,1,2,e) &
965 + dz(k,3) * u(i,1,3,e) &
966 + dz(k,4) * u(i,1,4,e)
967 end do
968 end do
969
970 do i = 1, lx * lx * lx
971 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
972 end do
973
974 do i = 1, lx * lx * lx
975 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
976 end do
977
978 end do
979
980 end subroutine cpu_dudxyz_lx4
981
982 subroutine cpu_dudxyz_lx3(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
983 integer, parameter :: lx = 3
984 integer, intent(in) :: nel
985 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
986 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
987 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
988 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
989 real(kind=rp), dimension(lx, lx, lx) :: drst
990 integer :: e, i, j, k
991
992 do e = 1, nel
993 do j = 1, lx * lx
994 do i = 1, lx
995 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
996 + dx(i,2) * u(2,j,1,e) &
997 + dx(i,3) * u(3,j,1,e)
998 end do
999 end do
1000
1001 do i = 1, lx * lx * lx
1002 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
1003 end do
1004
1005 do k = 1, lx
1006 do j = 1, lx
1007 do i = 1, lx
1008 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
1009 + dy(j,2) * u(i,2,k,e) &
1010 + dy(j,3) * u(i,3,k,e)
1011 end do
1012 end do
1013 end do
1014
1015 do i = 1, lx * lx * lx
1016 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
1017 end do
1018
1019 do k = 1, lx
1020 do i = 1, lx*lx
1021 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
1022 + dz(k,2) * u(i,1,2,e) &
1023 + dz(k,3) * u(i,1,3,e)
1024 end do
1025 end do
1026
1027 do i = 1, lx * lx * lx
1028 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
1029 end do
1030
1031 do i = 1, lx * lx * lx
1032 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
1033 end do
1034
1035 end do
1036
1037 end subroutine cpu_dudxyz_lx3
1038
1039 subroutine cpu_dudxyz_lx2(du, u, dr, ds, dt, dx, dy, dz, jacinv, nel)
1040 integer, parameter :: lx = 2
1041 integer, intent(in) :: nel
1042 real(kind=rp), dimension(lx, lx, lx, nel), intent(inout) :: du
1043 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: u, dr, ds, dt
1044 real(kind=rp), dimension(lx, lx, lx, nel), intent(in) :: jacinv
1045 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1046 real(kind=rp), dimension(lx, lx, lx) :: drst
1047 integer :: e, i, j, k
1048
1049 do e = 1, nel
1050 do j = 1, lx * lx
1051 do i = 1, lx
1052 du(i,j,1,e) = dx(i,1) * u(1,j,1,e) &
1053 + dx(i,2) * u(2,j,1,e)
1054 end do
1055 end do
1056
1057 do i = 1, lx * lx * lx
1058 du(i,1,1,e) = du(i,1,1,e) * dr(i,1,1,e)
1059 end do
1060
1061 do k = 1, lx
1062 do j = 1, lx
1063 do i = 1, lx
1064 drst(i,j,k) = dy(j,1) * u(i,1,k,e) &
1065 + dy(j,2) * u(i,2,k,e)
1066 end do
1067 end do
1068 end do
1069
1070 do i = 1, lx * lx * lx
1071 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * ds(i,1,1,e)
1072 end do
1073
1074 do k = 1, lx
1075 do i = 1, lx*lx
1076 drst(i,1,k) = dz(k,1) * u(i,1,1,e) &
1077 + dz(k,2) * u(i,1,2,e)
1078 end do
1079 end do
1080
1081 do i = 1, lx * lx * lx
1082 du(i,1,1,e) = du(i,1,1,e) + drst(i,1,1) * dt(i,1,1,e)
1083 end do
1084
1085 do i = 1, lx * lx * lx
1086 du(i,1,1,e) = du(i,1,1,e) * jacinv(i,1,1,e)
1087 end do
1088
1089 end do
1090
1091 end subroutine cpu_dudxyz_lx2
1092
1093end submodule cpu_dudxyz
1094
Operators CPU backend.
Definition opr_cpu.f90:34