Neko 1.99.2
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
cpu_opgrad.f90
Go to the documentation of this file.
1! Copyright (c) 2021-2024, 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_opgrad
35 implicit none
36
37contains
38
39 module subroutine opr_cpu_opgrad(ux, uy, uz, u, coef, e_start, e_end)
40 type(coef_t), intent(in) :: coef
41 integer, intent(in) :: e_start, e_end
42 real(kind=rp), dimension(coef%Xh%lxyz, e_end - e_start + 1), &
43 intent(inout) :: ux
44 real(kind=rp), dimension(coef%Xh%lxyz, e_end - e_start + 1), &
45 intent(inout) :: uy
46 real(kind=rp), dimension(coef%Xh%lxyz, e_end - e_start + 1), &
47 intent(inout) :: uz
48 real(kind=rp), dimension(coef%Xh%lxyz, e_end - e_start + 1), &
49 intent(in) :: u
50 integer :: e_len
51
52 e_len = e_end-e_start+1
53
54 if (e_len .eq. 1) then
55 call opr_cpu_opgrad_single(ux, uy, uz, u, coef, e_start)
56 else
57 call opr_cpu_opgrad_many(ux, uy, uz, u, coef, e_start, e_len)
58 end if
59 end subroutine opr_cpu_opgrad
60
61 subroutine opr_cpu_opgrad_many(ux, uy, uz, u, coef, e_start, e_len)
62 type(coef_t), intent(in) :: coef
63 integer, intent(in) :: e_start, e_len
64 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: ux
65 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: uy
66 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: uz
67 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(in) :: u
68
69 associate(xh => coef%Xh, msh => coef%msh, &
70 drdx => coef%drdx, drdy => coef%drdy, drdz => coef%drdz, &
71 dsdx => coef%dsdx, dsdy => coef%dsdy, dsdz => coef%dsdz, &
72 dtdx => coef%dtdx, dtdy => coef%dtdy, dtdz => coef%dtdz)
73
74 select case (xh%lx)
75 case (18)
76 call cpu_opgrad_lx18(ux, uy, uz, u, &
77 xh%dx, xh%dy, xh%dz, &
78 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
79 dtdx(1, 1, 1, e_start), &
80 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
81 dtdy(1, 1, 1, e_start), &
82 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
83 dtdz(1, 1, 1, e_start), &
84 xh%w3, e_len)
85 case (17)
86 call cpu_opgrad_lx17(ux, uy, uz, u, &
87 xh%dx, xh%dy, xh%dz, &
88 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
89 dtdx(1, 1, 1, e_start), &
90 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
91 dtdy(1, 1, 1, e_start), &
92 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
93 dtdz(1, 1, 1, e_start), &
94 xh%w3, e_len)
95 case (16)
96 call cpu_opgrad_lx16(ux, uy, uz, u, &
97 xh%dx, xh%dy, xh%dz, &
98 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
99 dtdx(1, 1, 1, e_start), &
100 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
101 dtdy(1, 1, 1, e_start), &
102 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
103 dtdz(1, 1, 1, e_start), &
104 xh%w3, e_len)
105 case (15)
106 call cpu_opgrad_lx15(ux, uy, uz, u, &
107 xh%dx, xh%dy, xh%dz, &
108 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
109 dtdx(1, 1, 1, e_start), &
110 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
111 dtdy(1, 1, 1, e_start), &
112 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
113 dtdz(1, 1, 1, e_start), &
114 xh%w3, e_len)
115 case (14)
116 call cpu_opgrad_lx14(ux, uy, uz, u, &
117 xh%dx, xh%dy, xh%dz, &
118 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
119 dtdx(1, 1, 1, e_start), &
120 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
121 dtdy(1, 1, 1, e_start), &
122 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
123 dtdz(1, 1, 1, e_start), &
124 xh%w3, e_len)
125 case (13)
126 call cpu_opgrad_lx13(ux, uy, uz, u, &
127 xh%dx, xh%dy, xh%dz, &
128 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
129 dtdx(1, 1, 1, e_start), &
130 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
131 dtdy(1, 1, 1, e_start), &
132 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
133 dtdz(1, 1, 1, e_start), &
134 xh%w3, e_len)
135 case (12)
136 call cpu_opgrad_lx12(ux, uy, uz, u, &
137 xh%dx, xh%dy, xh%dz, &
138 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
139 dtdx(1, 1, 1, e_start), &
140 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
141 dtdy(1, 1, 1, e_start), &
142 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
143 dtdz(1, 1, 1, e_start), &
144 xh%w3, e_len)
145 case (11)
146 call cpu_opgrad_lx11(ux, uy, uz, u, &
147 xh%dx, xh%dy, xh%dz, &
148 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
149 dtdx(1, 1, 1, e_start), &
150 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
151 dtdy(1, 1, 1, e_start), &
152 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
153 dtdz(1, 1, 1, e_start), &
154 xh%w3, e_len)
155 case (10)
156 call cpu_opgrad_lx10(ux, uy, uz, u, &
157 xh%dx, xh%dy, xh%dz, &
158 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
159 dtdx(1, 1, 1, e_start), &
160 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
161 dtdy(1, 1, 1, e_start), &
162 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
163 dtdz(1, 1, 1, e_start), &
164 xh%w3, e_len)
165
166 case (9)
167 call cpu_opgrad_lx9(ux, uy, uz, u, &
168 xh%dx, xh%dy, xh%dz, &
169 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
170 dtdx(1, 1, 1, e_start), &
171 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
172 dtdy(1, 1, 1, e_start), &
173 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
174 dtdz(1, 1, 1, e_start), &
175 xh%w3, e_len)
176 case (8)
177 call cpu_opgrad_lx8(ux, uy, uz, u, &
178 xh%dx, xh%dy, xh%dz, &
179 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
180 dtdx(1, 1, 1, e_start), &
181 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
182 dtdy(1, 1, 1, e_start), &
183 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
184 dtdz(1, 1, 1, e_start), &
185 xh%w3, e_len)
186 case (7)
187 call cpu_opgrad_lx7(ux, uy, uz, u, &
188 xh%dx, xh%dy, xh%dz, &
189 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
190 dtdx(1, 1, 1, e_start), &
191 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
192 dtdy(1, 1, 1, e_start), &
193 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
194 dtdz(1, 1, 1, e_start), &
195 xh%w3, e_len)
196 case (6)
197 call cpu_opgrad_lx6(ux, uy, uz, u, &
198 xh%dx, xh%dy, xh%dz, &
199 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
200 dtdx(1, 1, 1, e_start), &
201 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
202 dtdy(1, 1, 1, e_start), &
203 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
204 dtdz(1, 1, 1, e_start), &
205 xh%w3, e_len)
206 case (5)
207 call cpu_opgrad_lx5(ux, uy, uz, u, &
208 xh%dx, xh%dy, xh%dz, &
209 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
210 dtdx(1, 1, 1, e_start), &
211 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
212 dtdy(1, 1, 1, e_start), &
213 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
214 dtdz(1, 1, 1, e_start), &
215 xh%w3, e_len)
216 case (4)
217 call cpu_opgrad_lx4(ux, uy, uz, u, &
218 xh%dx, xh%dy, xh%dz, &
219 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
220 dtdx(1, 1, 1, e_start), &
221 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
222 dtdy(1, 1, 1, e_start), &
223 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
224 dtdz(1, 1, 1, e_start), &
225 xh%w3, e_len)
226 case (3)
227 call cpu_opgrad_lx3(ux, uy, uz, u, &
228 xh%dx, xh%dy, xh%dz, &
229 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
230 dtdx(1, 1, 1, e_start), &
231 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
232 dtdy(1, 1, 1, e_start), &
233 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
234 dtdz(1, 1, 1, e_start), &
235 xh%w3, e_len)
236 case (2)
237 call cpu_opgrad_lx2(ux, uy, uz, u, &
238 xh%dx, xh%dy, xh%dz, &
239 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
240 dtdx(1, 1, 1, e_start), &
241 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
242 dtdy(1, 1, 1, e_start), &
243 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
244 dtdz(1, 1, 1, e_start), &
245 xh%w3, e_len)
246 case default
247 call cpu_opgrad_lx(ux, uy, uz, u, &
248 xh%dx, xh%dy, xh%dz, &
249 drdx(1, 1, 1, e_start), dsdx(1, 1, 1, e_start), &
250 dtdx(1, 1, 1, e_start), &
251 drdy(1, 1, 1, e_start), dsdy(1, 1, 1, e_start), &
252 dtdy(1, 1, 1, e_start), &
253 drdz(1, 1, 1, e_start), dsdz(1, 1, 1, e_start), &
254 dtdz(1, 1, 1, e_start), &
255 xh%w3, e_len, xh%lx)
256 end select
257 end associate
258
259 end subroutine opr_cpu_opgrad_many
260
261 subroutine cpu_opgrad_lx(ux, uy, uz, u, dx, dy, dz, &
262 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n, lx)
263 integer, intent(in) :: n, lx
264 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
265 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
266 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
267 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
268 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
269 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
270 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
271 real(kind=rp) :: ur(lx, lx, lx)
272 real(kind=rp) :: us(lx, lx, lx)
273 real(kind=rp) :: ut(lx, lx, lx)
274 real(kind=rp) :: tmp
275 integer :: e, i, j, k, l
276
277 do e = 1, n
278 do j = 1, lx * lx
279 do i = 1, lx
280 tmp = 0.0_rp
281 do k = 1, lx
282 tmp = tmp + dx(i,k) * u(k,j,1,e)
283 end do
284 ur(i,j,1) = tmp
285 end do
286 end do
287
288 do k = 1, lx
289 do j = 1, lx
290 do i = 1, lx
291 tmp = 0.0_rp
292 do l = 1, lx
293 tmp = tmp + dy(j,l) * u(i,l,k,e)
294 end do
295 us(i,j,k) = tmp
296 end do
297 end do
298 end do
299
300 do k = 1, lx
301 do i = 1, lx*lx
302 tmp = 0.0_rp
303 do l = 1, lx
304 tmp = tmp + dz(k,l) * u(i,1,l,e)
305 end do
306 ut(i,1,k) = tmp
307 end do
308 end do
309
310 do i = 1, lx * lx * lx
311 ux(i,1,1,e) = w3(i,1,1) &
312 * ( drdx(i,1,1,e) * ur(i,1,1) &
313 + dsdx(i,1,1,e) * us(i,1,1) &
314 + dtdx(i,1,1,e) * ut(i,1,1) )
315 uy(i,1,1,e) = w3(i,1,1) &
316 * ( dsdy(i,1,1,e) * us(i,1,1) &
317 + drdy(i,1,1,e) * ur(i,1,1) &
318 + dtdy(i,1,1,e) * ut(i,1,1) )
319 uz(i,1,1,e) = w3(i,1,1) &
320 * ( dtdz(i,1,1,e) * ut(i,1,1) &
321 + drdz(i,1,1,e) * ur(i,1,1) &
322 + dsdz(i,1,1,e) * us(i,1,1) )
323 end do
324 end do
325 end subroutine cpu_opgrad_lx
326
327 subroutine cpu_opgrad_lx18(ux, uy, uz, u, dx, dy, dz, &
328 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
329 integer, parameter :: lx = 18
330 integer, intent(in) :: n
331 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
332 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
333 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
334 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
335 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
336 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
337 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
338 real(kind=rp) :: ur(lx, lx, lx)
339 real(kind=rp) :: us(lx, lx, lx)
340 real(kind=rp) :: ut(lx, lx, lx)
341 integer :: e, i, j, k
342
343 do e = 1, n
344 do j = 1, lx * lx
345 do i = 1, lx
346 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
347 + dx(i,2) * u(2,j,1,e) &
348 + dx(i,3) * u(3,j,1,e) &
349 + dx(i,4) * u(4,j,1,e) &
350 + dx(i,5) * u(5,j,1,e) &
351 + dx(i,6) * u(6,j,1,e) &
352 + dx(i,7) * u(7,j,1,e) &
353 + dx(i,8) * u(8,j,1,e) &
354 + dx(i,9) * u(9,j,1,e) &
355 + dx(i,10) * u(10,j,1,e) &
356 + dx(i,11) * u(11,j,1,e) &
357 + dx(i,12) * u(12,j,1,e) &
358 + dx(i,13) * u(13,j,1,e) &
359 + dx(i,14) * u(14,j,1,e) &
360 + dx(i,15) * u(15,j,1,e) &
361 + dx(i,16) * u(16,j,1,e) &
362 + dx(i,17) * u(17,j,1,e) &
363 + dx(i,18) * u(18,j,1,e)
364 end do
365 end do
366
367 do k = 1, lx
368 do j = 1, lx
369 do i = 1, lx
370 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
371 + dy(j,2) * u(i,2,k,e) &
372 + dy(j,3) * u(i,3,k,e) &
373 + dy(j,4) * u(i,4,k,e) &
374 + dy(j,5) * u(i,5,k,e) &
375 + dy(j,6) * u(i,6,k,e) &
376 + dy(j,7) * u(i,7,k,e) &
377 + dy(j,8) * u(i,8,k,e) &
378 + dy(j,9) * u(i,9,k,e) &
379 + dy(j,10) * u(i,10,k,e) &
380 + dy(j,11) * u(i,11,k,e) &
381 + dy(j,12) * u(i,12,k,e) &
382 + dy(j,13) * u(i,13,k,e) &
383 + dy(j,14) * u(i,14,k,e) &
384 + dy(j,15) * u(i,15,k,e) &
385 + dy(j,16) * u(i,16,k,e) &
386 + dy(j,17) * u(i,17,k,e) &
387 + dy(j,18) * u(i,18,k,e)
388 end do
389 end do
390 end do
391
392 do k = 1, lx
393 do i = 1, lx*lx
394 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
395 + dz(k,2) * u(i,1,2,e) &
396 + dz(k,3) * u(i,1,3,e) &
397 + dz(k,4) * u(i,1,4,e) &
398 + dz(k,5) * u(i,1,5,e) &
399 + dz(k,6) * u(i,1,6,e) &
400 + dz(k,7) * u(i,1,7,e) &
401 + dz(k,8) * u(i,1,8,e) &
402 + dz(k,9) * u(i,1,9,e) &
403 + dz(k,10) * u(i,1,10,e) &
404 + dz(k,11) * u(i,1,11,e) &
405 + dz(k,12) * u(i,1,12,e) &
406 + dz(k,13) * u(i,1,13,e) &
407 + dz(k,14) * u(i,1,14,e) &
408 + dz(k,15) * u(i,1,15,e) &
409 + dz(k,16) * u(i,1,16,e) &
410 + dz(k,17) * u(i,1,17,e) &
411 + dz(k,18) * u(i,1,18,e)
412 end do
413 end do
414
415 do i = 1, lx * lx * lx
416 ux(i,1,1,e) = w3(i,1,1) &
417 * ( drdx(i,1,1,e) * ur(i,1,1) &
418 + dsdx(i,1,1,e) * us(i,1,1) &
419 + dtdx(i,1,1,e) * ut(i,1,1) )
420 uy(i,1,1,e) = w3(i,1,1) &
421 * ( dsdy(i,1,1,e) * us(i,1,1) &
422 + drdy(i,1,1,e) * ur(i,1,1) &
423 + dtdy(i,1,1,e) * ut(i,1,1) )
424 uz(i,1,1,e) = w3(i,1,1) &
425 * ( dtdz(i,1,1,e) * ut(i,1,1) &
426 + drdz(i,1,1,e) * ur(i,1,1) &
427 + dsdz(i,1,1,e) * us(i,1,1) )
428 end do
429 end do
430 end subroutine cpu_opgrad_lx18
431
432 subroutine cpu_opgrad_lx17(ux, uy, uz, u, dx, dy, dz, &
433 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
434 integer, parameter :: lx = 17
435 integer, intent(in) :: n
436 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
437 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
438 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
439 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
440 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
441 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
442 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
443 real(kind=rp) :: ur(lx, lx, lx)
444 real(kind=rp) :: us(lx, lx, lx)
445 real(kind=rp) :: ut(lx, lx, lx)
446 integer :: e, i, j, k
447
448 do e = 1, n
449 do j = 1, lx * lx
450 do i = 1, lx
451 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
452 + dx(i,2) * u(2,j,1,e) &
453 + dx(i,3) * u(3,j,1,e) &
454 + dx(i,4) * u(4,j,1,e) &
455 + dx(i,5) * u(5,j,1,e) &
456 + dx(i,6) * u(6,j,1,e) &
457 + dx(i,7) * u(7,j,1,e) &
458 + dx(i,8) * u(8,j,1,e) &
459 + dx(i,9) * u(9,j,1,e) &
460 + dx(i,10) * u(10,j,1,e) &
461 + dx(i,11) * u(11,j,1,e) &
462 + dx(i,12) * u(12,j,1,e) &
463 + dx(i,13) * u(13,j,1,e) &
464 + dx(i,14) * u(14,j,1,e) &
465 + dx(i,15) * u(15,j,1,e) &
466 + dx(i,16) * u(16,j,1,e) &
467 + dx(i,17) * u(17,j,1,e)
468 end do
469 end do
470
471 do k = 1, lx
472 do j = 1, lx
473 do i = 1, lx
474 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
475 + dy(j,2) * u(i,2,k,e) &
476 + dy(j,3) * u(i,3,k,e) &
477 + dy(j,4) * u(i,4,k,e) &
478 + dy(j,5) * u(i,5,k,e) &
479 + dy(j,6) * u(i,6,k,e) &
480 + dy(j,7) * u(i,7,k,e) &
481 + dy(j,8) * u(i,8,k,e) &
482 + dy(j,9) * u(i,9,k,e) &
483 + dy(j,10) * u(i,10,k,e) &
484 + dy(j,11) * u(i,11,k,e) &
485 + dy(j,12) * u(i,12,k,e) &
486 + dy(j,13) * u(i,13,k,e) &
487 + dy(j,14) * u(i,14,k,e) &
488 + dy(j,15) * u(i,15,k,e) &
489 + dy(j,16) * u(i,16,k,e) &
490 + dy(j,17) * u(i,17,k,e)
491 end do
492 end do
493 end do
494
495 do k = 1, lx
496 do i = 1, lx*lx
497 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
498 + dz(k,2) * u(i,1,2,e) &
499 + dz(k,3) * u(i,1,3,e) &
500 + dz(k,4) * u(i,1,4,e) &
501 + dz(k,5) * u(i,1,5,e) &
502 + dz(k,6) * u(i,1,6,e) &
503 + dz(k,7) * u(i,1,7,e) &
504 + dz(k,8) * u(i,1,8,e) &
505 + dz(k,9) * u(i,1,9,e) &
506 + dz(k,10) * u(i,1,10,e) &
507 + dz(k,11) * u(i,1,11,e) &
508 + dz(k,12) * u(i,1,12,e) &
509 + dz(k,13) * u(i,1,13,e) &
510 + dz(k,14) * u(i,1,14,e) &
511 + dz(k,15) * u(i,1,15,e) &
512 + dz(k,16) * u(i,1,16,e) &
513 + dz(k,17) * u(i,1,17,e)
514 end do
515 end do
516
517 do i = 1, lx * lx * lx
518 ux(i,1,1,e) = w3(i,1,1) &
519 * ( drdx(i,1,1,e) * ur(i,1,1) &
520 + dsdx(i,1,1,e) * us(i,1,1) &
521 + dtdx(i,1,1,e) * ut(i,1,1) )
522 uy(i,1,1,e) = w3(i,1,1) &
523 * ( dsdy(i,1,1,e) * us(i,1,1) &
524 + drdy(i,1,1,e) * ur(i,1,1) &
525 + dtdy(i,1,1,e) * ut(i,1,1) )
526 uz(i,1,1,e) = w3(i,1,1) &
527 * ( dtdz(i,1,1,e) * ut(i,1,1) &
528 + drdz(i,1,1,e) * ur(i,1,1) &
529 + dsdz(i,1,1,e) * us(i,1,1) )
530 end do
531 end do
532 end subroutine cpu_opgrad_lx17
533
534 subroutine cpu_opgrad_lx16(ux, uy, uz, u, dx, dy, dz, &
535 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
536 integer, parameter :: lx = 16
537 integer, intent(in) :: n
538 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
539 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
540 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
541 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
542 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
543 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
544 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
545 real(kind=rp) :: ur(lx, lx, lx)
546 real(kind=rp) :: us(lx, lx, lx)
547 real(kind=rp) :: ut(lx, lx, lx)
548 integer :: e, i, j, k
549
550 do e = 1, n
551 do j = 1, lx * lx
552 do i = 1, lx
553 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
554 + dx(i,2) * u(2,j,1,e) &
555 + dx(i,3) * u(3,j,1,e) &
556 + dx(i,4) * u(4,j,1,e) &
557 + dx(i,5) * u(5,j,1,e) &
558 + dx(i,6) * u(6,j,1,e) &
559 + dx(i,7) * u(7,j,1,e) &
560 + dx(i,8) * u(8,j,1,e) &
561 + dx(i,9) * u(9,j,1,e) &
562 + dx(i,10) * u(10,j,1,e) &
563 + dx(i,11) * u(11,j,1,e) &
564 + dx(i,12) * u(12,j,1,e) &
565 + dx(i,13) * u(13,j,1,e) &
566 + dx(i,14) * u(14,j,1,e) &
567 + dx(i,15) * u(15,j,1,e) &
568 + dx(i,16) * u(16,j,1,e)
569 end do
570 end do
571
572 do k = 1, lx
573 do j = 1, lx
574 do i = 1, lx
575 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
576 + dy(j,2) * u(i,2,k,e) &
577 + dy(j,3) * u(i,3,k,e) &
578 + dy(j,4) * u(i,4,k,e) &
579 + dy(j,5) * u(i,5,k,e) &
580 + dy(j,6) * u(i,6,k,e) &
581 + dy(j,7) * u(i,7,k,e) &
582 + dy(j,8) * u(i,8,k,e) &
583 + dy(j,9) * u(i,9,k,e) &
584 + dy(j,10) * u(i,10,k,e) &
585 + dy(j,11) * u(i,11,k,e) &
586 + dy(j,12) * u(i,12,k,e) &
587 + dy(j,13) * u(i,13,k,e) &
588 + dy(j,14) * u(i,14,k,e) &
589 + dy(j,15) * u(i,15,k,e) &
590 + dy(j,16) * u(i,16,k,e)
591 end do
592 end do
593 end do
594
595 do k = 1, lx
596 do i = 1, lx*lx
597 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
598 + dz(k,2) * u(i,1,2,e) &
599 + dz(k,3) * u(i,1,3,e) &
600 + dz(k,4) * u(i,1,4,e) &
601 + dz(k,5) * u(i,1,5,e) &
602 + dz(k,6) * u(i,1,6,e) &
603 + dz(k,7) * u(i,1,7,e) &
604 + dz(k,8) * u(i,1,8,e) &
605 + dz(k,9) * u(i,1,9,e) &
606 + dz(k,10) * u(i,1,10,e) &
607 + dz(k,11) * u(i,1,11,e) &
608 + dz(k,12) * u(i,1,12,e) &
609 + dz(k,13) * u(i,1,13,e) &
610 + dz(k,14) * u(i,1,14,e) &
611 + dz(k,15) * u(i,1,15,e) &
612 + dz(k,16) * u(i,1,16,e)
613 end do
614 end do
615
616 do i = 1, lx * lx * lx
617 ux(i,1,1,e) = w3(i,1,1) &
618 * ( drdx(i,1,1,e) * ur(i,1,1) &
619 + dsdx(i,1,1,e) * us(i,1,1) &
620 + dtdx(i,1,1,e) * ut(i,1,1) )
621 uy(i,1,1,e) = w3(i,1,1) &
622 * ( dsdy(i,1,1,e) * us(i,1,1) &
623 + drdy(i,1,1,e) * ur(i,1,1) &
624 + dtdy(i,1,1,e) * ut(i,1,1) )
625 uz(i,1,1,e) = w3(i,1,1) &
626 * ( dtdz(i,1,1,e) * ut(i,1,1) &
627 + drdz(i,1,1,e) * ur(i,1,1) &
628 + dsdz(i,1,1,e) * us(i,1,1) )
629 end do
630 end do
631 end subroutine cpu_opgrad_lx16
632
633 subroutine cpu_opgrad_lx15(ux, uy, uz, u, dx, dy, dz, &
634 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
635 integer, parameter :: lx = 15
636 integer, intent(in) :: n
637 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
638 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
639 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
640 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
641 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
642 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
643 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
644 real(kind=rp) :: ur(lx, lx, lx)
645 real(kind=rp) :: us(lx, lx, lx)
646 real(kind=rp) :: ut(lx, lx, lx)
647 integer :: e, i, j, k
648
649 do e = 1, n
650 do j = 1, lx * lx
651 do i = 1, lx
652 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
653 + dx(i,2) * u(2,j,1,e) &
654 + dx(i,3) * u(3,j,1,e) &
655 + dx(i,4) * u(4,j,1,e) &
656 + dx(i,5) * u(5,j,1,e) &
657 + dx(i,6) * u(6,j,1,e) &
658 + dx(i,7) * u(7,j,1,e) &
659 + dx(i,8) * u(8,j,1,e) &
660 + dx(i,9) * u(9,j,1,e) &
661 + dx(i,10) * u(10,j,1,e) &
662 + dx(i,11) * u(11,j,1,e) &
663 + dx(i,12) * u(12,j,1,e) &
664 + dx(i,13) * u(13,j,1,e) &
665 + dx(i,14) * u(14,j,1,e) &
666 + dx(i,15) * u(15,j,1,e)
667 end do
668 end do
669
670 do k = 1, lx
671 do j = 1, lx
672 do i = 1, lx
673 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
674 + dy(j,2) * u(i,2,k,e) &
675 + dy(j,3) * u(i,3,k,e) &
676 + dy(j,4) * u(i,4,k,e) &
677 + dy(j,5) * u(i,5,k,e) &
678 + dy(j,6) * u(i,6,k,e) &
679 + dy(j,7) * u(i,7,k,e) &
680 + dy(j,8) * u(i,8,k,e) &
681 + dy(j,9) * u(i,9,k,e) &
682 + dy(j,10) * u(i,10,k,e) &
683 + dy(j,11) * u(i,11,k,e) &
684 + dy(j,12) * u(i,12,k,e) &
685 + dy(j,13) * u(i,13,k,e) &
686 + dy(j,14) * u(i,14,k,e) &
687 + dy(j,15) * u(i,15,k,e)
688 end do
689 end do
690 end do
691
692 do k = 1, lx
693 do i = 1, lx*lx
694 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
695 + dz(k,2) * u(i,1,2,e) &
696 + dz(k,3) * u(i,1,3,e) &
697 + dz(k,4) * u(i,1,4,e) &
698 + dz(k,5) * u(i,1,5,e) &
699 + dz(k,6) * u(i,1,6,e) &
700 + dz(k,7) * u(i,1,7,e) &
701 + dz(k,8) * u(i,1,8,e) &
702 + dz(k,9) * u(i,1,9,e) &
703 + dz(k,10) * u(i,1,10,e) &
704 + dz(k,11) * u(i,1,11,e) &
705 + dz(k,12) * u(i,1,12,e) &
706 + dz(k,13) * u(i,1,13,e) &
707 + dz(k,14) * u(i,1,14,e) &
708 + dz(k,15) * u(i,1,15,e)
709 end do
710 end do
711
712 do i = 1, lx * lx * lx
713 ux(i,1,1,e) = w3(i,1,1) &
714 * ( drdx(i,1,1,e) * ur(i,1,1) &
715 + dsdx(i,1,1,e) * us(i,1,1) &
716 + dtdx(i,1,1,e) * ut(i,1,1) )
717 uy(i,1,1,e) = w3(i,1,1) &
718 * ( dsdy(i,1,1,e) * us(i,1,1) &
719 + drdy(i,1,1,e) * ur(i,1,1) &
720 + dtdy(i,1,1,e) * ut(i,1,1) )
721 uz(i,1,1,e) = w3(i,1,1) &
722 * ( dtdz(i,1,1,e) * ut(i,1,1) &
723 + drdz(i,1,1,e) * ur(i,1,1) &
724 + dsdz(i,1,1,e) * us(i,1,1) )
725 end do
726 end do
727 end subroutine cpu_opgrad_lx15
728
729 subroutine cpu_opgrad_lx14(ux, uy, uz, u, dx, dy, dz, &
730 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
731 integer, parameter :: lx = 14
732 integer, intent(in) :: n
733 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
734 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
735 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
736 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
737 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
738 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
739 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
740 real(kind=rp) :: ur(lx, lx, lx)
741 real(kind=rp) :: us(lx, lx, lx)
742 real(kind=rp) :: ut(lx, lx, lx)
743 integer :: e, i, j, k
744
745 do e = 1, n
746 do j = 1, lx * lx
747 do i = 1, lx
748 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
749 + dx(i,2) * u(2,j,1,e) &
750 + dx(i,3) * u(3,j,1,e) &
751 + dx(i,4) * u(4,j,1,e) &
752 + dx(i,5) * u(5,j,1,e) &
753 + dx(i,6) * u(6,j,1,e) &
754 + dx(i,7) * u(7,j,1,e) &
755 + dx(i,8) * u(8,j,1,e) &
756 + dx(i,9) * u(9,j,1,e) &
757 + dx(i,10) * u(10,j,1,e) &
758 + dx(i,11) * u(11,j,1,e) &
759 + dx(i,12) * u(12,j,1,e) &
760 + dx(i,13) * u(13,j,1,e) &
761 + dx(i,14) * u(14,j,1,e)
762 end do
763 end do
764
765 do k = 1, lx
766 do j = 1, lx
767 do i = 1, lx
768 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
769 + dy(j,2) * u(i,2,k,e) &
770 + dy(j,3) * u(i,3,k,e) &
771 + dy(j,4) * u(i,4,k,e) &
772 + dy(j,5) * u(i,5,k,e) &
773 + dy(j,6) * u(i,6,k,e) &
774 + dy(j,7) * u(i,7,k,e) &
775 + dy(j,8) * u(i,8,k,e) &
776 + dy(j,9) * u(i,9,k,e) &
777 + dy(j,10) * u(i,10,k,e) &
778 + dy(j,11) * u(i,11,k,e) &
779 + dy(j,12) * u(i,12,k,e) &
780 + dy(j,13) * u(i,13,k,e) &
781 + dy(j,14) * u(i,14,k,e)
782 end do
783 end do
784 end do
785
786 do k = 1, lx
787 do i = 1, lx*lx
788 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
789 + dz(k,2) * u(i,1,2,e) &
790 + dz(k,3) * u(i,1,3,e) &
791 + dz(k,4) * u(i,1,4,e) &
792 + dz(k,5) * u(i,1,5,e) &
793 + dz(k,6) * u(i,1,6,e) &
794 + dz(k,7) * u(i,1,7,e) &
795 + dz(k,8) * u(i,1,8,e) &
796 + dz(k,9) * u(i,1,9,e) &
797 + dz(k,10) * u(i,1,10,e) &
798 + dz(k,11) * u(i,1,11,e) &
799 + dz(k,12) * u(i,1,12,e) &
800 + dz(k,13) * u(i,1,13,e) &
801 + dz(k,14) * u(i,1,14,e)
802 end do
803 end do
804
805 do i = 1, lx * lx * lx
806 ux(i,1,1,e) = w3(i,1,1) &
807 * ( drdx(i,1,1,e) * ur(i,1,1) &
808 + dsdx(i,1,1,e) * us(i,1,1) &
809 + dtdx(i,1,1,e) * ut(i,1,1) )
810 uy(i,1,1,e) = w3(i,1,1) &
811 * ( dsdy(i,1,1,e) * us(i,1,1) &
812 + drdy(i,1,1,e) * ur(i,1,1) &
813 + dtdy(i,1,1,e) * ut(i,1,1) )
814 uz(i,1,1,e) = w3(i,1,1) &
815 * ( dtdz(i,1,1,e) * ut(i,1,1) &
816 + drdz(i,1,1,e) * ur(i,1,1) &
817 + dsdz(i,1,1,e) * us(i,1,1) )
818 end do
819 end do
820 end subroutine cpu_opgrad_lx14
821
822 subroutine cpu_opgrad_lx13(ux, uy, uz, u, dx, dy, dz, &
823 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
824 integer, parameter :: lx = 13
825 integer, intent(in) :: n
826 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
827 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
828 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
829 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
830 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
831 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
832 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
833 real(kind=rp) :: ur(lx, lx, lx)
834 real(kind=rp) :: us(lx, lx, lx)
835 real(kind=rp) :: ut(lx, lx, lx)
836 integer :: e, i, j, k
837
838 do e = 1, n
839 do j = 1, lx * lx
840 do i = 1, lx
841 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
842 + dx(i,2) * u(2,j,1,e) &
843 + dx(i,3) * u(3,j,1,e) &
844 + dx(i,4) * u(4,j,1,e) &
845 + dx(i,5) * u(5,j,1,e) &
846 + dx(i,6) * u(6,j,1,e) &
847 + dx(i,7) * u(7,j,1,e) &
848 + dx(i,8) * u(8,j,1,e) &
849 + dx(i,9) * u(9,j,1,e) &
850 + dx(i,10) * u(10,j,1,e) &
851 + dx(i,11) * u(11,j,1,e) &
852 + dx(i,12) * u(12,j,1,e) &
853 + dx(i,13) * u(13,j,1,e)
854 end do
855 end do
856
857 do k = 1, lx
858 do j = 1, lx
859 do i = 1, lx
860 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
861 + dy(j,2) * u(i,2,k,e) &
862 + dy(j,3) * u(i,3,k,e) &
863 + dy(j,4) * u(i,4,k,e) &
864 + dy(j,5) * u(i,5,k,e) &
865 + dy(j,6) * u(i,6,k,e) &
866 + dy(j,7) * u(i,7,k,e) &
867 + dy(j,8) * u(i,8,k,e) &
868 + dy(j,9) * u(i,9,k,e) &
869 + dy(j,10) * u(i,10,k,e) &
870 + dy(j,11) * u(i,11,k,e) &
871 + dy(j,12) * u(i,12,k,e) &
872 + dy(j,13) * u(i,13,k,e)
873 end do
874 end do
875 end do
876
877 do k = 1, lx
878 do i = 1, lx*lx
879 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
880 + dz(k,2) * u(i,1,2,e) &
881 + dz(k,3) * u(i,1,3,e) &
882 + dz(k,4) * u(i,1,4,e) &
883 + dz(k,5) * u(i,1,5,e) &
884 + dz(k,6) * u(i,1,6,e) &
885 + dz(k,7) * u(i,1,7,e) &
886 + dz(k,8) * u(i,1,8,e) &
887 + dz(k,9) * u(i,1,9,e) &
888 + dz(k,10) * u(i,1,10,e) &
889 + dz(k,11) * u(i,1,11,e) &
890 + dz(k,12) * u(i,1,12,e) &
891 + dz(k,13) * u(i,1,13,e)
892 end do
893 end do
894
895 do i = 1, lx * lx * lx
896 ux(i,1,1,e) = w3(i,1,1) &
897 * ( drdx(i,1,1,e) * ur(i,1,1) &
898 + dsdx(i,1,1,e) * us(i,1,1) &
899 + dtdx(i,1,1,e) * ut(i,1,1) )
900 uy(i,1,1,e) = w3(i,1,1) &
901 * ( dsdy(i,1,1,e) * us(i,1,1) &
902 + drdy(i,1,1,e) * ur(i,1,1) &
903 + dtdy(i,1,1,e) * ut(i,1,1) )
904 uz(i,1,1,e) = w3(i,1,1) &
905 * ( dtdz(i,1,1,e) * ut(i,1,1) &
906 + drdz(i,1,1,e) * ur(i,1,1) &
907 + dsdz(i,1,1,e) * us(i,1,1) )
908 end do
909 end do
910 end subroutine cpu_opgrad_lx13
911
912 subroutine cpu_opgrad_lx12(ux, uy, uz, u, dx, dy, dz, &
913 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
914 integer, parameter :: lx = 12
915 integer, intent(in) :: n
916 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
917 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
918 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
919 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
920 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
921 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
922 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
923 real(kind=rp) :: ur(lx, lx, lx)
924 real(kind=rp) :: us(lx, lx, lx)
925 real(kind=rp) :: ut(lx, lx, lx)
926 integer :: e, i, j, k
927
928 do e = 1, n
929 do j = 1, lx * lx
930 do i = 1, lx
931 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
932 + dx(i,2) * u(2,j,1,e) &
933 + dx(i,3) * u(3,j,1,e) &
934 + dx(i,4) * u(4,j,1,e) &
935 + dx(i,5) * u(5,j,1,e) &
936 + dx(i,6) * u(6,j,1,e) &
937 + dx(i,7) * u(7,j,1,e) &
938 + dx(i,8) * u(8,j,1,e) &
939 + dx(i,9) * u(9,j,1,e) &
940 + dx(i,10) * u(10,j,1,e) &
941 + dx(i,11) * u(11,j,1,e) &
942 + dx(i,12) * u(12,j,1,e)
943 end do
944 end do
945
946 do k = 1, lx
947 do j = 1, lx
948 do i = 1, lx
949 us(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 + dy(j,5) * u(i,5,k,e) &
954 + dy(j,6) * u(i,6,k,e) &
955 + dy(j,7) * u(i,7,k,e) &
956 + dy(j,8) * u(i,8,k,e) &
957 + dy(j,9) * u(i,9,k,e) &
958 + dy(j,10) * u(i,10,k,e) &
959 + dy(j,11) * u(i,11,k,e) &
960 + dy(j,12) * u(i,12,k,e)
961 end do
962 end do
963 end do
964
965 do k = 1, lx
966 do i = 1, lx*lx
967 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
968 + dz(k,2) * u(i,1,2,e) &
969 + dz(k,3) * u(i,1,3,e) &
970 + dz(k,4) * u(i,1,4,e) &
971 + dz(k,5) * u(i,1,5,e) &
972 + dz(k,6) * u(i,1,6,e) &
973 + dz(k,7) * u(i,1,7,e) &
974 + dz(k,8) * u(i,1,8,e) &
975 + dz(k,9) * u(i,1,9,e) &
976 + dz(k,10) * u(i,1,10,e) &
977 + dz(k,11) * u(i,1,11,e) &
978 + dz(k,12) * u(i,1,12,e)
979 end do
980 end do
981
982 do i = 1, lx * lx * lx
983 ux(i,1,1,e) = w3(i,1,1) &
984 * ( drdx(i,1,1,e) * ur(i,1,1) &
985 + dsdx(i,1,1,e) * us(i,1,1) &
986 + dtdx(i,1,1,e) * ut(i,1,1) )
987 uy(i,1,1,e) = w3(i,1,1) &
988 * ( dsdy(i,1,1,e) * us(i,1,1) &
989 + drdy(i,1,1,e) * ur(i,1,1) &
990 + dtdy(i,1,1,e) * ut(i,1,1) )
991 uz(i,1,1,e) = w3(i,1,1) &
992 * ( dtdz(i,1,1,e) * ut(i,1,1) &
993 + drdz(i,1,1,e) * ur(i,1,1) &
994 + dsdz(i,1,1,e) * us(i,1,1) )
995 end do
996 end do
997 end subroutine cpu_opgrad_lx12
998
999 subroutine cpu_opgrad_lx11(ux, uy, uz, u, dx, dy, dz, &
1000 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1001 integer, parameter :: lx = 11
1002 integer, intent(in) :: n
1003 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1004 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1005 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1006 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1007 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1008 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1009 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1010 real(kind=rp) :: ur(lx, lx, lx)
1011 real(kind=rp) :: us(lx, lx, lx)
1012 real(kind=rp) :: ut(lx, lx, lx)
1013 integer :: e, i, j, k
1014
1015 do e = 1, n
1016 do j = 1, lx * lx
1017 do i = 1, lx
1018 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1019 + dx(i,2) * u(2,j,1,e) &
1020 + dx(i,3) * u(3,j,1,e) &
1021 + dx(i,4) * u(4,j,1,e) &
1022 + dx(i,5) * u(5,j,1,e) &
1023 + dx(i,6) * u(6,j,1,e) &
1024 + dx(i,7) * u(7,j,1,e) &
1025 + dx(i,8) * u(8,j,1,e) &
1026 + dx(i,9) * u(9,j,1,e) &
1027 + dx(i,10) * u(10,j,1,e) &
1028 + dx(i,11) * u(11,j,1,e)
1029 end do
1030 end do
1031
1032 do k = 1, lx
1033 do j = 1, lx
1034 do i = 1, lx
1035 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1036 + dy(j,2) * u(i,2,k,e) &
1037 + dy(j,3) * u(i,3,k,e) &
1038 + dy(j,4) * u(i,4,k,e) &
1039 + dy(j,5) * u(i,5,k,e) &
1040 + dy(j,6) * u(i,6,k,e) &
1041 + dy(j,7) * u(i,7,k,e) &
1042 + dy(j,8) * u(i,8,k,e) &
1043 + dy(j,9) * u(i,9,k,e) &
1044 + dy(j,10) * u(i,10,k,e) &
1045 + dy(j,11) * u(i,11,k,e)
1046 end do
1047 end do
1048 end do
1049
1050 do k = 1, lx
1051 do i = 1, lx*lx
1052 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1053 + dz(k,2) * u(i,1,2,e) &
1054 + dz(k,3) * u(i,1,3,e) &
1055 + dz(k,4) * u(i,1,4,e) &
1056 + dz(k,5) * u(i,1,5,e) &
1057 + dz(k,6) * u(i,1,6,e) &
1058 + dz(k,7) * u(i,1,7,e) &
1059 + dz(k,8) * u(i,1,8,e) &
1060 + dz(k,9) * u(i,1,9,e) &
1061 + dz(k,10) * u(i,1,10,e) &
1062 + dz(k,11) * u(i,1,11,e)
1063 end do
1064 end do
1065
1066 do i = 1, lx * lx * lx
1067 ux(i,1,1,e) = w3(i,1,1) &
1068 * ( drdx(i,1,1,e) * ur(i,1,1) &
1069 + dsdx(i,1,1,e) * us(i,1,1) &
1070 + dtdx(i,1,1,e) * ut(i,1,1) )
1071 uy(i,1,1,e) = w3(i,1,1) &
1072 * ( dsdy(i,1,1,e) * us(i,1,1) &
1073 + drdy(i,1,1,e) * ur(i,1,1) &
1074 + dtdy(i,1,1,e) * ut(i,1,1) )
1075 uz(i,1,1,e) = w3(i,1,1) &
1076 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1077 + drdz(i,1,1,e) * ur(i,1,1) &
1078 + dsdz(i,1,1,e) * us(i,1,1) )
1079 end do
1080 end do
1081 end subroutine cpu_opgrad_lx11
1082
1083 subroutine cpu_opgrad_lx10(ux, uy, uz, u, dx, dy, dz, &
1084 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1085 integer, parameter :: lx = 10
1086 integer, intent(in) :: n
1087 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1088 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1089 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1090 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1091 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1092 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1093 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1094 real(kind=rp) :: ur(lx, lx, lx)
1095 real(kind=rp) :: us(lx, lx, lx)
1096 real(kind=rp) :: ut(lx, lx, lx)
1097 integer :: e, i, j, k
1098
1099 do e = 1, n
1100 do j = 1, lx * lx
1101 do i = 1, lx
1102 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1103 + dx(i,2) * u(2,j,1,e) &
1104 + dx(i,3) * u(3,j,1,e) &
1105 + dx(i,4) * u(4,j,1,e) &
1106 + dx(i,5) * u(5,j,1,e) &
1107 + dx(i,6) * u(6,j,1,e) &
1108 + dx(i,7) * u(7,j,1,e) &
1109 + dx(i,8) * u(8,j,1,e) &
1110 + dx(i,9) * u(9,j,1,e) &
1111 + dx(i,10) * u(10,j,1,e)
1112 end do
1113 end do
1114
1115 do k = 1, lx
1116 do j = 1, lx
1117 do i = 1, lx
1118 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1119 + dy(j,2) * u(i,2,k,e) &
1120 + dy(j,3) * u(i,3,k,e) &
1121 + dy(j,4) * u(i,4,k,e) &
1122 + dy(j,5) * u(i,5,k,e) &
1123 + dy(j,6) * u(i,6,k,e) &
1124 + dy(j,7) * u(i,7,k,e) &
1125 + dy(j,8) * u(i,8,k,e) &
1126 + dy(j,9) * u(i,9,k,e) &
1127 + dy(j,10) * u(i,10,k,e)
1128 end do
1129 end do
1130 end do
1131
1132 do k = 1, lx
1133 do i = 1, lx*lx
1134 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1135 + dz(k,2) * u(i,1,2,e) &
1136 + dz(k,3) * u(i,1,3,e) &
1137 + dz(k,4) * u(i,1,4,e) &
1138 + dz(k,5) * u(i,1,5,e) &
1139 + dz(k,6) * u(i,1,6,e) &
1140 + dz(k,7) * u(i,1,7,e) &
1141 + dz(k,8) * u(i,1,8,e) &
1142 + dz(k,9) * u(i,1,9,e) &
1143 + dz(k,10) * u(i,1,10,e)
1144 end do
1145 end do
1146
1147 do i = 1, lx * lx * lx
1148 ux(i,1,1,e) = w3(i,1,1) &
1149 * ( drdx(i,1,1,e) * ur(i,1,1) &
1150 + dsdx(i,1,1,e) * us(i,1,1) &
1151 + dtdx(i,1,1,e) * ut(i,1,1) )
1152 uy(i,1,1,e) = w3(i,1,1) &
1153 * ( dsdy(i,1,1,e) * us(i,1,1) &
1154 + drdy(i,1,1,e) * ur(i,1,1) &
1155 + dtdy(i,1,1,e) * ut(i,1,1) )
1156 uz(i,1,1,e) = w3(i,1,1) &
1157 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1158 + drdz(i,1,1,e) * ur(i,1,1) &
1159 + dsdz(i,1,1,e) * us(i,1,1) )
1160 end do
1161 end do
1162 end subroutine cpu_opgrad_lx10
1163
1164 subroutine cpu_opgrad_lx9(ux, uy, uz, u, dx, dy, dz, &
1165 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1166 integer, parameter :: lx = 9
1167 integer, intent(in) :: n
1168 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1169 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1170 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1171 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1172 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1173 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1174 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1175 real(kind=rp) :: ur(lx, lx, lx)
1176 real(kind=rp) :: us(lx, lx, lx)
1177 real(kind=rp) :: ut(lx, lx, lx)
1178 integer :: e, i, j, k
1179
1180 do e = 1, n
1181 do j = 1, lx * lx
1182 do i = 1, lx
1183 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1184 + dx(i,2) * u(2,j,1,e) &
1185 + dx(i,3) * u(3,j,1,e) &
1186 + dx(i,4) * u(4,j,1,e) &
1187 + dx(i,5) * u(5,j,1,e) &
1188 + dx(i,6) * u(6,j,1,e) &
1189 + dx(i,7) * u(7,j,1,e) &
1190 + dx(i,8) * u(8,j,1,e) &
1191 + dx(i,9) * u(9,j,1,e)
1192 end do
1193 end do
1194
1195 do k = 1, lx
1196 do j = 1, lx
1197 do i = 1, lx
1198 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1199 + dy(j,2) * u(i,2,k,e) &
1200 + dy(j,3) * u(i,3,k,e) &
1201 + dy(j,4) * u(i,4,k,e) &
1202 + dy(j,5) * u(i,5,k,e) &
1203 + dy(j,6) * u(i,6,k,e) &
1204 + dy(j,7) * u(i,7,k,e) &
1205 + dy(j,8) * u(i,8,k,e) &
1206 + dy(j,9) * u(i,9,k,e)
1207 end do
1208 end do
1209 end do
1210
1211 do k = 1, lx
1212 do i = 1, lx*lx
1213 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1214 + dz(k,2) * u(i,1,2,e) &
1215 + dz(k,3) * u(i,1,3,e) &
1216 + dz(k,4) * u(i,1,4,e) &
1217 + dz(k,5) * u(i,1,5,e) &
1218 + dz(k,6) * u(i,1,6,e) &
1219 + dz(k,7) * u(i,1,7,e) &
1220 + dz(k,8) * u(i,1,8,e) &
1221 + dz(k,9) * u(i,1,9,e)
1222 end do
1223 end do
1224
1225 do i = 1, lx * lx * lx
1226 ux(i,1,1,e) = w3(i,1,1) &
1227 * ( drdx(i,1,1,e) * ur(i,1,1) &
1228 + dsdx(i,1,1,e) * us(i,1,1) &
1229 + dtdx(i,1,1,e) * ut(i,1,1) )
1230 uy(i,1,1,e) = w3(i,1,1) &
1231 * ( dsdy(i,1,1,e) * us(i,1,1) &
1232 + drdy(i,1,1,e) * ur(i,1,1) &
1233 + dtdy(i,1,1,e) * ut(i,1,1) )
1234 uz(i,1,1,e) = w3(i,1,1) &
1235 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1236 + drdz(i,1,1,e) * ur(i,1,1) &
1237 + dsdz(i,1,1,e) * us(i,1,1) )
1238 end do
1239 end do
1240 end subroutine cpu_opgrad_lx9
1241
1242 subroutine cpu_opgrad_lx8(ux, uy, uz, u, dx, dy, dz, &
1243 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1244 integer, parameter :: lx = 8
1245 integer, intent(in) :: n
1246 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1247 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1248 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1249 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1250 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1251 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1252 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1253 real(kind=rp) :: ur(lx, lx, lx)
1254 real(kind=rp) :: us(lx, lx, lx)
1255 real(kind=rp) :: ut(lx, lx, lx)
1256 integer :: e, i, j, k
1257
1258 do e = 1, n
1259 do j = 1, lx * lx
1260 do i = 1, lx
1261 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1262 + dx(i,2) * u(2,j,1,e) &
1263 + dx(i,3) * u(3,j,1,e) &
1264 + dx(i,4) * u(4,j,1,e) &
1265 + dx(i,5) * u(5,j,1,e) &
1266 + dx(i,6) * u(6,j,1,e) &
1267 + dx(i,7) * u(7,j,1,e) &
1268 + dx(i,8) * u(8,j,1,e)
1269 end do
1270 end do
1271
1272 do k = 1, lx
1273 do j = 1, lx
1274 do i = 1, lx
1275 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1276 + dy(j,2) * u(i,2,k,e) &
1277 + dy(j,3) * u(i,3,k,e) &
1278 + dy(j,4) * u(i,4,k,e) &
1279 + dy(j,5) * u(i,5,k,e) &
1280 + dy(j,6) * u(i,6,k,e) &
1281 + dy(j,7) * u(i,7,k,e) &
1282 + dy(j,8) * u(i,8,k,e)
1283 end do
1284 end do
1285 end do
1286
1287 do k = 1, lx
1288 do i = 1, lx*lx
1289 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1290 + dz(k,2) * u(i,1,2,e) &
1291 + dz(k,3) * u(i,1,3,e) &
1292 + dz(k,4) * u(i,1,4,e) &
1293 + dz(k,5) * u(i,1,5,e) &
1294 + dz(k,6) * u(i,1,6,e) &
1295 + dz(k,7) * u(i,1,7,e) &
1296 + dz(k,8) * u(i,1,8,e)
1297 end do
1298 end do
1299
1300 do i = 1, lx * lx * lx
1301 ux(i,1,1,e) = w3(i,1,1) &
1302 * ( drdx(i,1,1,e) * ur(i,1,1) &
1303 + dsdx(i,1,1,e) * us(i,1,1) &
1304 + dtdx(i,1,1,e) * ut(i,1,1) )
1305 uy(i,1,1,e) = w3(i,1,1) &
1306 * ( dsdy(i,1,1,e) * us(i,1,1) &
1307 + drdy(i,1,1,e) * ur(i,1,1) &
1308 + dtdy(i,1,1,e) * ut(i,1,1) )
1309 uz(i,1,1,e) = w3(i,1,1) &
1310 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1311 + drdz(i,1,1,e) * ur(i,1,1) &
1312 + dsdz(i,1,1,e) * us(i,1,1) )
1313 end do
1314 end do
1315 end subroutine cpu_opgrad_lx8
1316
1317 subroutine cpu_opgrad_lx7(ux, uy, uz, u, dx, dy, dz, &
1318 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1319 integer, parameter :: lx = 7
1320 integer, intent(in) :: n
1321 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1322 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1323 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1324 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1325 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1326 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1327 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1328 real(kind=rp) :: ur(lx, lx, lx)
1329 real(kind=rp) :: us(lx, lx, lx)
1330 real(kind=rp) :: ut(lx, lx, lx)
1331 integer :: e, i, j, k
1332
1333 do e = 1, n
1334 do j = 1, lx * lx
1335 do i = 1, lx
1336 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1337 + dx(i,2) * u(2,j,1,e) &
1338 + dx(i,3) * u(3,j,1,e) &
1339 + dx(i,4) * u(4,j,1,e) &
1340 + dx(i,5) * u(5,j,1,e) &
1341 + dx(i,6) * u(6,j,1,e) &
1342 + dx(i,7) * u(7,j,1,e)
1343 end do
1344 end do
1345
1346 do k = 1, lx
1347 do j = 1, lx
1348 do i = 1, lx
1349 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1350 + dy(j,2) * u(i,2,k,e) &
1351 + dy(j,3) * u(i,3,k,e) &
1352 + dy(j,4) * u(i,4,k,e) &
1353 + dy(j,5) * u(i,5,k,e) &
1354 + dy(j,6) * u(i,6,k,e) &
1355 + dy(j,7) * u(i,7,k,e)
1356 end do
1357 end do
1358 end do
1359
1360 do k = 1, lx
1361 do i = 1, lx*lx
1362 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1363 + dz(k,2) * u(i,1,2,e) &
1364 + dz(k,3) * u(i,1,3,e) &
1365 + dz(k,4) * u(i,1,4,e) &
1366 + dz(k,5) * u(i,1,5,e) &
1367 + dz(k,6) * u(i,1,6,e) &
1368 + dz(k,7) * u(i,1,7,e)
1369 end do
1370 end do
1371
1372 do i = 1, lx * lx * lx
1373 ux(i,1,1,e) = w3(i,1,1) &
1374 * ( drdx(i,1,1,e) * ur(i,1,1) &
1375 + dsdx(i,1,1,e) * us(i,1,1) &
1376 + dtdx(i,1,1,e) * ut(i,1,1) )
1377 uy(i,1,1,e) = w3(i,1,1) &
1378 * ( dsdy(i,1,1,e) * us(i,1,1) &
1379 + drdy(i,1,1,e) * ur(i,1,1) &
1380 + dtdy(i,1,1,e) * ut(i,1,1) )
1381 uz(i,1,1,e) = w3(i,1,1) &
1382 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1383 + drdz(i,1,1,e) * ur(i,1,1) &
1384 + dsdz(i,1,1,e) * us(i,1,1) )
1385 end do
1386 end do
1387 end subroutine cpu_opgrad_lx7
1388
1389 subroutine cpu_opgrad_lx6(ux, uy, uz, u, dx, dy, dz, &
1390 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1391 integer, parameter :: lx = 6
1392 integer, intent(in) :: n
1393 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1394 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1395 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1396 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1397 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1398 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1399 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1400 real(kind=rp) :: ur(lx, lx, lx)
1401 real(kind=rp) :: us(lx, lx, lx)
1402 real(kind=rp) :: ut(lx, lx, lx)
1403 integer :: e, i, j, k
1404
1405 do e = 1, n
1406 do j = 1, lx * lx
1407 do i = 1, lx
1408 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1409 + dx(i,2) * u(2,j,1,e) &
1410 + dx(i,3) * u(3,j,1,e) &
1411 + dx(i,4) * u(4,j,1,e) &
1412 + dx(i,5) * u(5,j,1,e) &
1413 + dx(i,6) * u(6,j,1,e)
1414 end do
1415 end do
1416
1417 do k = 1, lx
1418 do j = 1, lx
1419 do i = 1, lx
1420 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1421 + dy(j,2) * u(i,2,k,e) &
1422 + dy(j,3) * u(i,3,k,e) &
1423 + dy(j,4) * u(i,4,k,e) &
1424 + dy(j,5) * u(i,5,k,e) &
1425 + dy(j,6) * u(i,6,k,e)
1426 end do
1427 end do
1428 end do
1429
1430 do k = 1, lx
1431 do i = 1, lx*lx
1432 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1433 + dz(k,2) * u(i,1,2,e) &
1434 + dz(k,3) * u(i,1,3,e) &
1435 + dz(k,4) * u(i,1,4,e) &
1436 + dz(k,5) * u(i,1,5,e) &
1437 + dz(k,6) * u(i,1,6,e)
1438 end do
1439 end do
1440
1441 do i = 1, lx * lx * lx
1442 ux(i,1,1,e) = w3(i,1,1) &
1443 * ( drdx(i,1,1,e) * ur(i,1,1) &
1444 + dsdx(i,1,1,e) * us(i,1,1) &
1445 + dtdx(i,1,1,e) * ut(i,1,1) )
1446 uy(i,1,1,e) = w3(i,1,1) &
1447 * ( dsdy(i,1,1,e) * us(i,1,1) &
1448 + drdy(i,1,1,e) * ur(i,1,1) &
1449 + dtdy(i,1,1,e) * ut(i,1,1) )
1450 uz(i,1,1,e) = w3(i,1,1) &
1451 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1452 + drdz(i,1,1,e) * ur(i,1,1) &
1453 + dsdz(i,1,1,e) * us(i,1,1) )
1454 end do
1455 end do
1456 end subroutine cpu_opgrad_lx6
1457
1458 subroutine cpu_opgrad_lx5(ux, uy, uz, u, dx, dy, dz, &
1459 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1460 integer, parameter :: lx = 5
1461 integer, intent(in) :: n
1462 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1463 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1464 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1465 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1466 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1467 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1468 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1469 real(kind=rp) :: ur(lx, lx, lx)
1470 real(kind=rp) :: us(lx, lx, lx)
1471 real(kind=rp) :: ut(lx, lx, lx)
1472 integer :: e, i, j, k
1473
1474 do e = 1, n
1475 do j = 1, lx * lx
1476 do i = 1, lx
1477 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1478 + dx(i,2) * u(2,j,1,e) &
1479 + dx(i,3) * u(3,j,1,e) &
1480 + dx(i,4) * u(4,j,1,e) &
1481 + dx(i,5) * u(5,j,1,e)
1482 end do
1483 end do
1484
1485 do k = 1, lx
1486 do j = 1, lx
1487 do i = 1, lx
1488 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1489 + dy(j,2) * u(i,2,k,e) &
1490 + dy(j,3) * u(i,3,k,e) &
1491 + dy(j,4) * u(i,4,k,e) &
1492 + dy(j,5) * u(i,5,k,e)
1493 end do
1494 end do
1495 end do
1496
1497 do k = 1, lx
1498 do i = 1, lx*lx
1499 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1500 + dz(k,2) * u(i,1,2,e) &
1501 + dz(k,3) * u(i,1,3,e) &
1502 + dz(k,4) * u(i,1,4,e) &
1503 + dz(k,5) * u(i,1,5,e)
1504 end do
1505 end do
1506
1507 do i = 1, lx * lx * lx
1508 ux(i,1,1,e) = w3(i,1,1) &
1509 * ( drdx(i,1,1,e) * ur(i,1,1) &
1510 + dsdx(i,1,1,e) * us(i,1,1) &
1511 + dtdx(i,1,1,e) * ut(i,1,1) )
1512 uy(i,1,1,e) = w3(i,1,1) &
1513 * ( dsdy(i,1,1,e) * us(i,1,1) &
1514 + drdy(i,1,1,e) * ur(i,1,1) &
1515 + dtdy(i,1,1,e) * ut(i,1,1) )
1516 uz(i,1,1,e) = w3(i,1,1) &
1517 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1518 + drdz(i,1,1,e) * ur(i,1,1) &
1519 + dsdz(i,1,1,e) * us(i,1,1) )
1520 end do
1521 end do
1522 end subroutine cpu_opgrad_lx5
1523
1524 subroutine cpu_opgrad_lx4(ux, uy, uz, u, dx, dy, dz, &
1525 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1526 integer, parameter :: lx = 4
1527 integer, intent(in) :: n
1528 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1529 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1530 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1531 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1532 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1533 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1534 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1535 real(kind=rp) :: ur(lx, lx, lx)
1536 real(kind=rp) :: us(lx, lx, lx)
1537 real(kind=rp) :: ut(lx, lx, lx)
1538 integer :: e, i, j, k
1539
1540 do e = 1, n
1541 do j = 1, lx * lx
1542 do i = 1, lx
1543 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1544 + dx(i,2) * u(2,j,1,e) &
1545 + dx(i,3) * u(3,j,1,e) &
1546 + dx(i,4) * u(4,j,1,e)
1547 end do
1548 end do
1549
1550 do k = 1, lx
1551 do j = 1, lx
1552 do i = 1, lx
1553 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1554 + dy(j,2) * u(i,2,k,e) &
1555 + dy(j,3) * u(i,3,k,e) &
1556 + dy(j,4) * u(i,4,k,e)
1557 end do
1558 end do
1559 end do
1560
1561 do k = 1, lx
1562 do i = 1, lx*lx
1563 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1564 + dz(k,2) * u(i,1,2,e) &
1565 + dz(k,3) * u(i,1,3,e) &
1566 + dz(k,4) * u(i,1,4,e)
1567 end do
1568 end do
1569
1570 do i = 1, lx * lx * lx
1571 ux(i,1,1,e) = w3(i,1,1) &
1572 * ( drdx(i,1,1,e) * ur(i,1,1) &
1573 + dsdx(i,1,1,e) * us(i,1,1) &
1574 + dtdx(i,1,1,e) * ut(i,1,1) )
1575 uy(i,1,1,e) = w3(i,1,1) &
1576 * ( dsdy(i,1,1,e) * us(i,1,1) &
1577 + drdy(i,1,1,e) * ur(i,1,1) &
1578 + dtdy(i,1,1,e) * ut(i,1,1) )
1579 uz(i,1,1,e) = w3(i,1,1) &
1580 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1581 + drdz(i,1,1,e) * ur(i,1,1) &
1582 + dsdz(i,1,1,e) * us(i,1,1) )
1583 end do
1584 end do
1585 end subroutine cpu_opgrad_lx4
1586
1587 subroutine cpu_opgrad_lx3(ux, uy, uz, u, dx, dy, dz, &
1588 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1589 integer, parameter :: lx = 3
1590 integer, intent(in) :: n
1591 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1592 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1593 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1594 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1595 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1596 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1597 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1598 real(kind=rp) :: ur(lx, lx, lx)
1599 real(kind=rp) :: us(lx, lx, lx)
1600 real(kind=rp) :: ut(lx, lx, lx)
1601 integer :: e, i, j, k
1602
1603 do e = 1, n
1604 do j = 1, lx * lx
1605 do i = 1, lx
1606 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1607 + dx(i,2) * u(2,j,1,e) &
1608 + dx(i,3) * u(3,j,1,e)
1609 end do
1610 end do
1611
1612 do k = 1, lx
1613 do j = 1, lx
1614 do i = 1, lx
1615 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1616 + dy(j,2) * u(i,2,k,e) &
1617 + dy(j,3) * u(i,3,k,e)
1618 end do
1619 end do
1620 end do
1621
1622 do k = 1, lx
1623 do i = 1, lx*lx
1624 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1625 + dz(k,2) * u(i,1,2,e) &
1626 + dz(k,3) * u(i,1,3,e)
1627 end do
1628 end do
1629
1630 do i = 1, lx * lx * lx
1631 ux(i,1,1,e) = w3(i,1,1) &
1632 * ( drdx(i,1,1,e) * ur(i,1,1) &
1633 + dsdx(i,1,1,e) * us(i,1,1) &
1634 + dtdx(i,1,1,e) * ut(i,1,1) )
1635 uy(i,1,1,e) = w3(i,1,1) &
1636 * ( dsdy(i,1,1,e) * us(i,1,1) &
1637 + drdy(i,1,1,e) * ur(i,1,1) &
1638 + dtdy(i,1,1,e) * ut(i,1,1) )
1639 uz(i,1,1,e) = w3(i,1,1) &
1640 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1641 + drdz(i,1,1,e) * ur(i,1,1) &
1642 + dsdz(i,1,1,e) * us(i,1,1) )
1643 end do
1644 end do
1645 end subroutine cpu_opgrad_lx3
1646
1647 subroutine cpu_opgrad_lx2(ux, uy, uz, u, dx, dy, dz, &
1648 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, n)
1649 integer, parameter :: lx = 2
1650 integer, intent(in) :: n
1651 real(kind=rp), dimension(lx, lx, lx, n), intent(inout) :: ux, uy, uz
1652 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: u
1653 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1654 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdx, dsdx, dtdx
1655 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdy, dsdy, dtdy
1656 real(kind=rp), dimension(lx, lx, lx, n), intent(in) :: drdz, dsdz, dtdz
1657 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1658 real(kind=rp) :: ur(lx, lx, lx)
1659 real(kind=rp) :: us(lx, lx, lx)
1660 real(kind=rp) :: ut(lx, lx, lx)
1661 integer :: e, i, j, k
1662
1663 do e = 1, n
1664 do j = 1, lx * lx
1665 do i = 1, lx
1666 ur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1667 + dx(i,2) * u(2,j,1,e)
1668 end do
1669 end do
1670
1671 do k = 1, lx
1672 do j = 1, lx
1673 do i = 1, lx
1674 us(i,j,k) = dy(j,1) * u(i,1,k,e) &
1675 + dy(j,2) * u(i,2,k,e)
1676 end do
1677 end do
1678 end do
1679
1680 do k = 1, lx
1681 do i = 1, lx*lx
1682 ut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1683 + dz(k,2) * u(i,1,2,e)
1684 end do
1685 end do
1686
1687 do i = 1, lx * lx * lx
1688 ux(i,1,1,e) = w3(i,1,1) &
1689 * ( drdx(i,1,1,e) * ur(i,1,1) &
1690 + dsdx(i,1,1,e) * us(i,1,1) &
1691 + dtdx(i,1,1,e) * ut(i,1,1) )
1692 uy(i,1,1,e) = w3(i,1,1) &
1693 * ( dsdy(i,1,1,e) * us(i,1,1) &
1694 + drdy(i,1,1,e) * ur(i,1,1) &
1695 + dtdy(i,1,1,e) * ut(i,1,1) )
1696 uz(i,1,1,e) = w3(i,1,1) &
1697 * ( dtdz(i,1,1,e) * ut(i,1,1) &
1698 + drdz(i,1,1,e) * ur(i,1,1) &
1699 + dsdz(i,1,1,e) * us(i,1,1) )
1700 end do
1701 end do
1702 end subroutine cpu_opgrad_lx2
1703
1704 subroutine opr_cpu_opgrad_single(ux, uy, uz, u, coef, e)
1705 integer, parameter :: e_len = 1
1706 type(coef_t), intent(in) :: coef
1707 integer, intent(in) :: e
1708 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: ux
1709 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: uy
1710 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(inout) :: uz
1711 real(kind=rp), dimension(coef%Xh%lxyz, e_len), intent(in) :: u
1712
1713 associate(xh => coef%Xh, msh => coef%msh, &
1714 drdx => coef%drdx, drdy => coef%drdy, drdz => coef%drdz, &
1715 dsdx => coef%dsdx, dsdy => coef%dsdy, dsdz => coef%dsdz, &
1716 dtdx => coef%dtdx, dtdy => coef%dtdy, dtdz => coef%dtdz)
1717
1718 select case (xh%lx)
1719 case (18)
1720 call cpu_opgrad_lx18_single(ux, uy, uz, u, &
1721 xh%dx, xh%dy, xh%dz, &
1722 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1723 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1724 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1725 xh%w3)
1726 case (17)
1727 call cpu_opgrad_lx17_single(ux, uy, uz, u, &
1728 xh%dx, xh%dy, xh%dz, &
1729 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1730 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1731 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1732 xh%w3)
1733 case (16)
1734 call cpu_opgrad_lx16_single(ux, uy, uz, u, &
1735 xh%dx, xh%dy, xh%dz, &
1736 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1737 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1738 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1739 xh%w3)
1740 case (15)
1741 call cpu_opgrad_lx15_single(ux, uy, uz, u, &
1742 xh%dx, xh%dy, xh%dz, &
1743 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1744 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1745 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1746 xh%w3)
1747 case (14)
1748 call cpu_opgrad_lx14_single(ux, uy, uz, u, &
1749 xh%dx, xh%dy, xh%dz, &
1750 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1751 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1752 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1753 xh%w3)
1754 case (13)
1755 call cpu_opgrad_lx13_single(ux, uy, uz, u, &
1756 xh%dx, xh%dy, xh%dz, &
1757 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1758 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1759 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1760 xh%w3)
1761 case (12)
1762 call cpu_opgrad_lx12_single(ux, uy, uz, u, &
1763 xh%dx, xh%dy, xh%dz, &
1764 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1765 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1766 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1767 xh%w3)
1768 case (11)
1769 call cpu_opgrad_lx11_single(ux, uy, uz, u, &
1770 xh%dx, xh%dy, xh%dz, &
1771 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1772 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1773 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1774 xh%w3)
1775 case (10)
1776 call cpu_opgrad_lx10_single(ux, uy, uz, u, &
1777 xh%dx, xh%dy, xh%dz, &
1778 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1779 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1780 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1781 xh%w3)
1782
1783 case (9)
1784 call cpu_opgrad_lx9_single(ux, uy, uz, u, &
1785 xh%dx, xh%dy, xh%dz, &
1786 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1787 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1788 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1789 xh%w3)
1790 case (8)
1791 call cpu_opgrad_lx8_single(ux, uy, uz, u, &
1792 xh%dx, xh%dy, xh%dz, &
1793 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1794 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1795 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1796 xh%w3)
1797 case (7)
1798 call cpu_opgrad_lx7_single(ux, uy, uz, u, &
1799 xh%dx, xh%dy, xh%dz, &
1800 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1801 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1802 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1803 xh%w3)
1804 case (6)
1805 call cpu_opgrad_lx6_single(ux, uy, uz, u, &
1806 xh%dx, xh%dy, xh%dz, &
1807 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1808 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1809 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1810 xh%w3)
1811 case (5)
1812 call cpu_opgrad_lx5_single(ux, uy, uz, u, &
1813 xh%dx, xh%dy, xh%dz, &
1814 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1815 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1816 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1817 xh%w3)
1818 case (4)
1819 call cpu_opgrad_lx4_single(ux, uy, uz, u, &
1820 xh%dx, xh%dy, xh%dz, &
1821 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1822 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1823 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1824 xh%w3)
1825 case (3)
1826 call cpu_opgrad_lx3_single(ux, uy, uz, u, &
1827 xh%dx, xh%dy, xh%dz, &
1828 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1829 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1830 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1831 xh%w3)
1832 case (2)
1833 call cpu_opgrad_lx2_single(ux, uy, uz, u, &
1834 xh%dx, xh%dy, xh%dz, &
1835 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1836 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1837 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1838 xh%w3)
1839 case default
1840 call cpu_opgrad_lx_single(ux, uy, uz, u, &
1841 xh%dx, xh%dy, xh%dz, &
1842 drdx(1,1,1,e), dsdx(1,1,1,e), dtdx(1,1,1,e), &
1843 drdy(1,1,1,e), dsdy(1,1,1,e), dtdy(1,1,1,e), &
1844 drdz(1,1,1,e), dsdz(1,1,1,e), dtdz(1,1,1,e), &
1845 xh%w3, xh%lx)
1846 end select
1847 end associate
1848
1849 end subroutine opr_cpu_opgrad_single
1850
1851 subroutine cpu_opgrad_lx_single(ux, uy, uz, u, dx, dy, dz, &
1852 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3, lx)
1853 integer, intent(in) :: lx
1854 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
1855 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
1856 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1857 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
1858 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
1859 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
1860 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1861 real(kind=rp) :: ur(lx, lx, lx)
1862 real(kind=rp) :: us(lx, lx, lx)
1863 real(kind=rp) :: ut(lx, lx, lx)
1864 real(kind=rp) :: tmp
1865 integer :: i, j, k, l
1866
1867 do j = 1, lx * lx
1868 do i = 1, lx
1869 tmp = 0.0_rp
1870 do k = 1, lx
1871 tmp = tmp + dx(i,k) * u(k,j,1)
1872 end do
1873 ur(i,j,1) = tmp
1874 end do
1875 end do
1876
1877 do k = 1, lx
1878 do j = 1, lx
1879 do i = 1, lx
1880 tmp = 0.0_rp
1881 do l = 1, lx
1882 tmp = tmp + dy(j,l) * u(i,l,k)
1883 end do
1884 us(i,j,k) = tmp
1885 end do
1886 end do
1887 end do
1888
1889 do k = 1, lx
1890 do i = 1, lx*lx
1891 tmp = 0.0_rp
1892 do l = 1, lx
1893 tmp = tmp + dz(k,l) * u(i,1,l)
1894 end do
1895 ut(i,1,k) = tmp
1896 end do
1897 end do
1898
1899 do i = 1, lx * lx * lx
1900 ux(i,1,1) = w3(i,1,1) &
1901 * ( drdx(i,1,1) * ur(i,1,1) &
1902 + dsdx(i,1,1) * us(i,1,1) &
1903 + dtdx(i,1,1) * ut(i,1,1) )
1904 uy(i,1,1) = w3(i,1,1) &
1905 * ( dsdy(i,1,1) * us(i,1,1) &
1906 + drdy(i,1,1) * ur(i,1,1) &
1907 + dtdy(i,1,1) * ut(i,1,1) )
1908 uz(i,1,1) = w3(i,1,1) &
1909 * ( dtdz(i,1,1) * ut(i,1,1) &
1910 + drdz(i,1,1) * ur(i,1,1) &
1911 + dsdz(i,1,1) * us(i,1,1) )
1912 end do
1913
1914 end subroutine cpu_opgrad_lx_single
1915
1916 subroutine cpu_opgrad_lx18_single(ux, uy, uz, u, dx, dy, dz, &
1917 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
1918 integer, parameter :: lx = 18
1919 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
1920 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
1921 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1922 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
1923 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
1924 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
1925 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
1926 real(kind=rp) :: ur(lx, lx, lx)
1927 real(kind=rp) :: us(lx, lx, lx)
1928 real(kind=rp) :: ut(lx, lx, lx)
1929 integer :: i, j, k
1930
1931 do j = 1, lx * lx
1932 do i = 1, lx
1933 ur(i,j,1) = dx(i,1) * u(1,j,1) &
1934 + dx(i,2) * u(2,j,1) &
1935 + dx(i,3) * u(3,j,1) &
1936 + dx(i,4) * u(4,j,1) &
1937 + dx(i,5) * u(5,j,1) &
1938 + dx(i,6) * u(6,j,1) &
1939 + dx(i,7) * u(7,j,1) &
1940 + dx(i,8) * u(8,j,1) &
1941 + dx(i,9) * u(9,j,1) &
1942 + dx(i,10) * u(10,j,1) &
1943 + dx(i,11) * u(11,j,1) &
1944 + dx(i,12) * u(12,j,1) &
1945 + dx(i,13) * u(13,j,1) &
1946 + dx(i,14) * u(14,j,1) &
1947 + dx(i,15) * u(15,j,1) &
1948 + dx(i,16) * u(16,j,1) &
1949 + dx(i,17) * u(17,j,1) &
1950 + dx(i,18) * u(18,j,1)
1951 end do
1952 end do
1953
1954 do k = 1, lx
1955 do j = 1, lx
1956 do i = 1, lx
1957 us(i,j,k) = dy(j,1) * u(i,1,k) &
1958 + dy(j,2) * u(i,2,k) &
1959 + dy(j,3) * u(i,3,k) &
1960 + dy(j,4) * u(i,4,k) &
1961 + dy(j,5) * u(i,5,k) &
1962 + dy(j,6) * u(i,6,k) &
1963 + dy(j,7) * u(i,7,k) &
1964 + dy(j,8) * u(i,8,k) &
1965 + dy(j,9) * u(i,9,k) &
1966 + dy(j,10) * u(i,10,k) &
1967 + dy(j,11) * u(i,11,k) &
1968 + dy(j,12) * u(i,12,k) &
1969 + dy(j,13) * u(i,13,k) &
1970 + dy(j,14) * u(i,14,k) &
1971 + dy(j,15) * u(i,15,k) &
1972 + dy(j,16) * u(i,16,k) &
1973 + dy(j,17) * u(i,17,k) &
1974 + dy(j,18) * u(i,18,k)
1975 end do
1976 end do
1977 end do
1978
1979 do k = 1, lx
1980 do i = 1, lx*lx
1981 ut(i,1,k) = dz(k,1) * u(i,1,1) &
1982 + dz(k,2) * u(i,1,2) &
1983 + dz(k,3) * u(i,1,3) &
1984 + dz(k,4) * u(i,1,4) &
1985 + dz(k,5) * u(i,1,5) &
1986 + dz(k,6) * u(i,1,6) &
1987 + dz(k,7) * u(i,1,7) &
1988 + dz(k,8) * u(i,1,8) &
1989 + dz(k,9) * u(i,1,9) &
1990 + dz(k,10) * u(i,1,10) &
1991 + dz(k,11) * u(i,1,11) &
1992 + dz(k,12) * u(i,1,12) &
1993 + dz(k,13) * u(i,1,13) &
1994 + dz(k,14) * u(i,1,14) &
1995 + dz(k,15) * u(i,1,15) &
1996 + dz(k,16) * u(i,1,16) &
1997 + dz(k,17) * u(i,1,17) &
1998 + dz(k,18) * u(i,1,18)
1999 end do
2000 end do
2001
2002 do i = 1, lx * lx * lx
2003 ux(i,1,1) = w3(i,1,1) &
2004 * ( drdx(i,1,1) * ur(i,1,1) &
2005 + dsdx(i,1,1) * us(i,1,1) &
2006 + dtdx(i,1,1) * ut(i,1,1) )
2007 uy(i,1,1) = w3(i,1,1) &
2008 * ( dsdy(i,1,1) * us(i,1,1) &
2009 + drdy(i,1,1) * ur(i,1,1) &
2010 + dtdy(i,1,1) * ut(i,1,1) )
2011 uz(i,1,1) = w3(i,1,1) &
2012 * ( dtdz(i,1,1) * ut(i,1,1) &
2013 + drdz(i,1,1) * ur(i,1,1) &
2014 + dsdz(i,1,1) * us(i,1,1) )
2015 end do
2016
2017 end subroutine cpu_opgrad_lx18_single
2018
2019 subroutine cpu_opgrad_lx17_single(ux, uy, uz, u, dx, dy, dz, &
2020 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2021 integer, parameter :: lx = 17
2022 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2023 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2024 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2025 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2026 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2027 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2028 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2029 real(kind=rp) :: ur(lx, lx, lx)
2030 real(kind=rp) :: us(lx, lx, lx)
2031 real(kind=rp) :: ut(lx, lx, lx)
2032 integer :: i, j, k
2033
2034 do j = 1, lx * lx
2035 do i = 1, lx
2036 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2037 + dx(i,2) * u(2,j,1) &
2038 + dx(i,3) * u(3,j,1) &
2039 + dx(i,4) * u(4,j,1) &
2040 + dx(i,5) * u(5,j,1) &
2041 + dx(i,6) * u(6,j,1) &
2042 + dx(i,7) * u(7,j,1) &
2043 + dx(i,8) * u(8,j,1) &
2044 + dx(i,9) * u(9,j,1) &
2045 + dx(i,10) * u(10,j,1) &
2046 + dx(i,11) * u(11,j,1) &
2047 + dx(i,12) * u(12,j,1) &
2048 + dx(i,13) * u(13,j,1) &
2049 + dx(i,14) * u(14,j,1) &
2050 + dx(i,15) * u(15,j,1) &
2051 + dx(i,16) * u(16,j,1) &
2052 + dx(i,17) * u(17,j,1)
2053 end do
2054 end do
2055
2056 do k = 1, lx
2057 do j = 1, lx
2058 do i = 1, lx
2059 us(i,j,k) = dy(j,1) * u(i,1,k) &
2060 + dy(j,2) * u(i,2,k) &
2061 + dy(j,3) * u(i,3,k) &
2062 + dy(j,4) * u(i,4,k) &
2063 + dy(j,5) * u(i,5,k) &
2064 + dy(j,6) * u(i,6,k) &
2065 + dy(j,7) * u(i,7,k) &
2066 + dy(j,8) * u(i,8,k) &
2067 + dy(j,9) * u(i,9,k) &
2068 + dy(j,10) * u(i,10,k) &
2069 + dy(j,11) * u(i,11,k) &
2070 + dy(j,12) * u(i,12,k) &
2071 + dy(j,13) * u(i,13,k) &
2072 + dy(j,14) * u(i,14,k) &
2073 + dy(j,15) * u(i,15,k) &
2074 + dy(j,16) * u(i,16,k) &
2075 + dy(j,17) * u(i,17,k)
2076 end do
2077 end do
2078 end do
2079
2080 do k = 1, lx
2081 do i = 1, lx*lx
2082 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2083 + dz(k,2) * u(i,1,2) &
2084 + dz(k,3) * u(i,1,3) &
2085 + dz(k,4) * u(i,1,4) &
2086 + dz(k,5) * u(i,1,5) &
2087 + dz(k,6) * u(i,1,6) &
2088 + dz(k,7) * u(i,1,7) &
2089 + dz(k,8) * u(i,1,8) &
2090 + dz(k,9) * u(i,1,9) &
2091 + dz(k,10) * u(i,1,10) &
2092 + dz(k,11) * u(i,1,11) &
2093 + dz(k,12) * u(i,1,12) &
2094 + dz(k,13) * u(i,1,13) &
2095 + dz(k,14) * u(i,1,14) &
2096 + dz(k,15) * u(i,1,15) &
2097 + dz(k,16) * u(i,1,16) &
2098 + dz(k,17) * u(i,1,17)
2099 end do
2100 end do
2101
2102 do i = 1, lx * lx * lx
2103 ux(i,1,1) = w3(i,1,1) &
2104 * ( drdx(i,1,1) * ur(i,1,1) &
2105 + dsdx(i,1,1) * us(i,1,1) &
2106 + dtdx(i,1,1) * ut(i,1,1) )
2107 uy(i,1,1) = w3(i,1,1) &
2108 * ( dsdy(i,1,1) * us(i,1,1) &
2109 + drdy(i,1,1) * ur(i,1,1) &
2110 + dtdy(i,1,1) * ut(i,1,1) )
2111 uz(i,1,1) = w3(i,1,1) &
2112 * ( dtdz(i,1,1) * ut(i,1,1) &
2113 + drdz(i,1,1) * ur(i,1,1) &
2114 + dsdz(i,1,1) * us(i,1,1) )
2115 end do
2116 end subroutine cpu_opgrad_lx17_single
2117
2118 subroutine cpu_opgrad_lx16_single(ux, uy, uz, u, dx, dy, dz, &
2119 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2120 integer, parameter :: lx = 16
2121 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2122 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2123 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2124 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2125 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2126 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2127 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2128 real(kind=rp) :: ur(lx, lx, lx)
2129 real(kind=rp) :: us(lx, lx, lx)
2130 real(kind=rp) :: ut(lx, lx, lx)
2131 integer :: i, j, k
2132
2133 do j = 1, lx * lx
2134 do i = 1, lx
2135 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2136 + dx(i,2) * u(2,j,1) &
2137 + dx(i,3) * u(3,j,1) &
2138 + dx(i,4) * u(4,j,1) &
2139 + dx(i,5) * u(5,j,1) &
2140 + dx(i,6) * u(6,j,1) &
2141 + dx(i,7) * u(7,j,1) &
2142 + dx(i,8) * u(8,j,1) &
2143 + dx(i,9) * u(9,j,1) &
2144 + dx(i,10) * u(10,j,1) &
2145 + dx(i,11) * u(11,j,1) &
2146 + dx(i,12) * u(12,j,1) &
2147 + dx(i,13) * u(13,j,1) &
2148 + dx(i,14) * u(14,j,1) &
2149 + dx(i,15) * u(15,j,1) &
2150 + dx(i,16) * u(16,j,1)
2151 end do
2152 end do
2153
2154 do k = 1, lx
2155 do j = 1, lx
2156 do i = 1, lx
2157 us(i,j,k) = dy(j,1) * u(i,1,k) &
2158 + dy(j,2) * u(i,2,k) &
2159 + dy(j,3) * u(i,3,k) &
2160 + dy(j,4) * u(i,4,k) &
2161 + dy(j,5) * u(i,5,k) &
2162 + dy(j,6) * u(i,6,k) &
2163 + dy(j,7) * u(i,7,k) &
2164 + dy(j,8) * u(i,8,k) &
2165 + dy(j,9) * u(i,9,k) &
2166 + dy(j,10) * u(i,10,k) &
2167 + dy(j,11) * u(i,11,k) &
2168 + dy(j,12) * u(i,12,k) &
2169 + dy(j,13) * u(i,13,k) &
2170 + dy(j,14) * u(i,14,k) &
2171 + dy(j,15) * u(i,15,k) &
2172 + dy(j,16) * u(i,16,k)
2173 end do
2174 end do
2175 end do
2176
2177 do k = 1, lx
2178 do i = 1, lx*lx
2179 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2180 + dz(k,2) * u(i,1,2) &
2181 + dz(k,3) * u(i,1,3) &
2182 + dz(k,4) * u(i,1,4) &
2183 + dz(k,5) * u(i,1,5) &
2184 + dz(k,6) * u(i,1,6) &
2185 + dz(k,7) * u(i,1,7) &
2186 + dz(k,8) * u(i,1,8) &
2187 + dz(k,9) * u(i,1,9) &
2188 + dz(k,10) * u(i,1,10) &
2189 + dz(k,11) * u(i,1,11) &
2190 + dz(k,12) * u(i,1,12) &
2191 + dz(k,13) * u(i,1,13) &
2192 + dz(k,14) * u(i,1,14) &
2193 + dz(k,15) * u(i,1,15) &
2194 + dz(k,16) * u(i,1,16)
2195 end do
2196 end do
2197
2198 do i = 1, lx * lx * lx
2199 ux(i,1,1) = w3(i,1,1) &
2200 * ( drdx(i,1,1) * ur(i,1,1) &
2201 + dsdx(i,1,1) * us(i,1,1) &
2202 + dtdx(i,1,1) * ut(i,1,1) )
2203 uy(i,1,1) = w3(i,1,1) &
2204 * ( dsdy(i,1,1) * us(i,1,1) &
2205 + drdy(i,1,1) * ur(i,1,1) &
2206 + dtdy(i,1,1) * ut(i,1,1) )
2207 uz(i,1,1) = w3(i,1,1) &
2208 * ( dtdz(i,1,1) * ut(i,1,1) &
2209 + drdz(i,1,1) * ur(i,1,1) &
2210 + dsdz(i,1,1) * us(i,1,1) )
2211 end do
2212 end subroutine cpu_opgrad_lx16_single
2213
2214 subroutine cpu_opgrad_lx15_single(ux, uy, uz, u, dx, dy, dz, &
2215 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2216 integer, parameter :: lx = 15
2217 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2218 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2219 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2220 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2221 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2222 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2223 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2224 real(kind=rp) :: ur(lx, lx, lx)
2225 real(kind=rp) :: us(lx, lx, lx)
2226 real(kind=rp) :: ut(lx, lx, lx)
2227 integer :: i, j, k
2228
2229 do j = 1, lx * lx
2230 do i = 1, lx
2231 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2232 + dx(i,2) * u(2,j,1) &
2233 + dx(i,3) * u(3,j,1) &
2234 + dx(i,4) * u(4,j,1) &
2235 + dx(i,5) * u(5,j,1) &
2236 + dx(i,6) * u(6,j,1) &
2237 + dx(i,7) * u(7,j,1) &
2238 + dx(i,8) * u(8,j,1) &
2239 + dx(i,9) * u(9,j,1) &
2240 + dx(i,10) * u(10,j,1) &
2241 + dx(i,11) * u(11,j,1) &
2242 + dx(i,12) * u(12,j,1) &
2243 + dx(i,13) * u(13,j,1) &
2244 + dx(i,14) * u(14,j,1) &
2245 + dx(i,15) * u(15,j,1)
2246 end do
2247 end do
2248
2249 do k = 1, lx
2250 do j = 1, lx
2251 do i = 1, lx
2252 us(i,j,k) = dy(j,1) * u(i,1,k) &
2253 + dy(j,2) * u(i,2,k) &
2254 + dy(j,3) * u(i,3,k) &
2255 + dy(j,4) * u(i,4,k) &
2256 + dy(j,5) * u(i,5,k) &
2257 + dy(j,6) * u(i,6,k) &
2258 + dy(j,7) * u(i,7,k) &
2259 + dy(j,8) * u(i,8,k) &
2260 + dy(j,9) * u(i,9,k) &
2261 + dy(j,10) * u(i,10,k) &
2262 + dy(j,11) * u(i,11,k) &
2263 + dy(j,12) * u(i,12,k) &
2264 + dy(j,13) * u(i,13,k) &
2265 + dy(j,14) * u(i,14,k) &
2266 + dy(j,15) * u(i,15,k)
2267 end do
2268 end do
2269 end do
2270
2271 do k = 1, lx
2272 do i = 1, lx*lx
2273 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2274 + dz(k,2) * u(i,1,2) &
2275 + dz(k,3) * u(i,1,3) &
2276 + dz(k,4) * u(i,1,4) &
2277 + dz(k,5) * u(i,1,5) &
2278 + dz(k,6) * u(i,1,6) &
2279 + dz(k,7) * u(i,1,7) &
2280 + dz(k,8) * u(i,1,8) &
2281 + dz(k,9) * u(i,1,9) &
2282 + dz(k,10) * u(i,1,10) &
2283 + dz(k,11) * u(i,1,11) &
2284 + dz(k,12) * u(i,1,12) &
2285 + dz(k,13) * u(i,1,13) &
2286 + dz(k,14) * u(i,1,14) &
2287 + dz(k,15) * u(i,1,15)
2288 end do
2289 end do
2290
2291 do i = 1, lx * lx * lx
2292 ux(i,1,1) = w3(i,1,1) &
2293 * ( drdx(i,1,1) * ur(i,1,1) &
2294 + dsdx(i,1,1) * us(i,1,1) &
2295 + dtdx(i,1,1) * ut(i,1,1) )
2296 uy(i,1,1) = w3(i,1,1) &
2297 * ( dsdy(i,1,1) * us(i,1,1) &
2298 + drdy(i,1,1) * ur(i,1,1) &
2299 + dtdy(i,1,1) * ut(i,1,1) )
2300 uz(i,1,1) = w3(i,1,1) &
2301 * ( dtdz(i,1,1) * ut(i,1,1) &
2302 + drdz(i,1,1) * ur(i,1,1) &
2303 + dsdz(i,1,1) * us(i,1,1) )
2304 end do
2305 end subroutine cpu_opgrad_lx15_single
2306
2307 subroutine cpu_opgrad_lx14_single(ux, uy, uz, u, dx, dy, dz, &
2308 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2309 integer, parameter :: lx = 14
2310 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2311 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2312 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2313 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2314 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2315 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2316 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2317 real(kind=rp) :: ur(lx, lx, lx)
2318 real(kind=rp) :: us(lx, lx, lx)
2319 real(kind=rp) :: ut(lx, lx, lx)
2320 integer :: i, j, k
2321
2322 do j = 1, lx * lx
2323 do i = 1, lx
2324 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2325 + dx(i,2) * u(2,j,1) &
2326 + dx(i,3) * u(3,j,1) &
2327 + dx(i,4) * u(4,j,1) &
2328 + dx(i,5) * u(5,j,1) &
2329 + dx(i,6) * u(6,j,1) &
2330 + dx(i,7) * u(7,j,1) &
2331 + dx(i,8) * u(8,j,1) &
2332 + dx(i,9) * u(9,j,1) &
2333 + dx(i,10) * u(10,j,1) &
2334 + dx(i,11) * u(11,j,1) &
2335 + dx(i,12) * u(12,j,1) &
2336 + dx(i,13) * u(13,j,1) &
2337 + dx(i,14) * u(14,j,1)
2338 end do
2339 end do
2340
2341 do k = 1, lx
2342 do j = 1, lx
2343 do i = 1, lx
2344 us(i,j,k) = dy(j,1) * u(i,1,k) &
2345 + dy(j,2) * u(i,2,k) &
2346 + dy(j,3) * u(i,3,k) &
2347 + dy(j,4) * u(i,4,k) &
2348 + dy(j,5) * u(i,5,k) &
2349 + dy(j,6) * u(i,6,k) &
2350 + dy(j,7) * u(i,7,k) &
2351 + dy(j,8) * u(i,8,k) &
2352 + dy(j,9) * u(i,9,k) &
2353 + dy(j,10) * u(i,10,k) &
2354 + dy(j,11) * u(i,11,k) &
2355 + dy(j,12) * u(i,12,k) &
2356 + dy(j,13) * u(i,13,k) &
2357 + dy(j,14) * u(i,14,k)
2358 end do
2359 end do
2360 end do
2361
2362 do k = 1, lx
2363 do i = 1, lx*lx
2364 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2365 + dz(k,2) * u(i,1,2) &
2366 + dz(k,3) * u(i,1,3) &
2367 + dz(k,4) * u(i,1,4) &
2368 + dz(k,5) * u(i,1,5) &
2369 + dz(k,6) * u(i,1,6) &
2370 + dz(k,7) * u(i,1,7) &
2371 + dz(k,8) * u(i,1,8) &
2372 + dz(k,9) * u(i,1,9) &
2373 + dz(k,10) * u(i,1,10) &
2374 + dz(k,11) * u(i,1,11) &
2375 + dz(k,12) * u(i,1,12) &
2376 + dz(k,13) * u(i,1,13) &
2377 + dz(k,14) * u(i,1,14)
2378 end do
2379 end do
2380
2381 do i = 1, lx * lx * lx
2382 ux(i,1,1) = w3(i,1,1) &
2383 * ( drdx(i,1,1) * ur(i,1,1) &
2384 + dsdx(i,1,1) * us(i,1,1) &
2385 + dtdx(i,1,1) * ut(i,1,1) )
2386 uy(i,1,1) = w3(i,1,1) &
2387 * ( dsdy(i,1,1) * us(i,1,1) &
2388 + drdy(i,1,1) * ur(i,1,1) &
2389 + dtdy(i,1,1) * ut(i,1,1) )
2390 uz(i,1,1) = w3(i,1,1) &
2391 * ( dtdz(i,1,1) * ut(i,1,1) &
2392 + drdz(i,1,1) * ur(i,1,1) &
2393 + dsdz(i,1,1) * us(i,1,1) )
2394 end do
2395 end subroutine cpu_opgrad_lx14_single
2396
2397 subroutine cpu_opgrad_lx13_single(ux, uy, uz, u, dx, dy, dz, &
2398 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2399 integer, parameter :: lx = 13
2400 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2401 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2402 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2403 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2404 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2405 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2406 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2407 real(kind=rp) :: ur(lx, lx, lx)
2408 real(kind=rp) :: us(lx, lx, lx)
2409 real(kind=rp) :: ut(lx, lx, lx)
2410 integer :: i, j, k
2411
2412 do j = 1, lx * lx
2413 do i = 1, lx
2414 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2415 + dx(i,2) * u(2,j,1) &
2416 + dx(i,3) * u(3,j,1) &
2417 + dx(i,4) * u(4,j,1) &
2418 + dx(i,5) * u(5,j,1) &
2419 + dx(i,6) * u(6,j,1) &
2420 + dx(i,7) * u(7,j,1) &
2421 + dx(i,8) * u(8,j,1) &
2422 + dx(i,9) * u(9,j,1) &
2423 + dx(i,10) * u(10,j,1) &
2424 + dx(i,11) * u(11,j,1) &
2425 + dx(i,12) * u(12,j,1) &
2426 + dx(i,13) * u(13,j,1)
2427 end do
2428 end do
2429
2430 do k = 1, lx
2431 do j = 1, lx
2432 do i = 1, lx
2433 us(i,j,k) = dy(j,1) * u(i,1,k) &
2434 + dy(j,2) * u(i,2,k) &
2435 + dy(j,3) * u(i,3,k) &
2436 + dy(j,4) * u(i,4,k) &
2437 + dy(j,5) * u(i,5,k) &
2438 + dy(j,6) * u(i,6,k) &
2439 + dy(j,7) * u(i,7,k) &
2440 + dy(j,8) * u(i,8,k) &
2441 + dy(j,9) * u(i,9,k) &
2442 + dy(j,10) * u(i,10,k) &
2443 + dy(j,11) * u(i,11,k) &
2444 + dy(j,12) * u(i,12,k) &
2445 + dy(j,13) * u(i,13,k)
2446 end do
2447 end do
2448 end do
2449
2450 do k = 1, lx
2451 do i = 1, lx*lx
2452 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2453 + dz(k,2) * u(i,1,2) &
2454 + dz(k,3) * u(i,1,3) &
2455 + dz(k,4) * u(i,1,4) &
2456 + dz(k,5) * u(i,1,5) &
2457 + dz(k,6) * u(i,1,6) &
2458 + dz(k,7) * u(i,1,7) &
2459 + dz(k,8) * u(i,1,8) &
2460 + dz(k,9) * u(i,1,9) &
2461 + dz(k,10) * u(i,1,10) &
2462 + dz(k,11) * u(i,1,11) &
2463 + dz(k,12) * u(i,1,12) &
2464 + dz(k,13) * u(i,1,13)
2465 end do
2466 end do
2467
2468 do i = 1, lx * lx * lx
2469 ux(i,1,1) = w3(i,1,1) &
2470 * ( drdx(i,1,1) * ur(i,1,1) &
2471 + dsdx(i,1,1) * us(i,1,1) &
2472 + dtdx(i,1,1) * ut(i,1,1) )
2473 uy(i,1,1) = w3(i,1,1) &
2474 * ( dsdy(i,1,1) * us(i,1,1) &
2475 + drdy(i,1,1) * ur(i,1,1) &
2476 + dtdy(i,1,1) * ut(i,1,1) )
2477 uz(i,1,1) = w3(i,1,1) &
2478 * ( dtdz(i,1,1) * ut(i,1,1) &
2479 + drdz(i,1,1) * ur(i,1,1) &
2480 + dsdz(i,1,1) * us(i,1,1) )
2481 end do
2482 end subroutine cpu_opgrad_lx13_single
2483
2484 subroutine cpu_opgrad_lx12_single(ux, uy, uz, u, dx, dy, dz, &
2485 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2486 integer, parameter :: lx = 12
2487 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2488 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2489 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2490 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2491 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2492 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2493 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2494 real(kind=rp) :: ur(lx, lx, lx)
2495 real(kind=rp) :: us(lx, lx, lx)
2496 real(kind=rp) :: ut(lx, lx, lx)
2497 integer :: i, j, k
2498
2499 do j = 1, lx * lx
2500 do i = 1, lx
2501 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2502 + dx(i,2) * u(2,j,1) &
2503 + dx(i,3) * u(3,j,1) &
2504 + dx(i,4) * u(4,j,1) &
2505 + dx(i,5) * u(5,j,1) &
2506 + dx(i,6) * u(6,j,1) &
2507 + dx(i,7) * u(7,j,1) &
2508 + dx(i,8) * u(8,j,1) &
2509 + dx(i,9) * u(9,j,1) &
2510 + dx(i,10) * u(10,j,1) &
2511 + dx(i,11) * u(11,j,1) &
2512 + dx(i,12) * u(12,j,1)
2513 end do
2514 end do
2515
2516 do k = 1, lx
2517 do j = 1, lx
2518 do i = 1, lx
2519 us(i,j,k) = dy(j,1) * u(i,1,k) &
2520 + dy(j,2) * u(i,2,k) &
2521 + dy(j,3) * u(i,3,k) &
2522 + dy(j,4) * u(i,4,k) &
2523 + dy(j,5) * u(i,5,k) &
2524 + dy(j,6) * u(i,6,k) &
2525 + dy(j,7) * u(i,7,k) &
2526 + dy(j,8) * u(i,8,k) &
2527 + dy(j,9) * u(i,9,k) &
2528 + dy(j,10) * u(i,10,k) &
2529 + dy(j,11) * u(i,11,k) &
2530 + dy(j,12) * u(i,12,k)
2531 end do
2532 end do
2533 end do
2534
2535 do k = 1, lx
2536 do i = 1, lx*lx
2537 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2538 + dz(k,2) * u(i,1,2) &
2539 + dz(k,3) * u(i,1,3) &
2540 + dz(k,4) * u(i,1,4) &
2541 + dz(k,5) * u(i,1,5) &
2542 + dz(k,6) * u(i,1,6) &
2543 + dz(k,7) * u(i,1,7) &
2544 + dz(k,8) * u(i,1,8) &
2545 + dz(k,9) * u(i,1,9) &
2546 + dz(k,10) * u(i,1,10) &
2547 + dz(k,11) * u(i,1,11) &
2548 + dz(k,12) * u(i,1,12)
2549 end do
2550 end do
2551
2552 do i = 1, lx * lx * lx
2553 ux(i,1,1) = w3(i,1,1) &
2554 * ( drdx(i,1,1) * ur(i,1,1) &
2555 + dsdx(i,1,1) * us(i,1,1) &
2556 + dtdx(i,1,1) * ut(i,1,1) )
2557 uy(i,1,1) = w3(i,1,1) &
2558 * ( dsdy(i,1,1) * us(i,1,1) &
2559 + drdy(i,1,1) * ur(i,1,1) &
2560 + dtdy(i,1,1) * ut(i,1,1) )
2561 uz(i,1,1) = w3(i,1,1) &
2562 * ( dtdz(i,1,1) * ut(i,1,1) &
2563 + drdz(i,1,1) * ur(i,1,1) &
2564 + dsdz(i,1,1) * us(i,1,1) )
2565 end do
2566 end subroutine cpu_opgrad_lx12_single
2567
2568 subroutine cpu_opgrad_lx11_single(ux, uy, uz, u, dx, dy, dz, &
2569 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2570 integer, parameter :: lx = 11
2571 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2572 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2573 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2574 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2575 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2576 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2577 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2578 real(kind=rp) :: ur(lx, lx, lx)
2579 real(kind=rp) :: us(lx, lx, lx)
2580 real(kind=rp) :: ut(lx, lx, lx)
2581 integer :: i, j, k
2582
2583 do j = 1, lx * lx
2584 do i = 1, lx
2585 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2586 + dx(i,2) * u(2,j,1) &
2587 + dx(i,3) * u(3,j,1) &
2588 + dx(i,4) * u(4,j,1) &
2589 + dx(i,5) * u(5,j,1) &
2590 + dx(i,6) * u(6,j,1) &
2591 + dx(i,7) * u(7,j,1) &
2592 + dx(i,8) * u(8,j,1) &
2593 + dx(i,9) * u(9,j,1) &
2594 + dx(i,10) * u(10,j,1) &
2595 + dx(i,11) * u(11,j,1)
2596 end do
2597 end do
2598
2599 do k = 1, lx
2600 do j = 1, lx
2601 do i = 1, lx
2602 us(i,j,k) = dy(j,1) * u(i,1,k) &
2603 + dy(j,2) * u(i,2,k) &
2604 + dy(j,3) * u(i,3,k) &
2605 + dy(j,4) * u(i,4,k) &
2606 + dy(j,5) * u(i,5,k) &
2607 + dy(j,6) * u(i,6,k) &
2608 + dy(j,7) * u(i,7,k) &
2609 + dy(j,8) * u(i,8,k) &
2610 + dy(j,9) * u(i,9,k) &
2611 + dy(j,10) * u(i,10,k) &
2612 + dy(j,11) * u(i,11,k)
2613 end do
2614 end do
2615 end do
2616
2617 do k = 1, lx
2618 do i = 1, lx*lx
2619 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2620 + dz(k,2) * u(i,1,2) &
2621 + dz(k,3) * u(i,1,3) &
2622 + dz(k,4) * u(i,1,4) &
2623 + dz(k,5) * u(i,1,5) &
2624 + dz(k,6) * u(i,1,6) &
2625 + dz(k,7) * u(i,1,7) &
2626 + dz(k,8) * u(i,1,8) &
2627 + dz(k,9) * u(i,1,9) &
2628 + dz(k,10) * u(i,1,10) &
2629 + dz(k,11) * u(i,1,11)
2630 end do
2631 end do
2632
2633 do i = 1, lx * lx * lx
2634 ux(i,1,1) = w3(i,1,1) &
2635 * ( drdx(i,1,1) * ur(i,1,1) &
2636 + dsdx(i,1,1) * us(i,1,1) &
2637 + dtdx(i,1,1) * ut(i,1,1) )
2638 uy(i,1,1) = w3(i,1,1) &
2639 * ( dsdy(i,1,1) * us(i,1,1) &
2640 + drdy(i,1,1) * ur(i,1,1) &
2641 + dtdy(i,1,1) * ut(i,1,1) )
2642 uz(i,1,1) = w3(i,1,1) &
2643 * ( dtdz(i,1,1) * ut(i,1,1) &
2644 + drdz(i,1,1) * ur(i,1,1) &
2645 + dsdz(i,1,1) * us(i,1,1) )
2646 end do
2647 end subroutine cpu_opgrad_lx11_single
2648
2649 subroutine cpu_opgrad_lx10_single(ux, uy, uz, u, dx, dy, dz, &
2650 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2651 integer, parameter :: lx = 10
2652 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2653 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2654 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2655 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2656 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2657 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2658 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2659 real(kind=rp) :: ur(lx, lx, lx)
2660 real(kind=rp) :: us(lx, lx, lx)
2661 real(kind=rp) :: ut(lx, lx, lx)
2662 integer :: i, j, k
2663
2664 do j = 1, lx * lx
2665 do i = 1, lx
2666 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2667 + dx(i,2) * u(2,j,1) &
2668 + dx(i,3) * u(3,j,1) &
2669 + dx(i,4) * u(4,j,1) &
2670 + dx(i,5) * u(5,j,1) &
2671 + dx(i,6) * u(6,j,1) &
2672 + dx(i,7) * u(7,j,1) &
2673 + dx(i,8) * u(8,j,1) &
2674 + dx(i,9) * u(9,j,1) &
2675 + dx(i,10) * u(10,j,1)
2676 end do
2677 end do
2678
2679 do k = 1, lx
2680 do j = 1, lx
2681 do i = 1, lx
2682 us(i,j,k) = dy(j,1) * u(i,1,k) &
2683 + dy(j,2) * u(i,2,k) &
2684 + dy(j,3) * u(i,3,k) &
2685 + dy(j,4) * u(i,4,k) &
2686 + dy(j,5) * u(i,5,k) &
2687 + dy(j,6) * u(i,6,k) &
2688 + dy(j,7) * u(i,7,k) &
2689 + dy(j,8) * u(i,8,k) &
2690 + dy(j,9) * u(i,9,k) &
2691 + dy(j,10) * u(i,10,k)
2692 end do
2693 end do
2694 end do
2695
2696 do k = 1, lx
2697 do i = 1, lx*lx
2698 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2699 + dz(k,2) * u(i,1,2) &
2700 + dz(k,3) * u(i,1,3) &
2701 + dz(k,4) * u(i,1,4) &
2702 + dz(k,5) * u(i,1,5) &
2703 + dz(k,6) * u(i,1,6) &
2704 + dz(k,7) * u(i,1,7) &
2705 + dz(k,8) * u(i,1,8) &
2706 + dz(k,9) * u(i,1,9) &
2707 + dz(k,10) * u(i,1,10)
2708 end do
2709 end do
2710
2711 do i = 1, lx * lx * lx
2712 ux(i,1,1) = w3(i,1,1) &
2713 * ( drdx(i,1,1) * ur(i,1,1) &
2714 + dsdx(i,1,1) * us(i,1,1) &
2715 + dtdx(i,1,1) * ut(i,1,1) )
2716 uy(i,1,1) = w3(i,1,1) &
2717 * ( dsdy(i,1,1) * us(i,1,1) &
2718 + drdy(i,1,1) * ur(i,1,1) &
2719 + dtdy(i,1,1) * ut(i,1,1) )
2720 uz(i,1,1) = w3(i,1,1) &
2721 * ( dtdz(i,1,1) * ut(i,1,1) &
2722 + drdz(i,1,1) * ur(i,1,1) &
2723 + dsdz(i,1,1) * us(i,1,1) )
2724 end do
2725 end subroutine cpu_opgrad_lx10_single
2726
2727 subroutine cpu_opgrad_lx9_single(ux, uy, uz, u, dx, dy, dz, &
2728 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2729 integer, parameter :: lx = 9
2730 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2731 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2732 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2733 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2734 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2735 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2736 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2737 real(kind=rp) :: ur(lx, lx, lx)
2738 real(kind=rp) :: us(lx, lx, lx)
2739 real(kind=rp) :: ut(lx, lx, lx)
2740 integer :: i, j, k
2741
2742 do j = 1, lx * lx
2743 do i = 1, lx
2744 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2745 + dx(i,2) * u(2,j,1) &
2746 + dx(i,3) * u(3,j,1) &
2747 + dx(i,4) * u(4,j,1) &
2748 + dx(i,5) * u(5,j,1) &
2749 + dx(i,6) * u(6,j,1) &
2750 + dx(i,7) * u(7,j,1) &
2751 + dx(i,8) * u(8,j,1) &
2752 + dx(i,9) * u(9,j,1)
2753 end do
2754 end do
2755
2756 do k = 1, lx
2757 do j = 1, lx
2758 do i = 1, lx
2759 us(i,j,k) = dy(j,1) * u(i,1,k) &
2760 + dy(j,2) * u(i,2,k) &
2761 + dy(j,3) * u(i,3,k) &
2762 + dy(j,4) * u(i,4,k) &
2763 + dy(j,5) * u(i,5,k) &
2764 + dy(j,6) * u(i,6,k) &
2765 + dy(j,7) * u(i,7,k) &
2766 + dy(j,8) * u(i,8,k) &
2767 + dy(j,9) * u(i,9,k)
2768 end do
2769 end do
2770 end do
2771
2772 do k = 1, lx
2773 do i = 1, lx*lx
2774 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2775 + dz(k,2) * u(i,1,2) &
2776 + dz(k,3) * u(i,1,3) &
2777 + dz(k,4) * u(i,1,4) &
2778 + dz(k,5) * u(i,1,5) &
2779 + dz(k,6) * u(i,1,6) &
2780 + dz(k,7) * u(i,1,7) &
2781 + dz(k,8) * u(i,1,8) &
2782 + dz(k,9) * u(i,1,9)
2783 end do
2784 end do
2785
2786 do i = 1, lx * lx * lx
2787 ux(i,1,1) = w3(i,1,1) &
2788 * ( drdx(i,1,1) * ur(i,1,1) &
2789 + dsdx(i,1,1) * us(i,1,1) &
2790 + dtdx(i,1,1) * ut(i,1,1) )
2791 uy(i,1,1) = w3(i,1,1) &
2792 * ( dsdy(i,1,1) * us(i,1,1) &
2793 + drdy(i,1,1) * ur(i,1,1) &
2794 + dtdy(i,1,1) * ut(i,1,1) )
2795 uz(i,1,1) = w3(i,1,1) &
2796 * ( dtdz(i,1,1) * ut(i,1,1) &
2797 + drdz(i,1,1) * ur(i,1,1) &
2798 + dsdz(i,1,1) * us(i,1,1) )
2799 end do
2800 end subroutine cpu_opgrad_lx9_single
2801
2802 subroutine cpu_opgrad_lx8_single(ux, uy, uz, u, dx, dy, dz, &
2803 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2804 integer, parameter :: lx = 8
2805 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2806 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2807 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2808 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2809 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2810 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2811 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2812 real(kind=rp) :: ur(lx, lx, lx)
2813 real(kind=rp) :: us(lx, lx, lx)
2814 real(kind=rp) :: ut(lx, lx, lx)
2815 integer :: i, j, k
2816
2817 do j = 1, lx * lx
2818 do i = 1, lx
2819 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2820 + dx(i,2) * u(2,j,1) &
2821 + dx(i,3) * u(3,j,1) &
2822 + dx(i,4) * u(4,j,1) &
2823 + dx(i,5) * u(5,j,1) &
2824 + dx(i,6) * u(6,j,1) &
2825 + dx(i,7) * u(7,j,1) &
2826 + dx(i,8) * u(8,j,1)
2827 end do
2828 end do
2829
2830 do k = 1, lx
2831 do j = 1, lx
2832 do i = 1, lx
2833 us(i,j,k) = dy(j,1) * u(i,1,k) &
2834 + dy(j,2) * u(i,2,k) &
2835 + dy(j,3) * u(i,3,k) &
2836 + dy(j,4) * u(i,4,k) &
2837 + dy(j,5) * u(i,5,k) &
2838 + dy(j,6) * u(i,6,k) &
2839 + dy(j,7) * u(i,7,k) &
2840 + dy(j,8) * u(i,8,k)
2841 end do
2842 end do
2843 end do
2844
2845 do k = 1, lx
2846 do i = 1, lx*lx
2847 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2848 + dz(k,2) * u(i,1,2) &
2849 + dz(k,3) * u(i,1,3) &
2850 + dz(k,4) * u(i,1,4) &
2851 + dz(k,5) * u(i,1,5) &
2852 + dz(k,6) * u(i,1,6) &
2853 + dz(k,7) * u(i,1,7) &
2854 + dz(k,8) * u(i,1,8)
2855 end do
2856 end do
2857
2858 do i = 1, lx * lx * lx
2859 ux(i,1,1) = w3(i,1,1) &
2860 * ( drdx(i,1,1) * ur(i,1,1) &
2861 + dsdx(i,1,1) * us(i,1,1) &
2862 + dtdx(i,1,1) * ut(i,1,1) )
2863 uy(i,1,1) = w3(i,1,1) &
2864 * ( dsdy(i,1,1) * us(i,1,1) &
2865 + drdy(i,1,1) * ur(i,1,1) &
2866 + dtdy(i,1,1) * ut(i,1,1) )
2867 uz(i,1,1) = w3(i,1,1) &
2868 * ( dtdz(i,1,1) * ut(i,1,1) &
2869 + drdz(i,1,1) * ur(i,1,1) &
2870 + dsdz(i,1,1) * us(i,1,1) )
2871 end do
2872 end subroutine cpu_opgrad_lx8_single
2873
2874 subroutine cpu_opgrad_lx7_single(ux, uy, uz, u, dx, dy, dz, &
2875 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2876 integer, parameter :: lx = 7
2877 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2878 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2879 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2880 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2881 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2882 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2883 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2884 real(kind=rp) :: ur(lx, lx, lx)
2885 real(kind=rp) :: us(lx, lx, lx)
2886 real(kind=rp) :: ut(lx, lx, lx)
2887 integer :: i, j, k
2888
2889 do j = 1, lx * lx
2890 do i = 1, lx
2891 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2892 + dx(i,2) * u(2,j,1) &
2893 + dx(i,3) * u(3,j,1) &
2894 + dx(i,4) * u(4,j,1) &
2895 + dx(i,5) * u(5,j,1) &
2896 + dx(i,6) * u(6,j,1) &
2897 + dx(i,7) * u(7,j,1)
2898 end do
2899 end do
2900
2901 do k = 1, lx
2902 do j = 1, lx
2903 do i = 1, lx
2904 us(i,j,k) = dy(j,1) * u(i,1,k) &
2905 + dy(j,2) * u(i,2,k) &
2906 + dy(j,3) * u(i,3,k) &
2907 + dy(j,4) * u(i,4,k) &
2908 + dy(j,5) * u(i,5,k) &
2909 + dy(j,6) * u(i,6,k) &
2910 + dy(j,7) * u(i,7,k)
2911 end do
2912 end do
2913 end do
2914
2915 do k = 1, lx
2916 do i = 1, lx*lx
2917 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2918 + dz(k,2) * u(i,1,2) &
2919 + dz(k,3) * u(i,1,3) &
2920 + dz(k,4) * u(i,1,4) &
2921 + dz(k,5) * u(i,1,5) &
2922 + dz(k,6) * u(i,1,6) &
2923 + dz(k,7) * u(i,1,7)
2924 end do
2925 end do
2926
2927 do i = 1, lx * lx * lx
2928 ux(i,1,1) = w3(i,1,1) &
2929 * ( drdx(i,1,1) * ur(i,1,1) &
2930 + dsdx(i,1,1) * us(i,1,1) &
2931 + dtdx(i,1,1) * ut(i,1,1) )
2932 uy(i,1,1) = w3(i,1,1) &
2933 * ( dsdy(i,1,1) * us(i,1,1) &
2934 + drdy(i,1,1) * ur(i,1,1) &
2935 + dtdy(i,1,1) * ut(i,1,1) )
2936 uz(i,1,1) = w3(i,1,1) &
2937 * ( dtdz(i,1,1) * ut(i,1,1) &
2938 + drdz(i,1,1) * ur(i,1,1) &
2939 + dsdz(i,1,1) * us(i,1,1) )
2940 end do
2941 end subroutine cpu_opgrad_lx7_single
2942
2943 subroutine cpu_opgrad_lx6_single(ux, uy, uz, u, dx, dy, dz, &
2944 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
2945 integer, parameter :: lx = 6
2946 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
2947 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
2948 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
2949 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
2950 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
2951 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
2952 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
2953 real(kind=rp) :: ur(lx, lx, lx)
2954 real(kind=rp) :: us(lx, lx, lx)
2955 real(kind=rp) :: ut(lx, lx, lx)
2956 integer :: i, j, k
2957
2958 do j = 1, lx * lx
2959 do i = 1, lx
2960 ur(i,j,1) = dx(i,1) * u(1,j,1) &
2961 + dx(i,2) * u(2,j,1) &
2962 + dx(i,3) * u(3,j,1) &
2963 + dx(i,4) * u(4,j,1) &
2964 + dx(i,5) * u(5,j,1) &
2965 + dx(i,6) * u(6,j,1)
2966 end do
2967 end do
2968
2969 do k = 1, lx
2970 do j = 1, lx
2971 do i = 1, lx
2972 us(i,j,k) = dy(j,1) * u(i,1,k) &
2973 + dy(j,2) * u(i,2,k) &
2974 + dy(j,3) * u(i,3,k) &
2975 + dy(j,4) * u(i,4,k) &
2976 + dy(j,5) * u(i,5,k) &
2977 + dy(j,6) * u(i,6,k)
2978 end do
2979 end do
2980 end do
2981
2982 do k = 1, lx
2983 do i = 1, lx*lx
2984 ut(i,1,k) = dz(k,1) * u(i,1,1) &
2985 + dz(k,2) * u(i,1,2) &
2986 + dz(k,3) * u(i,1,3) &
2987 + dz(k,4) * u(i,1,4) &
2988 + dz(k,5) * u(i,1,5) &
2989 + dz(k,6) * u(i,1,6)
2990 end do
2991 end do
2992
2993 do i = 1, lx * lx * lx
2994 ux(i,1,1) = w3(i,1,1) &
2995 * ( drdx(i,1,1) * ur(i,1,1) &
2996 + dsdx(i,1,1) * us(i,1,1) &
2997 + dtdx(i,1,1) * ut(i,1,1) )
2998 uy(i,1,1) = w3(i,1,1) &
2999 * ( dsdy(i,1,1) * us(i,1,1) &
3000 + drdy(i,1,1) * ur(i,1,1) &
3001 + dtdy(i,1,1) * ut(i,1,1) )
3002 uz(i,1,1) = w3(i,1,1) &
3003 * ( dtdz(i,1,1) * ut(i,1,1) &
3004 + drdz(i,1,1) * ur(i,1,1) &
3005 + dsdz(i,1,1) * us(i,1,1) )
3006 end do
3007 end subroutine cpu_opgrad_lx6_single
3008
3009 subroutine cpu_opgrad_lx5_single(ux, uy, uz, u, dx, dy, dz, &
3010 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
3011 integer, parameter :: lx = 5
3012 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
3013 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
3014 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
3015 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
3016 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
3017 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
3018 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
3019 real(kind=rp) :: ur(lx, lx, lx)
3020 real(kind=rp) :: us(lx, lx, lx)
3021 real(kind=rp) :: ut(lx, lx, lx)
3022 integer :: i, j, k
3023
3024 do j = 1, lx * lx
3025 do i = 1, lx
3026 ur(i,j,1) = dx(i,1) * u(1,j,1) &
3027 + dx(i,2) * u(2,j,1) &
3028 + dx(i,3) * u(3,j,1) &
3029 + dx(i,4) * u(4,j,1) &
3030 + dx(i,5) * u(5,j,1)
3031 end do
3032 end do
3033
3034 do k = 1, lx
3035 do j = 1, lx
3036 do i = 1, lx
3037 us(i,j,k) = dy(j,1) * u(i,1,k) &
3038 + dy(j,2) * u(i,2,k) &
3039 + dy(j,3) * u(i,3,k) &
3040 + dy(j,4) * u(i,4,k) &
3041 + dy(j,5) * u(i,5,k)
3042 end do
3043 end do
3044 end do
3045
3046 do k = 1, lx
3047 do i = 1, lx*lx
3048 ut(i,1,k) = dz(k,1) * u(i,1,1) &
3049 + dz(k,2) * u(i,1,2) &
3050 + dz(k,3) * u(i,1,3) &
3051 + dz(k,4) * u(i,1,4) &
3052 + dz(k,5) * u(i,1,5)
3053 end do
3054 end do
3055
3056 do i = 1, lx * lx * lx
3057 ux(i,1,1) = w3(i,1,1) &
3058 * ( drdx(i,1,1) * ur(i,1,1) &
3059 + dsdx(i,1,1) * us(i,1,1) &
3060 + dtdx(i,1,1) * ut(i,1,1) )
3061 uy(i,1,1) = w3(i,1,1) &
3062 * ( dsdy(i,1,1) * us(i,1,1) &
3063 + drdy(i,1,1) * ur(i,1,1) &
3064 + dtdy(i,1,1) * ut(i,1,1) )
3065 uz(i,1,1) = w3(i,1,1) &
3066 * ( dtdz(i,1,1) * ut(i,1,1) &
3067 + drdz(i,1,1) * ur(i,1,1) &
3068 + dsdz(i,1,1) * us(i,1,1) )
3069 end do
3070 end subroutine cpu_opgrad_lx5_single
3071
3072 subroutine cpu_opgrad_lx4_single(ux, uy, uz, u, dx, dy, dz, &
3073 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
3074 integer, parameter :: lx = 4
3075 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
3076 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
3077 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
3078 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
3079 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
3080 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
3081 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
3082 real(kind=rp) :: ur(lx, lx, lx)
3083 real(kind=rp) :: us(lx, lx, lx)
3084 real(kind=rp) :: ut(lx, lx, lx)
3085 integer :: i, j, k
3086
3087 do j = 1, lx * lx
3088 do i = 1, lx
3089 ur(i,j,1) = dx(i,1) * u(1,j,1) &
3090 + dx(i,2) * u(2,j,1) &
3091 + dx(i,3) * u(3,j,1) &
3092 + dx(i,4) * u(4,j,1)
3093 end do
3094 end do
3095
3096 do k = 1, lx
3097 do j = 1, lx
3098 do i = 1, lx
3099 us(i,j,k) = dy(j,1) * u(i,1,k) &
3100 + dy(j,2) * u(i,2,k) &
3101 + dy(j,3) * u(i,3,k) &
3102 + dy(j,4) * u(i,4,k)
3103 end do
3104 end do
3105 end do
3106
3107 do k = 1, lx
3108 do i = 1, lx*lx
3109 ut(i,1,k) = dz(k,1) * u(i,1,1) &
3110 + dz(k,2) * u(i,1,2) &
3111 + dz(k,3) * u(i,1,3) &
3112 + dz(k,4) * u(i,1,4)
3113 end do
3114 end do
3115
3116 do i = 1, lx * lx * lx
3117 ux(i,1,1) = w3(i,1,1) &
3118 * ( drdx(i,1,1) * ur(i,1,1) &
3119 + dsdx(i,1,1) * us(i,1,1) &
3120 + dtdx(i,1,1) * ut(i,1,1) )
3121 uy(i,1,1) = w3(i,1,1) &
3122 * ( dsdy(i,1,1) * us(i,1,1) &
3123 + drdy(i,1,1) * ur(i,1,1) &
3124 + dtdy(i,1,1) * ut(i,1,1) )
3125 uz(i,1,1) = w3(i,1,1) &
3126 * ( dtdz(i,1,1) * ut(i,1,1) &
3127 + drdz(i,1,1) * ur(i,1,1) &
3128 + dsdz(i,1,1) * us(i,1,1) )
3129 end do
3130 end subroutine cpu_opgrad_lx4_single
3131
3132 subroutine cpu_opgrad_lx3_single(ux, uy, uz, u, dx, dy, dz, &
3133 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
3134 integer, parameter :: lx = 3
3135 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
3136 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
3137 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
3138 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
3139 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
3140 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
3141 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
3142 real(kind=rp) :: ur(lx, lx, lx)
3143 real(kind=rp) :: us(lx, lx, lx)
3144 real(kind=rp) :: ut(lx, lx, lx)
3145 integer :: i, j, k
3146
3147 do j = 1, lx * lx
3148 do i = 1, lx
3149 ur(i,j,1) = dx(i,1) * u(1,j,1) &
3150 + dx(i,2) * u(2,j,1) &
3151 + dx(i,3) * u(3,j,1)
3152 end do
3153 end do
3154
3155 do k = 1, lx
3156 do j = 1, lx
3157 do i = 1, lx
3158 us(i,j,k) = dy(j,1) * u(i,1,k) &
3159 + dy(j,2) * u(i,2,k) &
3160 + dy(j,3) * u(i,3,k)
3161 end do
3162 end do
3163 end do
3164
3165 do k = 1, lx
3166 do i = 1, lx*lx
3167 ut(i,1,k) = dz(k,1) * u(i,1,1) &
3168 + dz(k,2) * u(i,1,2) &
3169 + dz(k,3) * u(i,1,3)
3170 end do
3171 end do
3172
3173 do i = 1, lx * lx * lx
3174 ux(i,1,1) = w3(i,1,1) &
3175 * ( drdx(i,1,1) * ur(i,1,1) &
3176 + dsdx(i,1,1) * us(i,1,1) &
3177 + dtdx(i,1,1) * ut(i,1,1) )
3178 uy(i,1,1) = w3(i,1,1) &
3179 * ( dsdy(i,1,1) * us(i,1,1) &
3180 + drdy(i,1,1) * ur(i,1,1) &
3181 + dtdy(i,1,1) * ut(i,1,1) )
3182 uz(i,1,1) = w3(i,1,1) &
3183 * ( dtdz(i,1,1) * ut(i,1,1) &
3184 + drdz(i,1,1) * ur(i,1,1) &
3185 + dsdz(i,1,1) * us(i,1,1) )
3186 end do
3187 end subroutine cpu_opgrad_lx3_single
3188
3189 subroutine cpu_opgrad_lx2_single(ux, uy, uz, u, dx, dy, dz, &
3190 drdx, dsdx, dtdx, drdy, dsdy, dtdy, drdz, dsdz, dtdz, w3)
3191 integer, parameter :: lx = 2
3192 real(kind=rp), dimension(lx, lx, lx), intent(inout) :: ux, uy, uz
3193 real(kind=rp), dimension(lx, lx, lx), intent(in) :: u
3194 real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
3195 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdx, dsdx, dtdx
3196 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdy, dsdy, dtdy
3197 real(kind=rp), dimension(lx, lx, lx), intent(in) :: drdz, dsdz, dtdz
3198 real(kind=rp), dimension(lx, lx, lx), intent(in) :: w3
3199 real(kind=rp) :: ur(lx, lx, lx)
3200 real(kind=rp) :: us(lx, lx, lx)
3201 real(kind=rp) :: ut(lx, lx, lx)
3202 integer :: i, j, k
3203
3204 do j = 1, lx * lx
3205 do i = 1, lx
3206 ur(i,j,1) = dx(i,1) * u(1,j,1) &
3207 + dx(i,2) * u(2,j,1)
3208 end do
3209 end do
3210
3211 do k = 1, lx
3212 do j = 1, lx
3213 do i = 1, lx
3214 us(i,j,k) = dy(j,1) * u(i,1,k) &
3215 + dy(j,2) * u(i,2,k)
3216 end do
3217 end do
3218 end do
3219
3220 do k = 1, lx
3221 do i = 1, lx*lx
3222 ut(i,1,k) = dz(k,1) * u(i,1,1) &
3223 + dz(k,2) * u(i,1,2)
3224 end do
3225 end do
3226
3227 do i = 1, lx * lx * lx
3228 ux(i,1,1) = w3(i,1,1) &
3229 * ( drdx(i,1,1) * ur(i,1,1) &
3230 + dsdx(i,1,1) * us(i,1,1) &
3231 + dtdx(i,1,1) * ut(i,1,1) )
3232 uy(i,1,1) = w3(i,1,1) &
3233 * ( dsdy(i,1,1) * us(i,1,1) &
3234 + drdy(i,1,1) * ur(i,1,1) &
3235 + dtdy(i,1,1) * ut(i,1,1) )
3236 uz(i,1,1) = w3(i,1,1) &
3237 * ( dtdz(i,1,1) * ut(i,1,1) &
3238 + drdz(i,1,1) * ur(i,1,1) &
3239 + dsdz(i,1,1) * us(i,1,1) )
3240 end do
3241
3242 end subroutine cpu_opgrad_lx2_single
3243
3244end submodule cpu_opgrad
Operators CPU backend.
Definition opr_cpu.f90:34