Neko 1.99.9
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
ax_helm_compr_cpu.f90
Go to the documentation of this file.
1! Copyright (c) 2026, 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!
37 use ax_helm_cpu, only : ax_helm_cpu_t
38 use num_types, only : rp
39 use coefs, only : coef_t
40 use space, only : space_t
41 use mesh, only : mesh_t
42 use utils, only : neko_error
43 implicit none
44 private
45
55 type, public, extends(ax_helm_cpu_t) :: ax_helm_compr_cpu_t
56 contains
58 procedure, pass(this) :: compute => ax_helm_compr_compute
60 procedure, pass(this) :: compute_vector => &
61 ax_helm_compr_compute_vector
62 end type ax_helm_compr_cpu_t
63
64 interface
65
76 module subroutine ax_helm_compr_compute_vector(this, au, av, aw, &
77 u, v, w, coef, msh, xh)
78 class(ax_helm_compr_cpu_t), intent(in) :: this
79 type(mesh_t), intent(in) :: msh
80 type(space_t), intent(in) :: xh
81 type(coef_t), intent(in) :: coef
82 real(kind=rp), intent(inout) :: au(xh%lx, xh%ly, xh%lz, msh%nelv)
83 real(kind=rp), intent(inout) :: av(xh%lx, xh%ly, xh%lz, msh%nelv)
84 real(kind=rp), intent(inout) :: aw(xh%lx, xh%ly, xh%lz, msh%nelv)
85 real(kind=rp), intent(in) :: u(xh%lx, xh%ly, xh%lz, msh%nelv)
86 real(kind=rp), intent(in) :: v(xh%lx, xh%ly, xh%lz, msh%nelv)
87 real(kind=rp), intent(in) :: w(xh%lx, xh%ly, xh%lz, msh%nelv)
88 end subroutine ax_helm_compr_compute_vector
89 end interface
90
91contains
92
101 subroutine ax_helm_compr_compute(this, w, u, coef, msh, Xh)
102 class(ax_helm_compr_cpu_t), intent(in) :: this
103 type(mesh_t), intent(in) :: msh
104 type(space_t), intent(in) :: Xh
105 type(coef_t), intent(in) :: coef
106 real(kind=rp), intent(inout) :: w(xh%lx, xh%ly, xh%lz, msh%nelv)
107 real(kind=rp), intent(in) :: u(xh%lx, xh%ly, xh%lz, msh%nelv)
108 integer :: i, m
109
110 if (.not. allocated(coef%compression_inds)) then
111 call this%ax_helm_cpu_t%compute(w, u, coef, msh, xh)
112 return
113 end if
114
115 if (size(coef%compression_inds) .ne. msh%nelv) then
116 call neko_error('Compressed geometric factors do not match the mesh')
117 end if
118 m = size(coef%G11_compr, 4)
119
120 !$omp parallel
121 select case(xh%lx)
122 case (14)
123 call ax_helm_compr_lx14(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
124 xh%dzt, coef%h1, coef%G11_compr, &
125 coef%G22_compr, coef%G33_compr, &
126 coef%G12_compr, coef%G13_compr, &
127 coef%G23_compr, coef%compression_inds, &
128 msh%nelv, m)
129 case (13)
130 call ax_helm_compr_lx13(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
131 xh%dzt, coef%h1, coef%G11_compr, &
132 coef%G22_compr, coef%G33_compr, &
133 coef%G12_compr, coef%G13_compr, &
134 coef%G23_compr, coef%compression_inds, &
135 msh%nelv, m)
136 case (12)
137 call ax_helm_compr_lx12(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
138 xh%dzt, coef%h1, coef%G11_compr, &
139 coef%G22_compr, coef%G33_compr, &
140 coef%G12_compr, coef%G13_compr, &
141 coef%G23_compr, coef%compression_inds, &
142 msh%nelv, m)
143 case (11)
144 call ax_helm_compr_lx11(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
145 xh%dzt, coef%h1, coef%G11_compr, &
146 coef%G22_compr, coef%G33_compr, &
147 coef%G12_compr, coef%G13_compr, &
148 coef%G23_compr, coef%compression_inds, &
149 msh%nelv, m)
150 case (10)
151 call ax_helm_compr_lx10(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
152 xh%dzt, coef%h1, coef%G11_compr, &
153 coef%G22_compr, coef%G33_compr, &
154 coef%G12_compr, coef%G13_compr, &
155 coef%G23_compr, coef%compression_inds, &
156 msh%nelv, m)
157 case (9)
158 call ax_helm_compr_lx9(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
159 xh%dzt, coef%h1, coef%G11_compr, &
160 coef%G22_compr, coef%G33_compr, &
161 coef%G12_compr, coef%G13_compr, &
162 coef%G23_compr, coef%compression_inds, &
163 msh%nelv, m)
164 case (8)
165 call ax_helm_compr_lx8(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
166 xh%dzt, coef%h1, coef%G11_compr, &
167 coef%G22_compr, coef%G33_compr, &
168 coef%G12_compr, coef%G13_compr, &
169 coef%G23_compr, coef%compression_inds, &
170 msh%nelv, m)
171 case (7)
172 call ax_helm_compr_lx7(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
173 xh%dzt, coef%h1, coef%G11_compr, &
174 coef%G22_compr, coef%G33_compr, &
175 coef%G12_compr, coef%G13_compr, &
176 coef%G23_compr, coef%compression_inds, &
177 msh%nelv, m)
178 case (6)
179 call ax_helm_compr_lx6(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
180 xh%dzt, coef%h1, coef%G11_compr, &
181 coef%G22_compr, coef%G33_compr, &
182 coef%G12_compr, coef%G13_compr, &
183 coef%G23_compr, coef%compression_inds, &
184 msh%nelv, m)
185 case (5)
186 call ax_helm_compr_lx5(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
187 xh%dzt, coef%h1, coef%G11_compr, &
188 coef%G22_compr, coef%G33_compr, &
189 coef%G12_compr, coef%G13_compr, &
190 coef%G23_compr, coef%compression_inds, &
191 msh%nelv, m)
192 case (4)
193 call ax_helm_compr_lx4(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
194 xh%dzt, coef%h1, coef%G11_compr, &
195 coef%G22_compr, coef%G33_compr, &
196 coef%G12_compr, coef%G13_compr, &
197 coef%G23_compr, coef%compression_inds, &
198 msh%nelv, m)
199 case (3)
200 call ax_helm_compr_lx3(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
201 xh%dzt, coef%h1, coef%G11_compr, &
202 coef%G22_compr, coef%G33_compr, &
203 coef%G12_compr, coef%G13_compr, &
204 coef%G23_compr, coef%compression_inds, &
205 msh%nelv, m)
206 case (2)
207 call ax_helm_compr_lx2(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
208 xh%dzt, coef%h1, coef%G11_compr, &
209 coef%G22_compr, coef%G33_compr, &
210 coef%G12_compr, coef%G13_compr, &
211 coef%G23_compr, coef%compression_inds, &
212 msh%nelv, m)
213 case default
214 call ax_helm_compr_lx(w, u, xh%dx, xh%dy, xh%dz, xh%dxt, xh%dyt, &
215 xh%dzt, coef%h1, coef%G11_compr, &
216 coef%G22_compr, coef%G33_compr, &
217 coef%G12_compr, coef%G13_compr, &
218 coef%G23_compr, coef%compression_inds, &
219 msh%nelv, m, xh%lx)
220 end select
221
222 if (coef%ifh2) then
223 !$omp do private(i)
224 do i = 1, coef%dof%size()
225 w(i,1,1,1) = w(i,1,1,1) + &
226 coef%h2(i,1,1,1) * coef%B(i,1,1,1) * u(i,1,1,1)
227 end do
228 !$omp end do
229 end if
230 !$omp end parallel
231
232 end subroutine ax_helm_compr_compute
233
248 subroutine ax_helm_compr_lx(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
249 h1, G11, G22, G33, G12, G13, G23, ind, n, m, lx)
250 integer, intent(in) :: n, m, lx
251 integer, intent(in) :: ind(n)
252 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
253 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
254 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
255 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
256 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
257 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
258 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
259 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
260 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
261 real(kind=rp), intent(in) :: dx(lx, lx)
262 real(kind=rp), intent(in) :: dy(lx, lx)
263 real(kind=rp), intent(in) :: dz(lx, lx)
264 real(kind=rp), intent(in) :: dxt(lx, lx)
265 real(kind=rp), intent(in) :: dyt(lx, lx)
266 real(kind=rp), intent(in) :: dzt(lx, lx)
267 real(kind=rp) :: ur(lx, lx, lx)
268 real(kind=rp) :: us(lx, lx, lx)
269 real(kind=rp) :: ut(lx, lx, lx)
270 real(kind=rp) :: wur(lx, lx, lx)
271 real(kind=rp) :: wus(lx, lx, lx)
272 real(kind=rp) :: wut(lx, lx, lx)
273 real(kind=rp) :: tmp
274 integer :: e, c, i, j, k, l
275
276 !$omp do
277 do e = 1, n
278 c = ind(e)
279 do j = 1, lx * lx
280 do i = 1, lx
281 tmp = 0.0_rp
282 do k = 1, lx
283 tmp = tmp + dx(i,k) * u(k,j,1,e)
284 end do
285 wur(i,j,1) = tmp
286 end do
287 end do
288
289 do k = 1, lx
290 do j = 1, lx
291 do i = 1, lx
292 tmp = 0.0_rp
293 do l = 1, lx
294 tmp = tmp + dy(j,l) * u(i,l,k,e)
295 end do
296 wus(i,j,k) = tmp
297 end do
298 end do
299 end do
300
301 do k = 1, lx
302 do i = 1, lx*lx
303 tmp = 0.0_rp
304 do l = 1, lx
305 tmp = tmp + dz(k,l) * u(i,1,l,e)
306 end do
307 wut(i,1,k) = tmp
308 end do
309 end do
310
311 do i = 1, lx*lx*lx
312 ur(i,1,1) = h1(i,1,1,e) &
313 * ( g11(i,1,1,c) * wur(i,1,1) &
314 + g12(i,1,1,c) * wus(i,1,1) &
315 + g13(i,1,1,c) * wut(i,1,1) )
316 us(i,1,1) = h1(i,1,1,e) &
317 * ( g12(i,1,1,c) * wur(i,1,1) &
318 + g22(i,1,1,c) * wus(i,1,1) &
319 + g23(i,1,1,c) * wut(i,1,1) )
320 ut(i,1,1) = h1(i,1,1,e) &
321 * ( g13(i,1,1,c) * wur(i,1,1) &
322 + g23(i,1,1,c) * wus(i,1,1) &
323 + g33(i,1,1,c) * wut(i,1,1) )
324 end do
325
326 do j = 1, lx*lx
327 do i = 1, lx
328 tmp = 0.0_rp
329 do k = 1, lx
330 tmp = tmp + dxt(i,k) * ur(k,j,1)
331 end do
332 w(i,j,1,e) = tmp
333 end do
334 end do
335
336 do k = 1, lx
337 do j = 1, lx
338 do i = 1, lx
339 tmp = 0.0_rp
340 do l = 1, lx
341 tmp = tmp + dyt(j,l) * us(i,l,k)
342 end do
343 w(i,j,k,e) = w(i,j,k,e) + tmp
344 end do
345 end do
346 end do
347
348 do k = 1, lx
349 do i = 1, lx*lx
350 tmp = 0.0_rp
351 do l = 1, lx
352 tmp = tmp + dzt(k,l) * ut(i,1,l)
353 end do
354 w(i,1,k,e) = w(i,1,k,e) + tmp
355 end do
356 end do
357
358 end do
359 !$omp end do
360 end subroutine ax_helm_compr_lx
361
362 subroutine ax_helm_compr_lx14(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
363 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
364 integer, parameter :: lx = 14
365 integer, intent(in) :: n, m
366 integer, intent(in) :: ind(n)
367 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
368 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
369 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
370 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
371 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
372 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
373 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
374 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
375 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
376 real(kind=rp), intent(in) :: dx(lx, lx)
377 real(kind=rp), intent(in) :: dy(lx, lx)
378 real(kind=rp), intent(in) :: dz(lx, lx)
379 real(kind=rp), intent(in) :: dxt(lx, lx)
380 real(kind=rp), intent(in) :: dyt(lx, lx)
381 real(kind=rp), intent(in) :: dzt(lx, lx)
382 real(kind=rp) :: ur(lx, lx, lx)
383 real(kind=rp) :: us(lx, lx, lx)
384 real(kind=rp) :: ut(lx, lx, lx)
385 real(kind=rp) :: wur(lx, lx, lx)
386 real(kind=rp) :: wus(lx, lx, lx)
387 real(kind=rp) :: wut(lx, lx, lx)
388 integer :: e, c, i, j, k
389
390 !$omp do
391 do e = 1, n
392 c = ind(e)
393 do j = 1, lx * lx
394 do i = 1, lx
395 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
396 + dx(i,2) * u(2,j,1,e) &
397 + dx(i,3) * u(3,j,1,e) &
398 + dx(i,4) * u(4,j,1,e) &
399 + dx(i,5) * u(5,j,1,e) &
400 + dx(i,6) * u(6,j,1,e) &
401 + dx(i,7) * u(7,j,1,e) &
402 + dx(i,8) * u(8,j,1,e) &
403 + dx(i,9) * u(9,j,1,e) &
404 + dx(i,10) * u(10,j,1,e) &
405 + dx(i,11) * u(11,j,1,e) &
406 + dx(i,12) * u(12,j,1,e) &
407 + dx(i,13) * u(13,j,1,e) &
408 + dx(i,14) * u(14,j,1,e)
409 end do
410 end do
411
412 do k = 1, lx
413 do j = 1, lx
414 do i = 1, lx
415 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
416 + dy(j,2) * u(i,2,k,e) &
417 + dy(j,3) * u(i,3,k,e) &
418 + dy(j,4) * u(i,4,k,e) &
419 + dy(j,5) * u(i,5,k,e) &
420 + dy(j,6) * u(i,6,k,e) &
421 + dy(j,7) * u(i,7,k,e) &
422 + dy(j,8) * u(i,8,k,e) &
423 + dy(j,9) * u(i,9,k,e) &
424 + dy(j,10) * u(i,10,k,e) &
425 + dy(j,11) * u(i,11,k,e) &
426 + dy(j,12) * u(i,12,k,e) &
427 + dy(j,13) * u(i,13,k,e) &
428 + dy(j,14) * u(i,14,k,e)
429 end do
430 end do
431 end do
432
433 do k = 1, lx
434 do i = 1, lx*lx
435 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
436 + dz(k,2) * u(i,1,2,e) &
437 + dz(k,3) * u(i,1,3,e) &
438 + dz(k,4) * u(i,1,4,e) &
439 + dz(k,5) * u(i,1,5,e) &
440 + dz(k,6) * u(i,1,6,e) &
441 + dz(k,7) * u(i,1,7,e) &
442 + dz(k,8) * u(i,1,8,e) &
443 + dz(k,9) * u(i,1,9,e) &
444 + dz(k,10) * u(i,1,10,e) &
445 + dz(k,11) * u(i,1,11,e) &
446 + dz(k,12) * u(i,1,12,e) &
447 + dz(k,13) * u(i,1,13,e) &
448 + dz(k,14) * u(i,1,14,e)
449 end do
450 end do
451
452 do i = 1, lx*lx*lx
453 ur(i,1,1) = h1(i,1,1,e) &
454 * ( g11(i,1,1,c) * wur(i,1,1) &
455 + g12(i,1,1,c) * wus(i,1,1) &
456 + g13(i,1,1,c) * wut(i,1,1) )
457 us(i,1,1) = h1(i,1,1,e) &
458 * ( g12(i,1,1,c) * wur(i,1,1) &
459 + g22(i,1,1,c) * wus(i,1,1) &
460 + g23(i,1,1,c) * wut(i,1,1) )
461 ut(i,1,1) = h1(i,1,1,e) &
462 * ( g13(i,1,1,c) * wur(i,1,1) &
463 + g23(i,1,1,c) * wus(i,1,1) &
464 + g33(i,1,1,c) * wut(i,1,1) )
465 end do
466
467 do j = 1, lx*lx
468 do i = 1, lx
469 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
470 + dxt(i,2) * ur(2,j,1) &
471 + dxt(i,3) * ur(3,j,1) &
472 + dxt(i,4) * ur(4,j,1) &
473 + dxt(i,5) * ur(5,j,1) &
474 + dxt(i,6) * ur(6,j,1) &
475 + dxt(i,7) * ur(7,j,1) &
476 + dxt(i,8) * ur(8,j,1) &
477 + dxt(i,9) * ur(9,j,1) &
478 + dxt(i,10) * ur(10,j,1) &
479 + dxt(i,11) * ur(11,j,1) &
480 + dxt(i,12) * ur(12,j,1) &
481 + dxt(i,13) * ur(13,j,1) &
482 + dxt(i,14) * ur(14,j,1)
483 end do
484 end do
485
486 do k = 1, lx
487 do j = 1, lx
488 do i = 1, lx
489 w(i,j,k,e) = w(i,j,k,e) &
490 + dyt(j,1) * us(i,1,k) &
491 + dyt(j,2) * us(i,2,k) &
492 + dyt(j,3) * us(i,3,k) &
493 + dyt(j,4) * us(i,4,k) &
494 + dyt(j,5) * us(i,5,k) &
495 + dyt(j,6) * us(i,6,k) &
496 + dyt(j,7) * us(i,7,k) &
497 + dyt(j,8) * us(i,8,k) &
498 + dyt(j,9) * us(i,9,k) &
499 + dyt(j,10) * us(i,10,k) &
500 + dyt(j,11) * us(i,11,k) &
501 + dyt(j,12) * us(i,12,k) &
502 + dyt(j,13) * us(i,13,k) &
503 + dyt(j,14) * us(i,14,k)
504 end do
505 end do
506 end do
507
508 do k = 1, lx
509 do i = 1, lx*lx
510 w(i,1,k,e) = w(i,1,k,e) &
511 + dzt(k,1) * ut(i,1,1) &
512 + dzt(k,2) * ut(i,1,2) &
513 + dzt(k,3) * ut(i,1,3) &
514 + dzt(k,4) * ut(i,1,4) &
515 + dzt(k,5) * ut(i,1,5) &
516 + dzt(k,6) * ut(i,1,6) &
517 + dzt(k,7) * ut(i,1,7) &
518 + dzt(k,8) * ut(i,1,8) &
519 + dzt(k,9) * ut(i,1,9) &
520 + dzt(k,10) * ut(i,1,10) &
521 + dzt(k,11) * ut(i,1,11) &
522 + dzt(k,12) * ut(i,1,12) &
523 + dzt(k,13) * ut(i,1,13) &
524 + dzt(k,14) * ut(i,1,14)
525 end do
526 end do
527
528 end do
529 !$omp end do
530 end subroutine ax_helm_compr_lx14
531
532 subroutine ax_helm_compr_lx13(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
533 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
534 integer, parameter :: lx = 13
535 integer, intent(in) :: n, m
536 integer, intent(in) :: ind(n)
537 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
538 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
539 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
540 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
541 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
542 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
543 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
544 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
545 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
546 real(kind=rp), intent(in) :: dx(lx, lx)
547 real(kind=rp), intent(in) :: dy(lx, lx)
548 real(kind=rp), intent(in) :: dz(lx, lx)
549 real(kind=rp), intent(in) :: dxt(lx, lx)
550 real(kind=rp), intent(in) :: dyt(lx, lx)
551 real(kind=rp), intent(in) :: dzt(lx, lx)
552 real(kind=rp) :: ur(lx, lx, lx)
553 real(kind=rp) :: us(lx, lx, lx)
554 real(kind=rp) :: ut(lx, lx, lx)
555 real(kind=rp) :: wur(lx, lx, lx)
556 real(kind=rp) :: wus(lx, lx, lx)
557 real(kind=rp) :: wut(lx, lx, lx)
558 integer :: e, c, i, j, k
559
560 !$omp do
561 do e = 1, n
562 c = ind(e)
563 do j = 1, lx * lx
564 do i = 1, lx
565 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
566 + dx(i,2) * u(2,j,1,e) &
567 + dx(i,3) * u(3,j,1,e) &
568 + dx(i,4) * u(4,j,1,e) &
569 + dx(i,5) * u(5,j,1,e) &
570 + dx(i,6) * u(6,j,1,e) &
571 + dx(i,7) * u(7,j,1,e) &
572 + dx(i,8) * u(8,j,1,e) &
573 + dx(i,9) * u(9,j,1,e) &
574 + dx(i,10) * u(10,j,1,e) &
575 + dx(i,11) * u(11,j,1,e) &
576 + dx(i,12) * u(12,j,1,e) &
577 + dx(i,13) * u(13,j,1,e)
578
579 end do
580 end do
581
582 do k = 1, lx
583 do j = 1, lx
584 do i = 1, lx
585 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
586 + dy(j,2) * u(i,2,k,e) &
587 + dy(j,3) * u(i,3,k,e) &
588 + dy(j,4) * u(i,4,k,e) &
589 + dy(j,5) * u(i,5,k,e) &
590 + dy(j,6) * u(i,6,k,e) &
591 + dy(j,7) * u(i,7,k,e) &
592 + dy(j,8) * u(i,8,k,e) &
593 + dy(j,9) * u(i,9,k,e) &
594 + dy(j,10) * u(i,10,k,e) &
595 + dy(j,11) * u(i,11,k,e) &
596 + dy(j,12) * u(i,12,k,e) &
597 + dy(j,13) * u(i,13,k,e)
598 end do
599 end do
600 end do
601
602 do k = 1, lx
603 do i = 1, lx*lx
604 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
605 + dz(k,2) * u(i,1,2,e) &
606 + dz(k,3) * u(i,1,3,e) &
607 + dz(k,4) * u(i,1,4,e) &
608 + dz(k,5) * u(i,1,5,e) &
609 + dz(k,6) * u(i,1,6,e) &
610 + dz(k,7) * u(i,1,7,e) &
611 + dz(k,8) * u(i,1,8,e) &
612 + dz(k,9) * u(i,1,9,e) &
613 + dz(k,10) * u(i,1,10,e) &
614 + dz(k,11) * u(i,1,11,e) &
615 + dz(k,12) * u(i,1,12,e) &
616 + dz(k,13) * u(i,1,13,e)
617 end do
618 end do
619
620 do i = 1, lx*lx*lx
621 ur(i,1,1) = h1(i,1,1,e) &
622 * ( g11(i,1,1,c) * wur(i,1,1) &
623 + g12(i,1,1,c) * wus(i,1,1) &
624 + g13(i,1,1,c) * wut(i,1,1) )
625 us(i,1,1) = h1(i,1,1,e) &
626 * ( g12(i,1,1,c) * wur(i,1,1) &
627 + g22(i,1,1,c) * wus(i,1,1) &
628 + g23(i,1,1,c) * wut(i,1,1) )
629 ut(i,1,1) = h1(i,1,1,e) &
630 * ( g13(i,1,1,c) * wur(i,1,1) &
631 + g23(i,1,1,c) * wus(i,1,1) &
632 + g33(i,1,1,c) * wut(i,1,1) )
633 end do
634
635 do j = 1, lx*lx
636 do i = 1, lx
637 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
638 + dxt(i,2) * ur(2,j,1) &
639 + dxt(i,3) * ur(3,j,1) &
640 + dxt(i,4) * ur(4,j,1) &
641 + dxt(i,5) * ur(5,j,1) &
642 + dxt(i,6) * ur(6,j,1) &
643 + dxt(i,7) * ur(7,j,1) &
644 + dxt(i,8) * ur(8,j,1) &
645 + dxt(i,9) * ur(9,j,1) &
646 + dxt(i,10) * ur(10,j,1) &
647 + dxt(i,11) * ur(11,j,1) &
648 + dxt(i,12) * ur(12,j,1) &
649 + dxt(i,13) * ur(13,j,1)
650 end do
651 end do
652
653 do k = 1, lx
654 do j = 1, lx
655 do i = 1, lx
656 w(i,j,k,e) = w(i,j,k,e) &
657 + dyt(j,1) * us(i,1,k) &
658 + dyt(j,2) * us(i,2,k) &
659 + dyt(j,3) * us(i,3,k) &
660 + dyt(j,4) * us(i,4,k) &
661 + dyt(j,5) * us(i,5,k) &
662 + dyt(j,6) * us(i,6,k) &
663 + dyt(j,7) * us(i,7,k) &
664 + dyt(j,8) * us(i,8,k) &
665 + dyt(j,9) * us(i,9,k) &
666 + dyt(j,10) * us(i,10,k) &
667 + dyt(j,11) * us(i,11,k) &
668 + dyt(j,12) * us(i,12,k) &
669 + dyt(j,13) * us(i,13,k)
670 end do
671 end do
672 end do
673
674 do k = 1, lx
675 do i = 1, lx*lx
676 w(i,1,k,e) = w(i,1,k,e) &
677 + dzt(k,1) * ut(i,1,1) &
678 + dzt(k,2) * ut(i,1,2) &
679 + dzt(k,3) * ut(i,1,3) &
680 + dzt(k,4) * ut(i,1,4) &
681 + dzt(k,5) * ut(i,1,5) &
682 + dzt(k,6) * ut(i,1,6) &
683 + dzt(k,7) * ut(i,1,7) &
684 + dzt(k,8) * ut(i,1,8) &
685 + dzt(k,9) * ut(i,1,9) &
686 + dzt(k,10) * ut(i,1,10) &
687 + dzt(k,11) * ut(i,1,11) &
688 + dzt(k,12) * ut(i,1,12) &
689 + dzt(k,13) * ut(i,1,13)
690 end do
691 end do
692
693 end do
694 !$omp end do
695 end subroutine ax_helm_compr_lx13
696
697 subroutine ax_helm_compr_lx12(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
698 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
699 integer, parameter :: lx = 12
700 integer, intent(in) :: n, m
701 integer, intent(in) :: ind(n)
702 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
703 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
704 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
705 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
706 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
707 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
708 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
709 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
710 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
711 real(kind=rp), intent(in) :: dx(lx, lx)
712 real(kind=rp), intent(in) :: dy(lx, lx)
713 real(kind=rp), intent(in) :: dz(lx, lx)
714 real(kind=rp), intent(in) :: dxt(lx, lx)
715 real(kind=rp), intent(in) :: dyt(lx, lx)
716 real(kind=rp), intent(in) :: dzt(lx, lx)
717 real(kind=rp) :: ur(lx, lx, lx)
718 real(kind=rp) :: us(lx, lx, lx)
719 real(kind=rp) :: ut(lx, lx, lx)
720 real(kind=rp) :: wur(lx, lx, lx)
721 real(kind=rp) :: wus(lx, lx, lx)
722 real(kind=rp) :: wut(lx, lx, lx)
723 integer :: e, c, i, j, k
724
725 !$omp do
726 do e = 1, n
727 c = ind(e)
728 do j = 1, lx * lx
729 do i = 1, lx
730 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
731 + dx(i,2) * u(2,j,1,e) &
732 + dx(i,3) * u(3,j,1,e) &
733 + dx(i,4) * u(4,j,1,e) &
734 + dx(i,5) * u(5,j,1,e) &
735 + dx(i,6) * u(6,j,1,e) &
736 + dx(i,7) * u(7,j,1,e) &
737 + dx(i,8) * u(8,j,1,e) &
738 + dx(i,9) * u(9,j,1,e) &
739 + dx(i,10) * u(10,j,1,e) &
740 + dx(i,11) * u(11,j,1,e) &
741 + dx(i,12) * u(12,j,1,e)
742 end do
743 end do
744
745 do k = 1, lx
746 do j = 1, lx
747 do i = 1, lx
748 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
749 + dy(j,2) * u(i,2,k,e) &
750 + dy(j,3) * u(i,3,k,e) &
751 + dy(j,4) * u(i,4,k,e) &
752 + dy(j,5) * u(i,5,k,e) &
753 + dy(j,6) * u(i,6,k,e) &
754 + dy(j,7) * u(i,7,k,e) &
755 + dy(j,8) * u(i,8,k,e) &
756 + dy(j,9) * u(i,9,k,e) &
757 + dy(j,10) * u(i,10,k,e) &
758 + dy(j,11) * u(i,11,k,e) &
759 + dy(j,12) * u(i,12,k,e)
760 end do
761 end do
762 end do
763
764 do k = 1, lx
765 do i = 1, lx*lx
766 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
767 + dz(k,2) * u(i,1,2,e) &
768 + dz(k,3) * u(i,1,3,e) &
769 + dz(k,4) * u(i,1,4,e) &
770 + dz(k,5) * u(i,1,5,e) &
771 + dz(k,6) * u(i,1,6,e) &
772 + dz(k,7) * u(i,1,7,e) &
773 + dz(k,8) * u(i,1,8,e) &
774 + dz(k,9) * u(i,1,9,e) &
775 + dz(k,10) * u(i,1,10,e) &
776 + dz(k,11) * u(i,1,11,e) &
777 + dz(k,12) * u(i,1,12,e)
778 end do
779 end do
780
781 do i = 1, lx*lx*lx
782 ur(i,1,1) = h1(i,1,1,e) &
783 * ( g11(i,1,1,c) * wur(i,1,1) &
784 + g12(i,1,1,c) * wus(i,1,1) &
785 + g13(i,1,1,c) * wut(i,1,1) )
786 us(i,1,1) = h1(i,1,1,e) &
787 * ( g12(i,1,1,c) * wur(i,1,1) &
788 + g22(i,1,1,c) * wus(i,1,1) &
789 + g23(i,1,1,c) * wut(i,1,1) )
790 ut(i,1,1) = h1(i,1,1,e) &
791 * ( g13(i,1,1,c) * wur(i,1,1) &
792 + g23(i,1,1,c) * wus(i,1,1) &
793 + g33(i,1,1,c) * wut(i,1,1) )
794 end do
795
796 do j = 1, lx*lx
797 do i = 1, lx
798 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
799 + dxt(i,2) * ur(2,j,1) &
800 + dxt(i,3) * ur(3,j,1) &
801 + dxt(i,4) * ur(4,j,1) &
802 + dxt(i,5) * ur(5,j,1) &
803 + dxt(i,6) * ur(6,j,1) &
804 + dxt(i,7) * ur(7,j,1) &
805 + dxt(i,8) * ur(8,j,1) &
806 + dxt(i,9) * ur(9,j,1) &
807 + dxt(i,10) * ur(10,j,1) &
808 + dxt(i,11) * ur(11,j,1) &
809 + dxt(i,12) * ur(12,j,1)
810 end do
811 end do
812
813 do k = 1, lx
814 do j = 1, lx
815 do i = 1, lx
816 w(i,j,k,e) = w(i,j,k,e) &
817 + dyt(j,1) * us(i,1,k) &
818 + dyt(j,2) * us(i,2,k) &
819 + dyt(j,3) * us(i,3,k) &
820 + dyt(j,4) * us(i,4,k) &
821 + dyt(j,5) * us(i,5,k) &
822 + dyt(j,6) * us(i,6,k) &
823 + dyt(j,7) * us(i,7,k) &
824 + dyt(j,8) * us(i,8,k) &
825 + dyt(j,9) * us(i,9,k) &
826 + dyt(j,10) * us(i,10,k) &
827 + dyt(j,11) * us(i,11,k) &
828 + dyt(j,12) * us(i,12,k)
829 end do
830 end do
831 end do
832
833 do k = 1, lx
834 do i = 1, lx*lx
835 w(i,1,k,e) = w(i,1,k,e) &
836 + dzt(k,1) * ut(i,1,1) &
837 + dzt(k,2) * ut(i,1,2) &
838 + dzt(k,3) * ut(i,1,3) &
839 + dzt(k,4) * ut(i,1,4) &
840 + dzt(k,5) * ut(i,1,5) &
841 + dzt(k,6) * ut(i,1,6) &
842 + dzt(k,7) * ut(i,1,7) &
843 + dzt(k,8) * ut(i,1,8) &
844 + dzt(k,9) * ut(i,1,9) &
845 + dzt(k,10) * ut(i,1,10) &
846 + dzt(k,11) * ut(i,1,11) &
847 + dzt(k,12) * ut(i,1,12)
848 end do
849 end do
850
851 end do
852 !$omp end do
853 end subroutine ax_helm_compr_lx12
854
855 subroutine ax_helm_compr_lx11(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
856 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
857 integer, parameter :: lx = 11
858 integer, intent(in) :: n, m
859 integer, intent(in) :: ind(n)
860 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
861 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
862 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
863 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
864 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
865 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
866 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
867 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
868 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
869 real(kind=rp), intent(in) :: dx(lx, lx)
870 real(kind=rp), intent(in) :: dy(lx, lx)
871 real(kind=rp), intent(in) :: dz(lx, lx)
872 real(kind=rp), intent(in) :: dxt(lx, lx)
873 real(kind=rp), intent(in) :: dyt(lx, lx)
874 real(kind=rp), intent(in) :: dzt(lx, lx)
875 real(kind=rp) :: ur(lx, lx, lx)
876 real(kind=rp) :: us(lx, lx, lx)
877 real(kind=rp) :: ut(lx, lx, lx)
878 real(kind=rp) :: wur(lx, lx, lx)
879 real(kind=rp) :: wus(lx, lx, lx)
880 real(kind=rp) :: wut(lx, lx, lx)
881 integer :: e, c, i, j, k
882
883 !$omp do
884 do e = 1, n
885 c = ind(e)
886 do j = 1, lx * lx
887 do i = 1, lx
888 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
889 + dx(i,2) * u(2,j,1,e) &
890 + dx(i,3) * u(3,j,1,e) &
891 + dx(i,4) * u(4,j,1,e) &
892 + dx(i,5) * u(5,j,1,e) &
893 + dx(i,6) * u(6,j,1,e) &
894 + dx(i,7) * u(7,j,1,e) &
895 + dx(i,8) * u(8,j,1,e) &
896 + dx(i,9) * u(9,j,1,e) &
897 + dx(i,10) * u(10,j,1,e) &
898 + dx(i,11) * u(11,j,1,e)
899 end do
900 end do
901
902 do k = 1, lx
903 do j = 1, lx
904 do i = 1, lx
905 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
906 + dy(j,2) * u(i,2,k,e) &
907 + dy(j,3) * u(i,3,k,e) &
908 + dy(j,4) * u(i,4,k,e) &
909 + dy(j,5) * u(i,5,k,e) &
910 + dy(j,6) * u(i,6,k,e) &
911 + dy(j,7) * u(i,7,k,e) &
912 + dy(j,8) * u(i,8,k,e) &
913 + dy(j,9) * u(i,9,k,e) &
914 + dy(j,10) * u(i,10,k,e) &
915 + dy(j,11) * u(i,11,k,e)
916 end do
917 end do
918 end do
919
920 do k = 1, lx
921 do i = 1, lx*lx
922 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
923 + dz(k,2) * u(i,1,2,e) &
924 + dz(k,3) * u(i,1,3,e) &
925 + dz(k,4) * u(i,1,4,e) &
926 + dz(k,5) * u(i,1,5,e) &
927 + dz(k,6) * u(i,1,6,e) &
928 + dz(k,7) * u(i,1,7,e) &
929 + dz(k,8) * u(i,1,8,e) &
930 + dz(k,9) * u(i,1,9,e) &
931 + dz(k,10) * u(i,1,10,e) &
932 + dz(k,11) * u(i,1,11,e)
933 end do
934 end do
935
936 do i = 1, lx*lx*lx
937 ur(i,1,1) = h1(i,1,1,e) &
938 * ( g11(i,1,1,c) * wur(i,1,1) &
939 + g12(i,1,1,c) * wus(i,1,1) &
940 + g13(i,1,1,c) * wut(i,1,1) )
941 us(i,1,1) = h1(i,1,1,e) &
942 * ( g12(i,1,1,c) * wur(i,1,1) &
943 + g22(i,1,1,c) * wus(i,1,1) &
944 + g23(i,1,1,c) * wut(i,1,1) )
945 ut(i,1,1) = h1(i,1,1,e) &
946 * ( g13(i,1,1,c) * wur(i,1,1) &
947 + g23(i,1,1,c) * wus(i,1,1) &
948 + g33(i,1,1,c) * wut(i,1,1) )
949 end do
950
951 do j = 1, lx*lx
952 do i = 1, lx
953 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
954 + dxt(i,2) * ur(2,j,1) &
955 + dxt(i,3) * ur(3,j,1) &
956 + dxt(i,4) * ur(4,j,1) &
957 + dxt(i,5) * ur(5,j,1) &
958 + dxt(i,6) * ur(6,j,1) &
959 + dxt(i,7) * ur(7,j,1) &
960 + dxt(i,8) * ur(8,j,1) &
961 + dxt(i,9) * ur(9,j,1) &
962 + dxt(i,10) * ur(10,j,1) &
963 + dxt(i,11) * ur(11,j,1)
964 end do
965 end do
966
967 do k = 1, lx
968 do j = 1, lx
969 do i = 1, lx
970 w(i,j,k,e) = w(i,j,k,e) &
971 + dyt(j,1) * us(i,1,k) &
972 + dyt(j,2) * us(i,2,k) &
973 + dyt(j,3) * us(i,3,k) &
974 + dyt(j,4) * us(i,4,k) &
975 + dyt(j,5) * us(i,5,k) &
976 + dyt(j,6) * us(i,6,k) &
977 + dyt(j,7) * us(i,7,k) &
978 + dyt(j,8) * us(i,8,k) &
979 + dyt(j,9) * us(i,9,k) &
980 + dyt(j,10) * us(i,10,k) &
981 + dyt(j,11) * us(i,11,k)
982 end do
983 end do
984 end do
985
986 do k = 1, lx
987 do i = 1, lx*lx
988 w(i,1,k,e) = w(i,1,k,e) &
989 + dzt(k,1) * ut(i,1,1) &
990 + dzt(k,2) * ut(i,1,2) &
991 + dzt(k,3) * ut(i,1,3) &
992 + dzt(k,4) * ut(i,1,4) &
993 + dzt(k,5) * ut(i,1,5) &
994 + dzt(k,6) * ut(i,1,6) &
995 + dzt(k,7) * ut(i,1,7) &
996 + dzt(k,8) * ut(i,1,8) &
997 + dzt(k,9) * ut(i,1,9) &
998 + dzt(k,10) * ut(i,1,10) &
999 + dzt(k,11) * ut(i,1,11)
1000 end do
1001 end do
1002
1003 end do
1004 !$omp end do
1005 end subroutine ax_helm_compr_lx11
1006
1007 subroutine ax_helm_compr_lx10(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1008 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1009 integer, parameter :: lx = 10
1010 integer, intent(in) :: n, m
1011 integer, intent(in) :: ind(n)
1012 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1013 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1014 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1015 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1016 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1017 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1018 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1019 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1020 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1021 real(kind=rp), intent(in) :: dx(lx, lx)
1022 real(kind=rp), intent(in) :: dy(lx, lx)
1023 real(kind=rp), intent(in) :: dz(lx, lx)
1024 real(kind=rp), intent(in) :: dxt(lx, lx)
1025 real(kind=rp), intent(in) :: dyt(lx, lx)
1026 real(kind=rp), intent(in) :: dzt(lx, lx)
1027 real(kind=rp) :: ur(lx, lx, lx)
1028 real(kind=rp) :: us(lx, lx, lx)
1029 real(kind=rp) :: ut(lx, lx, lx)
1030 real(kind=rp) :: wur(lx, lx, lx)
1031 real(kind=rp) :: wus(lx, lx, lx)
1032 real(kind=rp) :: wut(lx, lx, lx)
1033 integer :: e, c, i, j, k
1034
1035 !$omp do
1036 do e = 1, n
1037 c = ind(e)
1038 do j = 1, lx * lx
1039 do i = 1, lx
1040 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1041 + dx(i,2) * u(2,j,1,e) &
1042 + dx(i,3) * u(3,j,1,e) &
1043 + dx(i,4) * u(4,j,1,e) &
1044 + dx(i,5) * u(5,j,1,e) &
1045 + dx(i,6) * u(6,j,1,e) &
1046 + dx(i,7) * u(7,j,1,e) &
1047 + dx(i,8) * u(8,j,1,e) &
1048 + dx(i,9) * u(9,j,1,e) &
1049 + dx(i,10) * u(10,j,1,e)
1050 end do
1051 end do
1052
1053 do k = 1, lx
1054 do j = 1, lx
1055 do i = 1, lx
1056 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1057 + dy(j,2) * u(i,2,k,e) &
1058 + dy(j,3) * u(i,3,k,e) &
1059 + dy(j,4) * u(i,4,k,e) &
1060 + dy(j,5) * u(i,5,k,e) &
1061 + dy(j,6) * u(i,6,k,e) &
1062 + dy(j,7) * u(i,7,k,e) &
1063 + dy(j,8) * u(i,8,k,e) &
1064 + dy(j,9) * u(i,9,k,e) &
1065 + dy(j,10) * u(i,10,k,e)
1066 end do
1067 end do
1068 end do
1069
1070 do k = 1, lx
1071 do i = 1, lx*lx
1072 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1073 + dz(k,2) * u(i,1,2,e) &
1074 + dz(k,3) * u(i,1,3,e) &
1075 + dz(k,4) * u(i,1,4,e) &
1076 + dz(k,5) * u(i,1,5,e) &
1077 + dz(k,6) * u(i,1,6,e) &
1078 + dz(k,7) * u(i,1,7,e) &
1079 + dz(k,8) * u(i,1,8,e) &
1080 + dz(k,9) * u(i,1,9,e) &
1081 + dz(k,10) * u(i,1,10,e)
1082 end do
1083 end do
1084
1085 do i = 1, lx*lx*lx
1086 ur(i,1,1) = h1(i,1,1,e) &
1087 * ( g11(i,1,1,c) * wur(i,1,1) &
1088 + g12(i,1,1,c) * wus(i,1,1) &
1089 + g13(i,1,1,c) * wut(i,1,1) )
1090 us(i,1,1) = h1(i,1,1,e) &
1091 * ( g12(i,1,1,c) * wur(i,1,1) &
1092 + g22(i,1,1,c) * wus(i,1,1) &
1093 + g23(i,1,1,c) * wut(i,1,1) )
1094 ut(i,1,1) = h1(i,1,1,e) &
1095 * ( g13(i,1,1,c) * wur(i,1,1) &
1096 + g23(i,1,1,c) * wus(i,1,1) &
1097 + g33(i,1,1,c) * wut(i,1,1) )
1098 end do
1099
1100 do j = 1, lx*lx
1101 do i = 1, lx
1102 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1103 + dxt(i,2) * ur(2,j,1) &
1104 + dxt(i,3) * ur(3,j,1) &
1105 + dxt(i,4) * ur(4,j,1) &
1106 + dxt(i,5) * ur(5,j,1) &
1107 + dxt(i,6) * ur(6,j,1) &
1108 + dxt(i,7) * ur(7,j,1) &
1109 + dxt(i,8) * ur(8,j,1) &
1110 + dxt(i,9) * ur(9,j,1) &
1111 + dxt(i,10) * ur(10,j,1)
1112 end do
1113 end do
1114
1115 do k = 1, lx
1116 do j = 1, lx
1117 do i = 1, lx
1118 w(i,j,k,e) = w(i,j,k,e) &
1119 + dyt(j,1) * us(i,1,k) &
1120 + dyt(j,2) * us(i,2,k) &
1121 + dyt(j,3) * us(i,3,k) &
1122 + dyt(j,4) * us(i,4,k) &
1123 + dyt(j,5) * us(i,5,k) &
1124 + dyt(j,6) * us(i,6,k) &
1125 + dyt(j,7) * us(i,7,k) &
1126 + dyt(j,8) * us(i,8,k) &
1127 + dyt(j,9) * us(i,9,k) &
1128 + dyt(j,10) * us(i,10,k)
1129 end do
1130 end do
1131 end do
1132
1133 do k = 1, lx
1134 do i = 1, lx*lx
1135 w(i,1,k,e) = w(i,1,k,e) &
1136 + dzt(k,1) * ut(i,1,1) &
1137 + dzt(k,2) * ut(i,1,2) &
1138 + dzt(k,3) * ut(i,1,3) &
1139 + dzt(k,4) * ut(i,1,4) &
1140 + dzt(k,5) * ut(i,1,5) &
1141 + dzt(k,6) * ut(i,1,6) &
1142 + dzt(k,7) * ut(i,1,7) &
1143 + dzt(k,8) * ut(i,1,8) &
1144 + dzt(k,9) * ut(i,1,9) &
1145 + dzt(k,10) * ut(i,1,10)
1146 end do
1147 end do
1148
1149 end do
1150 !$omp end do
1151 end subroutine ax_helm_compr_lx10
1152
1153 subroutine ax_helm_compr_lx9(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1154 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1155 integer, parameter :: lx = 9
1156 integer, intent(in) :: n, m
1157 integer, intent(in) :: ind(n)
1158 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1159 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1160 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1161 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1162 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1163 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1164 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1165 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1166 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1167 real(kind=rp), intent(in) :: dx(lx, lx)
1168 real(kind=rp), intent(in) :: dy(lx, lx)
1169 real(kind=rp), intent(in) :: dz(lx, lx)
1170 real(kind=rp), intent(in) :: dxt(lx, lx)
1171 real(kind=rp), intent(in) :: dyt(lx, lx)
1172 real(kind=rp), intent(in) :: dzt(lx, lx)
1173 real(kind=rp) :: ur(lx, lx, lx)
1174 real(kind=rp) :: us(lx, lx, lx)
1175 real(kind=rp) :: ut(lx, lx, lx)
1176 real(kind=rp) :: wur(lx, lx, lx)
1177 real(kind=rp) :: wus(lx, lx, lx)
1178 real(kind=rp) :: wut(lx, lx, lx)
1179 integer :: e, c, i, j, k
1180
1181 !$omp do
1182 do e = 1, n
1183 c = ind(e)
1184 do j = 1, lx * lx
1185 do i = 1, lx
1186 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1187 + dx(i,2) * u(2,j,1,e) &
1188 + dx(i,3) * u(3,j,1,e) &
1189 + dx(i,4) * u(4,j,1,e) &
1190 + dx(i,5) * u(5,j,1,e) &
1191 + dx(i,6) * u(6,j,1,e) &
1192 + dx(i,7) * u(7,j,1,e) &
1193 + dx(i,8) * u(8,j,1,e) &
1194 + dx(i,9) * u(9,j,1,e)
1195 end do
1196 end do
1197
1198 do k = 1, lx
1199 do j = 1, lx
1200 do i = 1, lx
1201 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1202 + dy(j,2) * u(i,2,k,e) &
1203 + dy(j,3) * u(i,3,k,e) &
1204 + dy(j,4) * u(i,4,k,e) &
1205 + dy(j,5) * u(i,5,k,e) &
1206 + dy(j,6) * u(i,6,k,e) &
1207 + dy(j,7) * u(i,7,k,e) &
1208 + dy(j,8) * u(i,8,k,e) &
1209 + dy(j,9) * u(i,9,k,e)
1210 end do
1211 end do
1212 end do
1213
1214 do k = 1, lx
1215 do i = 1, lx*lx
1216 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1217 + dz(k,2) * u(i,1,2,e) &
1218 + dz(k,3) * u(i,1,3,e) &
1219 + dz(k,4) * u(i,1,4,e) &
1220 + dz(k,5) * u(i,1,5,e) &
1221 + dz(k,6) * u(i,1,6,e) &
1222 + dz(k,7) * u(i,1,7,e) &
1223 + dz(k,8) * u(i,1,8,e) &
1224 + dz(k,9) * u(i,1,9,e)
1225 end do
1226 end do
1227
1228 do i = 1, lx*lx*lx
1229 ur(i,1,1) = h1(i,1,1,e) &
1230 * ( g11(i,1,1,c) * wur(i,1,1) &
1231 + g12(i,1,1,c) * wus(i,1,1) &
1232 + g13(i,1,1,c) * wut(i,1,1) )
1233 us(i,1,1) = h1(i,1,1,e) &
1234 * ( g12(i,1,1,c) * wur(i,1,1) &
1235 + g22(i,1,1,c) * wus(i,1,1) &
1236 + g23(i,1,1,c) * wut(i,1,1) )
1237 ut(i,1,1) = h1(i,1,1,e) &
1238 * ( g13(i,1,1,c) * wur(i,1,1) &
1239 + g23(i,1,1,c) * wus(i,1,1) &
1240 + g33(i,1,1,c) * wut(i,1,1) )
1241 end do
1242
1243 do j = 1, lx*lx
1244 do i = 1, lx
1245 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1246 + dxt(i,2) * ur(2,j,1) &
1247 + dxt(i,3) * ur(3,j,1) &
1248 + dxt(i,4) * ur(4,j,1) &
1249 + dxt(i,5) * ur(5,j,1) &
1250 + dxt(i,6) * ur(6,j,1) &
1251 + dxt(i,7) * ur(7,j,1) &
1252 + dxt(i,8) * ur(8,j,1) &
1253 + dxt(i,9) * ur(9,j,1)
1254 end do
1255 end do
1256
1257 do k = 1, lx
1258 do j = 1, lx
1259 do i = 1, lx
1260 w(i,j,k,e) = w(i,j,k,e) &
1261 + dyt(j,1) * us(i,1,k) &
1262 + dyt(j,2) * us(i,2,k) &
1263 + dyt(j,3) * us(i,3,k) &
1264 + dyt(j,4) * us(i,4,k) &
1265 + dyt(j,5) * us(i,5,k) &
1266 + dyt(j,6) * us(i,6,k) &
1267 + dyt(j,7) * us(i,7,k) &
1268 + dyt(j,8) * us(i,8,k) &
1269 + dyt(j,9) * us(i,9,k)
1270 end do
1271 end do
1272 end do
1273
1274 do k = 1, lx
1275 do i = 1, lx*lx
1276 w(i,1,k,e) = w(i,1,k,e) &
1277 + dzt(k,1) * ut(i,1,1) &
1278 + dzt(k,2) * ut(i,1,2) &
1279 + dzt(k,3) * ut(i,1,3) &
1280 + dzt(k,4) * ut(i,1,4) &
1281 + dzt(k,5) * ut(i,1,5) &
1282 + dzt(k,6) * ut(i,1,6) &
1283 + dzt(k,7) * ut(i,1,7) &
1284 + dzt(k,8) * ut(i,1,8) &
1285 + dzt(k,9) * ut(i,1,9)
1286 end do
1287 end do
1288
1289 end do
1290 !$omp end do
1291 end subroutine ax_helm_compr_lx9
1292
1293 subroutine ax_helm_compr_lx8(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1294 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1295 integer, parameter :: lx = 8
1296 integer, intent(in) :: n, m
1297 integer, intent(in) :: ind(n)
1298 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1299 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1300 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1301 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1302 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1303 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1304 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1305 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1306 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1307 real(kind=rp), intent(in) :: dx(lx, lx)
1308 real(kind=rp), intent(in) :: dy(lx, lx)
1309 real(kind=rp), intent(in) :: dz(lx, lx)
1310 real(kind=rp), intent(in) :: dxt(lx, lx)
1311 real(kind=rp), intent(in) :: dyt(lx, lx)
1312 real(kind=rp), intent(in) :: dzt(lx, lx)
1313 real(kind=rp) :: ur(lx, lx, lx)
1314 real(kind=rp) :: us(lx, lx, lx)
1315 real(kind=rp) :: ut(lx, lx, lx)
1316 real(kind=rp) :: wur(lx, lx, lx)
1317 real(kind=rp) :: wus(lx, lx, lx)
1318 real(kind=rp) :: wut(lx, lx, lx)
1319 integer :: e, c, i, j, k
1320
1321 !$omp do
1322 do e = 1, n
1323 c = ind(e)
1324 do j = 1, lx * lx
1325 do i = 1, lx
1326 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1327 + dx(i,2) * u(2,j,1,e) &
1328 + dx(i,3) * u(3,j,1,e) &
1329 + dx(i,4) * u(4,j,1,e) &
1330 + dx(i,5) * u(5,j,1,e) &
1331 + dx(i,6) * u(6,j,1,e) &
1332 + dx(i,7) * u(7,j,1,e) &
1333 + dx(i,8) * u(8,j,1,e)
1334 end do
1335 end do
1336
1337 do k = 1, lx
1338 do j = 1, lx
1339 do i = 1, lx
1340 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1341 + dy(j,2) * u(i,2,k,e) &
1342 + dy(j,3) * u(i,3,k,e) &
1343 + dy(j,4) * u(i,4,k,e) &
1344 + dy(j,5) * u(i,5,k,e) &
1345 + dy(j,6) * u(i,6,k,e) &
1346 + dy(j,7) * u(i,7,k,e) &
1347 + dy(j,8) * u(i,8,k,e)
1348 end do
1349 end do
1350 end do
1351
1352 do k = 1, lx
1353 do i = 1, lx*lx
1354 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1355 + dz(k,2) * u(i,1,2,e) &
1356 + dz(k,3) * u(i,1,3,e) &
1357 + dz(k,4) * u(i,1,4,e) &
1358 + dz(k,5) * u(i,1,5,e) &
1359 + dz(k,6) * u(i,1,6,e) &
1360 + dz(k,7) * u(i,1,7,e) &
1361 + dz(k,8) * u(i,1,8,e)
1362 end do
1363 end do
1364
1365 do i = 1, lx*lx*lx
1366 ur(i,1,1) = h1(i,1,1,e) &
1367 * ( g11(i,1,1,c) * wur(i,1,1) &
1368 + g12(i,1,1,c) * wus(i,1,1) &
1369 + g13(i,1,1,c) * wut(i,1,1) )
1370 us(i,1,1) = h1(i,1,1,e) &
1371 * ( g12(i,1,1,c) * wur(i,1,1) &
1372 + g22(i,1,1,c) * wus(i,1,1) &
1373 + g23(i,1,1,c) * wut(i,1,1) )
1374 ut(i,1,1) = h1(i,1,1,e) &
1375 * ( g13(i,1,1,c) * wur(i,1,1) &
1376 + g23(i,1,1,c) * wus(i,1,1) &
1377 + g33(i,1,1,c) * wut(i,1,1) )
1378 end do
1379
1380 do j = 1, lx*lx
1381 do i = 1, lx
1382 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1383 + dxt(i,2) * ur(2,j,1) &
1384 + dxt(i,3) * ur(3,j,1) &
1385 + dxt(i,4) * ur(4,j,1) &
1386 + dxt(i,5) * ur(5,j,1) &
1387 + dxt(i,6) * ur(6,j,1) &
1388 + dxt(i,7) * ur(7,j,1) &
1389 + dxt(i,8) * ur(8,j,1)
1390 end do
1391 end do
1392
1393 do k = 1, lx
1394 do j = 1, lx
1395 do i = 1, lx
1396 w(i,j,k,e) = w(i,j,k,e) &
1397 + dyt(j,1) * us(i,1,k) &
1398 + dyt(j,2) * us(i,2,k) &
1399 + dyt(j,3) * us(i,3,k) &
1400 + dyt(j,4) * us(i,4,k) &
1401 + dyt(j,5) * us(i,5,k) &
1402 + dyt(j,6) * us(i,6,k) &
1403 + dyt(j,7) * us(i,7,k) &
1404 + dyt(j,8) * us(i,8,k)
1405 end do
1406 end do
1407 end do
1408
1409 do k = 1, lx
1410 do i = 1, lx*lx
1411 w(i,1,k,e) = w(i,1,k,e) &
1412 + dzt(k,1) * ut(i,1,1) &
1413 + dzt(k,2) * ut(i,1,2) &
1414 + dzt(k,3) * ut(i,1,3) &
1415 + dzt(k,4) * ut(i,1,4) &
1416 + dzt(k,5) * ut(i,1,5) &
1417 + dzt(k,6) * ut(i,1,6) &
1418 + dzt(k,7) * ut(i,1,7) &
1419 + dzt(k,8) * ut(i,1,8)
1420 end do
1421 end do
1422
1423 end do
1424 !$omp end do
1425 end subroutine ax_helm_compr_lx8
1426
1427 subroutine ax_helm_compr_lx7(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1428 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1429 integer, parameter :: lx = 7
1430 integer, intent(in) :: n, m
1431 integer, intent(in) :: ind(n)
1432 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1433 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1434 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1435 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1436 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1437 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1438 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1439 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1440 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1441 real(kind=rp), intent(in) :: dx(lx, lx)
1442 real(kind=rp), intent(in) :: dy(lx, lx)
1443 real(kind=rp), intent(in) :: dz(lx, lx)
1444 real(kind=rp), intent(in) :: dxt(lx, lx)
1445 real(kind=rp), intent(in) :: dyt(lx, lx)
1446 real(kind=rp), intent(in) :: dzt(lx, lx)
1447 real(kind=rp) :: ur(lx, lx, lx)
1448 real(kind=rp) :: us(lx, lx, lx)
1449 real(kind=rp) :: ut(lx, lx, lx)
1450 real(kind=rp) :: wur(lx, lx, lx)
1451 real(kind=rp) :: wus(lx, lx, lx)
1452 real(kind=rp) :: wut(lx, lx, lx)
1453 integer :: e, c, i, j, k
1454
1455 !$omp do
1456 do e = 1, n
1457 c = ind(e)
1458 do j = 1, lx * lx
1459 do i = 1, lx
1460 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1461 + dx(i,2) * u(2,j,1,e) &
1462 + dx(i,3) * u(3,j,1,e) &
1463 + dx(i,4) * u(4,j,1,e) &
1464 + dx(i,5) * u(5,j,1,e) &
1465 + dx(i,6) * u(6,j,1,e) &
1466 + dx(i,7) * u(7,j,1,e)
1467 end do
1468 end do
1469
1470 do k = 1, lx
1471 do j = 1, lx
1472 do i = 1, lx
1473 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1474 + dy(j,2) * u(i,2,k,e) &
1475 + dy(j,3) * u(i,3,k,e) &
1476 + dy(j,4) * u(i,4,k,e) &
1477 + dy(j,5) * u(i,5,k,e) &
1478 + dy(j,6) * u(i,6,k,e) &
1479 + dy(j,7) * u(i,7,k,e)
1480 end do
1481 end do
1482 end do
1483
1484 do k = 1, lx
1485 do i = 1, lx*lx
1486 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1487 + dz(k,2) * u(i,1,2,e) &
1488 + dz(k,3) * u(i,1,3,e) &
1489 + dz(k,4) * u(i,1,4,e) &
1490 + dz(k,5) * u(i,1,5,e) &
1491 + dz(k,6) * u(i,1,6,e) &
1492 + dz(k,7) * u(i,1,7,e)
1493 end do
1494 end do
1495
1496 do i = 1, lx*lx*lx
1497 ur(i,1,1) = h1(i,1,1,e) &
1498 * ( g11(i,1,1,c) * wur(i,1,1) &
1499 + g12(i,1,1,c) * wus(i,1,1) &
1500 + g13(i,1,1,c) * wut(i,1,1) )
1501 us(i,1,1) = h1(i,1,1,e) &
1502 * ( g12(i,1,1,c) * wur(i,1,1) &
1503 + g22(i,1,1,c) * wus(i,1,1) &
1504 + g23(i,1,1,c) * wut(i,1,1) )
1505 ut(i,1,1) = h1(i,1,1,e) &
1506 * ( g13(i,1,1,c) * wur(i,1,1) &
1507 + g23(i,1,1,c) * wus(i,1,1) &
1508 + g33(i,1,1,c) * wut(i,1,1) )
1509 end do
1510
1511 do j = 1, lx*lx
1512 do i = 1, lx
1513 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1514 + dxt(i,2) * ur(2,j,1) &
1515 + dxt(i,3) * ur(3,j,1) &
1516 + dxt(i,4) * ur(4,j,1) &
1517 + dxt(i,5) * ur(5,j,1) &
1518 + dxt(i,6) * ur(6,j,1) &
1519 + dxt(i,7) * ur(7,j,1)
1520 end do
1521 end do
1522
1523 do k = 1, lx
1524 do j = 1, lx
1525 do i = 1, lx
1526 w(i,j,k,e) = w(i,j,k,e) &
1527 + dyt(j,1) * us(i,1,k) &
1528 + dyt(j,2) * us(i,2,k) &
1529 + dyt(j,3) * us(i,3,k) &
1530 + dyt(j,4) * us(i,4,k) &
1531 + dyt(j,5) * us(i,5,k) &
1532 + dyt(j,6) * us(i,6,k) &
1533 + dyt(j,7) * us(i,7,k)
1534 end do
1535 end do
1536 end do
1537
1538 do k = 1, lx
1539 do i = 1, lx*lx
1540 w(i,1,k,e) = w(i,1,k,e) &
1541 + dzt(k,1) * ut(i,1,1) &
1542 + dzt(k,2) * ut(i,1,2) &
1543 + dzt(k,3) * ut(i,1,3) &
1544 + dzt(k,4) * ut(i,1,4) &
1545 + dzt(k,5) * ut(i,1,5) &
1546 + dzt(k,6) * ut(i,1,6) &
1547 + dzt(k,7) * ut(i,1,7)
1548 end do
1549 end do
1550
1551 end do
1552 !$omp end do
1553 end subroutine ax_helm_compr_lx7
1554
1555 subroutine ax_helm_compr_lx6(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1556 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1557 integer, parameter :: lx = 6
1558 integer, intent(in) :: n, m
1559 integer, intent(in) :: ind(n)
1560 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1561 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1562 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1563 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1564 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1565 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1566 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1567 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1568 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1569 real(kind=rp), intent(in) :: dx(lx, lx)
1570 real(kind=rp), intent(in) :: dy(lx, lx)
1571 real(kind=rp), intent(in) :: dz(lx, lx)
1572 real(kind=rp), intent(in) :: dxt(lx, lx)
1573 real(kind=rp), intent(in) :: dyt(lx, lx)
1574 real(kind=rp), intent(in) :: dzt(lx, lx)
1575 real(kind=rp) :: ur(lx, lx, lx)
1576 real(kind=rp) :: us(lx, lx, lx)
1577 real(kind=rp) :: ut(lx, lx, lx)
1578 real(kind=rp) :: wur(lx, lx, lx)
1579 real(kind=rp) :: wus(lx, lx, lx)
1580 real(kind=rp) :: wut(lx, lx, lx)
1581 integer :: e, c, i, j, k
1582
1583 !$omp do
1584 do e = 1, n
1585 c = ind(e)
1586 do j = 1, lx * lx
1587 do i = 1, lx
1588 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1589 + dx(i,2) * u(2,j,1,e) &
1590 + dx(i,3) * u(3,j,1,e) &
1591 + dx(i,4) * u(4,j,1,e) &
1592 + dx(i,5) * u(5,j,1,e) &
1593 + dx(i,6) * u(6,j,1,e)
1594 end do
1595 end do
1596
1597 do k = 1, lx
1598 do j = 1, lx
1599 do i = 1, lx
1600 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1601 + dy(j,2) * u(i,2,k,e) &
1602 + dy(j,3) * u(i,3,k,e) &
1603 + dy(j,4) * u(i,4,k,e) &
1604 + dy(j,5) * u(i,5,k,e) &
1605 + dy(j,6) * u(i,6,k,e)
1606 end do
1607 end do
1608 end do
1609
1610 do k = 1, lx
1611 do i = 1, lx*lx
1612 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1613 + dz(k,2) * u(i,1,2,e) &
1614 + dz(k,3) * u(i,1,3,e) &
1615 + dz(k,4) * u(i,1,4,e) &
1616 + dz(k,5) * u(i,1,5,e) &
1617 + dz(k,6) * u(i,1,6,e)
1618 end do
1619 end do
1620
1621 do i = 1, lx*lx*lx
1622 ur(i,1,1) = h1(i,1,1,e) &
1623 * ( g11(i,1,1,c) * wur(i,1,1) &
1624 + g12(i,1,1,c) * wus(i,1,1) &
1625 + g13(i,1,1,c) * wut(i,1,1) )
1626 us(i,1,1) = h1(i,1,1,e) &
1627 * ( g12(i,1,1,c) * wur(i,1,1) &
1628 + g22(i,1,1,c) * wus(i,1,1) &
1629 + g23(i,1,1,c) * wut(i,1,1) )
1630 ut(i,1,1) = h1(i,1,1,e) &
1631 * ( g13(i,1,1,c) * wur(i,1,1) &
1632 + g23(i,1,1,c) * wus(i,1,1) &
1633 + g33(i,1,1,c) * wut(i,1,1) )
1634 end do
1635
1636 do j = 1, lx*lx
1637 do i = 1, lx
1638 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1639 + dxt(i,2) * ur(2,j,1) &
1640 + dxt(i,3) * ur(3,j,1) &
1641 + dxt(i,4) * ur(4,j,1) &
1642 + dxt(i,5) * ur(5,j,1) &
1643 + dxt(i,6) * ur(6,j,1)
1644 end do
1645 end do
1646
1647 do k = 1, lx
1648 do j = 1, lx
1649 do i = 1, lx
1650 w(i,j,k,e) = w(i,j,k,e) &
1651 + dyt(j,1) * us(i,1,k) &
1652 + dyt(j,2) * us(i,2,k) &
1653 + dyt(j,3) * us(i,3,k) &
1654 + dyt(j,4) * us(i,4,k) &
1655 + dyt(j,5) * us(i,5,k) &
1656 + dyt(j,6) * us(i,6,k)
1657 end do
1658 end do
1659 end do
1660
1661 do k = 1, lx
1662 do i = 1, lx*lx
1663 w(i,1,k,e) = w(i,1,k,e) &
1664 + dzt(k,1) * ut(i,1,1) &
1665 + dzt(k,2) * ut(i,1,2) &
1666 + dzt(k,3) * ut(i,1,3) &
1667 + dzt(k,4) * ut(i,1,4) &
1668 + dzt(k,5) * ut(i,1,5) &
1669 + dzt(k,6) * ut(i,1,6)
1670 end do
1671 end do
1672
1673 end do
1674 !$omp end do
1675 end subroutine ax_helm_compr_lx6
1676
1677 subroutine ax_helm_compr_lx5(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1678 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1679 integer, parameter :: lx = 5
1680 integer, intent(in) :: n, m
1681 integer, intent(in) :: ind(n)
1682 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1683 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1684 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1685 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1686 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1687 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1688 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1689 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1690 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1691 real(kind=rp), intent(in) :: dx(lx, lx)
1692 real(kind=rp), intent(in) :: dy(lx, lx)
1693 real(kind=rp), intent(in) :: dz(lx, lx)
1694 real(kind=rp), intent(in) :: dxt(lx, lx)
1695 real(kind=rp), intent(in) :: dyt(lx, lx)
1696 real(kind=rp), intent(in) :: dzt(lx, lx)
1697 real(kind=rp) :: ur(lx, lx, lx)
1698 real(kind=rp) :: us(lx, lx, lx)
1699 real(kind=rp) :: ut(lx, lx, lx)
1700 real(kind=rp) :: wur(lx, lx, lx)
1701 real(kind=rp) :: wus(lx, lx, lx)
1702 real(kind=rp) :: wut(lx, lx, lx)
1703 integer :: e, c, i, j, k
1704
1705 !$omp do
1706 do e = 1, n
1707 c = ind(e)
1708 do j = 1, lx * lx
1709 do i = 1, lx
1710 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1711 + dx(i,2) * u(2,j,1,e) &
1712 + dx(i,3) * u(3,j,1,e) &
1713 + dx(i,4) * u(4,j,1,e) &
1714 + dx(i,5) * u(5,j,1,e)
1715 end do
1716 end do
1717
1718 do k = 1, lx
1719 do j = 1, lx
1720 do i = 1, lx
1721 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1722 + dy(j,2) * u(i,2,k,e) &
1723 + dy(j,3) * u(i,3,k,e) &
1724 + dy(j,4) * u(i,4,k,e) &
1725 + dy(j,5) * u(i,5,k,e)
1726 end do
1727 end do
1728 end do
1729
1730 do k = 1, lx
1731 do i = 1, lx*lx
1732 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1733 + dz(k,2) * u(i,1,2,e) &
1734 + dz(k,3) * u(i,1,3,e) &
1735 + dz(k,4) * u(i,1,4,e) &
1736 + dz(k,5) * u(i,1,5,e)
1737 end do
1738 end do
1739
1740 do i = 1, lx*lx*lx
1741 ur(i,1,1) = h1(i,1,1,e) &
1742 * ( g11(i,1,1,c) * wur(i,1,1) &
1743 + g12(i,1,1,c) * wus(i,1,1) &
1744 + g13(i,1,1,c) * wut(i,1,1) )
1745 us(i,1,1) = h1(i,1,1,e) &
1746 * ( g12(i,1,1,c) * wur(i,1,1) &
1747 + g22(i,1,1,c) * wus(i,1,1) &
1748 + g23(i,1,1,c) * wut(i,1,1) )
1749 ut(i,1,1) = h1(i,1,1,e) &
1750 * ( g13(i,1,1,c) * wur(i,1,1) &
1751 + g23(i,1,1,c) * wus(i,1,1) &
1752 + g33(i,1,1,c) * wut(i,1,1) )
1753 end do
1754
1755 do j = 1, lx*lx
1756 do i = 1, lx
1757 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1758 + dxt(i,2) * ur(2,j,1) &
1759 + dxt(i,3) * ur(3,j,1) &
1760 + dxt(i,4) * ur(4,j,1) &
1761 + dxt(i,5) * ur(5,j,1)
1762 end do
1763 end do
1764
1765 do k = 1, lx
1766 do j = 1, lx
1767 do i = 1, lx
1768 w(i,j,k,e) = w(i,j,k,e) &
1769 + dyt(j,1) * us(i,1,k) &
1770 + dyt(j,2) * us(i,2,k) &
1771 + dyt(j,3) * us(i,3,k) &
1772 + dyt(j,4) * us(i,4,k) &
1773 + dyt(j,5) * us(i,5,k)
1774 end do
1775 end do
1776 end do
1777
1778 do k = 1, lx
1779 do i = 1, lx*lx
1780 w(i,1,k,e) = w(i,1,k,e) &
1781 + dzt(k,1) * ut(i,1,1) &
1782 + dzt(k,2) * ut(i,1,2) &
1783 + dzt(k,3) * ut(i,1,3) &
1784 + dzt(k,4) * ut(i,1,4) &
1785 + dzt(k,5) * ut(i,1,5)
1786 end do
1787 end do
1788
1789 end do
1790 !$omp end do
1791 end subroutine ax_helm_compr_lx5
1792
1793 subroutine ax_helm_compr_lx4(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1794 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1795 integer, parameter :: lx = 4
1796 integer, intent(in) :: n, m
1797 integer, intent(in) :: ind(n)
1798 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1799 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1800 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1801 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1802 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1803 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1804 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1805 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1806 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1807 real(kind=rp), intent(in) :: dx(lx, lx)
1808 real(kind=rp), intent(in) :: dy(lx, lx)
1809 real(kind=rp), intent(in) :: dz(lx, lx)
1810 real(kind=rp), intent(in) :: dxt(lx, lx)
1811 real(kind=rp), intent(in) :: dyt(lx, lx)
1812 real(kind=rp), intent(in) :: dzt(lx, lx)
1813 real(kind=rp) :: ur(lx, lx, lx)
1814 real(kind=rp) :: us(lx, lx, lx)
1815 real(kind=rp) :: ut(lx, lx, lx)
1816 real(kind=rp) :: wur(lx, lx, lx)
1817 real(kind=rp) :: wus(lx, lx, lx)
1818 real(kind=rp) :: wut(lx, lx, lx)
1819 integer :: e, c, i, j, k
1820
1821 !$omp do
1822 do e = 1, n
1823 c = ind(e)
1824 do j = 1, lx * lx
1825 do i = 1, lx
1826 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1827 + dx(i,2) * u(2,j,1,e) &
1828 + dx(i,3) * u(3,j,1,e) &
1829 + dx(i,4) * u(4,j,1,e)
1830 end do
1831 end do
1832
1833 do k = 1, lx
1834 do j = 1, lx
1835 do i = 1, lx
1836 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1837 + dy(j,2) * u(i,2,k,e) &
1838 + dy(j,3) * u(i,3,k,e) &
1839 + dy(j,4) * u(i,4,k,e)
1840 end do
1841 end do
1842 end do
1843
1844 do k = 1, lx
1845 do i = 1, lx*lx
1846 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1847 + dz(k,2) * u(i,1,2,e) &
1848 + dz(k,3) * u(i,1,3,e) &
1849 + dz(k,4) * u(i,1,4,e)
1850 end do
1851 end do
1852
1853 do i = 1, lx*lx*lx
1854 ur(i,1,1) = h1(i,1,1,e) &
1855 * ( g11(i,1,1,c) * wur(i,1,1) &
1856 + g12(i,1,1,c) * wus(i,1,1) &
1857 + g13(i,1,1,c) * wut(i,1,1) )
1858 us(i,1,1) = h1(i,1,1,e) &
1859 * ( g12(i,1,1,c) * wur(i,1,1) &
1860 + g22(i,1,1,c) * wus(i,1,1) &
1861 + g23(i,1,1,c) * wut(i,1,1) )
1862 ut(i,1,1) = h1(i,1,1,e) &
1863 * ( g13(i,1,1,c) * wur(i,1,1) &
1864 + g23(i,1,1,c) * wus(i,1,1) &
1865 + g33(i,1,1,c) * wut(i,1,1) )
1866 end do
1867
1868 do j = 1, lx*lx
1869 do i = 1, lx
1870 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1871 + dxt(i,2) * ur(2,j,1) &
1872 + dxt(i,3) * ur(3,j,1) &
1873 + dxt(i,4) * ur(4,j,1)
1874 end do
1875 end do
1876
1877 do k = 1, lx
1878 do j = 1, lx
1879 do i = 1, lx
1880 w(i,j,k,e) = w(i,j,k,e) &
1881 + dyt(j,1) * us(i,1,k) &
1882 + dyt(j,2) * us(i,2,k) &
1883 + dyt(j,3) * us(i,3,k) &
1884 + dyt(j,4) * us(i,4,k)
1885 end do
1886 end do
1887 end do
1888
1889 do k = 1, lx
1890 do i = 1, lx*lx
1891 w(i,1,k,e) = w(i,1,k,e) &
1892 + dzt(k,1) * ut(i,1,1) &
1893 + dzt(k,2) * ut(i,1,2) &
1894 + dzt(k,3) * ut(i,1,3) &
1895 + dzt(k,4) * ut(i,1,4)
1896 end do
1897 end do
1898
1899 end do
1900 !$omp end do
1901 end subroutine ax_helm_compr_lx4
1902
1903 subroutine ax_helm_compr_lx3(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
1904 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
1905 integer, parameter :: lx = 3
1906 integer, intent(in) :: n, m
1907 integer, intent(in) :: ind(n)
1908 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
1909 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
1910 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
1911 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
1912 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
1913 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
1914 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
1915 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
1916 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
1917 real(kind=rp), intent(in) :: dx(lx, lx)
1918 real(kind=rp), intent(in) :: dy(lx, lx)
1919 real(kind=rp), intent(in) :: dz(lx, lx)
1920 real(kind=rp), intent(in) :: dxt(lx, lx)
1921 real(kind=rp), intent(in) :: dyt(lx, lx)
1922 real(kind=rp), intent(in) :: dzt(lx, lx)
1923 real(kind=rp) :: ur(lx, lx, lx)
1924 real(kind=rp) :: us(lx, lx, lx)
1925 real(kind=rp) :: ut(lx, lx, lx)
1926 real(kind=rp) :: wur(lx, lx, lx)
1927 real(kind=rp) :: wus(lx, lx, lx)
1928 real(kind=rp) :: wut(lx, lx, lx)
1929 integer :: e, c, i, j, k
1930
1931 !$omp do
1932 do e = 1, n
1933 c = ind(e)
1934 do j = 1, lx * lx
1935 do i = 1, lx
1936 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
1937 + dx(i,2) * u(2,j,1,e) &
1938 + dx(i,3) * u(3,j,1,e)
1939 end do
1940 end do
1941
1942 do k = 1, lx
1943 do j = 1, lx
1944 do i = 1, lx
1945 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
1946 + dy(j,2) * u(i,2,k,e) &
1947 + dy(j,3) * u(i,3,k,e)
1948 end do
1949 end do
1950 end do
1951
1952 do k = 1, lx
1953 do i = 1, lx*lx
1954 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
1955 + dz(k,2) * u(i,1,2,e) &
1956 + dz(k,3) * u(i,1,3,e)
1957 end do
1958 end do
1959
1960 do i = 1, lx*lx*lx
1961 ur(i,1,1) = h1(i,1,1,e) &
1962 * ( g11(i,1,1,c) * wur(i,1,1) &
1963 + g12(i,1,1,c) * wus(i,1,1) &
1964 + g13(i,1,1,c) * wut(i,1,1) )
1965 us(i,1,1) = h1(i,1,1,e) &
1966 * ( g12(i,1,1,c) * wur(i,1,1) &
1967 + g22(i,1,1,c) * wus(i,1,1) &
1968 + g23(i,1,1,c) * wut(i,1,1) )
1969 ut(i,1,1) = h1(i,1,1,e) &
1970 * ( g13(i,1,1,c) * wur(i,1,1) &
1971 + g23(i,1,1,c) * wus(i,1,1) &
1972 + g33(i,1,1,c) * wut(i,1,1) )
1973 end do
1974
1975 do j = 1, lx*lx
1976 do i = 1, lx
1977 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
1978 + dxt(i,2) * ur(2,j,1) &
1979 + dxt(i,3) * ur(3,j,1)
1980 end do
1981 end do
1982
1983 do k = 1, lx
1984 do j = 1, lx
1985 do i = 1, lx
1986 w(i,j,k,e) = w(i,j,k,e) &
1987 + dyt(j,1) * us(i,1,k) &
1988 + dyt(j,2) * us(i,2,k) &
1989 + dyt(j,3) * us(i,3,k)
1990 end do
1991 end do
1992 end do
1993
1994 do k = 1, lx
1995 do i = 1, lx*lx
1996 w(i,1,k,e) = w(i,1,k,e) &
1997 + dzt(k,1) * ut(i,1,1) &
1998 + dzt(k,2) * ut(i,1,2) &
1999 + dzt(k,3) * ut(i,1,3)
2000 end do
2001 end do
2002
2003 end do
2004 !$omp end do
2005 end subroutine ax_helm_compr_lx3
2006
2007 subroutine ax_helm_compr_lx2(w, u, Dx, Dy, Dz, Dxt, Dyt, Dzt, &
2008 h1, G11, G22, G33, G12, G13, G23, ind, n, m)
2009 integer, parameter :: lx = 2
2010 integer, intent(in) :: n, m
2011 integer, intent(in) :: ind(n)
2012 real(kind=rp), intent(inout) :: w(lx, lx, lx, n)
2013 real(kind=rp), intent(in) :: u(lx, lx, lx, n)
2014 real(kind=rp), intent(in) :: h1(lx, lx, lx, n)
2015 real(kind=rp), intent(in) :: g11(lx, lx, lx, m)
2016 real(kind=rp), intent(in) :: g22(lx, lx, lx, m)
2017 real(kind=rp), intent(in) :: g33(lx, lx, lx, m)
2018 real(kind=rp), intent(in) :: g12(lx, lx, lx, m)
2019 real(kind=rp), intent(in) :: g13(lx, lx, lx, m)
2020 real(kind=rp), intent(in) :: g23(lx, lx, lx, m)
2021 real(kind=rp), intent(in) :: dx(lx, lx)
2022 real(kind=rp), intent(in) :: dy(lx, lx)
2023 real(kind=rp), intent(in) :: dz(lx, lx)
2024 real(kind=rp), intent(in) :: dxt(lx, lx)
2025 real(kind=rp), intent(in) :: dyt(lx, lx)
2026 real(kind=rp), intent(in) :: dzt(lx, lx)
2027 real(kind=rp) :: ur(lx, lx, lx)
2028 real(kind=rp) :: us(lx, lx, lx)
2029 real(kind=rp) :: ut(lx, lx, lx)
2030 real(kind=rp) :: wur(lx, lx, lx)
2031 real(kind=rp) :: wus(lx, lx, lx)
2032 real(kind=rp) :: wut(lx, lx, lx)
2033 integer :: e, c, i, j, k
2034
2035 !$omp do
2036 do e = 1, n
2037 c = ind(e)
2038 do j = 1, lx * lx
2039 do i = 1, lx
2040 wur(i,j,1) = dx(i,1) * u(1,j,1,e) &
2041 + dx(i,2) * u(2,j,1,e)
2042 end do
2043 end do
2044
2045 do k = 1, lx
2046 do j = 1, lx
2047 do i = 1, lx
2048 wus(i,j,k) = dy(j,1) * u(i,1,k,e) &
2049 + dy(j,2) * u(i,2,k,e)
2050 end do
2051 end do
2052 end do
2053
2054 do k = 1, lx
2055 do i = 1, lx*lx
2056 wut(i,1,k) = dz(k,1) * u(i,1,1,e) &
2057 + dz(k,2) * u(i,1,2,e)
2058 end do
2059 end do
2060
2061 do i = 1, lx*lx*lx
2062 ur(i,1,1) = h1(i,1,1,e) &
2063 * ( g11(i,1,1,c) * wur(i,1,1) &
2064 + g12(i,1,1,c) * wus(i,1,1) &
2065 + g13(i,1,1,c) * wut(i,1,1) )
2066 us(i,1,1) = h1(i,1,1,e) &
2067 * ( g12(i,1,1,c) * wur(i,1,1) &
2068 + g22(i,1,1,c) * wus(i,1,1) &
2069 + g23(i,1,1,c) * wut(i,1,1) )
2070 ut(i,1,1) = h1(i,1,1,e) &
2071 * ( g13(i,1,1,c) * wur(i,1,1) &
2072 + g23(i,1,1,c) * wus(i,1,1) &
2073 + g33(i,1,1,c) * wut(i,1,1) )
2074 end do
2075
2076 do j = 1, lx*lx
2077 do i = 1, lx
2078 w(i,j,1,e) = dxt(i,1) * ur(1,j,1) &
2079 + dxt(i,2) * ur(2,j,1)
2080 end do
2081 end do
2082
2083 do k = 1, lx
2084 do j = 1, lx
2085 do i = 1, lx
2086 w(i,j,k,e) = w(i,j,k,e) &
2087 + dyt(j,1) * us(i,1,k) &
2088 + dyt(j,2) * us(i,2,k)
2089 end do
2090 end do
2091 end do
2092
2093 do k = 1, lx
2094 do i = 1, lx*lx
2095 w(i,1,k,e) = w(i,1,k,e) &
2096 + dzt(k,1) * ut(i,1,1) &
2097 + dzt(k,2) * ut(i,1,2)
2098 end do
2099 end do
2100
2101 end do
2102 !$omp end do
2103 end subroutine ax_helm_compr_lx2
2104
2105end module ax_helm_compr_cpu
CPU Helmholtz matrix-vector product on compressed geometric factors.
subroutine ax_helm_compr_compute(this, w, u, coef, msh, xh)
Compute product, taking three components of a vector field in an uncoupled manner.
Coefficients.
Definition coef.f90:34
Defines a mesh.
Definition mesh.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:14
Defines a function space.
Definition space.f90:34
Utilities.
Definition utils.f90:35
CPU matrix-vector product for a Helmholtz problem, reading the geometric factors from their compress...
CPU matrix-vector product for a Helmholtz problem.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:136
The function space for the SEM solution fields.
Definition space.f90:64