70 use,
intrinsic :: iso_c_binding, only : c_ptr
90 subroutine tensr3(v, nv, u, nu, A, Bt, Ct, w)
93 real(kind=
rp),
intent(inout) :: v(nv, nv, nv)
94 real(kind=
rp),
intent(inout) :: u(nu, nu, nu)
95 real(kind=
rp),
intent(inout) :: w(nu*nu*nv)
96 real(kind=
rp),
intent(inout) :: a(nv, nu)
97 real(kind=
rp),
intent(inout) :: bt(nu, nv)
98 real(kind=
rp),
intent(inout) :: ct(nu, nv)
99 integer :: j, k, l, nunu, nvnv, nunv
107 call mxm(a, nv, u, nu, v, nunu)
111 call mxm(v(k, 1, 1), nv, bt, nu, w(l), nv)
115 call mxm(w, nvnv, ct, nu, v, nv)
120 subroutine trsp(a, lda, b, ldb)
121 integer,
intent(in) :: lda
122 integer,
intent(in) :: ldb
123 real(kind=
rp),
intent(inout) :: a(lda, ldb)
124 real(kind=
rp),
intent(in) :: b(ldb, lda)
137 integer,
intent(in) :: n
138 real(kind=
rp),
intent(inout) :: a(n, n)
154 integer,
intent(in) :: nv, nu
155 real(kind=
rp),
intent(inout) :: v(nv*nv), u(nu*nu)
156 real(kind=
rp),
intent(inout) :: a(nv,nu), bt(nu,nv)
171 integer,
intent(in) :: nv, nu
172 real(kind=
rp),
intent(inout) :: v(nv*nv*nv), u(nu*nu*nu)
173 real(kind=
rp),
intent(inout) :: a(nv,nu),bt(nu, nv),ct(nu,nv)
188 integer,
intent(in) :: nv, nu, n_pt, el_list(n_pt)
189 real(kind=
rp),
intent(inout) :: v(nv*nv*nv, n_pt), u(nu*nu*nu,1)
190 real(kind=
rp),
intent(inout) :: a(nv,nu,n_pt),bt(nu, nv,n_pt),ct(nu,nv,n_pt)
191 type(c_ptr) :: v_d, u_d, a_d, bt_d, ct_d, el_list_d
196 call tnsr3d_el_sx(v(1,i), nv, u(1,el_list(i)), nu, a(1,1,i), bt(1,1,i), ct(1,1,i))
200 call tnsr3d_el_xsmm(v(1,i), nv, u(1,el_list(i)), nu, a(1,1,i), bt(1,1,i), ct(1,1,i))
213 call tnsr3d_el_cpu(v(1,i), nv, u(1,el_list(i)+1), nu, a(1,1,i), bt(1,1,i), ct(1,1,i))
222 subroutine tnsr3d(v, nv, u, nu, A, Bt, Ct, nelv)
223 integer,
intent(inout) :: nv, nu, nelv
224 real(kind=
rp),
intent(inout) :: v(nv*nv*nv,nelv), u(nu*nu*nu,nelv)
225 real(kind=
rp),
intent(inout) :: a(nv,nu), bt(nu, nv), ct(nu,nv)
226 type(c_ptr) :: v_d, u_d, a_d, bt_d, ct_d
230 call tnsr3d_sx(v, nv, u, nu, a, bt, ct, nelv)
241 call tnsr3d_cpu(v, nv, u, nu, a, bt, ct, nelv)
248 integer,
intent(inout) :: nv, nu, nelv
249 real(kind=
rp),
intent(inout) :: v(nv*nv*nv*nelv)
250 real(kind=
rp),
intent(inout) :: a(nv,nu), bt(nu, nv), ct(nu,nv)
266 integer,
intent(in) :: nx, ny, nz
267 real(kind=
rp),
intent(in) :: h1(nx), h2(ny), h3(nz)
268 real(kind=
rp),
intent(inout) :: s(nx, ny, nz)
270 integer :: ix, iy, iz
276 s(ix,iy,iz) = s(ix,iy,iz)+hh*h1(ix)
298 real(kind=
rp),
intent(inout) :: v
299 integer,
intent(in) :: nu
300 real(kind=
rp),
intent(inout) :: u(nu,nu,nu)
301 real(kind=
rp),
intent(inout) :: hr(nu)
302 real(kind=
rp),
intent(inout) :: hs(nu)
303 real(kind=
rp),
intent(inout) :: ht(nu)
307 real(kind=
rp) :: vv(1)
334 real(kind=
rp),
intent(inout) :: v(3)
335 integer,
intent(in) :: nu
336 real(kind=
rp),
intent(inout) :: u1(nu,nu,nu)
337 real(kind=
rp),
intent(inout) :: u2(nu,nu,nu)
338 real(kind=
rp),
intent(inout) :: u3(nu,nu,nu)
339 real(kind=
rp),
intent(inout) :: hr(nu)
340 real(kind=
rp),
intent(inout) :: hs(nu)
341 real(kind=
rp),
intent(inout) :: ht(nu)
Return the device pointer for an associated Fortran array.
Device abstraction, common interface for various accelerators.
Wrapper for all matrix-matrix product implementations.
subroutine, public mxm(a, n1, b, n2, c, n3)
Compute matrix-matrix product for contiguously packed matrices A,B, and C.
integer, parameter neko_bcknd_sx
integer, parameter neko_bcknd_device
integer, parameter neko_bcknd_xsmm
integer, parameter, public rp
Global precision used in computations.
subroutine, public tnsr1_3d_cpu(v, nv, nu, A, Bt, Ct, nelv)
subroutine, public tnsr2d_el_cpu(v, nv, u, nu, A, Bt)
subroutine, public tnsr3d_el_cpu(v, nv, u, nu, A, Bt, Ct)
subroutine, public tnsr3d_cpu(v, nv, u, nu, A, Bt, Ct, nelv)
subroutine, public tnsr3d_el_list_device(v_d, nv, u_d, nu, A_d, Bt_d, Ct_d, elements, n_points)
subroutine, public tnsr3d_device(v_d, nv, u_d, nu, A_d, Bt_d, Ct_d, nelv)
Tensor operations SX-Aurora backend.
subroutine, public tnsr2d_el_sx(v, nv, u, nu, A, Bt)
subroutine, public tnsr1_3d_sx(v, nv, nu, A, Bt, Ct, nelv)
subroutine, public tnsr3d_el_sx(v, nv, u, nu, A, Bt, Ct)
subroutine, public tnsr3d_sx(v, nv, u, nu, A, Bt, Ct, nelv)
Tensor operations libxsmm backend.
subroutine, public tnsr3d_el_xsmm(v, nv, u, nu, A, Bt, Ct)
subroutine, public tnsr2d_el_xsmm(v, nv, u, nu, A, Bt)
subroutine, public tnsr1_3d_xsmm(v, nv, nu, A, Bt, Ct, nelv)
subroutine, public tnsr3d_xsmm(v, nv, u, nu, A, Bt, Ct, nelv)
subroutine, public tensr3(v, nv, u, nu, A, Bt, Ct, w)
Tensor product .
subroutine, public tnsr3d(v, nv, u, nu, A, Bt, Ct, nelv)
Tensor product performed on nelv elements.
subroutine, public tnsr1_3d(v, nv, nu, A, Bt, Ct, nelv)
In place tensor product .
subroutine, public trsp1(a, n)
In-place transpose of a square tensor.
subroutine, public addtnsr(s, h1, h2, h3, nx, ny, nz)
Maps and adds to S a tensor product form of the three functions H1,H2,H3. This is a single element ro...
subroutine, public trsp(a, lda, b, ldb)
Transpose of a rectangular tensor .
subroutine triple_tensor_product_vector(v, u1, u2, u3, nu, Hr, Hs, Ht)
Computes the tensor product on a vector field . This operation is usually performed for spectral inte...
subroutine, public tnsr3d_el_list(v, nv, u, nu, A, Bt, Ct, el_list, n_pt)
Tensor product performed on a subset of the elements.
subroutine triple_tensor_product_scalar(v, u, nu, Hr, Hs, Ht)
Computes the tensor product . This operation is usually performed for spectral interpolation of a sca...
subroutine, public tnsr2d_el(v, nv, u, nu, A, Bt)
Computes .
subroutine, public tnsr3d_el(v, nv, u, nu, A, Bt, Ct)
Tensor product performed on a single element.