Line data Source code
1 : !{\src2tex{textfont=tt}}
2 : !!****m* ABINIT/m_abi_gpu_linalg
3 : !! NAME
4 : !! m_abi_gpu_linalg
5 : !!
6 : !! FUNCTION
7 : !! Interfaces of GPU subroutines wrapper
8 : !!
9 : !! COPYRIGHT
10 : !! Copyright (C) 2011-2026 ABINIT group (FDahm)
11 : !! This file is distributed under the terms of the
12 : !! GNU General Public License, see ~abinit/COPYING
13 : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : !!
15 : !! NOTES
16 : !! This is an include file.
17 : !!
18 : !! SOURCE
19 :
20 : !!***
21 :
22 : #ifndef HAVE_GPU
23 :
24 : !!****f* m_abi_gpu_linalg/gpu_linalg_init
25 : !! NAME
26 : !! gpu_linalg_init
27 : !!
28 : !! FUNCTION
29 : !! initialisation of linalg environnement on GPU
30 : !!
31 : !! INPUTS
32 : !!
33 : !! OUTPUT
34 : !!
35 : !! SIDE EFFECTS
36 : !! WARNING! : this routine is a dummy one when HAVE_GPU is not enabled
37 : !! the correct one is in 17_gpu_toolbox/gpu_linalg.cu
38 : !!
39 : !! SOURCE
40 :
41 0 : subroutine gpu_linalg_init()
42 :
43 :
44 0 : end subroutine gpu_linalg_init
45 : !!***
46 :
47 : !!****f* m_abi_gpu_linalg/gpu_linalg_shutdown
48 : !! NAME
49 : !! gpu_linalg_shutdown
50 : !!
51 : !! FUNCTION
52 : !! close linalg environnement on GPU
53 : !!
54 : !! INPUTS
55 : !!
56 : !! OUTPUT
57 : !!
58 : !! SIDE EFFECTS
59 : !! WARNING! : this routine is a dummy one when HAVE_GPU is not enabled
60 : !! the correct one is in 17_gpu_toolbox/gpu_linalg.cu
61 : !!
62 : !! SOURCE
63 0 : subroutine gpu_linalg_shutdown()
64 :
65 0 : end subroutine gpu_linalg_shutdown
66 : !!***
67 :
68 : !!****f* m_abi_gpu_linalg/gpu_xgemm
69 : !! NAME
70 : !! gpu_xgemm
71 : !!
72 : !! FUNCTION
73 : !! Compute a scalar-matrix-matrix product and return a scalar-matrix product on GPU
74 : !! c = alpha * op(a) * op(b) + beta * c
75 : !!
76 : !! INPUTS
77 : !! cplx = 1 if real 2 if complex
78 : !! transa= from of op(a) to be used in the matrix multiplication
79 : !! transb= from of op(b) to be used in the matrix multiplication
80 : !! m = number of rows of the matrix op(a) and of the matrix c
81 : !! n = number of rows of the matrix op(b) and the number of columns of the matrix c
82 : !! k = number of columns of the matrix op(a) and the number of rows of the matrix op(b)
83 : !! alpha = alpha scalar coefficient for matrix op(a)
84 : !! a_gpu = pointer to gpu memory location of matrix a
85 : !! lda = first dimension of a
86 : !! b_gpu = pointer to gpu memory location of matrix b
87 : !! ldb = first dimension of b
88 : !! beta = beta scalar coefficient for matrix c
89 : !! c_gpu = pointer to gpu memory location of matrix c
90 : !! ldc = first dimension of c
91 : !!
92 : !! OUTPUT
93 : !! c = c matrix
94 : !!
95 : !! SIDE EFFECTS
96 : !! WARNING! : this routine is a dummy one when HAVE_GPU is not enabled
97 : !! the correct one is in 17_gpu_toolbox/gpu_linalg.cu
98 : !!
99 : !! SOURCE
100 :
101 0 : subroutine gpu_xgemm(cplx,transa,transb,m,n,k,alpha,a_gpu,lda,b_gpu,ldb,beta,c_gpu,ldc)
102 :
103 : !Arguments ------------------------------------
104 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
105 : complex(dp),intent(in) :: alpha,beta
106 : character(len=1),intent(in) :: transa,transb
107 : type(c_ptr),intent(in) :: a_gpu,b_gpu
108 : type(c_ptr),intent(inout) :: c_gpu
109 : !Local variables ------------------------------
110 : type(c_ptr) :: cptr
111 : ! *********************************************************************
112 :
113 : if (.false.) then
114 : cptr=a_gpu;cptr=b_gpu;cptr=c_gpu
115 : write(std_out,*) transa,transb,cplx,lda,ldb,ldc,m,n,k,alpha,beta
116 : end if
117 :
118 0 : end subroutine gpu_xgemm
119 : !!***
120 :
121 : !!****f* m_abi_gpu_linalg/gpu_xtrsm
122 : !! NAME
123 : !! gpu_xtrsm
124 : !!
125 : !! FUNCTION
126 : !! Solves a matrix equation (one matrix operand is triangular) on GPU.
127 : !! The xtrsm routines solve one of the following matrix equations
128 : !! op(a)*x = alpha*b
129 : !! or
130 : !! x*op(a) = alpha*b,
131 : !!
132 : !! INPUTS
133 : !! cplx= 1 if real 2 if complex
134 : !! side= Specifies whether op(a) appears on the left or right of x for
135 : !! the operation to be performed as follows:
136 : !! L or l op(a)*x = alpha*b
137 : !! R or r x*op(a) = alpha*b
138 : !! uplo= Specifies whether the matrix a is an upper or lower triangular
139 : !! matrix as follows:
140 : !! U or u Matrix a is an upper triangular matrix.
141 : !! L or l Matrix a is a lower triangular matrix
142 : !! transa= Specifies the form of op(a) to be used in the matrix
143 : !! multiplication as follows:
144 : !! N or n op(a) = a
145 : !! T or t op(a) = a'
146 : !! C or c op(a) = conjg(a')
147 : !! diag= Specifies whether or not a is unit triangular as follows:
148 : !! U or u Matrix a is assumed to be unit triangular.
149 : !! N or n Matrix a is not assumed to be unit triangular.
150 : !! m= Specifies the number of rows of b. The value of m must be at least zero
151 : !! n= Specifies the number of columns of b. The value of n must be at least zero
152 : !! alpha= Specifies the scalar alpha. When alpha is zero, then a is not referenced and b
153 : !! need not be set before entry.
154 : !! a_gpu = pointer to gpu memory location of array a, DIMENSION (lda, k), where k is m when side = 'L' or 'l' and is n
155 : !! when side = 'R' or 'r'.
156 : !! lda= Specifies the first dimension of a as declared in the calling
157 : !! (sub)program. When side = 'L' or 'l', then lda must be at least max(1,
158 : !! m), when side = 'R' or 'r', then lda must be at least max(1, n).
159 : !! b_gpu = pointer to gpu memory location of b Array, DIMENSION (ldb,n). Before entry, the leading m-by-n part of the array
160 : !! b must contain the right-hand side matrix b.
161 : !! ldb= Specifies the first dimension of b as declared in the calling
162 : !! (sub)program. The value of ldb must be at least max(1, m).
163 : !!
164 : !! OUTPUT
165 : !! b_gpu
166 : !!
167 : !! SIDE EFFECTS
168 : !! WARNING! : this routine is a dummy one when HAVE_GPU is not enabled
169 : !! the correct one is in 17_gpu_toolbox/gpu_linalg.cu
170 : !!
171 : !! SOURCE
172 :
173 0 : subroutine gpu_xtrsm(cplx,side,uplo,transa,diag,m,n,alpha,a_gpu,lda,b_gpu,ldb)
174 :
175 : ! !Arguments ------------------------------------
176 : integer, intent(in) :: cplx,lda,ldb,m,n
177 : complex(dp), intent(in) :: alpha
178 : character(len=1), intent(in) :: side,uplo,transa,diag
179 : type(c_ptr),intent(in) :: a_gpu
180 : type(c_ptr),intent(inout) :: b_gpu
181 : !Local variables ------------------------------
182 : type(c_ptr) :: cptr
183 : ! *********************************************************************
184 :
185 : if (.false.) then
186 : cptr=a_gpu;cptr=b_gpu
187 : write(std_out,*) side,uplo,transa,diag,cplx,lda,ldb,m,n,alpha
188 : end if
189 :
190 0 : end subroutine gpu_xtrsm
191 : !!***
192 :
193 : !!****f* m_abi_gpu_linalg/gpu_xaxpy
194 : !! NAME
195 : !! gpu_xaxpy
196 : !!
197 : !! FUNCTION
198 : !! Compute a BLAS-1 AXPY operation on GPU
199 : !! y = alpha * x + y
200 : !!
201 : !! INPUTS
202 : !! cplx = 1 if real 2 if complex
203 : !! size = vector size
204 : !! alpha = scalar complex value
205 : !! x_gpu = pointer to gpu memory location of array x
206 : !! incrx = stride between consecutive elements of x
207 : !! y_gpu = pointer to gpu memory location of array y
208 : !! incry = stride between consecutive elements of y
209 : !!
210 : !! SOURCE
211 0 : subroutine gpu_xaxpy(cplx, size, alpha, x_gpu, incrx, y_gpu, incry)
212 :
213 : ! !Arguments ------------------------------------
214 : integer, intent(in) :: cplx
215 : integer, intent(in) :: size
216 : complex(dp), intent(in) :: alpha
217 : type(c_ptr), intent(in) :: x_gpu
218 : integer, intent(in) :: incrx
219 : type(c_ptr), intent(inout) :: y_gpu
220 : integer, intent(in) :: incry
221 :
222 : ABI_UNUSED((/cplx,size,incrx,incry/))
223 : ABI_UNUSED(alpha)
224 0 : ABI_UNUSED_A(x_gpu)
225 0 : ABI_UNUSED_A(y_gpu)
226 0 : end subroutine gpu_xaxpy
227 : !!***
228 :
229 : !!****f* m_abi_gpu_linalg/gpu_xcopy
230 : !! NAME
231 : !! gpu_xcopy
232 : !!
233 : !! FUNCTION
234 : !! Compute a BLAS-1 COPY operation on GPU
235 : !! y = x (copy x into y)
236 : !!
237 : !! INPUTS
238 : !! cplx = 1 if real 2 if complex
239 : !! size = input vector size
240 : !! x_gpu = pointer to gpu memory location of array x
241 : !! incrx = stride between consecutive elements of x
242 : !! y_gpu = pointer to gpu memory location of array y
243 : !! incry = stride between consecutive elements of y
244 : !!
245 : !! SOURCE
246 0 : subroutine gpu_xcopy(cplx, size, x_gpu, incrx, y_gpu, incry)
247 :
248 : ! !Arguments ------------------------------------
249 : integer, intent(in) :: cplx
250 : integer, intent(in) :: size
251 : type(c_ptr), intent(in) :: x_gpu
252 : integer, intent(in) :: incrx
253 : type(c_ptr), intent(inout) :: y_gpu
254 : integer, intent(in) :: incry
255 :
256 : ABI_UNUSED((/cplx,size,incrx,incry/))
257 0 : ABI_UNUSED_A(x_gpu)
258 0 : ABI_UNUSED_A(y_gpu)
259 0 : end subroutine gpu_xcopy
260 : !!***
261 :
262 : !!****f* m_abi_gpu_linalg/gpu_xscal
263 : !! NAME
264 : !! gpu_xscal
265 : !!
266 : !! FUNCTION
267 : !! Compute a BLAS-1 SCAL operation on GPU
268 : !! x = alpha * x
269 : !!
270 : !! INPUTS
271 : !! cplx = 1 if real 2 if complex
272 : !! size = vector size
273 : !! alpha = scalar complex value
274 : !! x_gpu = pointer to gpu memory location of array x
275 : !! incrx = stride between consecutive elements of x
276 : !!
277 : !! SOURCE
278 0 : subroutine gpu_xscal(cplx, size, alpha, x_gpu, incrx)
279 :
280 : ! !Arguments ------------------------------------
281 : integer, intent(in) :: cplx
282 : integer, intent(in) :: size
283 : complex(dp), intent(in) :: alpha
284 : type(c_ptr), intent(in) :: x_gpu
285 : integer, intent(in) :: incrx
286 :
287 : ABI_UNUSED((/cplx,size,incrx/))
288 : ABI_UNUSED(alpha)
289 0 : ABI_UNUSED_A(x_gpu)
290 0 : end subroutine gpu_xscal
291 : !!***
292 :
293 : !!****f* m_abi_gpu_linalg/gpu_xdot
294 : !! NAME
295 : !! gpu_xdot
296 : !!
297 : !! FUNCTION
298 : !! Compute blas-3 DOT on GPU
299 : !! trace(X^H*Y)
300 : !!
301 : !! INPUTS
302 : !! cplx = 1 if real 2 if complex
303 : !! N = number of elements n*s
304 : !! X = array (n,s)
305 : !! incrx = stride for X
306 : !! Y = array (n,s)
307 : !! incry = stride for Y
308 : !! alpha = output pointer
309 : !!
310 : !! SOURCE
311 0 : subroutine gpu_xdot(cplx, size, alpha, x_gpu, incrx, y_gpu, incry)
312 :
313 : ! !Arguments ------------------------------------
314 : integer, intent(in) :: cplx
315 : integer, intent(in) :: size
316 : complex(dp), intent(in) :: alpha
317 : type(c_ptr), intent(in) :: x_gpu
318 : integer, intent(in) :: incrx
319 : type(c_ptr), intent(inout) :: y_gpu
320 : integer, intent(in) :: incry
321 :
322 : ABI_UNUSED((/cplx,size,incrx,incry/))
323 : ABI_UNUSED(alpha)
324 0 : ABI_UNUSED_A(x_gpu)
325 0 : ABI_UNUSED_A(y_gpu)
326 0 : end subroutine gpu_xdot
327 : !!***
328 :
329 : !!****f* m_abi_gpu_linalg/gpu_xsygvd
330 : !! NAME
331 : !! gpu_xsygvd
332 : !!
333 : !! FUNCTION
334 : !! Compute a LAPACK SYGVD operation on GPU
335 : !! compute eigen values/vectors of a real generalized
336 : !! symmetric-definite eigenproblem
337 : !!
338 : !! See cusolver documentation
339 : !! https://docs.nvidia.com/cuda/cusolver/index.html#cuSolverDN-lt-t-gt-sygvd
340 : !!
341 : !! See also LAPACK doc in reference implementation:
342 : !! https://github.com/Reference-LAPACK/lapack/blob/master/SRC/dsygvd.f
343 : !!
344 : !! INPUTS
345 : !! cplx = 1 if real 2 if complex
346 : !! itype = integer, type of problem
347 : !! jobz = character, 'n'(eigenvalues only) or 'v' (eigenvalues + eigenvectors)
348 : !! uplo = character, 'u' or 'l'
349 : !! A_nrows = matrix size
350 : !! A_ptr = pointer to gpu memory location of matrix A
351 : !! lda = leading dimension of matrix A
352 : !! B_ptr = pointer to gpu memory location of matrix B
353 : !! ldb = leading dimension of matrix B
354 : !! W_ptr = pointer to gpu memory location of matrix W (output eigen values)
355 : !! work_ptr =
356 : !! lwork =
357 : !! devInfo =
358 : !!
359 : !! SOURCE
360 0 : subroutine gpu_xsygvd(cplx, itype, jobz, uplo, A_nrows, &
361 : & A_ptr, lda, &
362 : & B_ptr, ldb, &
363 : & W_ptr, &
364 : & work_ptr, lwork, &
365 : & devInfo)
366 :
367 : ! Arguments ------------------------------------
368 : integer, intent(in ) :: cplx
369 : integer, intent(in ) :: itype
370 : character(len=1),intent(in ) :: jobz
371 : character(len=1),intent(in ) :: uplo
372 : integer, intent(in ) :: A_nrows
373 : type(c_ptr), intent(in ) :: A_ptr
374 : integer, intent(in ) :: lda
375 : type(c_ptr), intent(in ) :: B_ptr
376 : integer, intent(in ) :: ldb
377 : type(c_ptr), intent(inout) :: W_ptr
378 : type(c_ptr), intent(inout) :: work_ptr
379 : integer, intent(in ) :: lwork
380 : integer, intent(inout) :: devInfo
381 :
382 : ABI_UNUSED((/cplx,itype,A_nrows,lda,ldb,lwork,devInfo/))
383 : ABI_UNUSED((/jobz,uplo/))
384 : ABI_UNUSED((/A_ptr,B_ptr,W_ptr,work_ptr/))
385 0 : end subroutine gpu_xsygvd
386 : !!***
387 :
388 : !!****f* m_abi_gpu_linalg/gpu_xsygvd_bufferSize
389 : !! NAME
390 : !! gpu_xsygvd_bufferSize
391 : !!
392 : !! FUNCTION
393 : !! Compute required size for auxiliary work buffer used internally by cusolver
394 : !! in cusolverDnDsygvd / cusolverDnZhegvd
395 : !!
396 : !! See cusolver documentation
397 : !! https://docs.nvidia.com/cuda/cusolver/index.html#cuSolverDN-lt-t-gt-sygvd
398 : !!
399 : !! INPUTS
400 : !! cplx = 1 if real 2 if complex
401 : !! itype = integer, type of problem
402 : !! jobz = character, 'n'(eigenvalues only) or 'v' (eigenvalues + eigenvectors)
403 : !! uplo = character, 'u' or 'l'
404 : !! A_nrows = matrix size
405 : !! A_ptr = pointer to gpu memory location of matrix A
406 : !! lda = leading dimension of matrix A
407 : !! B_ptr = pointer to gpu memory location of matrix B
408 : !! ldb = leading dimension of matrix B
409 : !! W_ptr = pointer to gpu memory location of matrix W (output eigen values)
410 : !! work_ptr =
411 : !! lwork =
412 : !! devInfo =
413 : !!
414 : !! SOURCE
415 0 : subroutine gpu_xsygvd_buffersize(cplx, itype, jobz, uplo, A_nrows, &
416 : & A_ptr, lda, &
417 : & B_ptr, ldb, &
418 : & W_ptr, &
419 : & lwork)
420 :
421 : ! Arguments ------------------------------------
422 : integer, intent(in ) :: cplx
423 : integer, intent(in ) :: itype
424 : character(len=1),intent(in ) :: jobz
425 : character(len=1),intent(in ) :: uplo
426 : integer, intent(in ) :: A_nrows
427 : type(c_ptr), intent(in ) :: A_ptr
428 : integer, intent(in ) :: lda
429 : type(c_ptr), intent(in ) :: B_ptr
430 : integer, intent(in ) :: ldb
431 : type(c_ptr), intent(inout) :: W_ptr
432 : integer, intent(in ) :: lwork
433 :
434 : ABI_UNUSED((/cplx,itype,A_nrows,lda,ldb,lwork/))
435 : ABI_UNUSED((/jobz,uplo/))
436 : ABI_UNUSED((/A_ptr,B_ptr,W_ptr/))
437 0 : end subroutine gpu_xsygvd_bufferSize
438 : !!***
439 : #endif
440 :
441 : !!***
442 :
443 : !!****f* m_abi_gpu_linalg/abi_gpu_xgemm
444 : !! NAME
445 : !! abi_gpu_xgemm
446 : !!
447 : !! FUNCTION
448 : !! Compute a scalar-matrix-matrix product and return a scalar-matrix product on GPU
449 : !! c = alpha * op(a) * op(b) + beta * c
450 : !!
451 : !! INPUTS
452 : !! cplx = 1 if real 2 if complex
453 : !! transa= from of op(a) to be used in the matrix multiplication
454 : !! transb= from of op(b) to be used in the matrix multiplication
455 : !! m = number of rows of the matrix op(a) and of the matrix c
456 : !! n = number of rows of the matrix op(b) and the number of columns of the matrix c
457 : !! k = number of columns of the matrix op(a) and the number of rows of the matrix op(b)
458 : !! alpha = alpha scalar coefficient for matrix op(a)
459 : !! a = pointer to gpu memory location of matrix a
460 : !! lda = first dimension of a
461 : !! b = pointer to gpu memory location of matrix b
462 : !! ldb = first dimension of b
463 : !! beta = beta scalar coefficient for matrix c
464 : !! c = pointer to gpu memory location of matrix c
465 : !! ldc = first dimension of c
466 : !!
467 : !! OUTPUT
468 : !! c = c matrix
469 : !!
470 : !! SOURCE
471 :
472 0 : subroutine abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,a,lda,b,ldb,beta,c,ldc)
473 :
474 : !Arguments ------------------------------------
475 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
476 : complex(dp),intent(in) :: alpha,beta
477 : character(len=1),intent(in) :: transa,transb
478 : type(c_ptr),intent(in) :: a,b
479 : type(c_ptr),intent(in) :: c
480 : ! *********************************************************************
481 :
482 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
483 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
484 : end if
485 :
486 : #ifdef HAVE_GPU
487 : call gpu_xgemm(cplx,transa,transb,m,n,k,alpha,&
488 : a,lda,b,ldb,beta,c,ldc)
489 :
490 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
491 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
492 : ! Therefore, we issue a stream sync here to avoid
493 : !potential mistakes in calling context.
494 : call gpu_linalg_stream_synchronize()
495 : end if
496 :
497 : #else
498 : ! Unused if GPU code disabled
499 : ABI_UNUSED((/cplx,lda,ldb,ldc,m,n,k/))
500 : ABI_UNUSED((/alpha,beta/))
501 : ABI_UNUSED((/transa,transb/))
502 : ABI_UNUSED((/a,b,c/))
503 : #endif
504 :
505 0 : end subroutine abi_gpu_xgemm_cptr
506 : !!***
507 :
508 0 : subroutine abi_gpu_xgemm_d(cplx,transa,transb,m,n,k,alpha,a,lda,b,ldb,beta,c,ldc)
509 :
510 : !Arguments ------------------------------------
511 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
512 : complex(dp),intent(in) :: alpha,beta
513 : character(len=1),intent(in) :: transa,transb
514 : real(dp), intent(in),target :: a(*),b(*)
515 : real(dp), intent(inout),target :: c(*)
516 : ! *********************************************************************
517 :
518 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
519 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
520 : end if
521 :
522 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
523 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
524 : c_loc(a),lda,&
525 : c_loc(b),ldb,&
526 : beta,&
527 0 : c_loc(c),ldc)
528 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
529 : #ifdef HAVE_OPENMP_OFFLOAD
530 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
531 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
532 : c_loc(a),lda,&
533 : c_loc(b),ldb,&
534 : beta,&
535 : c_loc(c),ldc)
536 : !$OMP END TARGET DATA
537 : #endif
538 : else
539 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
540 : end if
541 :
542 0 : end subroutine abi_gpu_xgemm_d
543 : !!***
544 :
545 0 : subroutine abi_gpu_xgemm_z(cplx,transa,transb,m,n,k,alpha,a,lda,b,ldb,beta,c,ldc)
546 :
547 : !Arguments ------------------------------------
548 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
549 : complex(dp),intent(in) :: alpha,beta
550 : character(len=1),intent(in) :: transa,transb
551 : complex(dp),intent(in),target :: a(*),b(*)
552 : complex(dp),intent(inout),target :: c(*)
553 : ! *********************************************************************
554 :
555 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
556 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
557 : end if
558 :
559 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
560 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
561 : c_loc(a),lda,&
562 : c_loc(b),ldb,&
563 : beta,&
564 0 : c_loc(c),ldc)
565 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
566 : #ifdef HAVE_OPENMP_OFFLOAD
567 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
568 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
569 : c_loc(a),lda,&
570 : c_loc(b),ldb,&
571 : beta,&
572 : c_loc(c),ldc)
573 : !$OMP END TARGET DATA
574 : #endif
575 : else
576 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
577 : end if
578 :
579 0 : end subroutine abi_gpu_xgemm_z
580 : !!***
581 :
582 0 : subroutine abi_gpu_xgemm_2d(cplx,transa,transb,m,n,k,alpha,a,lda,b,ldb,beta,c,ldc)
583 :
584 : !Arguments ------------------------------------
585 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
586 : complex(dp),intent(in) :: alpha,beta
587 : character(len=1),intent(in) :: transa,transb
588 : real(dp), intent(in),target :: a(lda,*),b(ldb,*)
589 : real(dp), intent(inout),target :: c(ldc,*)
590 : ! *********************************************************************
591 :
592 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
593 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
594 : end if
595 :
596 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
597 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
598 : c_loc(a),lda,&
599 : c_loc(b),ldb,&
600 : beta,&
601 0 : c_loc(c),ldc)
602 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
603 : #ifdef HAVE_OPENMP_OFFLOAD
604 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
605 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
606 : c_loc(a),lda,&
607 : c_loc(b),ldb,&
608 : beta,&
609 : c_loc(c),ldc)
610 : !$OMP END TARGET DATA
611 : #endif
612 : else
613 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
614 : end if
615 :
616 0 : end subroutine abi_gpu_xgemm_2d
617 : !!***
618 :
619 0 : subroutine abi_gpu_xgemm_2z(cplx,transa,transb,m,n,k,alpha,a,lda,b,ldb,beta,c,ldc)
620 :
621 : !Arguments ------------------------------------
622 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
623 : complex(dp),intent(in) :: alpha,beta
624 : character(len=1),intent(in) :: transa,transb
625 : complex(dp),intent(in),target :: a(lda,*),b(ldb,*)
626 : complex(dp),intent(inout),target :: c(ldc,*)
627 : ! *********************************************************************
628 :
629 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
630 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
631 : end if
632 :
633 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
634 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
635 : c_loc(a),lda,&
636 : c_loc(b),ldb,&
637 : beta,&
638 0 : c_loc(c),ldc)
639 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
640 : #ifdef HAVE_OPENMP_OFFLOAD
641 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
642 : call abi_gpu_xgemm_cptr(cplx,transa,transb,m,n,k,alpha,&
643 : c_loc(a),lda,&
644 : c_loc(b),ldb,&
645 : beta,&
646 : c_loc(c),ldc)
647 : !$OMP END TARGET DATA
648 : #endif
649 : else
650 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
651 : end if
652 :
653 0 : end subroutine abi_gpu_xgemm_2z
654 : !!***
655 :
656 : !------------------------------------------------------------------------------
657 :
658 : !!****f* m_abi_gpu_linalg/abi_gpu_xgemm_strided
659 : !! NAME
660 : !! abi_gpu_xgemm_strided
661 : !!
662 : !! FUNCTION
663 : !! Compute a batched scalar-matrix-matrix product and return a scalar-matrix product on GPU.
664 : !! Meant to be used on non-contiguous matrixes with data is uniformly split in the same number of batches in each matrix.
665 : !! c = alpha * op(a) * op(b) + beta * c
666 : !!
667 : !! INPUTS
668 : !! cplx = 1 if real 2 if complex
669 : !! transa = from of op(a) to be used in the matrix multiplication
670 : !! transb = from of op(b) to be used in the matrix multiplication
671 : !! m = number of rows of the matrix op(a) and of the matrix c
672 : !! n = number of rows of the matrix op(b) and the number of columns of the matrix c
673 : !! k = number of columns of the matrix op(a) and the number of rows of the matrix op(b)
674 : !! alpha = alpha scalar coefficient for matrix op(a)
675 : !! a = pointer to gpu memory location of matrix a
676 : !! lda = first dimension of a
677 : !! strideC = stride between each batch in matrix a
678 : !! b = pointer to gpu memory location of matrix b
679 : !! ldb = first dimension of b
680 : !! strideC = stride between each batch in matrix b
681 : !! beta = beta scalar coefficient for matrix c
682 : !! c = pointer to gpu memory location of matrix c
683 : !! ldc = first dimension of c
684 : !! strideC = stride between each batch in matrix c
685 : !! batchCount = number of batches in any matrix
686 : !!
687 : !! OUTPUT
688 : !! c = c matrix
689 : !!
690 : !! SOURCE
691 :
692 0 : subroutine abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
693 : a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount,async,stream_id)
694 :
695 : !Arguments ------------------------------------
696 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
697 : integer,intent(in) :: strideA,strideB,strideC,batchCount
698 : logical,intent(in),optional :: async
699 : integer,intent(in),optional :: stream_id
700 : complex(dp),intent(in) :: alpha,beta
701 : character(len=1),intent(in) :: transa,transb
702 : type(c_ptr),intent(in) :: a,b
703 : type(c_ptr),intent(in) :: c
704 : !Locals ---------------------------------------
705 : logical :: async_
706 : integer :: stream_id_
707 : ! *********************************************************************
708 :
709 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
710 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
711 : end if
712 : async_=.false.
713 0 : if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) async_=.false.
714 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
715 : ! Therefore, we issue a stream sync here to avoid
716 : ! potential mistakes in calling context.
717 : if(present(async)) async_=async
718 0 : stream_id_=-1;
719 : if(present(stream_id)) then
720 : stream_id_=modulo(stream_id,32);
721 : async_=.true.
722 : end if
723 :
724 : #ifdef HAVE_GPU
725 : call gpu_xgemm_strided_batched(cplx,transa,transb,m,n,k,alpha,&
726 : a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount,stream_id_)
727 :
728 : if (.not. async_) then
729 : call gpu_linalg_stream_synchronize()
730 : end if
731 :
732 : #else
733 : ! Unused if GPU code disabled
734 : ABI_UNUSED((/cplx,lda,ldb,ldc,m,n,k/))
735 : ABI_UNUSED((/strideA,strideB,strideC,batchCount/))
736 : ABI_UNUSED((/alpha,beta/))
737 : ABI_UNUSED((/transa,transb/))
738 : ABI_UNUSED((/a,b,c/))
739 : #endif
740 :
741 0 : end subroutine abi_gpu_xgemm_strided_cptr
742 : !!***
743 :
744 0 : subroutine abi_gpu_xgemm_strided_d(cplx,transa,transb,m,n,k,alpha,a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount)
745 :
746 : !Arguments ------------------------------------
747 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
748 : integer,intent(in) :: strideA,strideB,strideC,batchCount
749 : complex(dp),intent(in) :: alpha,beta
750 : character(len=1),intent(in) :: transa,transb
751 : real(dp), intent(in),target :: a(*),b(*)
752 : real(dp), intent(inout),target :: c(*)
753 : ! *********************************************************************
754 :
755 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
756 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
757 : end if
758 :
759 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
760 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
761 : c_loc(a),lda,strideA,&
762 : c_loc(b),ldb,strideB,&
763 : beta,&
764 0 : c_loc(c),ldc,strideC,batchCount)
765 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
766 : #ifdef HAVE_OPENMP_OFFLOAD
767 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
768 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
769 : c_loc(a),lda,strideA,&
770 : c_loc(b),ldb,strideB,&
771 : beta,&
772 : c_loc(c),ldc,strideC,batchCount)
773 : !$OMP END TARGET DATA
774 : #endif
775 : else
776 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
777 : end if
778 :
779 0 : end subroutine abi_gpu_xgemm_strided_d
780 : !!***
781 :
782 0 : subroutine abi_gpu_xgemm_strided_z(cplx,transa,transb,m,n,k,alpha,a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount)
783 :
784 : !Arguments ------------------------------------
785 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
786 : integer,intent(in) :: strideA,strideB,strideC,batchCount
787 : complex(dp),intent(in) :: alpha,beta
788 : character(len=1),intent(in) :: transa,transb
789 : complex(dp),intent(in),target :: a(*),b(*)
790 : complex(dp),intent(inout),target :: c(*)
791 : ! *********************************************************************
792 :
793 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
794 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
795 : end if
796 :
797 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
798 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
799 : c_loc(a),lda,strideA,&
800 : c_loc(b),ldb,strideB,&
801 : beta,&
802 0 : c_loc(c),ldc,strideC,batchCount)
803 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
804 : #ifdef HAVE_OPENMP_OFFLOAD
805 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
806 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
807 : c_loc(a),lda,strideA,&
808 : c_loc(b),ldb,strideB,&
809 : beta,&
810 : c_loc(c),ldc,strideC,batchCount)
811 : !$OMP END TARGET DATA
812 : #endif
813 : else
814 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
815 : end if
816 :
817 0 : end subroutine abi_gpu_xgemm_strided_z
818 : !!***
819 :
820 0 : subroutine abi_gpu_xgemm_strided_2d(cplx,transa,transb,m,n,k,alpha,a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount)
821 :
822 : !Arguments ------------------------------------
823 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
824 : integer,intent(in) :: strideA,strideB,strideC,batchCount
825 : complex(dp),intent(in) :: alpha,beta
826 : character(len=1),intent(in) :: transa,transb
827 : real(dp), intent(in),target :: a(lda,*),b(ldb,*)
828 : real(dp), intent(inout),target :: c(ldc,*)
829 : ! *********************************************************************
830 :
831 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
832 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
833 : end if
834 :
835 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
836 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
837 : c_loc(a),lda,strideA,&
838 : c_loc(b),ldb,strideB,&
839 : beta,&
840 0 : c_loc(c),ldc,strideC,batchCount)
841 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
842 : #ifdef HAVE_OPENMP_OFFLOAD
843 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
844 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
845 : c_loc(a),lda,strideA,&
846 : c_loc(b),ldb,strideB,&
847 : beta,&
848 : c_loc(c),ldc,strideC,batchCount)
849 : !$OMP END TARGET DATA
850 : #endif
851 : else
852 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
853 : end if
854 :
855 0 : end subroutine abi_gpu_xgemm_strided_2d
856 : !!***
857 :
858 0 : subroutine abi_gpu_xgemm_strided_2z(cplx,transa,transb,m,n,k,alpha,a,lda,strideA,b,ldb,strideB,beta,c,ldc,strideC,batchCount)
859 :
860 : !Arguments ------------------------------------
861 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n,k
862 : integer,intent(in) :: strideA,strideB,strideC,batchCount
863 : complex(dp),intent(in) :: alpha,beta
864 : character(len=1),intent(in) :: transa,transb
865 : complex(dp),intent(in),target :: a(lda,*),b(ldb,*)
866 : complex(dp),intent(inout),target :: c(ldc,*)
867 : ! *********************************************************************
868 :
869 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
870 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
871 : end if
872 :
873 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
874 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
875 : c_loc(a),lda,strideA,&
876 : c_loc(b),ldb,strideB,&
877 : beta,&
878 0 : c_loc(c),ldc,strideC,batchCount)
879 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
880 : #ifdef HAVE_OPENMP_OFFLOAD
881 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
882 : call abi_gpu_xgemm_strided_cptr(cplx,transa,transb,m,n,k,alpha,&
883 : c_loc(a),lda,strideA,&
884 : c_loc(b),ldb,strideB,&
885 : beta,&
886 : c_loc(c),ldc,strideC,batchCount)
887 : !$OMP END TARGET DATA
888 : #endif
889 : else
890 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
891 : end if
892 :
893 0 : end subroutine abi_gpu_xgemm_strided_2z
894 : !!***
895 :
896 : !------------------------------------------------------------------------------
897 : !!****f* m_abi_gpu_linalg/abi_gpu_xsymm
898 : !! NAME
899 : !! abi_gpu_xsymm
900 : !!
901 : !! FUNCTION
902 : !! Compute a symmetric scalar-matrix-matrix product and return a scalar-matrix product on GPU
903 : !! c = alpha * op(a) * op(b) + beta * c , if side == L
904 : !! c = alpha * op(b) * op(a) + beta * c , if side == R
905 : !!
906 : !! INPUTS
907 : !! cplx = 1 if real 2 if complex
908 : !! side = Specifies whether op(a) appears on the left or right of x for
909 : !! the operation to be performed as follows:
910 : !! L or l op(a)*x = alpha*b
911 : !! R or r x*op(a) = alpha*b
912 : !! uplo = Specifies whether the matrix a is an upper or lower triangular
913 : !! matrix as follows:
914 : !! U or u Matrix a is an upper triangular matrix.
915 : !! L or l Matrix a is a lower triangular matrix
916 : !! m = number of rows of the matrix op(a) and of the matrix c
917 : !! n = number of rows of the matrix op(b) and the number of columns of the matrix c
918 : !! alpha = alpha scalar coefficient for matrix op(a)
919 : !! a = pointer to gpu memory location of matrix a
920 : !! lda = first dimension of a
921 : !! b = pointer to gpu memory location of matrix b
922 : !! ldb = first dimension of b
923 : !! beta = beta scalar coefficient for matrix c
924 : !! c = pointer to gpu memory location of matrix c
925 : !! ldc = first dimension of c
926 : !!
927 : !! OUTPUT
928 : !! c = c matrix
929 : !!
930 : !! SOURCE
931 :
932 0 : subroutine abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
933 :
934 : !Arguments ------------------------------------
935 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n
936 : complex(dp),intent(in) :: alpha,beta
937 : character(len=1),intent(in) :: side,uplo
938 : type(c_ptr),intent(in) :: a,b
939 : type(c_ptr),intent(in) :: c
940 : ! *********************************************************************
941 :
942 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
943 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
944 : end if
945 :
946 : #ifdef HAVE_GPU
947 :
948 : call gpu_xsymm(cplx,side,uplo,m,n,alpha,&
949 : a,lda,b,ldb,beta,c,ldc)
950 :
951 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
952 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
953 : ! Therefore, we issue a stream sync here to avoid
954 : !potential mistakes in calling context.
955 : call gpu_linalg_stream_synchronize()
956 : end if
957 :
958 : #else
959 : ! Unused if GPU code disabled
960 : ABI_UNUSED((/cplx,lda,ldb,ldc,m,n/))
961 : ABI_UNUSED((/alpha,beta/))
962 : ABI_UNUSED((/side,uplo/))
963 : ABI_UNUSED((/a,b,c/))
964 : #endif
965 :
966 0 : end subroutine abi_gpu_xsymm_cptr
967 : !!***
968 :
969 0 : subroutine abi_gpu_xsymm_d(cplx,side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
970 :
971 : !Arguments ------------------------------------
972 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n
973 : complex(dp),intent(in) :: alpha,beta
974 : character(len=1),intent(in) :: side,uplo
975 : real(dp), intent(in),target :: a(*),b(*)
976 : real(dp), intent(inout),target :: c(*)
977 : ! *********************************************************************
978 :
979 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
980 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
981 : end if
982 :
983 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
984 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
985 : c_loc(a),lda,&
986 : c_loc(b),ldb,&
987 : beta,&
988 0 : c_loc(c),ldc)
989 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
990 : #ifdef HAVE_OPENMP_OFFLOAD
991 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
992 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
993 : c_loc(a),lda,&
994 : c_loc(b),ldb,&
995 : beta,&
996 : c_loc(c),ldc)
997 : !$OMP END TARGET DATA
998 : #endif
999 : else
1000 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1001 : end if
1002 :
1003 0 : end subroutine abi_gpu_xsymm_d
1004 : !!***
1005 :
1006 0 : subroutine abi_gpu_xsymm_z(cplx,side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1007 :
1008 : !Arguments ------------------------------------
1009 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n
1010 : complex(dp),intent(in) :: alpha,beta
1011 : character(len=1),intent(in) :: side,uplo
1012 : complex(dp),intent(in),target :: a(*),b(*)
1013 : complex(dp),intent(inout),target :: c(*)
1014 : ! *********************************************************************
1015 :
1016 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1017 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1018 : end if
1019 :
1020 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1021 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1022 : c_loc(a),lda,&
1023 : c_loc(b),ldb,&
1024 : beta,&
1025 0 : c_loc(c),ldc)
1026 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1027 : #ifdef HAVE_OPENMP_OFFLOAD
1028 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1029 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1030 : c_loc(a),lda,&
1031 : c_loc(b),ldb,&
1032 : beta,&
1033 : c_loc(c),ldc)
1034 : !$OMP END TARGET DATA
1035 : #endif
1036 : else
1037 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1038 : end if
1039 :
1040 0 : end subroutine abi_gpu_xsymm_z
1041 : !!***
1042 :
1043 0 : subroutine abi_gpu_xsymm_2d(cplx,side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1044 :
1045 : !Arguments ------------------------------------
1046 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n
1047 : complex(dp),intent(in) :: alpha,beta
1048 : character(len=1),intent(in) :: side,uplo
1049 : real(dp), intent(in),target :: a(lda,*),b(ldb,*)
1050 : real(dp), intent(inout),target :: c(ldc,*)
1051 : ! *********************************************************************
1052 :
1053 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1054 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1055 : end if
1056 :
1057 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1058 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1059 : c_loc(a),lda,&
1060 : c_loc(b),ldb,&
1061 : beta,&
1062 0 : c_loc(c),ldc)
1063 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1064 : #ifdef HAVE_OPENMP_OFFLOAD
1065 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1066 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1067 : c_loc(a),lda,&
1068 : c_loc(b),ldb,&
1069 : beta,&
1070 : c_loc(c),ldc)
1071 : !$OMP END TARGET DATA
1072 : #endif
1073 : else
1074 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1075 : end if
1076 :
1077 0 : end subroutine abi_gpu_xsymm_2d
1078 : !!***
1079 :
1080 0 : subroutine abi_gpu_xsymm_2z(cplx,side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1081 :
1082 : !Arguments ------------------------------------
1083 : integer,intent(in) :: cplx,lda,ldb,ldc,m,n
1084 : complex(dp),intent(in) :: alpha,beta
1085 : character(len=1),intent(in) :: side,uplo
1086 : complex(dp),intent(in),target :: a(lda,*),b(ldb,*)
1087 : complex(dp),intent(inout),target :: c(ldc,*)
1088 : ! *********************************************************************
1089 :
1090 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1091 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1092 : end if
1093 :
1094 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1095 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1096 : c_loc(a),lda,&
1097 : c_loc(b),ldb,&
1098 : beta,&
1099 0 : c_loc(c),ldc)
1100 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1101 : #ifdef HAVE_OPENMP_OFFLOAD
1102 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1103 : call abi_gpu_xsymm_cptr(cplx,side,uplo,m,n,alpha,&
1104 : c_loc(a),lda,&
1105 : c_loc(b),ldb,&
1106 : beta,&
1107 : c_loc(c),ldc)
1108 : !$OMP END TARGET DATA
1109 : #endif
1110 : else
1111 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1112 : end if
1113 :
1114 0 : end subroutine abi_gpu_xsymm_2z
1115 : !!***
1116 :
1117 : !------------------------------------------------------------------------------
1118 : !!****f* m_abi_gpu_linalg/abi_gpu_zhemm
1119 : !! NAME
1120 : !! abi_gpu_xhemm
1121 : !!
1122 : !! FUNCTION
1123 : !! Compute a Hermitian scalar-matrix-matrix product and return a scalar-matrix product on GPU
1124 : !! c = alpha * op(a) * op(b) + beta * c , if side == L
1125 : !! c = alpha * op(b) * op(a) + beta * c , if side == R
1126 : !!
1127 : !! INPUTS
1128 : !! side = Specifies whether op(a) appears on the left or right of x for
1129 : !! the operation to be performed as follows:
1130 : !! L or l op(a)*x = alpha*b
1131 : !! R or r x*op(a) = alpha*b
1132 : !! uplo = Specifies whether the matrix a is an upper or lower triangular
1133 : !! matrix as follows:
1134 : !! U or u Matrix a is an upper triangular matrix.
1135 : !! L or l Matrix a is a lower triangular matrix
1136 : !! m = number of rows of the matrix op(a) and of the matrix c
1137 : !! n = number of rows of the matrix op(b) and the number of columns of the matrix c
1138 : !! alpha = alpha scalar coefficient for matrix op(a)
1139 : !! a_gpu = pointer to gpu memory location of matrix a
1140 : !! lda = first dimension of a
1141 : !! b_gpu = pointer to gpu memory location of matrix b
1142 : !! ldb = first dimension of b
1143 : !! beta = beta scalar coefficient for matrix c
1144 : !! c_gpu = pointer to gpu memory location of matrix c
1145 : !! ldc = first dimension of c
1146 : !!
1147 : !! OUTPUT
1148 : !! c = c matrix
1149 : !!
1150 : !! SOURCE
1151 :
1152 0 : subroutine abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1153 :
1154 : !Arguments ------------------------------------
1155 : integer,intent(in) :: lda,ldb,ldc,m,n
1156 : complex(dp),intent(in) :: alpha,beta
1157 : character(len=1),intent(in) :: side,uplo
1158 : type(c_ptr),intent(in) :: a,b
1159 : type(c_ptr),intent(in) :: c
1160 : ! *********************************************************************
1161 :
1162 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1163 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1164 : end if
1165 :
1166 : #ifdef HAVE_GPU
1167 :
1168 : call gpu_zhemm(side,uplo,m,n,alpha,&
1169 : a,lda,b,ldb,beta,c,ldc)
1170 :
1171 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1172 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
1173 : ! Therefore, we issue a stream sync here to avoid
1174 : !potential mistakes in calling context.
1175 : call gpu_linalg_stream_synchronize()
1176 : end if
1177 :
1178 : #else
1179 : ! Unused if GPU code disabled
1180 : ABI_UNUSED((/lda,ldb,ldc,m,n/))
1181 : ABI_UNUSED((/alpha,beta/))
1182 : ABI_UNUSED((/side,uplo/))
1183 : ABI_UNUSED((/a,b,c/))
1184 : #endif
1185 :
1186 0 : end subroutine abi_gpu_zhemm_cptr
1187 : !!***
1188 :
1189 0 : subroutine abi_gpu_zhemm_d(side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1190 :
1191 : !Arguments ------------------------------------
1192 : integer,intent(in) :: lda,ldb,ldc,m,n
1193 : complex(dp),intent(in) :: alpha,beta
1194 : character(len=1),intent(in) :: side,uplo
1195 : real(dp),intent(in),target :: a(*),b(*)
1196 : real(dp),intent(inout),target :: c(*)
1197 : ! *********************************************************************
1198 :
1199 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1200 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1201 : end if
1202 :
1203 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1204 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1205 : c_loc(a),lda,&
1206 : c_loc(b),ldb,&
1207 : beta,&
1208 0 : c_loc(c),ldc)
1209 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1210 : #ifdef HAVE_OPENMP_OFFLOAD
1211 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1212 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1213 : c_loc(a),lda,&
1214 : c_loc(b),ldb,&
1215 : beta,&
1216 : c_loc(c),ldc)
1217 : !$OMP END TARGET DATA
1218 : #endif
1219 : else
1220 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1221 : end if
1222 :
1223 0 : end subroutine abi_gpu_zhemm_d
1224 : !!***
1225 :
1226 0 : subroutine abi_gpu_zhemm_z(side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1227 :
1228 : !Arguments ------------------------------------
1229 : integer,intent(in) :: lda,ldb,ldc,m,n
1230 : complex(dp),intent(in) :: alpha,beta
1231 : character(len=1),intent(in) :: side,uplo
1232 : complex(dp),intent(in),target :: a(*),b(*)
1233 : complex(dp),intent(inout),target :: c(*)
1234 : ! *********************************************************************
1235 :
1236 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1237 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1238 : end if
1239 :
1240 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1241 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1242 : c_loc(a),lda,&
1243 : c_loc(b),ldb,&
1244 : beta,&
1245 0 : c_loc(c),ldc)
1246 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1247 : #ifdef HAVE_OPENMP_OFFLOAD
1248 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1249 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1250 : c_loc(a),lda,&
1251 : c_loc(b),ldb,&
1252 : beta,&
1253 : c_loc(c),ldc)
1254 : !$OMP END TARGET DATA
1255 : #endif
1256 : else
1257 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1258 : end if
1259 :
1260 0 : end subroutine abi_gpu_zhemm_z
1261 : !!***
1262 :
1263 0 : subroutine abi_gpu_zhemm_2d(side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1264 :
1265 : !Arguments ------------------------------------
1266 : integer,intent(in) :: lda,ldb,ldc,m,n
1267 : complex(dp),intent(in) :: alpha,beta
1268 : character(len=1),intent(in) :: side,uplo
1269 : real(dp),intent(in),target :: a(lda,*),b(ldb,*)
1270 : real(dp),intent(inout),target :: c(ldc,*)
1271 : ! *********************************************************************
1272 :
1273 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1274 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1275 : end if
1276 :
1277 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1278 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1279 : c_loc(a),lda,&
1280 : c_loc(b),ldb,&
1281 : beta,&
1282 0 : c_loc(c),ldc)
1283 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1284 : #ifdef HAVE_OPENMP_OFFLOAD
1285 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1286 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1287 : c_loc(a),lda,&
1288 : c_loc(b),ldb,&
1289 : beta,&
1290 : c_loc(c),ldc)
1291 : !$OMP END TARGET DATA
1292 : #endif
1293 : else
1294 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1295 : end if
1296 :
1297 0 : end subroutine abi_gpu_zhemm_2d
1298 : !!***
1299 :
1300 0 : subroutine abi_gpu_zhemm_2z(side,uplo,m,n,alpha,a,lda,b,ldb,beta,c,ldc)
1301 :
1302 : !Arguments ------------------------------------
1303 : integer,intent(in) :: lda,ldb,ldc,m,n
1304 : complex(dp),intent(in) :: alpha,beta
1305 : character(len=1),intent(in) :: side,uplo
1306 : complex(dp),intent(in),target :: a(lda,*),b(ldb,*)
1307 : complex(dp),intent(inout),target :: c(ldc,*)
1308 : ! *********************************************************************
1309 :
1310 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1311 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1312 : end if
1313 :
1314 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1315 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1316 : c_loc(a),lda,&
1317 : c_loc(b),ldb,&
1318 : beta,&
1319 0 : c_loc(c),ldc)
1320 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1321 : #ifdef HAVE_OPENMP_OFFLOAD
1322 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b,c)
1323 : call abi_gpu_zhemm_cptr(side,uplo,m,n,alpha,&
1324 : c_loc(a),lda,&
1325 : c_loc(b),ldb,&
1326 : beta,&
1327 : c_loc(c),ldc)
1328 : !$OMP END TARGET DATA
1329 : #endif
1330 : else
1331 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1332 : end if
1333 :
1334 0 : end subroutine abi_gpu_zhemm_2z
1335 : !!***
1336 :
1337 : !------------------------------------------------------------------------------
1338 : !!****f* m_abi_gpu_linalg/abi_gpu_xscal
1339 : !! NAME
1340 : !! abi_gpu_xscal
1341 : !!
1342 : !! FUNCTION
1343 : !! Compute a BLAS-1 SCAL operation on GPU
1344 : !! x = alpha * x
1345 : !!
1346 : !! INPUTS
1347 : !! cplx = 1 if real 2 if complex
1348 : !! size = vector size
1349 : !! alpha = scalar complex value
1350 : !! x_gpu = pointer to gpu memory location of array x
1351 : !! incrx = stride between consecutive elements of x
1352 : !!
1353 : !! SOURCE
1354 0 : subroutine abi_gpu_xscal_cptr(cplx, size, alpha, x, incrx)
1355 :
1356 : ! !Arguments ------------------------------------
1357 : integer, intent(in) :: cplx
1358 : integer, intent(in) :: size
1359 : complex(dp), intent(in) :: alpha
1360 : type(c_ptr), intent(in) :: x
1361 : integer, intent(in) :: incrx
1362 : ! *************************************************************************
1363 :
1364 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1365 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1366 : end if
1367 :
1368 : #ifdef HAVE_GPU
1369 :
1370 : call gpu_xscal(cplx, size, alpha, x, incrx)
1371 :
1372 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1373 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
1374 : ! Therefore, we issue a stream sync here to avoid
1375 : !potential mistakes in calling context.
1376 : call gpu_linalg_stream_synchronize()
1377 : end if
1378 :
1379 : #else
1380 : ! Unused if GPU code disabled
1381 : ABI_UNUSED((/cplx,incrx,size/))
1382 : ABI_UNUSED((/alpha/))
1383 : ABI_UNUSED((/x/))
1384 : #endif
1385 :
1386 0 : end subroutine abi_gpu_xscal_cptr
1387 : !!***
1388 :
1389 0 : subroutine abi_gpu_xscal_d(cplx, size, alpha, x, incrx)
1390 :
1391 : ! !Arguments ------------------------------------
1392 : integer, intent(in) :: cplx
1393 : integer, intent(in) :: size
1394 : complex(dp), intent(in) :: alpha
1395 : real(dp), intent(inout), target :: x(*)
1396 : integer, intent(in) :: incrx
1397 : ! *************************************************************************
1398 :
1399 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1400 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1401 : end if
1402 :
1403 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1404 0 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1405 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1406 : #ifdef HAVE_OPENMP_OFFLOAD
1407 : !$OMP TARGET DATA USE_DEVICE_ADDR(x)
1408 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1409 : !$OMP END TARGET DATA
1410 : #endif
1411 : else
1412 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1413 : end if
1414 :
1415 0 : end subroutine abi_gpu_xscal_d
1416 : !!***
1417 :
1418 0 : subroutine abi_gpu_xscal_z(cplx, size, alpha, x, incrx)
1419 :
1420 : ! !Arguments ------------------------------------
1421 : integer, intent(in) :: cplx
1422 : integer, intent(in) :: size
1423 : complex(dp), intent(in) :: alpha
1424 : complex(dp), intent(inout), target :: x(*)
1425 : integer, intent(in) :: incrx
1426 : ! *************************************************************************
1427 :
1428 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1429 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1430 : end if
1431 :
1432 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1433 0 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1434 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1435 : #ifdef HAVE_OPENMP_OFFLOAD
1436 : !$OMP TARGET DATA USE_DEVICE_ADDR(x)
1437 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1438 : !$OMP END TARGET DATA
1439 : #endif
1440 : else
1441 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1442 : end if
1443 :
1444 0 : end subroutine abi_gpu_xscal_z
1445 : !!***
1446 :
1447 0 : subroutine abi_gpu_xscal_2d(cplx, size, alpha, x, incrx)
1448 :
1449 : ! !Arguments ------------------------------------
1450 : integer, intent(in) :: cplx
1451 : integer, intent(in) :: size
1452 : complex(dp), intent(in) :: alpha
1453 : real(dp), intent(inout), target :: x(size,*)
1454 : integer, intent(in) :: incrx
1455 : ! *************************************************************************
1456 :
1457 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1458 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1459 : end if
1460 :
1461 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1462 0 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1463 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1464 : #ifdef HAVE_OPENMP_OFFLOAD
1465 : !$OMP TARGET DATA USE_DEVICE_ADDR(x)
1466 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1467 : !$OMP END TARGET DATA
1468 : #endif
1469 : else
1470 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1471 : end if
1472 :
1473 0 : end subroutine abi_gpu_xscal_2d
1474 : !!***
1475 :
1476 0 : subroutine abi_gpu_xscal_2z(cplx, size, alpha, x, incrx)
1477 :
1478 : ! !Arguments ------------------------------------
1479 : integer, intent(in) :: cplx
1480 : integer, intent(in) :: size
1481 : complex(dp), intent(in) :: alpha
1482 : complex(dp), intent(inout), target :: x(size,*)
1483 : integer, intent(in) :: incrx
1484 : ! *************************************************************************
1485 :
1486 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1487 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1488 : end if
1489 :
1490 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1491 0 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1492 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1493 : #ifdef HAVE_OPENMP_OFFLOAD
1494 : !$OMP TARGET DATA USE_DEVICE_ADDR(x)
1495 : call abi_gpu_xscal_cptr(cplx, size, alpha, c_loc(x), incrx)
1496 : !$OMP END TARGET DATA
1497 : #endif
1498 : else
1499 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1500 : end if
1501 :
1502 0 : end subroutine abi_gpu_xscal_2z
1503 : !!***
1504 :
1505 : !!****f* m_abi_gpu_linalg/abi_gpu_xdot
1506 : !! NAME
1507 : !! abi_gpu_xdot
1508 : !!
1509 : !! FUNCTION
1510 : !! Compute blas-3 DOT on GPU
1511 : !! trace(X^H*Y)
1512 : !!
1513 : !! INPUTS
1514 : !! cplx = 1 if real 2 if complex
1515 : !! size = vector size
1516 : !! alpha = scalar complex value
1517 : !! x_gpu = pointer to gpu memory location of array x
1518 : !! incrx = stride between consecutive elements of x
1519 : !! y_gpu = pointer to gpu memory location of array y
1520 : !! incry = stride between consecutive elements of y
1521 : !!
1522 : !! SOURCE
1523 0 : subroutine abi_gpu_xdot_cptr(cplx, size, alpha, x, incrx, y, incry)
1524 :
1525 : ! !Arguments ------------------------------------
1526 : integer, intent(in) :: cplx
1527 : integer, intent(in) :: size
1528 : complex(dp), intent(in) :: alpha
1529 : type(c_ptr), intent(in) :: x
1530 : integer, intent(in) :: incrx
1531 : type(c_ptr), intent(in) :: y
1532 : integer, intent(in) :: incry
1533 : ! *************************************************************************
1534 :
1535 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1536 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1537 : end if
1538 :
1539 : #ifdef HAVE_GPU
1540 :
1541 : call gpu_xdot(cplx, size, alpha, x, incrx, y, incry)
1542 :
1543 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1544 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
1545 : ! Therefore, we issue a stream sync here to avoid
1546 : !potential mistakes in calling context.
1547 : call gpu_linalg_stream_synchronize()
1548 : end if
1549 :
1550 : #else
1551 : ! Unused if GPU code disabled
1552 : ABI_UNUSED((/cplx,incrx,incry,size/))
1553 : ABI_UNUSED((/alpha/))
1554 : ABI_UNUSED((/x,y/))
1555 : #endif
1556 :
1557 0 : end subroutine abi_gpu_xdot_cptr
1558 : !!***
1559 :
1560 0 : subroutine abi_gpu_xdot_d(cplx, size, alpha, x, incrx, y, incry)
1561 :
1562 : !Arguments ------------------------------------
1563 : integer, intent(in) :: cplx
1564 : integer, intent(in) :: size
1565 : complex(dp), intent(in) :: alpha
1566 : real(dp), intent(in), target :: x(*)
1567 : integer, intent(in) :: incrx
1568 : real(dp), intent(inout), target :: y(*)
1569 : integer, intent(in) :: incry
1570 : ! *************************************************************************
1571 :
1572 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1573 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1574 : end if
1575 :
1576 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1577 0 : call abi_gpu_xdot_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1578 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1579 : #ifdef HAVE_OPENMP_OFFLOAD
1580 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1581 : call abi_gpu_xdot_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1582 : !$OMP END TARGET DATA
1583 : #endif
1584 : else
1585 0 : ABI_BUG("Unhandled GPU mode !")
1586 : end if
1587 :
1588 0 : end subroutine abi_gpu_xdot_d
1589 : !!***
1590 :
1591 0 : subroutine abi_gpu_xdot_z(cplx, size, alpha, x, incrx, y, incry)
1592 :
1593 : !Arguments ------------------------------------
1594 : integer, intent(in) :: cplx
1595 : integer, intent(in) :: size
1596 : complex(dp), intent(in) :: alpha
1597 : complex(dp), intent(in), target :: x(*)
1598 : integer, intent(in) :: incrx
1599 : complex(dp), intent(inout), target :: y(*)
1600 : integer, intent(in) :: incry
1601 : ! *************************************************************************
1602 :
1603 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1604 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1605 : end if
1606 :
1607 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1608 0 : call abi_gpu_xdot_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1609 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1610 : #ifdef HAVE_OPENMP_OFFLOAD
1611 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1612 : call abi_gpu_xdot_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1613 : !$OMP END TARGET DATA
1614 : #endif
1615 : else
1616 0 : ABI_BUG("Unhandled GPU mode !")
1617 : end if
1618 :
1619 0 : end subroutine abi_gpu_xdot_z
1620 : !!***
1621 :
1622 : !!****f* m_abi_gpu_linalg/abi_gpu_xaxpy
1623 : !! NAME
1624 : !! abi_gpu_xaxpy
1625 : !!
1626 : !! FUNCTION
1627 : !! Compute a BLAS-1 AXPY operation on GPU
1628 : !! y = alpha * x + y
1629 : !!
1630 : !! INPUTS
1631 : !! cplx = 1 if real 2 if complex
1632 : !! size = vector size
1633 : !! alpha = scalar complex value
1634 : !! x_gpu = pointer to gpu memory location of array x
1635 : !! incrx = stride between consecutive elements of x
1636 : !! y_gpu = pointer to gpu memory location of array y
1637 : !! incry = stride between consecutive elements of y
1638 : !!
1639 : !! SOURCE
1640 0 : subroutine abi_gpu_xaxpy_cptr(cplx, size, alpha, x, incrx, y, incry)
1641 :
1642 : ! !Arguments ------------------------------------
1643 : integer, intent(in) :: cplx
1644 : integer, intent(in) :: size
1645 : complex(dp), intent(in) :: alpha
1646 : type(c_ptr), intent(in) :: x
1647 : integer, intent(in) :: incrx
1648 : type(c_ptr), intent(in) :: y
1649 : integer, intent(in) :: incry
1650 : ! *************************************************************************
1651 :
1652 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1653 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1654 : end if
1655 :
1656 : #ifdef HAVE_GPU
1657 :
1658 : call gpu_xaxpy(cplx, size, alpha, x, incrx, y, incry)
1659 :
1660 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1661 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
1662 : ! Therefore, we issue a stream sync here to avoid
1663 : !potential mistakes in calling context.
1664 : call gpu_linalg_stream_synchronize()
1665 : end if
1666 :
1667 : #else
1668 : ! Unused if GPU code disabled
1669 : ABI_UNUSED((/cplx,incrx,incry,size/))
1670 : ABI_UNUSED((/alpha/))
1671 : ABI_UNUSED((/x,y/))
1672 : #endif
1673 :
1674 0 : end subroutine abi_gpu_xaxpy_cptr
1675 : !!***
1676 :
1677 0 : subroutine abi_gpu_xaxpy_d(cplx, size, alpha, x, incrx, y, incry)
1678 :
1679 : !Arguments ------------------------------------
1680 : integer, intent(in) :: cplx
1681 : integer, intent(in) :: size
1682 : complex(dp), intent(in) :: alpha
1683 : real(dp), intent(in), target :: x(*)
1684 : integer, intent(in) :: incrx
1685 : real(dp), intent(inout), target :: y(*)
1686 : integer, intent(in) :: incry
1687 : ! *************************************************************************
1688 :
1689 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1690 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1691 : end if
1692 :
1693 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1694 0 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1695 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1696 : #ifdef HAVE_OPENMP_OFFLOAD
1697 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1698 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1699 : !$OMP END TARGET DATA
1700 : #endif
1701 : else
1702 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1703 : end if
1704 :
1705 0 : end subroutine abi_gpu_xaxpy_d
1706 : !!***
1707 :
1708 0 : subroutine abi_gpu_xaxpy_z(cplx, size, alpha, x, incrx, y, incry)
1709 :
1710 : !Arguments ------------------------------------
1711 : integer, intent(in) :: cplx
1712 : integer, intent(in) :: size
1713 : complex(dp), intent(in) :: alpha
1714 : complex(dp), intent(in), target :: x(*)
1715 : integer, intent(in) :: incrx
1716 : complex(dp), intent(inout), target :: y(*)
1717 : integer, intent(in) :: incry
1718 : ! *************************************************************************
1719 :
1720 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1721 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1722 : end if
1723 :
1724 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1725 0 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1726 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1727 : #ifdef HAVE_OPENMP_OFFLOAD
1728 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1729 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1730 : !$OMP END TARGET DATA
1731 : #endif
1732 : else
1733 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1734 : end if
1735 :
1736 0 : end subroutine abi_gpu_xaxpy_z
1737 : !!***
1738 :
1739 0 : subroutine abi_gpu_xaxpy_2d(cplx, size, alpha, x, incrx, y, incry)
1740 :
1741 : !Arguments ------------------------------------
1742 : integer, intent(in) :: cplx
1743 : integer, intent(in) :: size
1744 : complex(dp), intent(in) :: alpha
1745 : real(dp), intent(in), target :: x(size,*)
1746 : integer, intent(in) :: incrx
1747 : real(dp), intent(inout), target :: y(size,*)
1748 : integer, intent(in) :: incry
1749 : ! *************************************************************************
1750 :
1751 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1752 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1753 : end if
1754 :
1755 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1756 0 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1757 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1758 : #ifdef HAVE_OPENMP_OFFLOAD
1759 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1760 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1761 : !$OMP END TARGET DATA
1762 : #endif
1763 : else
1764 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1765 : end if
1766 :
1767 0 : end subroutine abi_gpu_xaxpy_2d
1768 : !!***
1769 :
1770 0 : subroutine abi_gpu_xaxpy_2z(cplx, size, alpha, x, incrx, y, incry)
1771 :
1772 : !Arguments ------------------------------------
1773 : integer, intent(in) :: cplx
1774 : integer, intent(in) :: size
1775 : complex(dp), intent(in) :: alpha
1776 : complex(dp), intent(in), target :: x(size,*)
1777 : integer, intent(in) :: incrx
1778 : complex(dp), intent(inout), target :: y(size,*)
1779 : integer, intent(in) :: incry
1780 : ! *************************************************************************
1781 :
1782 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1783 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1784 : end if
1785 :
1786 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1787 0 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1788 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1789 : #ifdef HAVE_OPENMP_OFFLOAD
1790 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1791 : call abi_gpu_xaxpy_cptr(cplx, size, alpha, c_loc(x), incrx, c_loc(y), incry)
1792 : !$OMP END TARGET DATA
1793 : #endif
1794 : else
1795 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1796 : end if
1797 :
1798 0 : end subroutine abi_gpu_xaxpy_2z
1799 : !!***
1800 :
1801 : !!****f* m_abi_gpu_linalg/abi_gpu_xcopy
1802 : !! NAME
1803 : !! abi_gpu_xcopy
1804 : !!
1805 : !! FUNCTION
1806 : !! Compute a BLAS-1 COPY operation on GPU
1807 : !! y = x (copy x into y)
1808 : !!
1809 : !! INPUTS
1810 : !! cplx = 1 if real 2 if complex
1811 : !! size = input vector size
1812 : !! x = pointer to gpu memory location of array x
1813 : !! incrx = stride between consecutive elements of x
1814 : !! y = pointer to gpu memory location of array y
1815 : !! incry = stride between consecutive elements of y
1816 : !!
1817 : !! SOURCE
1818 :
1819 0 : subroutine abi_gpu_xcopy_cptr(cplx, size, x, incrx, y, incry)
1820 :
1821 : ! !Arguments ------------------------------------
1822 : integer, intent(in) :: cplx
1823 : integer, intent(in) :: size
1824 : type(c_ptr), intent(in) :: x
1825 : integer, intent(in) :: incrx
1826 : type(c_ptr), intent(in) :: y
1827 : integer, intent(in) :: incry
1828 :
1829 : ! *************************************************************************
1830 :
1831 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1832 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1833 : end if
1834 :
1835 : #ifdef HAVE_GPU
1836 :
1837 : call gpu_xcopy(cplx, size, x, incrx, y, incry)
1838 :
1839 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1840 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
1841 : ! Therefore, we issue a stream sync here to avoid
1842 : !potential mistakes in calling context.
1843 : call gpu_linalg_stream_synchronize()
1844 : end if
1845 :
1846 : #else
1847 : ! Unused if GPU code disabled
1848 : ABI_UNUSED((/cplx,incrx,incry,size/))
1849 : ABI_UNUSED((/x,y/))
1850 : #endif
1851 :
1852 0 : end subroutine abi_gpu_xcopy_cptr
1853 : !!***
1854 :
1855 0 : subroutine abi_gpu_xcopy_d(cplx, size, x, incrx, y, incry)
1856 :
1857 : ! !Arguments ------------------------------------
1858 : integer, intent(in) :: cplx
1859 : integer, intent(in) :: size
1860 : real(dp), intent(in),target :: x(*)
1861 : integer, intent(in) :: incrx
1862 : real(dp), intent(inout),target :: y(*)
1863 : integer, intent(in) :: incry
1864 :
1865 : ! *************************************************************************
1866 :
1867 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1868 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1869 : end if
1870 :
1871 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1872 0 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1873 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1874 : #ifdef HAVE_OPENMP_OFFLOAD
1875 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1876 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1877 : !$OMP END TARGET DATA
1878 : #endif
1879 : else
1880 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1881 : end if
1882 :
1883 0 : end subroutine abi_gpu_xcopy_d
1884 : !!***
1885 :
1886 0 : subroutine abi_gpu_xcopy_z(cplx, size, x, incrx, y, incry)
1887 :
1888 : ! !Arguments ------------------------------------
1889 : integer, intent(in) :: cplx
1890 : integer, intent(in) :: size
1891 : complex(dp), intent(in),target :: x(*)
1892 : integer, intent(in) :: incrx
1893 : complex(dp), intent(inout),target :: y(*)
1894 : integer, intent(in) :: incry
1895 : ! *************************************************************************
1896 :
1897 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1898 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1899 : end if
1900 :
1901 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1902 0 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1903 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1904 : #ifdef HAVE_OPENMP_OFFLOAD
1905 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1906 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1907 : !$OMP END TARGET DATA
1908 : #endif
1909 : else
1910 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1911 : end if
1912 :
1913 0 : end subroutine abi_gpu_xcopy_z
1914 : !!***
1915 :
1916 0 : subroutine abi_gpu_xcopy_2d(cplx, size, x, incrx, y, incry)
1917 :
1918 : ! !Arguments ------------------------------------
1919 : integer, intent(in) :: cplx
1920 : integer, intent(in) :: size
1921 : real(dp), intent(in),target :: x(size,*)
1922 : integer, intent(in) :: incrx
1923 : real(dp), intent(inout),target :: y(size,*)
1924 : integer, intent(in) :: incry
1925 :
1926 : ! *************************************************************************
1927 :
1928 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1929 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1930 : end if
1931 :
1932 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1933 0 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1934 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1935 : #ifdef HAVE_OPENMP_OFFLOAD
1936 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1937 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1938 : !$OMP END TARGET DATA
1939 : #endif
1940 : else
1941 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1942 : end if
1943 :
1944 0 : end subroutine abi_gpu_xcopy_2d
1945 : !!***
1946 :
1947 0 : subroutine abi_gpu_xcopy_2z(cplx, size, x, incrx, y, incry)
1948 :
1949 : !Arguments ------------------------------------
1950 : integer, intent(in) :: cplx
1951 : integer, intent(in) :: size
1952 : complex(dp), intent(in),target :: x(size,*)
1953 : integer, intent(in) :: incrx
1954 : complex(dp), intent(inout),target :: y(size,*)
1955 : integer, intent(in) :: incry
1956 : ! *************************************************************************
1957 :
1958 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
1959 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
1960 : end if
1961 :
1962 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
1963 0 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1964 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
1965 : #ifdef HAVE_OPENMP_OFFLOAD
1966 : !$OMP TARGET DATA USE_DEVICE_ADDR(x,y)
1967 : call abi_gpu_xcopy_cptr(cplx, size, c_loc(x), incrx, c_loc(y), incry)
1968 : !$OMP END TARGET DATA
1969 : #endif
1970 : else
1971 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
1972 : end if
1973 :
1974 0 : end subroutine abi_gpu_xcopy_2z
1975 : !!***
1976 :
1977 : !!****f* m_abi_gpu_linalg/abi_gpu_xtrsm
1978 : !! NAME
1979 : !! abi_gpu_xtrsm
1980 : !!
1981 : !! FUNCTION
1982 : !! Solves a matrix equation (one matrix operand is triangular) on GPU.
1983 : !! The xtrsm routines solve one of the following matrix equations
1984 : !! op(a)*x = alpha*b
1985 : !! or
1986 : !! x*op(a) = alpha*b,
1987 : !!
1988 : !! INPUTS
1989 : !! cplx= 1 if real 2 if complex
1990 : !! side= Specifies whether op(a) appears on the left or right of x for
1991 : !! the operation to be performed as follows:
1992 : !! L or l op(a)*x = alpha*b
1993 : !! R or r x*op(a) = alpha*b
1994 : !! uplo= Specifies whether the matrix a is an upper or lower triangular
1995 : !! matrix as follows:
1996 : !! U or u Matrix a is an upper triangular matrix.
1997 : !! L or l Matrix a is a lower triangular matrix
1998 : !! transa= Specifies the form of op(a) to be used in the matrix
1999 : !! multiplication as follows:
2000 : !! N or n op(a) = a
2001 : !! T or t op(a) = a'
2002 : !! C or c op(a) = conjg(a')
2003 : !! diag= Specifies whether or not a is unit triangular as follows:
2004 : !! U or u Matrix a is assumed to be unit triangular.
2005 : !! N or n Matrix a is not assumed to be unit triangular.
2006 : !! m= Specifies the number of rows of b. The value of m must be at least zero
2007 : !! n= Specifies the number of columns of b. The value of n must be at least zero
2008 : !! alpha= Specifies the scalar alpha. When alpha is zero, then a is not referenced and b
2009 : !! need not be set before entry.
2010 : !! a = pointer to gpu memory location of array a, DIMENSION (lda, k), where k is m when side = 'L' or 'l' and is n
2011 : !! when side = 'R' or 'r'.
2012 : !! lda= Specifies the first dimension of a as declared in the calling
2013 : !! (sub)program. When side = 'L' or 'l', then lda must be at least max(1,
2014 : !! m), when side = 'R' or 'r', then lda must be at least max(1, n).
2015 : !! b = pointer to gpu memory location of b Array, DIMENSION (ldb,n). Before entry, the leading m-by-n part of the array
2016 : !! b must contain the right-hand side matrix b.
2017 : !! ldb= Specifies the first dimension of b as declared in the calling
2018 : !! (sub)program. The value of ldb must be at least max(1, m).
2019 : !!
2020 : !! OUTPUT
2021 : !! b
2022 : !!
2023 : !! SIDE EFFECTS
2024 : !! WARNING! : this routine is a dummy one when HAVE_GPU_CUDA is not enabled
2025 : !! the correct one is in 17_toolbox/gpu_linalg.cu
2026 : !!
2027 : !! SOURCE
2028 :
2029 0 : subroutine abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2030 :
2031 : ! !Arguments ------------------------------------
2032 : integer, intent(in) :: cplx,lda,ldb,m,n
2033 : complex(dp), intent(in) :: alpha
2034 : character(len=1), intent(in) :: side,uplo,transa,diag
2035 : type(c_ptr),intent(in) :: a
2036 : type(c_ptr),intent(in) :: b
2037 : ! *********************************************************************
2038 :
2039 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2040 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2041 : end if
2042 :
2043 : #ifdef HAVE_GPU
2044 :
2045 : call gpu_xtrsm(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2046 :
2047 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2048 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
2049 : ! Therefore, we issue a stream sync here to avoid
2050 : !potential mistakes in calling context.
2051 : call gpu_linalg_stream_synchronize()
2052 : end if
2053 :
2054 : #else
2055 : ! Unused if GPU code disabled
2056 : ABI_UNUSED((/cplx,lda,ldb,m,n/))
2057 : ABI_UNUSED((/alpha/))
2058 : ABI_UNUSED((/side,uplo,transa,diag/))
2059 : ABI_UNUSED((/a,b/))
2060 : #endif
2061 :
2062 0 : end subroutine abi_gpu_xtrsm_cptr
2063 : !!***
2064 :
2065 0 : subroutine abi_gpu_xtrsm_d(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2066 :
2067 : ! !Arguments ------------------------------------
2068 : integer, intent(in) :: cplx,lda,ldb,m,n
2069 : complex(dp), intent(in) :: alpha
2070 : character(len=1), intent(in) :: side,uplo,transa,diag
2071 : real(dp), intent(in),target :: a(*)
2072 : real(dp), intent(inout),target :: b(*)
2073 : ! *********************************************************************
2074 :
2075 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2076 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2077 : end if
2078 :
2079 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2080 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2081 0 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2082 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2083 : #ifdef HAVE_OPENMP_OFFLOAD
2084 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b)
2085 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2086 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2087 : !$OMP END TARGET DATA
2088 : #endif
2089 : else
2090 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2091 : end if
2092 :
2093 0 : end subroutine abi_gpu_xtrsm_d
2094 : !!***
2095 :
2096 0 : subroutine abi_gpu_xtrsm_z(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2097 :
2098 : ! !Arguments ------------------------------------
2099 : integer, intent(in) :: cplx,lda,ldb,m,n
2100 : complex(dp), intent(in) :: alpha
2101 : character(len=1), intent(in) :: side,uplo,transa,diag
2102 : complex(dp),intent(in),target :: a(*)
2103 : complex(dp),intent(inout),target :: b(*)
2104 : ! *********************************************************************
2105 :
2106 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2107 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2108 : end if
2109 :
2110 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2111 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2112 0 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2113 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2114 : #ifdef HAVE_OPENMP_OFFLOAD
2115 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b)
2116 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2117 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2118 : !$OMP END TARGET DATA
2119 : #endif
2120 : else
2121 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2122 : end if
2123 :
2124 0 : end subroutine abi_gpu_xtrsm_z
2125 : !!***
2126 :
2127 0 : subroutine abi_gpu_xtrsm_2d(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2128 :
2129 : !Arguments ------------------------------------
2130 : integer, intent(in) :: cplx,lda,ldb,m,n
2131 : complex(dp), intent(in) :: alpha
2132 : character(len=1), intent(in) :: side,uplo,transa,diag
2133 : real(dp), intent(in),target :: a(lda,*)
2134 : real(dp), intent(inout),target :: b(ldb,*)
2135 : ! *********************************************************************
2136 :
2137 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2138 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2139 : end if
2140 :
2141 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2142 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2143 0 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2144 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2145 : #ifdef HAVE_OPENMP_OFFLOAD
2146 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b)
2147 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2148 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2149 : !$OMP END TARGET DATA
2150 : #endif
2151 : else
2152 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2153 : end if
2154 :
2155 0 : end subroutine abi_gpu_xtrsm_2d
2156 : !!***
2157 :
2158 0 : subroutine abi_gpu_xtrsm_2z(cplx,side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb)
2159 :
2160 : !Arguments ------------------------------------
2161 : integer, intent(in) :: cplx,lda,ldb,m,n
2162 : complex(dp), intent(in) :: alpha
2163 : character(len=1), intent(in) :: side,uplo,transa,diag
2164 : complex(dp),intent(in),target :: a(lda,*)
2165 : complex(dp),intent(inout),target :: b(ldb,*)
2166 : ! *********************************************************************
2167 :
2168 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2169 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2170 : end if
2171 :
2172 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2173 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2174 0 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2175 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2176 : #ifdef HAVE_OPENMP_OFFLOAD
2177 : !$OMP TARGET DATA USE_DEVICE_ADDR(a,b)
2178 : call abi_gpu_xtrsm_cptr(cplx,side,uplo,transa,&
2179 : diag,m,n,alpha,c_loc(a),lda,c_loc(b),ldb)
2180 : !$OMP END TARGET DATA
2181 : #endif
2182 : else
2183 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2184 : end if
2185 :
2186 0 : end subroutine abi_gpu_xtrsm_2z
2187 : !!***
2188 :
2189 :
2190 : !!****f* m_abi_gpu_linalg/abi_gpu_work_resizeI
2191 : !!
2192 : !! NAME
2193 : !! abi_gpu_work_resizeI
2194 :
2195 : subroutine abi_gpu_work_resizeI(array,array_managed,current_dim,asked_dim)
2196 :
2197 : integer, allocatable, intent(inout) :: array(:)
2198 : integer(kind=c_int32_t), contiguous, pointer, intent(inout) :: array_managed(:)
2199 : integer, intent(inout) :: current_dim
2200 : integer, intent(in ) :: asked_dim
2201 : ! *************************************************************************
2202 :
2203 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2204 :
2205 : #ifdef HAVE_YAKL
2206 : if ( current_dim < asked_dim ) then
2207 : current_dim = asked_dim
2208 : if ( associated(array_managed) ) then
2209 : ABI_FREE_MANAGED(array_managed)
2210 : end if
2211 : ABI_MALLOC_MANAGED(array_managed,(/asked_dim/))
2212 : end if
2213 : #endif
2214 :
2215 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2216 :
2217 : #ifdef HAVE_OPENMP_OFFLOAD
2218 : if ( current_dim < asked_dim ) then
2219 : current_dim = asked_dim
2220 : if ( allocated(array) ) then
2221 : !$OMP TARGET EXIT DATA MAP(release:array)
2222 : ABI_FREE(array)
2223 : end if
2224 : ABI_MALLOC(array,(asked_dim))
2225 : !$OMP TARGET ENTER DATA MAP(alloc:array)
2226 : end if
2227 : #endif
2228 :
2229 : else
2230 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2231 : end if
2232 :
2233 : #ifndef HAVE_GPU
2234 : ! Unused if GPU code disabled
2235 : ABI_UNUSED(array)
2236 : ABI_UNUSED(array_managed)
2237 : ABI_UNUSED((/current_dim,asked_dim/))
2238 : #endif
2239 :
2240 : end subroutine abi_gpu_work_resizeI
2241 : !!***
2242 :
2243 : !!****f* m_abi_gpu_linalg/abi_gpu_work_resizeR
2244 : !!
2245 : !! NAME
2246 : !! abi_gpu_work_resizeR
2247 :
2248 : subroutine abi_gpu_work_resizeR(array,array_managed,current_dim,asked_dim)
2249 :
2250 : double precision, allocatable, intent(inout) :: array(:)
2251 : real(kind=c_double), contiguous, pointer, intent(inout) :: array_managed(:)
2252 : integer, intent(inout) :: current_dim
2253 : integer, intent(in ) :: asked_dim
2254 :
2255 : ! *************************************************************************
2256 :
2257 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2258 :
2259 : #ifdef HAVE_YAKL
2260 : if ( current_dim < asked_dim ) then
2261 : current_dim = asked_dim
2262 : if ( associated(array_managed) ) then
2263 : ABI_FREE_MANAGED(array_managed)
2264 : end if
2265 : ABI_MALLOC_MANAGED(array_managed,(/asked_dim/))
2266 : end if
2267 : #endif
2268 :
2269 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2270 :
2271 : #ifdef HAVE_OPENMP_OFFLOAD
2272 : if ( current_dim < asked_dim ) then
2273 : current_dim = asked_dim
2274 : if ( allocated(array) ) then
2275 : !$OMP TARGET EXIT DATA MAP(release:array)
2276 : ABI_FREE(array)
2277 : end if
2278 : ABI_MALLOC(array,(asked_dim))
2279 : !$OMP TARGET ENTER DATA MAP(alloc:array)
2280 : end if
2281 : #endif
2282 :
2283 : else
2284 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2285 : end if
2286 :
2287 : #ifndef HAVE_GPU
2288 : ! Unused if GPU code disabled
2289 : ABI_UNUSED(array)
2290 : ABI_UNUSED(array_managed)
2291 : ABI_UNUSED((/current_dim,asked_dim/))
2292 : #endif
2293 :
2294 : end subroutine abi_gpu_work_resizeR
2295 : !!***
2296 :
2297 :
2298 : !!****f* m_abi_gpu_linalg/abi_gpu_work_resizeC
2299 : !!
2300 : !! NAME
2301 : !! abi_gpu_work_resizeC
2302 :
2303 : subroutine abi_gpu_work_resizeC(array,array_managed,current_dim,asked_dim)
2304 :
2305 : complex(kind=8), allocatable, intent(inout) :: array(:)
2306 : complex(kind=c_double_complex), contiguous, pointer, intent(inout) :: array_managed(:)
2307 : integer, intent(inout) :: current_dim
2308 : integer, intent(in ) :: asked_dim
2309 : ! *************************************************************************
2310 :
2311 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2312 :
2313 : #ifdef HAVE_YAKL
2314 : if ( current_dim < asked_dim ) then
2315 : current_dim = asked_dim
2316 : if ( associated(array_managed) ) then
2317 : ABI_FREE_MANAGED(array_managed)
2318 : end if
2319 : ABI_MALLOC_MANAGED(array_managed,(/asked_dim/))
2320 : end if
2321 : #endif
2322 :
2323 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2324 :
2325 : #ifdef HAVE_OPENMP_OFFLOAD
2326 : if ( current_dim < asked_dim ) then
2327 : current_dim = asked_dim
2328 : if ( allocated(array) ) then
2329 : !$OMP TARGET EXIT DATA MAP(release:array)
2330 : ABI_FREE(array)
2331 : end if
2332 : ABI_MALLOC(array,(asked_dim))
2333 : !$OMP TARGET ENTER DATA MAP(alloc:array)
2334 : end if
2335 : #endif
2336 :
2337 : else
2338 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2339 : end if
2340 :
2341 : #ifndef HAVE_GPU
2342 : ! Unused if GPU code disabled
2343 : ABI_UNUSED(array)
2344 : ABI_UNUSED(array_managed)
2345 : ABI_UNUSED((/current_dim,asked_dim/))
2346 : #endif
2347 :
2348 : end subroutine abi_gpu_work_resizeC
2349 : !!***
2350 :
2351 : !!****f* m_abi_gpu_linalg/abi_gpu_work_resizeCptr
2352 : !!
2353 : !! NAME
2354 : !! abi_gpu_work_resizeCptr
2355 :
2356 : subroutine abi_gpu_work_resizeCptr(array,current_dim,asked_dim)
2357 : use, intrinsic :: iso_c_binding
2358 :
2359 : type(c_ptr), intent(inout) :: array
2360 : integer(c_size_t), intent(inout) :: current_dim
2361 : integer(c_size_t), intent(in ) :: asked_dim
2362 : ! *************************************************************************
2363 :
2364 : if ( current_dim < asked_dim ) then
2365 : if(current_dim /= 0) then
2366 : call dealloc_on_gpu(array)
2367 : end if
2368 : current_dim = asked_dim
2369 : call alloc_on_gpu(array, asked_dim)
2370 : end if
2371 :
2372 : #ifndef HAVE_GPU
2373 : ! Unused if GPU code disabled
2374 : ABI_UNUSED(array)
2375 : ABI_UNUSED((/current_dim,asked_dim/))
2376 : #endif
2377 :
2378 : end subroutine abi_gpu_work_resizeCptr
2379 : !!***
2380 :
2381 : subroutine abi_gpu_work_finalize()
2382 : use, intrinsic :: iso_c_binding
2383 :
2384 : #ifdef HAVE_GPU
2385 : !FIXME Assuming managed here ?
2386 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2387 :
2388 : #ifdef HAVE_YAKL
2389 : if ( associated(i_work_managed) ) then
2390 : ABI_FREE_MANAGED(i_work_managed)
2391 : end if
2392 : if ( associated(r_work_managed) ) then
2393 : ABI_FREE_MANAGED(r_work_managed)
2394 : end if
2395 : if ( associated(c_work_managed) ) then
2396 : ABI_FREE_MANAGED(c_work_managed)
2397 : end if
2398 : #endif
2399 :
2400 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2401 :
2402 : #ifdef HAVE_OPENMP_OFFLOAD
2403 : if ( allocated(i_work) ) then
2404 : !$OMP TARGET EXIT DATA MAP(release:i_work)
2405 : ABI_FREE(i_work)
2406 : end if
2407 : if ( allocated(r_work) ) then
2408 : !$OMP TARGET EXIT DATA MAP(release:r_work)
2409 : ABI_FREE(r_work)
2410 : end if
2411 : if ( allocated(c_work) ) then
2412 : !$OMP TARGET EXIT DATA MAP(release:c_work)
2413 : ABI_FREE(c_work)
2414 : end if
2415 : #endif
2416 :
2417 : if(gpu_work_len > 0) then
2418 : call dealloc_on_gpu(gpu_work)
2419 : end if
2420 :
2421 : end if
2422 :
2423 : i_work_len = 0
2424 : r_work_len = 0
2425 : c_work_len = 0
2426 : gpu_work_len = 0
2427 :
2428 : #endif
2429 :
2430 : end subroutine abi_gpu_work_finalize
2431 : !!***
2432 :
2433 : !!****f* m_abi_gpu_linalg/abi_gpu_xhegvd
2434 : !! NAME
2435 : !! abi_gpu_xhegvd
2436 : !!
2437 : !! FUNCTION
2438 : !! Compute a LAPACK SYGVD operation on GPU
2439 : !! compute eigen values/vectors of a real generalized
2440 : !! symmetric-definite eigenproblem
2441 : !!
2442 : !! See cusolver documentation
2443 : !! https://docs.nvidia.com/cuda/cusolver/index.html#cuSolverDN-lt-t-gt-hegvd
2444 : !!
2445 : !! See also LAPACK doc in reference implementation:
2446 : !! https://github.com/Reference-LAPACK/lapack/blob/master/SRC/dhegvd.f
2447 : !!
2448 : !! INPUTS
2449 : !! cplx = 1 if real 2 if complex
2450 : !! itype = integer, type of problem
2451 : !! jobz = character, 'n'(eigenvalues only) or 'v' (eigenvalues + eigenvectors)
2452 : !! uplo = character, 'u' or 'l'
2453 : !! A_nrows = matrix size
2454 : !! A = pointer to gpu memory location of matrix A
2455 : !! lda = leading dimension of matrix A
2456 : !! B = pointer to gpu memory location of matrix B
2457 : !! ldb = leading dimension of matrix B
2458 : !! W = pointer to gpu memory location of matrix W (output eigen values)
2459 : !! devInfo =
2460 : !!
2461 : !! SOURCE
2462 :
2463 0 : subroutine abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, A_nrows, &
2464 : A, lda, &
2465 : B, ldb, &
2466 : W, &
2467 : devInfo)
2468 :
2469 : ! Arguments ------------------------------------
2470 : integer, intent(in ) :: cplx
2471 : integer, intent(in ) :: itype
2472 : character(len=1),intent(in ) :: jobz
2473 : character(len=1),intent(in ) :: uplo
2474 : integer, intent(in ) :: A_nrows
2475 : type(c_ptr), intent(in ) :: A
2476 : integer, intent(in ) :: lda
2477 : type(c_ptr), intent(in ) :: B
2478 : integer, intent(in ) :: ldb
2479 : type(c_ptr), intent(in ) :: W
2480 : integer, intent(inout) :: devInfo
2481 :
2482 : ! Local variables ------------------------------
2483 : integer :: bufferSize
2484 : type(c_ptr) :: gpu_ptr
2485 : ! *************************************************************************
2486 :
2487 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2488 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2489 : end if
2490 :
2491 : #ifdef HAVE_GPU
2492 :
2493 : ! probe needed bufferSize
2494 : call gpu_xsygvd_buffersize(cplx, itype, jobz, uplo, &
2495 : A_nrows, &
2496 : A, lda, &
2497 : B, ldb, &
2498 : W, &
2499 : bufferSize)
2500 :
2501 : select case(cplx)
2502 :
2503 : case (1)
2504 : ! resize work array if needed and retrieve work pointer to use
2505 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
2506 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2507 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
2508 : gpu_ptr = c_loc(r_work_managed)
2509 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2510 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(1,c_size_t)*bufferSize*dp)
2511 : !gpu_ptr = gpu_work
2512 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
2513 : gpu_ptr = xomp_get_mapped_ptr(c_loc(r_work))
2514 : end if
2515 :
2516 : case (2)
2517 : ! resize work array if needed and retrieve work pointer to use
2518 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
2519 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2520 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
2521 : gpu_ptr = c_loc(c_work_managed)
2522 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2523 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(2,c_size_t)*bufferSize*dp)
2524 : !gpu_ptr = gpu_work
2525 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
2526 : gpu_ptr = xomp_get_mapped_ptr(c_loc(c_work))
2527 : end if
2528 :
2529 : end select
2530 :
2531 : ! and compute (finally)
2532 : call gpu_xsygvd(cplx, itype, jobz, uplo, &
2533 : A_nrows, &
2534 : A, lda, &
2535 : B, ldb, &
2536 : W, &
2537 : gpu_ptr, bufferSize, devInfo)
2538 :
2539 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2540 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
2541 : ! Therefore, we issue a stream sync here to avoid
2542 : !potential mistakes in calling context.
2543 : call gpu_linalg_stream_synchronize()
2544 : !FIXME Free memory to help GPU memory constraint (temporary?)
2545 : call abi_gpu_work_finalize()
2546 : end if
2547 :
2548 : #else
2549 : ! Unused if GPU code disabled
2550 : ABI_UNUSED((/cplx,itype,A_nrows,lda,ldb,devInfo,bufferSize/))
2551 : ABI_UNUSED((/jobz,uplo/))
2552 : ABI_UNUSED((/A,B,W,gpu_ptr/))
2553 : #endif
2554 :
2555 0 : end subroutine abi_gpu_xhegvd_cptr
2556 : !!***
2557 :
2558 0 : subroutine abi_gpu_xhegvd_d(cplx, itype, jobz, uplo, A_nrows, &
2559 : A, lda, &
2560 : B, ldb, &
2561 : W, &
2562 : devInfo)
2563 :
2564 : ! Arguments ------------------------------------
2565 : integer, intent(in ) :: cplx
2566 : integer, intent(in ) :: itype
2567 : character(len=1),intent(in ) :: jobz
2568 : character(len=1),intent(in ) :: uplo
2569 : integer, intent(in ) :: A_nrows
2570 : real(dp), intent(in ),target :: A(*)
2571 : integer, intent(in ) :: lda
2572 : real(dp), intent(in ),target :: B(*)
2573 : integer, intent(in ) :: ldb
2574 : real(dp), intent(inout),target :: W(*)
2575 : integer, intent(inout) :: devInfo
2576 : ! *************************************************************************
2577 :
2578 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2579 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2580 : A_nrows, &
2581 : c_loc(A), lda, &
2582 : c_loc(B), ldb, &
2583 : c_loc(W), &
2584 0 : devInfo)
2585 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2586 : #ifdef HAVE_OPENMP_OFFLOAD
2587 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,B,W)
2588 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2589 : A_nrows, &
2590 : c_loc(A), lda, &
2591 : c_loc(B), ldb, &
2592 : c_loc(W), &
2593 : devInfo)
2594 : !$OMP END TARGET DATA
2595 : #endif
2596 : else
2597 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2598 : end if
2599 :
2600 0 : end subroutine abi_gpu_xhegvd_d
2601 : !!***
2602 :
2603 0 : subroutine abi_gpu_xhegvd_z(cplx, itype, jobz, uplo, A_nrows, &
2604 : A, lda, &
2605 : B, ldb, &
2606 : W, &
2607 : devInfo)
2608 :
2609 : ! Arguments ------------------------------------
2610 : integer, intent(in ) :: cplx
2611 : integer, intent(in ) :: itype
2612 : character(len=1),intent(in ) :: jobz
2613 : character(len=1),intent(in ) :: uplo
2614 : integer, intent(in ) :: A_nrows
2615 : complex(dp), intent(in ),target :: A(*)
2616 : integer, intent(in ) :: lda
2617 : complex(dp), intent(in ),target :: B(*)
2618 : integer, intent(in ) :: ldb
2619 : real(dp), intent(inout),target :: W(*)
2620 : integer, intent(inout) :: devInfo
2621 : ! *************************************************************************
2622 :
2623 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2624 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2625 : A_nrows, &
2626 : c_loc(A), lda, &
2627 : c_loc(B), ldb, &
2628 : c_loc(W), &
2629 0 : devInfo)
2630 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2631 : #ifdef HAVE_OPENMP_OFFLOAD
2632 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,B,W)
2633 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2634 : A_nrows, &
2635 : c_loc(A), lda, &
2636 : c_loc(B), ldb, &
2637 : c_loc(W), &
2638 : devInfo)
2639 : !$OMP END TARGET DATA
2640 : #endif
2641 : else
2642 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2643 : end if
2644 :
2645 0 : end subroutine abi_gpu_xhegvd_z
2646 : !!***
2647 :
2648 0 : subroutine abi_gpu_xhegvd_2d(cplx, itype, jobz, uplo, A_nrows, &
2649 0 : A, lda, &
2650 0 : B, ldb, &
2651 0 : W, &
2652 : devInfo)
2653 :
2654 : ! Arguments ------------------------------------
2655 : integer, intent(in ) :: cplx
2656 : integer, intent(in ) :: itype
2657 : character(len=1),intent(in ) :: jobz
2658 : character(len=1),intent(in ) :: uplo
2659 : integer, intent(in ) :: A_nrows,lda,ldb
2660 : real(dp), intent(in ),target :: A(lda,*)
2661 : real(dp), intent(in ),target :: B(ldb,*)
2662 : real(dp), intent(inout),target :: W(A_nrows,*)
2663 : integer, intent(inout) :: devInfo
2664 : ! *************************************************************************
2665 :
2666 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2667 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2668 : A_nrows, &
2669 : c_loc(A), lda, &
2670 : c_loc(B), ldb, &
2671 : c_loc(W), &
2672 0 : devInfo)
2673 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2674 : #ifdef HAVE_OPENMP_OFFLOAD
2675 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,B,W)
2676 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2677 : A_nrows, &
2678 : c_loc(A), lda, &
2679 : c_loc(B), ldb, &
2680 : c_loc(W), &
2681 : devInfo)
2682 : !$OMP END TARGET DATA
2683 : #endif
2684 : else
2685 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2686 : end if
2687 :
2688 0 : end subroutine abi_gpu_xhegvd_2d
2689 : !!***
2690 :
2691 0 : subroutine abi_gpu_xhegvd_2z(cplx, itype, jobz, uplo, A_nrows, &
2692 0 : A, lda, &
2693 0 : B, ldb, &
2694 0 : W, &
2695 : devInfo)
2696 :
2697 : ! Arguments ------------------------------------
2698 : integer, intent(in ) :: cplx
2699 : integer, intent(in ) :: itype
2700 : character(len=1),intent(in ) :: jobz
2701 : character(len=1),intent(in ) :: uplo
2702 : integer, intent(in ) :: A_nrows,lda,ldb
2703 : complex(dp), intent(in ),target :: A(lda,*)
2704 : complex(dp), intent(in ),target :: B(ldb,*)
2705 : real(dp), intent(inout),target :: W(A_nrows,*)
2706 : integer, intent(inout) :: devInfo
2707 : ! *************************************************************************
2708 :
2709 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2710 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2711 : A_nrows, &
2712 : c_loc(A), lda, &
2713 : c_loc(B), ldb, &
2714 : c_loc(W), &
2715 0 : devInfo)
2716 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2717 : #ifdef HAVE_OPENMP_OFFLOAD
2718 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,B,W)
2719 : call abi_gpu_xhegvd_cptr(cplx, itype, jobz, uplo, &
2720 : A_nrows, &
2721 : c_loc(A), lda, &
2722 : c_loc(B), ldb, &
2723 : c_loc(W), &
2724 : devInfo)
2725 : !$OMP END TARGET DATA
2726 : #endif
2727 : else
2728 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2729 : end if
2730 :
2731 0 : end subroutine abi_gpu_xhegvd_2z
2732 : !!***
2733 :
2734 : !!****f* m_abi_gpu_linalg/abi_gpu_xheevd
2735 : !! NAME
2736 : !! abi_gpu_xheevd
2737 : !!
2738 : !! FUNCTION
2739 : !! Compute a LAPACK SYEVD operation on GPU
2740 : !! compute eigen values/vectors of a real
2741 : !! symmetric-definite eigenproblem
2742 : !!
2743 : !! See cusolver documentation
2744 : !! https://docs.nvidia.com/cuda/cusolver/index.html#cuSolverDN-lt-t-gt-heevd
2745 : !!
2746 : !! See also LAPACK doc in reference implementation:
2747 : !! https://github.com/Reference-LAPACK/lapack/blob/master/SRC/dheevd.f
2748 : !!
2749 : !! INPUTS
2750 : !! cplx = 1 if real 2 if complex
2751 : !! jobz = character, 'n'(eigenvalues only) or 'v' (eigenvalues + eigenvectors)
2752 : !! uplo = character, 'u' or 'l'
2753 : !! A_nrows = matrix size
2754 : !! A = pointer to gpu memory location of matrix A
2755 : !! lda = leading dimension of matrix A
2756 : !! W = pointer to gpu memory location of matrix W (output eigen values)
2757 : !! devInfo =
2758 : !!
2759 : !! SOURCE
2760 :
2761 0 : subroutine abi_gpu_xheevd_cptr(cplx, jobz, uplo, A_nrows, &
2762 : A, lda, &
2763 : W, &
2764 : devInfo)
2765 :
2766 : ! Arguments ------------------------------------
2767 : integer, intent(in ) :: cplx
2768 : character(len=1),intent(in ) :: jobz
2769 : character(len=1),intent(in ) :: uplo
2770 : integer, intent(in ) :: A_nrows
2771 : type(c_ptr), intent(in ) :: A
2772 : integer, intent(in ) :: lda
2773 : type(c_ptr), intent(in) :: W
2774 : integer, intent(inout) :: devInfo
2775 :
2776 : ! Local variables ------------------------------
2777 : integer :: bufferSize
2778 : type(c_ptr) :: gpu_ptr
2779 : ! *************************************************************************
2780 :
2781 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
2782 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
2783 : end if
2784 :
2785 : #ifdef HAVE_GPU
2786 :
2787 : ! probe needed bufferSize
2788 : call gpu_xsyevd_buffersize(cplx, jobz, uplo, &
2789 : A_nrows, &
2790 : A, lda, &
2791 : W, &
2792 : bufferSize)
2793 :
2794 : select case(cplx)
2795 :
2796 : case (1)
2797 : ! resize work array if needed and retrieve work pointer to use
2798 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
2799 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2800 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
2801 : gpu_ptr = c_loc(r_work_managed)
2802 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2803 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(1,c_size_t)*bufferSize*dp)
2804 : !gpu_ptr = gpu_work
2805 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
2806 : gpu_ptr = xomp_get_mapped_ptr(c_loc(r_work))
2807 : end if
2808 :
2809 : case (2)
2810 : ! resize work array if needed and retrieve work pointer to use
2811 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
2812 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2813 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
2814 : gpu_ptr = c_loc(c_work_managed)
2815 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2816 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(2,c_size_t)*bufferSize*dp)
2817 : !gpu_ptr = gpu_work
2818 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
2819 : gpu_ptr = xomp_get_mapped_ptr(c_loc(c_work))
2820 : end if
2821 :
2822 : end select
2823 :
2824 : ! and compute (finally)
2825 : call gpu_xsyevd(cplx, jobz, uplo, &
2826 : A_nrows, &
2827 : A, lda, &
2828 : W, &
2829 : gpu_ptr, bufferSize, devInfo)
2830 :
2831 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2832 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
2833 : ! Therefore, we issue a stream sync here to avoid
2834 : !potential mistakes in calling context.
2835 : call gpu_linalg_stream_synchronize()
2836 : end if
2837 :
2838 : #else
2839 : ! Unused if GPU code disabled
2840 : ABI_UNUSED((/cplx,A_nrows,lda,devInfo,bufferSize/))
2841 : ABI_UNUSED((/jobz,uplo/))
2842 : ABI_UNUSED((/A,W,gpu_ptr/))
2843 : #endif
2844 :
2845 0 : end subroutine abi_gpu_xheevd_cptr
2846 : !!***
2847 :
2848 0 : subroutine abi_gpu_xheevd_d(cplx, jobz, uplo, A_nrows, &
2849 : A, lda, &
2850 : W, &
2851 : devInfo)
2852 :
2853 : ! Arguments ------------------------------------
2854 : integer, intent(in ) :: cplx
2855 : character(len=1),intent(in ) :: jobz
2856 : character(len=1),intent(in ) :: uplo
2857 : integer, intent(in ) :: A_nrows
2858 : real(dp), intent(in ),target :: A(*)
2859 : integer, intent(in ) :: lda
2860 : real(dp), intent(inout),target :: W(*)
2861 : integer, intent(inout) :: devInfo
2862 : ! *************************************************************************
2863 :
2864 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2865 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2866 : A_nrows, &
2867 : c_loc(A), lda, &
2868 : c_loc(W), &
2869 0 : devInfo)
2870 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2871 : #ifdef HAVE_OPENMP_OFFLOAD
2872 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,W)
2873 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2874 : A_nrows, &
2875 : c_loc(A), lda, &
2876 : c_loc(W), &
2877 : devInfo)
2878 : !$OMP END TARGET DATA
2879 : #endif
2880 : else
2881 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2882 : end if
2883 :
2884 0 : end subroutine abi_gpu_xheevd_d
2885 : !!***
2886 :
2887 0 : subroutine abi_gpu_xheevd_z(cplx, jobz, uplo, A_nrows, &
2888 : A, lda, &
2889 : W, &
2890 : devInfo)
2891 :
2892 : ! Arguments ------------------------------------
2893 : integer, intent(in ) :: cplx
2894 : character(len=1),intent(in ) :: jobz
2895 : character(len=1),intent(in ) :: uplo
2896 : integer, intent(in ) :: A_nrows
2897 : complex(dp), intent(in ),target :: A(*)
2898 : integer, intent(in ) :: lda
2899 : real(dp), intent(inout),target :: W(*)
2900 : integer, intent(inout) :: devInfo
2901 : ! *************************************************************************
2902 :
2903 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2904 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2905 : A_nrows, &
2906 : c_loc(A), lda, &
2907 : c_loc(W), &
2908 0 : devInfo)
2909 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2910 : #ifdef HAVE_OPENMP_OFFLOAD
2911 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,W)
2912 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2913 : A_nrows, &
2914 : c_loc(A), lda, &
2915 : c_loc(W), &
2916 : devInfo)
2917 : !$OMP END TARGET DATA
2918 : #endif
2919 : else
2920 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2921 : end if
2922 :
2923 0 : end subroutine abi_gpu_xheevd_z
2924 : !!***
2925 :
2926 0 : subroutine abi_gpu_xheevd_2d(cplx, jobz, uplo, A_nrows, &
2927 0 : A, lda, &
2928 0 : W, &
2929 : devInfo)
2930 :
2931 : ! Arguments ------------------------------------
2932 : integer, intent(in ) :: cplx
2933 : character(len=1),intent(in ) :: jobz
2934 : character(len=1),intent(in ) :: uplo
2935 : integer, intent(in ) :: A_nrows,lda
2936 : real(dp), intent(in ),target :: A(lda,*)
2937 : real(dp), intent(inout),target :: W(A_nrows,*)
2938 : integer, intent(inout) :: devInfo
2939 : ! *************************************************************************
2940 :
2941 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2942 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2943 : A_nrows, &
2944 : c_loc(A), lda, &
2945 : c_loc(W), &
2946 0 : devInfo)
2947 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2948 : #ifdef HAVE_OPENMP_OFFLOAD
2949 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,W)
2950 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2951 : A_nrows, &
2952 : c_loc(A), lda, &
2953 : c_loc(W), &
2954 : devInfo)
2955 : !$OMP END TARGET DATA
2956 : #endif
2957 : else
2958 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2959 : end if
2960 :
2961 0 : end subroutine abi_gpu_xheevd_2d
2962 : !!***
2963 :
2964 0 : subroutine abi_gpu_xheevd_2z(cplx, jobz, uplo, A_nrows, &
2965 0 : A, lda, &
2966 0 : W, &
2967 : devInfo)
2968 :
2969 : ! Arguments ------------------------------------
2970 : integer, intent(in ) :: cplx
2971 : character(len=1),intent(in ) :: jobz
2972 : character(len=1),intent(in ) :: uplo
2973 : integer, intent(in ) :: A_nrows,lda
2974 : complex(dp), intent(in ),target :: A(lda,*)
2975 : real(dp), intent(inout),target :: W(A_nrows,*)
2976 : integer, intent(inout) :: devInfo
2977 : ! *************************************************************************
2978 :
2979 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
2980 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2981 : A_nrows, &
2982 : c_loc(A), lda, &
2983 : c_loc(W), &
2984 0 : devInfo)
2985 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
2986 : #ifdef HAVE_OPENMP_OFFLOAD
2987 : !$OMP TARGET DATA USE_DEVICE_ADDR(A,W)
2988 : call abi_gpu_xheevd_cptr(cplx, jobz, uplo, &
2989 : A_nrows, &
2990 : c_loc(A), lda, &
2991 : c_loc(W), &
2992 : devInfo)
2993 : !$OMP END TARGET DATA
2994 : #endif
2995 : else
2996 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
2997 : end if
2998 :
2999 0 : end subroutine abi_gpu_xheevd_2z
3000 : !!***
3001 :
3002 : !!****f* m_abi_gpu_linalg/abi_gpu_xpotrf
3003 : !! NAME
3004 : !! abi_gpu_xpotrf
3005 : !!
3006 : !! FUNCTION
3007 : !! Compute a LAPACK SYGVD operation on GPU
3008 : !! compute eigen values/vectors of a real
3009 : !! symmetric-definite eigenproblem
3010 : !!
3011 : !! See cusolver documentation
3012 : !! https://docs.nvidia.com/cuda/cusolver/index.html#cuSolverDN-lt-t-gt-potrf
3013 : !!
3014 : !! See also LAPACK doc in reference implementation:
3015 : !! https://github.com/Reference-LAPACK/lapack/blob/master/SRC/dpotrf.f
3016 : !!
3017 : !! INPUTS
3018 : !! cplx = 1 if real 2 if complex
3019 : !! uplo = character, 'u' or 'l'
3020 : !! A_nrows = matrix size
3021 : !! A = pointer to gpu memory location of matrix A
3022 : !! lda = leading dimension of matrix A
3023 : !! devInfo =
3024 : !!
3025 : !! SOURCE
3026 0 : subroutine abi_gpu_xpotrf_cptr(cplx, uplo, A_nrows, &
3027 : A, lda, &
3028 : devInfo)
3029 :
3030 : ! Arguments ------------------------------------
3031 : integer, intent(in ) :: cplx
3032 : character(len=1),intent(in ) :: uplo
3033 : integer, intent(in ) :: A_nrows
3034 : type(c_ptr), intent(in ) :: A
3035 : integer, intent(in ) :: lda
3036 : integer, intent(inout) :: devInfo
3037 :
3038 : ! Local variables ------------------------------
3039 : integer :: bufferSize
3040 : type(c_ptr) :: gpu_ptr
3041 : ! *************************************************************************
3042 :
3043 0 : if (abi_linalg_gpu_mode == ABI_GPU_DISABLED) then
3044 0 : ABI_BUG("You requested to run on CPU to a GPU wrapper :/")
3045 : end if
3046 :
3047 : #ifdef HAVE_GPU
3048 :
3049 : ! probe needed bufferSize
3050 : call gpu_xpotrf_buffersize(cplx, uplo, &
3051 : A_nrows, &
3052 : A, lda, &
3053 : bufferSize)
3054 :
3055 : select case(cplx)
3056 :
3057 : case (1)
3058 : ! resize work array if needed and retrieve work pointer to use
3059 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
3060 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3061 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
3062 : gpu_ptr = c_loc(r_work_managed)
3063 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3064 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(1,c_size_t)*bufferSize*dp)
3065 : !gpu_ptr = gpu_work
3066 : call abi_gpu_work_resize(r_work,r_work_managed,r_work_len,bufferSize)
3067 : gpu_ptr = xomp_get_mapped_ptr(c_loc(r_work))
3068 : end if
3069 :
3070 : case (2)
3071 : ! resize work array if needed and retrieve work pointer to use
3072 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY &
3073 : .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3074 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
3075 : gpu_ptr = c_loc(c_work_managed)
3076 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3077 : !call abi_gpu_work_resizeCptr(gpu_work,gpu_work_len,INT(1,c_size_t)*bufferSize*dp)
3078 : !gpu_ptr = gpu_work
3079 : call abi_gpu_work_resize(c_work,c_work_managed,c_work_len,bufferSize)
3080 : gpu_ptr = xomp_get_mapped_ptr(c_loc(c_work))
3081 : end if
3082 :
3083 : end select
3084 :
3085 : ! and compute (finally)
3086 : call gpu_xpotrf(cplx, uplo, &
3087 : A_nrows, &
3088 : A, lda, &
3089 : gpu_ptr, bufferSize, devInfo)
3090 :
3091 : if (abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3092 : ! CUDA/HIP linalg calls are run asynchronously and OpenMP is unaware of them.
3093 : ! Therefore, we issue a stream sync here to avoid
3094 : !potential mistakes in calling context.
3095 : call gpu_linalg_stream_synchronize()
3096 : end if
3097 :
3098 : #else
3099 : ! Unused if GPU code disabled
3100 : ABI_UNUSED((/cplx,A_nrows,lda,devInfo,bufferSize/))
3101 : ABI_UNUSED((/uplo/))
3102 : ABI_UNUSED((/A,gpu_ptr/))
3103 : #endif
3104 :
3105 0 : end subroutine abi_gpu_xpotrf_cptr
3106 : !!***
3107 :
3108 0 : subroutine abi_gpu_xpotrf_d(cplx, uplo, A_nrows, &
3109 : A, lda, &
3110 : devInfo)
3111 :
3112 : ! Arguments ------------------------------------
3113 : integer, intent(in ) :: cplx
3114 : character(len=1),intent(in ) :: uplo
3115 : integer, intent(in ) :: A_nrows
3116 : real(dp), intent(in ),target :: A(*)
3117 : integer, intent(in ) :: lda
3118 : integer, intent(inout) :: devInfo
3119 : ! *************************************************************************
3120 :
3121 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3122 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3123 : A_nrows, &
3124 : c_loc(A), lda, &
3125 0 : devInfo)
3126 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3127 : #ifdef HAVE_OPENMP_OFFLOAD
3128 : !$OMP TARGET DATA USE_DEVICE_ADDR(A)
3129 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3130 : A_nrows, &
3131 : c_loc(A), lda, &
3132 : devInfo)
3133 : !$OMP END TARGET DATA
3134 : #endif
3135 : else
3136 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
3137 : end if
3138 :
3139 0 : end subroutine abi_gpu_xpotrf_d
3140 : !!***
3141 :
3142 0 : subroutine abi_gpu_xpotrf_z(cplx, uplo, A_nrows, &
3143 : A, lda, &
3144 : devInfo)
3145 :
3146 : ! Arguments ------------------------------------
3147 : integer, intent(in ) :: cplx
3148 : character(len=1),intent(in ) :: uplo
3149 : integer, intent(in ) :: A_nrows
3150 : complex(dp), intent(in ),target :: A(*)
3151 : integer, intent(in ) :: lda
3152 : integer, intent(inout) :: devInfo
3153 : ! *************************************************************************
3154 :
3155 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3156 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3157 : A_nrows, &
3158 : c_loc(A), lda, &
3159 0 : devInfo)
3160 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3161 : #ifdef HAVE_OPENMP_OFFLOAD
3162 : !$OMP TARGET DATA USE_DEVICE_ADDR(A)
3163 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3164 : A_nrows, &
3165 : c_loc(A), lda, &
3166 : devInfo)
3167 : !$OMP END TARGET DATA
3168 : #endif
3169 : else
3170 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
3171 : end if
3172 :
3173 0 : end subroutine abi_gpu_xpotrf_z
3174 : !!***
3175 :
3176 0 : subroutine abi_gpu_xpotrf_2d(cplx, uplo, A_nrows, &
3177 0 : A, lda, &
3178 : devInfo)
3179 :
3180 : ! Arguments ------------------------------------
3181 : integer, intent(in ) :: cplx
3182 : character(len=1),intent(in ) :: uplo
3183 : integer, intent(in ) :: A_nrows,lda
3184 : real(dp), intent(in ),target :: A(lda,*)
3185 : integer, intent(inout) :: devInfo
3186 : ! *************************************************************************
3187 :
3188 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3189 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3190 : A_nrows, &
3191 : c_loc(A), lda, &
3192 0 : devInfo)
3193 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3194 : #ifdef HAVE_OPENMP_OFFLOAD
3195 : !$OMP TARGET DATA USE_DEVICE_ADDR(A)
3196 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3197 : A_nrows, &
3198 : c_loc(A), lda, &
3199 : devInfo)
3200 : !$OMP END TARGET DATA
3201 : #endif
3202 : else
3203 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
3204 : end if
3205 :
3206 0 : end subroutine abi_gpu_xpotrf_2d
3207 : !!***
3208 :
3209 0 : subroutine abi_gpu_xpotrf_2z(cplx, uplo, A_nrows, &
3210 0 : A, lda, &
3211 : devInfo)
3212 :
3213 : ! Arguments ------------------------------------
3214 : integer, intent(in ) :: cplx
3215 : character(len=1),intent(in ) :: uplo
3216 : integer, intent(in ) :: A_nrows,lda
3217 : complex(dp), intent(in ),target :: A(lda,*)
3218 : integer, intent(inout) :: devInfo
3219 : ! *************************************************************************
3220 :
3221 0 : if(abi_linalg_gpu_mode == ABI_GPU_LEGACY .or. abi_linalg_gpu_mode == ABI_GPU_KOKKOS) then
3222 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3223 : A_nrows, &
3224 : c_loc(A), lda, &
3225 0 : devInfo)
3226 0 : else if(abi_linalg_gpu_mode == ABI_GPU_OPENMP) then
3227 : #ifdef HAVE_OPENMP_OFFLOAD
3228 : !$OMP TARGET DATA USE_DEVICE_ADDR(A)
3229 : call abi_gpu_xpotrf_cptr(cplx, uplo, &
3230 : A_nrows, &
3231 : c_loc(A), lda, &
3232 : devInfo)
3233 : !$OMP END TARGET DATA
3234 : #endif
3235 : else
3236 0 : ABI_BUG(sjoin("Unhandled GPU mode:", itoa(abi_linalg_gpu_mode)))
3237 : end if
3238 :
3239 0 : end subroutine abi_gpu_xpotrf_2z
3240 : !!***
3241 :
3242 : !!****f* m_abi_linalg/gpu_xorthonormalize
3243 : !! NAME
3244 : !! gpu_xorthonormalize
3245 : !!
3246 : !! FUNCTION
3247 : !! This routine computes the overlap of two complex wavefunctions (for a given number of bands)
3248 : !! and orthonormalizes it using gpu:
3249 : !! - Computes the products of two rectangular matrices
3250 : !! containing the wavefunctions psi and S.psi (where S is the
3251 : !! overlap (with the PAW terms if necessary)).
3252 : !! - Does a Cholesky decomposition of this overlap
3253 : !! - rotates the initial matrix blockvectorx by the triangular matrix to
3254 : !! have an orthonormal set of wavefunctions
3255 : !!
3256 : !! INPUTS
3257 : !! blockvectorbx = matrix of dimension (blocksize,vectsize) as a GPU ptr
3258 : !! (e.g. block of overlap*wavefunction)
3259 : !! blocksize = dimension of matrices (e.g number of bands)
3260 : !! spaceComm = communicator used for MPI parallelization
3261 : !! vectsize = dimension of matrices (e.g number of G vector)
3262 : !!
3263 : !! OUTPUT
3264 : !! sqgram = Choleski decomposition of transpose(blockvector)*blockvectorx as a GPU ptr
3265 : !!
3266 : !! SIDE EFFECTS
3267 : !! blockvectorx = on input, matrix of dimension (vectsize,blocksize) as a GPU ptr
3268 : !! (e.g block of wavefunction)
3269 : !! blockvectorx = on output, orthonormalized wavefunction. as a GPU ptr
3270 : !!
3271 : !! SOURCE
3272 :
3273 0 : subroutine gpu_xorthonormalize(blockvectorx_gpu,blockvectorbx_gpu,blocksize,spaceComm,&
3274 : & sqgram_gpu,vectsize,&
3275 : & x_cplx,timopt,tim_xortho) ! optional arguments
3276 :
3277 : !Arguments ------------------------------------
3278 : !scalars
3279 : integer,intent(in) :: blocksize,spaceComm,vectsize,x_cplx
3280 : integer, intent(in), optional :: timopt,tim_xortho
3281 : !arrays
3282 : type(c_ptr),intent(inout) :: blockvectorbx_gpu, blockvectorx_gpu, sqgram_gpu
3283 : !Local variables-------------------------------
3284 : #if defined HAVE_GPU
3285 : integer :: ierr,info
3286 : real(dp), dimension(:,:),allocatable, target :: d_sqgram
3287 : complex(dp),dimension(:,:),allocatable, target :: z_sqgram
3288 : character :: tr
3289 : real(dp) :: tsec(2)
3290 : integer(c_size_t) :: size
3291 : #else
3292 : type(c_ptr) :: cptr_a
3293 : #endif
3294 : character(len=500) :: message
3295 : ! *********************************************************************
3296 :
3297 : #if defined HAVE_GPU
3298 : if (present(tim_xortho).and.present(timopt)) then
3299 : if(abs(timopt)==3) then
3300 : call timab(tim_xortho,1,tsec)
3301 : end if
3302 : end if
3303 :
3304 : if ( x_cplx == 1 ) then
3305 : tr='t'
3306 : ABI_MALLOC(d_sqgram,(blocksize,blocksize))
3307 : else
3308 : tr='c'
3309 : ABI_MALLOC(z_sqgram,(blocksize,blocksize))
3310 : end if
3311 :
3312 : call gpu_xgemm(x_cplx,tr,'n',blocksize,blocksize,vectsize, &
3313 : & cone,blockvectorx_gpu,vectsize,blockvectorbx_gpu,vectsize,czero,sqgram_gpu,blocksize)
3314 : call gpu_device_synchronize()
3315 : size=x_cplx*dp*blocksize*blocksize
3316 :
3317 : if ( x_cplx == 1 ) then
3318 : call copy_from_gpu(d_sqgram, sqgram_gpu, INT(x_cplx, c_size_t)*dp*blocksize*blocksize)
3319 : call xmpi_sum(d_sqgram,spaceComm,ierr)
3320 : call abi_xpotrf('u',blocksize,d_sqgram,blocksize,info)
3321 : call copy_on_gpu(d_sqgram, sqgram_gpu, INT(x_cplx, c_size_t)*dp*blocksize*blocksize)
3322 : else
3323 : call copy_from_gpu(z_sqgram, sqgram_gpu, size)
3324 : call xmpi_sum(z_sqgram,spaceComm,ierr)
3325 : call abi_xpotrf('u',blocksize,z_sqgram,blocksize,info)
3326 : call copy_on_gpu(z_sqgram, sqgram_gpu, size)
3327 : end if
3328 :
3329 : if (info /= 0 ) then
3330 : write(message,'(a,i3)') ' xpotrf, info=',info
3331 : ABI_WARNING(message)
3332 : end if
3333 :
3334 : call gpu_xtrsm(x_cplx,'r','u','n','n',vectsize,blocksize,cone,sqgram_gpu,blocksize,&
3335 : & blockvectorx_gpu,vectsize)
3336 :
3337 : if(x_cplx==1) then
3338 : ABI_FREE(d_sqgram)
3339 : else
3340 : ABI_FREE(z_sqgram)
3341 : end if
3342 : if (present(tim_xortho).and.present(timopt)) then
3343 : if(abs(timopt)==3) then
3344 : call timab(tim_xortho,2,tsec)
3345 : end if
3346 : end if
3347 : return
3348 :
3349 : #else
3350 0 : message=' This routine is not allowed when running on GPU is disabled !'
3351 0 : ABI_BUG(message)
3352 : if (.false.) then
3353 : write(std_out,*) blocksize,vectsize,spaceComm,x_cplx
3354 : if (present(timopt)) write(std_out,*) timopt
3355 : if (present(tim_xortho)) write(std_out,*) tim_xortho
3356 : cptr_a=blockvectorbx_gpu;cptr_a=blockvectorx_gpu;cptr_a=sqgram_gpu
3357 : end if
3358 : #endif
3359 :
3360 0 : end subroutine gpu_xorthonormalize
3361 : !!***
|
