Neko  0.9.0
A portable framework for high-order spectral element flow simulations
sx_convect_scalar.f90
Go to the documentation of this file.
1 ! Copyright (c) 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 !
33 submodule(opr_sx) sx_convect_scalar
34  use num_types, only : rp
35  use math, only : col2
36  implicit none
37 
38 contains
39 
40  module subroutine opr_sx_convect_scalar(du, u, c, xh_gll, xh_gl, coef_gll, &
41  coef_gl, gll_to_gl)
42  type(space_t), intent(in) :: Xh_GL
43  type(space_t), intent(in) :: Xh_GLL
44  type(coef_t), intent(in) :: coef_GLL
45  type(coef_t), intent(in) :: coef_GL
46  type(interpolator_t), intent(inout) :: GLL_to_GL
47  real(kind=rp), intent(inout) :: &
48  du(xh_gll%lx, xh_gll%ly, xh_gll%lz, coef_gl%msh%nelv)
49  real(kind=rp), intent(inout) :: &
50  u(xh_gl%lx, xh_gl%lx, xh_gl%lx, coef_gl%msh%nelv)
51  real(kind=rp), intent(inout) :: c(xh_gl%lxyz, coef_gl%msh%nelv, 3)
52  associate(dx => xh_gl%dx, dy => xh_gl%dy, dz => xh_gl%dz, &
53  lx => xh_gl%lx, nelv => coef_gl%msh%nelv)
54 
55  select case (lx)
56  case (18)
57  call sx_convect_scalar_lx18(du, u, c, dx, dy, dz, &
58  xh_gll, coef_gll, gll_to_gl, nelv)
59  case (17)
60  call sx_convect_scalar_lx17(du, u, c, dx, dy, dz, &
61  xh_gll, coef_gll, gll_to_gl, nelv)
62  case (16)
63  call sx_convect_scalar_lx16(du, u, c, dx, dy, dz, &
64  xh_gll, coef_gll, gll_to_gl, nelv)
65  case (15)
66  call sx_convect_scalar_lx15(du, u, c, dx, dy, dz, &
67  xh_gll, coef_gll, gll_to_gl, nelv)
68  case (14)
69  call sx_convect_scalar_lx14(du, u, c, dx, dy, dz, &
70  xh_gll, coef_gll, gll_to_gl, nelv)
71  case (13)
72  call sx_convect_scalar_lx13(du, u, c, dx, dy, dz, &
73  xh_gll, coef_gll, gll_to_gl, nelv)
74  case (12)
75  call sx_convect_scalar_lx12(du, u, c, dx, dy, dz, &
76  xh_gll, coef_gll, gll_to_gl, nelv)
77  case (11)
78  call sx_convect_scalar_lx11(du, u, c, dx, dy, dz, &
79  xh_gll, coef_gll, gll_to_gl, nelv)
80  case (10)
81  call sx_convect_scalar_lx10(du, u, c, dx, dy, dz, &
82  xh_gll, coef_gll, gll_to_gl, nelv)
83  case (9)
84  call sx_convect_scalar_lx9(du, u, c, dx, dy, dz, &
85  xh_gll, coef_gll, gll_to_gl, nelv)
86  case (8)
87  call sx_convect_scalar_lx8(du, u, c, dx, dy, dz, &
88  xh_gll, coef_gll, gll_to_gl, nelv)
89  case (7)
90  call sx_convect_scalar_lx7(du, u, c, dx, dy, dz, &
91  xh_gll, coef_gll, gll_to_gl, nelv)
92  case (6)
93  call sx_convect_scalar_lx6(du, u, c, dx, dy, dz, &
94  xh_gll, coef_gll, gll_to_gl, nelv)
95  case (5)
96  call sx_convect_scalar_lx5(du, u, c, dx, dy, dz, &
97  xh_gll, coef_gll, gll_to_gl, nelv)
98  case (4)
99  call sx_convect_scalar_lx4(du, u, c, dx, dy, dz, &
100  xh_gll, coef_gll, gll_to_gl, nelv)
101  case (3)
102  call sx_convect_scalar_lx3(du, u, c, dx, dy, dz, &
103  xh_gll, coef_gll, gll_to_gl, nelv)
104  case (2)
105  call sx_convect_scalar_lx2(du, u, c, dx, dy, dz, &
106  xh_gll, coef_gll, gll_to_gl, nelv)
107  case default
108  call sx_convect_scalar_lx(du, u, c, dx, dy, dz, &
109  xh_gll, coef_gll, gll_to_gl, nelv, lx)
110  end select
111  end associate
112 
113  end subroutine opr_sx_convect_scalar
114 
115  subroutine sx_convect_scalar_lx(du, u, c, dx, dy, dz, Xh_GLL, &
116  coef_GLL, GLL_to_GL, nelv, lx)
117  integer, intent(in) :: nelv, lx
118  type(space_t), intent(in) :: Xh_GLL
119  type(coef_t), intent(in) :: coef_GLL
120  type(interpolator_t), intent(inout) :: GLL_to_GL
121  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
122  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
123  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
124  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
125  real(kind=rp) :: ur(lx, lx, lx, nelv)
126  real(kind=rp) :: us(lx, lx, lx, nelv)
127  real(kind=rp) :: ut(lx, lx, lx, nelv)
128  real(kind=rp) :: ud(lx, lx, lx, nelv)
129  real(kind=rp) :: wr, ws, wt
130  integer :: e, i, j, k, jj, kk, n_GLL
131  n_gll = nelv * xh_gll%lxyz
132  do i = 1, lx
133  do jj = 1, lx * lx * nelv
134  wr = 0d0
135  do kk = 1, lx
136  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
137  end do
138  ur(i, jj, 1, 1) = wr
139  end do
140  end do
141 
142  do k = 1, lx
143  do i = 1, lx
144  do j = 1, lx
145  do e = 1, nelv
146  ws = 0d0
147  !NEC$ unroll_completely
148  do kk = 1, lx
149  ws = ws + dy(j, kk) * u(i, kk, k, e)
150  end do
151  us(i,j,k,e) = ws
152  end do
153  end do
154  end do
155  end do
156 
157  do j = 1, lx
158  do i = 1, lx
159  do k = 1, lx
160  do e = 1, nelv
161  wt = 0d0
162  !NEC$ unroll_completely
163  do kk = 1, lx
164  wt = wt + dz(k, kk) * u(i, j, kk, e)
165  end do
166  ut(i,j,k,e) = wt
167  end do
168  end do
169  end do
170  end do
171 
172  do i = 1, lx * lx * lx
173  do e = 1, nelv
174  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
175  + c(i,e,2) * us(i,1,1,e) &
176  + c(i,e,3) * ut(i,1,1,e) )
177  end do
178  end do
179  call gll_to_gl%map(du, ud, nelv, xh_gll)
180  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
181  call col2(du, coef_gll%Binv, n_gll)
182 
183  end subroutine sx_convect_scalar_lx
184 
185  subroutine sx_convect_scalar_lx18(du, u, c, dx, dy, dz, Xh_GLL, &
186  coef_GLL, GLL_to_GL, nelv)
187  integer, parameter :: lx = 18
188  integer, intent(in) :: nelv
189  type(space_t), intent(in) :: Xh_GLL
190  type(coef_t), intent(in) :: coef_GLL
191  type(interpolator_t), intent(inout) :: GLL_to_GL
192  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
193  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
194  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
195  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
196  real(kind=rp) :: ur(lx, lx, lx, nelv)
197  real(kind=rp) :: us(lx, lx, lx, nelv)
198  real(kind=rp) :: ut(lx, lx, lx, nelv)
199  real(kind=rp) :: ud(lx, lx, lx, nelv)
200  real(kind=rp) :: wr, ws, wt
201  integer :: e, i, j, k, jj, kk, n_GLL
202  n_gll = nelv * xh_gll%lxyz
203  do i = 1, lx
204  do jj = 1, lx * lx * nelv
205  wr = 0d0
206  do kk = 1, lx
207  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
208  end do
209  ur(i, jj, 1, 1) = wr
210  end do
211  end do
212 
213  do k = 1, lx
214  do i = 1, lx
215  do j = 1, lx
216  do e = 1, nelv
217  ws = 0d0
218  !NEC$ unroll_completely
219  do kk = 1, lx
220  ws = ws + dy(j, kk) * u(i, kk, k, e)
221  end do
222  us(i,j,k,e) = ws
223  end do
224  end do
225  end do
226  end do
227 
228  do j = 1, lx
229  do i = 1, lx
230  do k = 1, lx
231  do e = 1, nelv
232  wt = 0d0
233  !NEC$ unroll_completely
234  do kk = 1, lx
235  wt = wt + dz(k, kk) * u(i, j, kk, e)
236  end do
237  ut(i,j,k,e) = wt
238  end do
239  end do
240  end do
241  end do
242 
243  do i = 1, lx * lx * lx
244  do e = 1, nelv
245  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
246  + c(i,e,2) * us(i,1,1,e) &
247  + c(i,e,3) * ut(i,1,1,e) )
248  end do
249  end do
250  call gll_to_gl%map(du, ud, nelv, xh_gll)
251  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
252  call col2(du, coef_gll%Binv, n_gll)
253 
254  end subroutine sx_convect_scalar_lx18
255 
256  subroutine sx_convect_scalar_lx17(du, u, c, dx, dy, dz, Xh_GLL, &
257  coef_GLL, GLL_to_GL, nelv)
258  integer, parameter :: lx = 17
259  integer, intent(in) :: nelv
260  type(space_t), intent(in) :: Xh_GLL
261  type(coef_t), intent(in) :: coef_GLL
262  type(interpolator_t), intent(inout) :: GLL_to_GL
263  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
264  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
265  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
266  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
267  real(kind=rp) :: ur(lx, lx, lx, nelv)
268  real(kind=rp) :: us(lx, lx, lx, nelv)
269  real(kind=rp) :: ut(lx, lx, lx, nelv)
270  real(kind=rp) :: ud(lx, lx, lx, nelv)
271  real(kind=rp) :: wr, ws, wt
272  integer :: e, i, j, k, jj, kk, n_GLL
273  n_gll = nelv * xh_gll%lxyz
274  do i = 1, lx
275  do jj = 1, lx * lx * nelv
276  wr = 0d0
277  do kk = 1, lx
278  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
279  end do
280  ur(i, jj, 1, 1) = wr
281  end do
282  end do
283 
284  do k = 1, lx
285  do i = 1, lx
286  do j = 1, lx
287  do e = 1, nelv
288  ws = 0d0
289  !NEC$ unroll_completely
290  do kk = 1, lx
291  ws = ws + dy(j, kk) * u(i, kk, k, e)
292  end do
293  us(i,j,k,e) = ws
294  end do
295  end do
296  end do
297  end do
298 
299  do j = 1, lx
300  do i = 1, lx
301  do k = 1, lx
302  do e = 1, nelv
303  wt = 0d0
304  !NEC$ unroll_completely
305  do kk = 1, lx
306  wt = wt + dz(k, kk) * u(i, j, kk, e)
307  end do
308  ut(i,j,k,e) = wt
309  end do
310  end do
311  end do
312  end do
313 
314  do i = 1, lx * lx * lx
315  do e = 1, nelv
316  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
317  + c(i,e,2) * us(i,1,1,e) &
318  + c(i,e,3) * ut(i,1,1,e) )
319  end do
320  end do
321  call gll_to_gl%map(du, ud, nelv, xh_gll)
322  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
323  call col2(du, coef_gll%Binv, n_gll)
324 
325  end subroutine sx_convect_scalar_lx17
326 
327  subroutine sx_convect_scalar_lx16(du, u, c, dx, dy, dz, Xh_GLL, &
328  coef_GLL, GLL_to_GL, nelv)
329  integer, parameter :: lx = 16
330  integer, intent(in) :: nelv
331  type(space_t), intent(in) :: Xh_GLL
332  type(coef_t), intent(in) :: coef_GLL
333  type(interpolator_t), intent(inout) :: GLL_to_GL
334  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
335  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
336  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
337  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
338  real(kind=rp) :: ur(lx, lx, lx, nelv)
339  real(kind=rp) :: us(lx, lx, lx, nelv)
340  real(kind=rp) :: ut(lx, lx, lx, nelv)
341  real(kind=rp) :: ud(lx, lx, lx, nelv)
342  real(kind=rp) :: wr, ws, wt
343  integer :: e, i, j, k, jj, kk, n_GLL
344  n_gll = nelv * xh_gll%lxyz
345  do i = 1, lx
346  do jj = 1, lx * lx * nelv
347  wr = 0d0
348  do kk = 1, lx
349  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
350  end do
351  ur(i, jj, 1, 1) = wr
352  end do
353  end do
354 
355  do k = 1, lx
356  do i = 1, lx
357  do j = 1, lx
358  do e = 1, nelv
359  ws = 0d0
360  !NEC$ unroll_completely
361  do kk = 1, lx
362  ws = ws + dy(j, kk) * u(i, kk, k, e)
363  end do
364  us(i,j,k,e) = ws
365  end do
366  end do
367  end do
368  end do
369 
370  do j = 1, lx
371  do i = 1, lx
372  do k = 1, lx
373  do e = 1, nelv
374  wt = 0d0
375  !NEC$ unroll_completely
376  do kk = 1, lx
377  wt = wt + dz(k, kk) * u(i, j, kk, e)
378  end do
379  ut(i,j,k,e) = wt
380  end do
381  end do
382  end do
383  end do
384 
385  do i = 1, lx * lx * lx
386  do e = 1, nelv
387  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
388  + c(i,e,2) * us(i,1,1,e) &
389  + c(i,e,3) * ut(i,1,1,e) )
390  end do
391  end do
392  call gll_to_gl%map(du, ud, nelv, xh_gll)
393  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
394  call col2(du, coef_gll%Binv, n_gll)
395 
396  end subroutine sx_convect_scalar_lx16
397 
398  subroutine sx_convect_scalar_lx15(du, u, c, dx, dy, dz, Xh_GLL, &
399  coef_GLL, GLL_to_GL, nelv)
400  integer, parameter :: lx = 15
401  integer, intent(in) :: nelv
402  type(space_t), intent(in) :: Xh_GLL
403  type(coef_t), intent(in) :: coef_GLL
404  type(interpolator_t), intent(inout) :: GLL_to_GL
405  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
406  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
407  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
408  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
409  real(kind=rp) :: ur(lx, lx, lx, nelv)
410  real(kind=rp) :: us(lx, lx, lx, nelv)
411  real(kind=rp) :: ut(lx, lx, lx, nelv)
412  real(kind=rp) :: ud(lx, lx, lx, nelv)
413  real(kind=rp) :: wr, ws, wt
414  integer :: e, i, j, k, jj, kk, n_GLL
415  n_gll = nelv * xh_gll%lxyz
416  do i = 1, lx
417  do jj = 1, lx * lx * nelv
418  wr = 0d0
419  do kk = 1, lx
420  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
421  end do
422  ur(i, jj, 1, 1) = wr
423  end do
424  end do
425 
426  do k = 1, lx
427  do i = 1, lx
428  do j = 1, lx
429  do e = 1, nelv
430  ws = 0d0
431  !NEC$ unroll_completely
432  do kk = 1, lx
433  ws = ws + dy(j, kk) * u(i, kk, k, e)
434  end do
435  us(i,j,k,e) = ws
436  end do
437  end do
438  end do
439  end do
440 
441  do j = 1, lx
442  do i = 1, lx
443  do k = 1, lx
444  do e = 1, nelv
445  wt = 0d0
446  !NEC$ unroll_completely
447  do kk = 1, lx
448  wt = wt + dz(k, kk) * u(i, j, kk, e)
449  end do
450  ut(i,j,k,e) = wt
451  end do
452  end do
453  end do
454  end do
455 
456  do i = 1, lx * lx * lx
457  do e = 1, nelv
458  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
459  + c(i,e,2) * us(i,1,1,e) &
460  + c(i,e,3) * ut(i,1,1,e) )
461  end do
462  end do
463  call gll_to_gl%map(du, ud, nelv, xh_gll)
464  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
465  call col2(du, coef_gll%Binv, n_gll)
466 
467  end subroutine sx_convect_scalar_lx15
468 
469  subroutine sx_convect_scalar_lx14(du, u, c, dx, dy, dz, Xh_GLL, &
470  coef_GLL, GLL_to_GL, nelv)
471  integer, parameter :: lx = 14
472  integer, intent(in) :: nelv
473  type(space_t), intent(in) :: Xh_GLL
474  type(coef_t), intent(in) :: coef_GLL
475  type(interpolator_t), intent(inout) :: GLL_to_GL
476  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
477  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
478  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
479  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
480  real(kind=rp) :: ur(lx, lx, lx, nelv)
481  real(kind=rp) :: us(lx, lx, lx, nelv)
482  real(kind=rp) :: ut(lx, lx, lx, nelv)
483  real(kind=rp) :: ud(lx, lx, lx, nelv)
484  real(kind=rp) :: wr, ws, wt
485  integer :: e, i, j, k, jj, kk, n_GLL
486  n_gll = nelv * xh_gll%lxyz
487  do i = 1, lx
488  do jj = 1, lx * lx * nelv
489  wr = 0d0
490  do kk = 1, lx
491  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
492  end do
493  ur(i, jj, 1, 1) = wr
494  end do
495  end do
496 
497  do k = 1, lx
498  do i = 1, lx
499  do j = 1, lx
500  do e = 1, nelv
501  ws = 0d0
502  !NEC$ unroll_completely
503  do kk = 1, lx
504  ws = ws + dy(j, kk) * u(i, kk, k, e)
505  end do
506  us(i,j,k,e) = ws
507  end do
508  end do
509  end do
510  end do
511 
512  do j = 1, lx
513  do i = 1, lx
514  do k = 1, lx
515  do e = 1, nelv
516  wt = 0d0
517  !NEC$ unroll_completely
518  do kk = 1, lx
519  wt = wt + dz(k, kk) * u(i, j, kk, e)
520  end do
521  ut(i,j,k,e) = wt
522  end do
523  end do
524  end do
525  end do
526 
527  do i = 1, lx * lx * lx
528  do e = 1, nelv
529  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
530  + c(i,e,2) * us(i,1,1,e) &
531  + c(i,e,3) * ut(i,1,1,e) )
532  end do
533  end do
534  call gll_to_gl%map(du, ud, nelv, xh_gll)
535  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
536  call col2(du, coef_gll%Binv, n_gll)
537 
538  end subroutine sx_convect_scalar_lx14
539 
540  subroutine sx_convect_scalar_lx13(du, u, c, dx, dy, dz, Xh_GLL, &
541  coef_GLL, GLL_to_GL, nelv)
542  integer, parameter :: lx = 13
543  integer, intent(in) :: nelv
544  type(space_t), intent(in) :: Xh_GLL
545  type(coef_t), intent(in) :: coef_GLL
546  type(interpolator_t), intent(inout) :: GLL_to_GL
547  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
548  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
549  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
550  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
551  real(kind=rp) :: ur(lx, lx, lx, nelv)
552  real(kind=rp) :: us(lx, lx, lx, nelv)
553  real(kind=rp) :: ut(lx, lx, lx, nelv)
554  real(kind=rp) :: ud(lx, lx, lx, nelv)
555  real(kind=rp) :: wr, ws, wt
556  integer :: e, i, j, k, jj, kk, n_GLL
557  n_gll = nelv * xh_gll%lxyz
558  do i = 1, lx
559  do jj = 1, lx * lx * nelv
560  wr = 0d0
561  do kk = 1, lx
562  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
563  end do
564  ur(i, jj, 1, 1) = wr
565  end do
566  end do
567 
568  do k = 1, lx
569  do i = 1, lx
570  do j = 1, lx
571  do e = 1, nelv
572  ws = 0d0
573  !NEC$ unroll_completely
574  do kk = 1, lx
575  ws = ws + dy(j, kk) * u(i, kk, k, e)
576  end do
577  us(i,j,k,e) = ws
578  end do
579  end do
580  end do
581  end do
582 
583  do j = 1, lx
584  do i = 1, lx
585  do k = 1, lx
586  do e = 1, nelv
587  wt = 0d0
588  !NEC$ unroll_completely
589  do kk = 1, lx
590  wt = wt + dz(k, kk) * u(i, j, kk, e)
591  end do
592  ut(i,j,k,e) = wt
593  end do
594  end do
595  end do
596  end do
597 
598  do i = 1, lx * lx * lx
599  do e = 1, nelv
600  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
601  + c(i,e,2) * us(i,1,1,e) &
602  + c(i,e,3) * ut(i,1,1,e) )
603  end do
604  end do
605  call gll_to_gl%map(du, ud, nelv, xh_gll)
606  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
607  call col2(du, coef_gll%Binv, n_gll)
608 
609  end subroutine sx_convect_scalar_lx13
610 
611  subroutine sx_convect_scalar_lx12(du, u, c, dx, dy, dz, Xh_GLL, &
612  coef_GLL, GLL_to_GL, nelv)
613  integer, parameter :: lx = 12
614  integer, intent(in) :: nelv
615  type(space_t), intent(in) :: Xh_GLL
616  type(coef_t), intent(in) :: coef_GLL
617  type(interpolator_t), intent(inout) :: GLL_to_GL
618  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
619  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
620  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
621  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
622  real(kind=rp) :: ur(lx, lx, lx, nelv)
623  real(kind=rp) :: us(lx, lx, lx, nelv)
624  real(kind=rp) :: ut(lx, lx, lx, nelv)
625  real(kind=rp) :: ud(lx, lx, lx, nelv)
626  real(kind=rp) :: wr, ws, wt
627  integer :: e, i, j, k, jj, kk, n_GLL
628  n_gll = nelv * xh_gll%lxyz
629  do i = 1, lx
630  do jj = 1, lx * lx * nelv
631  wr = 0d0
632  do kk = 1, lx
633  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
634  end do
635  ur(i, jj, 1, 1) = wr
636  end do
637  end do
638 
639  do k = 1, lx
640  do i = 1, lx
641  do j = 1, lx
642  do e = 1, nelv
643  ws = 0d0
644  !NEC$ unroll_completely
645  do kk = 1, lx
646  ws = ws + dy(j, kk) * u(i, kk, k, e)
647  end do
648  us(i,j,k,e) = ws
649  end do
650  end do
651  end do
652  end do
653 
654  do j = 1, lx
655  do i = 1, lx
656  do k = 1, lx
657  do e = 1, nelv
658  wt = 0d0
659  !NEC$ unroll_completely
660  do kk = 1, lx
661  wt = wt + dz(k, kk) * u(i, j, kk, e)
662  end do
663  ut(i,j,k,e) = wt
664  end do
665  end do
666  end do
667  end do
668 
669  do i = 1, lx * lx * lx
670  do e = 1, nelv
671  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
672  + c(i,e,2) * us(i,1,1,e) &
673  + c(i,e,3) * ut(i,1,1,e) )
674  end do
675  end do
676  call gll_to_gl%map(du, ud, nelv, xh_gll)
677  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
678  call col2(du, coef_gll%Binv, n_gll)
679 
680  end subroutine sx_convect_scalar_lx12
681 
682  subroutine sx_convect_scalar_lx11(du, u, c, dx, dy, dz, Xh_GLL, &
683  coef_GLL, GLL_to_GL, nelv)
684  integer, parameter :: lx = 11
685  integer, intent(in) :: nelv
686  type(space_t), intent(in) :: Xh_GLL
687  type(coef_t), intent(in) :: coef_GLL
688  type(interpolator_t), intent(inout) :: GLL_to_GL
689  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
690  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
691  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
692  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
693  real(kind=rp) :: ur(lx, lx, lx, nelv)
694  real(kind=rp) :: us(lx, lx, lx, nelv)
695  real(kind=rp) :: ut(lx, lx, lx, nelv)
696  real(kind=rp) :: ud(lx, lx, lx, nelv)
697  real(kind=rp) :: wr, ws, wt
698  integer :: e, i, j, k, jj, kk, n_GLL
699  n_gll = nelv * xh_gll%lxyz
700  do i = 1, lx
701  do jj = 1, lx * lx * nelv
702  wr = 0d0
703  do kk = 1, lx
704  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
705  end do
706  ur(i, jj, 1, 1) = wr
707  end do
708  end do
709 
710  do k = 1, lx
711  do i = 1, lx
712  do j = 1, lx
713  do e = 1, nelv
714  ws = 0d0
715  !NEC$ unroll_completely
716  do kk = 1, lx
717  ws = ws + dy(j, kk) * u(i, kk, k, e)
718  end do
719  us(i,j,k,e) = ws
720  end do
721  end do
722  end do
723  end do
724 
725  do j = 1, lx
726  do i = 1, lx
727  do k = 1, lx
728  do e = 1, nelv
729  wt = 0d0
730  !NEC$ unroll_completely
731  do kk = 1, lx
732  wt = wt + dz(k, kk) * u(i, j, kk, e)
733  end do
734  ut(i,j,k,e) = wt
735  end do
736  end do
737  end do
738  end do
739 
740  do i = 1, lx * lx * lx
741  do e = 1, nelv
742  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
743  + c(i,e,2) * us(i,1,1,e) &
744  + c(i,e,3) * ut(i,1,1,e) )
745  end do
746  end do
747  call gll_to_gl%map(du, ud, nelv, xh_gll)
748  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
749  call col2(du, coef_gll%Binv, n_gll)
750 
751  end subroutine sx_convect_scalar_lx11
752 
753  subroutine sx_convect_scalar_lx10(du, u, c, dx, dy, dz, Xh_GLL, &
754  coef_GLL, GLL_to_GL, nelv)
755  integer, parameter :: lx = 10
756  integer, intent(in) :: nelv
757  type(space_t), intent(in) :: Xh_GLL
758  type(coef_t), intent(in) :: coef_GLL
759  type(interpolator_t), intent(inout) :: GLL_to_GL
760  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
761  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
762  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
763  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
764  real(kind=rp) :: ur(lx, lx, lx, nelv)
765  real(kind=rp) :: us(lx, lx, lx, nelv)
766  real(kind=rp) :: ut(lx, lx, lx, nelv)
767  real(kind=rp) :: ud(lx, lx, lx, nelv)
768  real(kind=rp) :: wr, ws, wt
769  integer :: e, i, j, k, jj, kk, n_GLL
770  n_gll = nelv * xh_gll%lxyz
771  do i = 1, lx
772  do jj = 1, lx * lx * nelv
773  wr = 0d0
774  do kk = 1, lx
775  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
776  end do
777  ur(i, jj, 1, 1) = wr
778  end do
779  end do
780 
781  do k = 1, lx
782  do i = 1, lx
783  do j = 1, lx
784  do e = 1, nelv
785  ws = 0d0
786  !NEC$ unroll_completely
787  do kk = 1, lx
788  ws = ws + dy(j, kk) * u(i, kk, k, e)
789  end do
790  us(i,j,k,e) = ws
791  end do
792  end do
793  end do
794  end do
795 
796  do j = 1, lx
797  do i = 1, lx
798  do k = 1, lx
799  do e = 1, nelv
800  wt = 0d0
801  !NEC$ unroll_completely
802  do kk = 1, lx
803  wt = wt + dz(k, kk) * u(i, j, kk, e)
804  end do
805  ut(i,j,k,e) = wt
806  end do
807  end do
808  end do
809  end do
810 
811  do i = 1, lx * lx * lx
812  do e = 1, nelv
813  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
814  + c(i,e,2) * us(i,1,1,e) &
815  + c(i,e,3) * ut(i,1,1,e) )
816  end do
817  end do
818  call gll_to_gl%map(du, ud, nelv, xh_gll)
819  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
820  call col2(du, coef_gll%Binv, n_gll)
821 
822  end subroutine sx_convect_scalar_lx10
823 
824  subroutine sx_convect_scalar_lx9(du, u, c, dx, dy, dz, Xh_GLL, &
825  coef_GLL, GLL_to_GL, nelv)
826  integer, parameter :: lx = 9
827  integer, intent(in) :: nelv
828  type(space_t), intent(in) :: Xh_GLL
829  type(coef_t), intent(in) :: coef_GLL
830  type(interpolator_t), intent(inout) :: GLL_to_GL
831  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
832  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
833  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
834  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
835  real(kind=rp) :: ur(lx, lx, lx, nelv)
836  real(kind=rp) :: us(lx, lx, lx, nelv)
837  real(kind=rp) :: ut(lx, lx, lx, nelv)
838  real(kind=rp) :: ud(lx, lx, lx, nelv)
839  real(kind=rp) :: wr, ws, wt
840  integer :: e, i, j, k, jj, kk, n_GLL
841  n_gll = nelv * xh_gll%lxyz
842  do i = 1, lx
843  do jj = 1, lx * lx * nelv
844  wr = 0d0
845  do kk = 1, lx
846  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
847  end do
848  ur(i, jj, 1, 1) = wr
849  end do
850  end do
851 
852  do k = 1, lx
853  do i = 1, lx
854  do j = 1, lx
855  do e = 1, nelv
856  ws = 0d0
857  !NEC$ unroll_completely
858  do kk = 1, lx
859  ws = ws + dy(j, kk) * u(i, kk, k, e)
860  end do
861  us(i,j,k,e) = ws
862  end do
863  end do
864  end do
865  end do
866 
867  do j = 1, lx
868  do i = 1, lx
869  do k = 1, lx
870  do e = 1, nelv
871  wt = 0d0
872  !NEC$ unroll_completely
873  do kk = 1, lx
874  wt = wt + dz(k, kk) * u(i, j, kk, e)
875  end do
876  ut(i,j,k,e) = wt
877  end do
878  end do
879  end do
880  end do
881 
882  do i = 1, lx * lx * lx
883  do e = 1, nelv
884  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
885  + c(i,e,2) * us(i,1,1,e) &
886  + c(i,e,3) * ut(i,1,1,e) )
887  end do
888  end do
889  call gll_to_gl%map(du, ud, nelv, xh_gll)
890  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
891  call col2(du, coef_gll%Binv, n_gll)
892 
893  end subroutine sx_convect_scalar_lx9
894 
895  subroutine sx_convect_scalar_lx8(du, u, c, dx, dy, dz, Xh_GLL, &
896  coef_GLL, GLL_to_GL, nelv)
897  integer, parameter :: lx = 8
898  integer, intent(in) :: nelv
899  type(space_t), intent(in) :: Xh_GLL
900  type(coef_t), intent(in) :: coef_GLL
901  type(interpolator_t), intent(inout) :: GLL_to_GL
902  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
903  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
904  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
905  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
906  real(kind=rp) :: ur(lx, lx, lx, nelv)
907  real(kind=rp) :: us(lx, lx, lx, nelv)
908  real(kind=rp) :: ut(lx, lx, lx, nelv)
909  real(kind=rp) :: ud(lx, lx, lx, nelv)
910  real(kind=rp) :: wr, ws, wt
911  integer :: e, i, j, k, jj, kk, n_GLL
912  n_gll = nelv * xh_gll%lxyz
913  do i = 1, lx
914  do jj = 1, lx * lx * nelv
915  wr = 0d0
916  do kk = 1, lx
917  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
918  end do
919  ur(i, jj, 1, 1) = wr
920  end do
921  end do
922 
923  do k = 1, lx
924  do i = 1, lx
925  do j = 1, lx
926  do e = 1, nelv
927  ws = 0d0
928  !NEC$ unroll_completely
929  do kk = 1, lx
930  ws = ws + dy(j, kk) * u(i, kk, k, e)
931  end do
932  us(i,j,k,e) = ws
933  end do
934  end do
935  end do
936  end do
937 
938  do j = 1, lx
939  do i = 1, lx
940  do k = 1, lx
941  do e = 1, nelv
942  wt = 0d0
943  !NEC$ unroll_completely
944  do kk = 1, lx
945  wt = wt + dz(k, kk) * u(i, j, kk, e)
946  end do
947  ut(i,j,k,e) = wt
948  end do
949  end do
950  end do
951  end do
952 
953  do i = 1, lx * lx * lx
954  do e = 1, nelv
955  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
956  + c(i,e,2) * us(i,1,1,e) &
957  + c(i,e,3) * ut(i,1,1,e) )
958  end do
959  end do
960  call gll_to_gl%map(du, ud, nelv, xh_gll)
961  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
962  call col2(du, coef_gll%Binv, n_gll)
963 
964  end subroutine sx_convect_scalar_lx8
965 
966  subroutine sx_convect_scalar_lx7(du, u, c, dx, dy, dz, Xh_GLL, &
967  coef_GLL, GLL_to_GL, nelv)
968  integer, parameter :: lx = 7
969  integer, intent(in) :: nelv
970  type(space_t), intent(in) :: Xh_GLL
971  type(coef_t), intent(in) :: coef_GLL
972  type(interpolator_t), intent(inout) :: GLL_to_GL
973  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
974  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
975  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
976  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
977  real(kind=rp) :: ur(lx, lx, lx, nelv)
978  real(kind=rp) :: us(lx, lx, lx, nelv)
979  real(kind=rp) :: ut(lx, lx, lx, nelv)
980  real(kind=rp) :: ud(lx, lx, lx, nelv)
981  real(kind=rp) :: wr, ws, wt
982  integer :: e, i, j, k, jj, kk, n_GLL
983  n_gll = nelv * xh_gll%lxyz
984  do i = 1, lx
985  do jj = 1, lx * lx * nelv
986  wr = 0d0
987  do kk = 1, lx
988  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
989  end do
990  ur(i, jj, 1, 1) = wr
991  end do
992  end do
993 
994  do k = 1, lx
995  do i = 1, lx
996  do j = 1, lx
997  do e = 1, nelv
998  ws = 0d0
999  !NEC$ unroll_completely
1000  do kk = 1, lx
1001  ws = ws + dy(j, kk) * u(i, kk, k, e)
1002  end do
1003  us(i,j,k,e) = ws
1004  end do
1005  end do
1006  end do
1007  end do
1008 
1009  do j = 1, lx
1010  do i = 1, lx
1011  do k = 1, lx
1012  do e = 1, nelv
1013  wt = 0d0
1014  !NEC$ unroll_completely
1015  do kk = 1, lx
1016  wt = wt + dz(k, kk) * u(i, j, kk, e)
1017  end do
1018  ut(i,j,k,e) = wt
1019  end do
1020  end do
1021  end do
1022  end do
1023 
1024  do i = 1, lx * lx * lx
1025  do e = 1, nelv
1026  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1027  + c(i,e,2) * us(i,1,1,e) &
1028  + c(i,e,3) * ut(i,1,1,e) )
1029  end do
1030  end do
1031  call gll_to_gl%map(du, ud, nelv, xh_gll)
1032  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1033  call col2(du, coef_gll%Binv, n_gll)
1034 
1035  end subroutine sx_convect_scalar_lx7
1036 
1037  subroutine sx_convect_scalar_lx6(du, u, c, dx, dy, dz, Xh_GLL, &
1038  coef_GLL, GLL_to_GL, nelv)
1039  integer, parameter :: lx = 6
1040  integer, intent(in) :: nelv
1041  type(space_t), intent(in) :: Xh_GLL
1042  type(coef_t), intent(in) :: coef_GLL
1043  type(interpolator_t), intent(inout) :: GLL_to_GL
1044  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
1045  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
1046  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
1047  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1048  real(kind=rp) :: ur(lx, lx, lx, nelv)
1049  real(kind=rp) :: us(lx, lx, lx, nelv)
1050  real(kind=rp) :: ut(lx, lx, lx, nelv)
1051  real(kind=rp) :: ud(lx, lx, lx, nelv)
1052  real(kind=rp) :: wr, ws, wt
1053  integer :: e, i, j, k, jj, kk, n_GLL
1054  n_gll = nelv * xh_gll%lxyz
1055  do i = 1, lx
1056  do jj = 1, lx * lx * nelv
1057  wr = 0d0
1058  do kk = 1, lx
1059  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
1060  end do
1061  ur(i, jj, 1, 1) = wr
1062  end do
1063  end do
1064 
1065  do k = 1, lx
1066  do i = 1, lx
1067  do j = 1, lx
1068  do e = 1, nelv
1069  ws = 0d0
1070  !NEC$ unroll_completely
1071  do kk = 1, lx
1072  ws = ws + dy(j, kk) * u(i, kk, k, e)
1073  end do
1074  us(i,j,k,e) = ws
1075  end do
1076  end do
1077  end do
1078  end do
1079 
1080  do j = 1, lx
1081  do i = 1, lx
1082  do k = 1, lx
1083  do e = 1, nelv
1084  wt = 0d0
1085  !NEC$ unroll_completely
1086  do kk = 1, lx
1087  wt = wt + dz(k, kk) * u(i, j, kk, e)
1088  end do
1089  ut(i,j,k,e) = wt
1090  end do
1091  end do
1092  end do
1093  end do
1094 
1095  do i = 1, lx * lx * lx
1096  do e = 1, nelv
1097  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1098  + c(i,e,2) * us(i,1,1,e) &
1099  + c(i,e,3) * ut(i,1,1,e) )
1100  end do
1101  end do
1102  call gll_to_gl%map(du, ud, nelv, xh_gll)
1103  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1104  call col2(du, coef_gll%Binv, n_gll)
1105 
1106  end subroutine sx_convect_scalar_lx6
1107 
1108  subroutine sx_convect_scalar_lx5(du, u, c, dx, dy, dz, Xh_GLL, &
1109  coef_GLL, GLL_to_GL, nelv)
1110  integer, parameter :: lx = 5
1111  integer, intent(in) :: nelv
1112  type(space_t), intent(in) :: Xh_GLL
1113  type(coef_t), intent(in) :: coef_GLL
1114  type(interpolator_t), intent(inout) :: GLL_to_GL
1115  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
1116  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
1117  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
1118  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1119  real(kind=rp) :: ur(lx, lx, lx, nelv)
1120  real(kind=rp) :: us(lx, lx, lx, nelv)
1121  real(kind=rp) :: ut(lx, lx, lx, nelv)
1122  real(kind=rp) :: ud(lx, lx, lx, nelv)
1123  real(kind=rp) :: wr, ws, wt
1124  integer :: e, i, j, k, jj, kk, n_GLL
1125  n_gll = nelv * xh_gll%lxyz
1126  do i = 1, lx
1127  do jj = 1, lx * lx * nelv
1128  wr = 0d0
1129  do kk = 1, lx
1130  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
1131  end do
1132  ur(i, jj, 1, 1) = wr
1133  end do
1134  end do
1135 
1136  do k = 1, lx
1137  do i = 1, lx
1138  do j = 1, lx
1139  do e = 1, nelv
1140  ws = 0d0
1141  !NEC$ unroll_completely
1142  do kk = 1, lx
1143  ws = ws + dy(j, kk) * u(i, kk, k, e)
1144  end do
1145  us(i,j,k,e) = ws
1146  end do
1147  end do
1148  end do
1149  end do
1150 
1151  do j = 1, lx
1152  do i = 1, lx
1153  do k = 1, lx
1154  do e = 1, nelv
1155  wt = 0d0
1156  !NEC$ unroll_completely
1157  do kk = 1, lx
1158  wt = wt + dz(k, kk) * u(i, j, kk, e)
1159  end do
1160  ut(i,j,k,e) = wt
1161  end do
1162  end do
1163  end do
1164  end do
1165 
1166  do i = 1, lx * lx * lx
1167  do e = 1, nelv
1168  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1169  + c(i,e,2) * us(i,1,1,e) &
1170  + c(i,e,3) * ut(i,1,1,e) )
1171  end do
1172  end do
1173  call gll_to_gl%map(du, ud, nelv, xh_gll)
1174  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1175  call col2(du, coef_gll%Binv, n_gll)
1176 
1177  end subroutine sx_convect_scalar_lx5
1178 
1179  subroutine sx_convect_scalar_lx4(du, u, c, dx, dy, dz, Xh_GLL, &
1180  coef_GLL, GLL_to_GL, nelv)
1181  integer, parameter :: lx = 4
1182  integer, intent(in) :: nelv
1183  type(space_t), intent(in) :: Xh_GLL
1184  type(coef_t), intent(in) :: coef_GLL
1185  type(interpolator_t), intent(inout) :: GLL_to_GL
1186  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
1187  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
1188  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
1189  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1190  real(kind=rp) :: ur(lx, lx, lx, nelv)
1191  real(kind=rp) :: us(lx, lx, lx, nelv)
1192  real(kind=rp) :: ut(lx, lx, lx, nelv)
1193  real(kind=rp) :: ud(lx, lx, lx, nelv)
1194  real(kind=rp) :: wr, ws, wt
1195  integer :: e, i, j, k, jj, kk, n_GLL
1196  n_gll = nelv * xh_gll%lxyz
1197  do i = 1, lx
1198  do jj = 1, lx * lx * nelv
1199  wr = 0d0
1200  do kk = 1, lx
1201  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
1202  end do
1203  ur(i, jj, 1, 1) = wr
1204  end do
1205  end do
1206 
1207  do k = 1, lx
1208  do i = 1, lx
1209  do j = 1, lx
1210  do e = 1, nelv
1211  ws = 0d0
1212  !NEC$ unroll_completely
1213  do kk = 1, lx
1214  ws = ws + dy(j, kk) * u(i, kk, k, e)
1215  end do
1216  us(i,j,k,e) = ws
1217  end do
1218  end do
1219  end do
1220  end do
1221 
1222  do j = 1, lx
1223  do i = 1, lx
1224  do k = 1, lx
1225  do e = 1, nelv
1226  wt = 0d0
1227  !NEC$ unroll_completely
1228  do kk = 1, lx
1229  wt = wt + dz(k, kk) * u(i, j, kk, e)
1230  end do
1231  ut(i,j,k,e) = wt
1232  end do
1233  end do
1234  end do
1235  end do
1236 
1237  do i = 1, lx * lx * lx
1238  do e = 1, nelv
1239  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1240  + c(i,e,2) * us(i,1,1,e) &
1241  + c(i,e,3) * ut(i,1,1,e) )
1242  end do
1243  end do
1244  call gll_to_gl%map(du, ud, nelv, xh_gll)
1245  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1246  call col2(du, coef_gll%Binv, n_gll)
1247 
1248  end subroutine sx_convect_scalar_lx4
1249 
1250  subroutine sx_convect_scalar_lx3(du, u, c, dx, dy, dz, Xh_GLL, &
1251  coef_GLL, GLL_to_GL, nelv)
1252  integer, parameter :: lx = 3
1253  integer, intent(in) :: nelv
1254  type(space_t), intent(in) :: Xh_GLL
1255  type(coef_t), intent(in) :: coef_GLL
1256  type(interpolator_t), intent(inout) :: GLL_to_GL
1257  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
1258  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
1259  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
1260  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1261  real(kind=rp) :: ur(lx, lx, lx, nelv)
1262  real(kind=rp) :: us(lx, lx, lx, nelv)
1263  real(kind=rp) :: ut(lx, lx, lx, nelv)
1264  real(kind=rp) :: ud(lx, lx, lx, nelv)
1265  real(kind=rp) :: wr, ws, wt
1266  integer :: e, i, j, k, jj, kk, n_GLL
1267  n_gll = nelv * xh_gll%lxyz
1268  do i = 1, lx
1269  do jj = 1, lx * lx * nelv
1270  wr = 0d0
1271  do kk = 1, lx
1272  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
1273  end do
1274  ur(i, jj, 1, 1) = wr
1275  end do
1276  end do
1277 
1278  do k = 1, lx
1279  do i = 1, lx
1280  do j = 1, lx
1281  do e = 1, nelv
1282  ws = 0d0
1283  !NEC$ unroll_completely
1284  do kk = 1, lx
1285  ws = ws + dy(j, kk) * u(i, kk, k, e)
1286  end do
1287  us(i,j,k,e) = ws
1288  end do
1289  end do
1290  end do
1291  end do
1292 
1293  do j = 1, lx
1294  do i = 1, lx
1295  do k = 1, lx
1296  do e = 1, nelv
1297  wt = 0d0
1298  !NEC$ unroll_completely
1299  do kk = 1, lx
1300  wt = wt + dz(k, kk) * u(i, j, kk, e)
1301  end do
1302  ut(i,j,k,e) = wt
1303  end do
1304  end do
1305  end do
1306  end do
1307 
1308  do i = 1, lx * lx * lx
1309  do e = 1, nelv
1310  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1311  + c(i,e,2) * us(i,1,1,e) &
1312  + c(i,e,3) * ut(i,1,1,e) )
1313  end do
1314  end do
1315  call gll_to_gl%map(du, ud, nelv, xh_gll)
1316  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1317  call col2(du, coef_gll%Binv, n_gll)
1318 
1319  end subroutine sx_convect_scalar_lx3
1320 
1321  subroutine sx_convect_scalar_lx2(du, u, c, dx, dy, dz, Xh_GLL, &
1322  coef_GLL, GLL_to_GL, nelv)
1323  integer, parameter :: lx = 2
1324  integer, intent(in) :: nelv
1325  type(space_t), intent(in) :: Xh_GLL
1326  type(coef_t), intent(in) :: coef_GLL
1327  type(interpolator_t), intent(inout) :: GLL_to_GL
1328  real(kind=rp), intent(inout) :: du(xh_gll%lx, xh_gll%lx, xh_gll%lx, nelv)
1329  real(kind=rp), intent(in) :: u(lx, lx, lx, nelv)
1330  real(kind=rp), intent(in) :: c(lx*lx*lx, nelv, 3)
1331  real(kind=rp), dimension(lx, lx), intent(in) :: dx, dy, dz
1332  real(kind=rp) :: ur(lx, lx, lx, nelv)
1333  real(kind=rp) :: us(lx, lx, lx, nelv)
1334  real(kind=rp) :: ut(lx, lx, lx, nelv)
1335  real(kind=rp) :: ud(lx, lx, lx, nelv)
1336  real(kind=rp) :: wr, ws, wt
1337  integer :: e, i, j, k, jj, kk, n_GLL
1338  n_gll = nelv * xh_gll%lxyz
1339  do i = 1, lx
1340  do jj = 1, lx * lx * nelv
1341  wr = 0d0
1342  do kk = 1, lx
1343  wr = wr + dx(i, kk) * u(kk, jj, 1, 1)
1344  end do
1345  ur(i, jj, 1, 1) = wr
1346  end do
1347  end do
1348 
1349  do k = 1, lx
1350  do i = 1, lx
1351  do j = 1, lx
1352  do e = 1, nelv
1353  ws = 0d0
1354  !NEC$ unroll_completely
1355  do kk = 1, lx
1356  ws = ws + dy(j, kk) * u(i, kk, k, e)
1357  end do
1358  us(i,j,k,e) = ws
1359  end do
1360  end do
1361  end do
1362  end do
1363 
1364  do j = 1, lx
1365  do i = 1, lx
1366  do k = 1, lx
1367  do e = 1, nelv
1368  wt = 0d0
1369  !NEC$ unroll_completely
1370  do kk = 1, lx
1371  wt = wt + dz(k, kk) * u(i, j, kk, e)
1372  end do
1373  ut(i,j,k,e) = wt
1374  end do
1375  end do
1376  end do
1377  end do
1378 
1379  do i = 1, lx * lx * lx
1380  do e = 1, nelv
1381  ud(i,1,1,e) = ( c(i,e,1) * ur(i,1,1,e) &
1382  + c(i,e,2) * us(i,1,1,e) &
1383  + c(i,e,3) * ut(i,1,1,e) )
1384  end do
1385  end do
1386  call gll_to_gl%map(du, ud, nelv, xh_gll)
1387  call coef_gll%gs_h%op(du, n_gll, gs_op_add)
1388  call col2(du, coef_gll%Binv, n_gll)
1389 
1390  end subroutine sx_convect_scalar_lx2
1391 
1392 end submodule sx_convect_scalar
math
Definition: math.f90:60
math::col2
subroutine, public col2(a, b, n)
Vector multiplication .
Definition: math.f90:729
num_types
Definition: num_types.f90:1
num_types::rp
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
opr_sx
Operators SX-Aurora backend.
Definition: opr_sx.f90:2