Actual source code: petsclegacycupmblas.h
1: #ifndef PETSCLEGACYCUPMBLAS_H
2: #define PETSCLEGACYCUPMBLAS_H
4: #include <petscdevice_cupm.h>
6: #if defined(PETSC_HAVE_CUDA)
7: /* complex single */
8: #if defined(PETSC_USE_COMPLEX)
9: #if defined(PETSC_USE_REAL_SINGLE)
10: #define cublasXaxpy(a, b, c, d, e, f, g) cublasCaxpy((a), (b), (cuComplex *)(c), (cuComplex *)(d), (e), (cuComplex *)(f), (g))
11: #define cublasXscal(a, b, c, d, e) cublasCscal((a), (b), (cuComplex *)(c), (cuComplex *)(d), (e))
12: #define cublasXdotu(a, b, c, d, e, f, g) cublasCdotu((a), (b), (cuComplex *)(c), (d), (cuComplex *)(e), (f), (cuComplex *)(g))
13: #define cublasXdot(a, b, c, d, e, f, g) cublasCdotc((a), (b), (cuComplex *)(c), (d), (cuComplex *)(e), (f), (cuComplex *)(g))
14: #define cublasXswap(a, b, c, d, e, f) cublasCswap((a), (b), (cuComplex *)(c), (d), (cuComplex *)(e), (f))
15: #define cublasXnrm2(a, b, c, d, e) cublasScnrm2((a), (b), (cuComplex *)(c), (d), (e))
16: #define cublasIXamax(a, b, c, d, e) cublasIcamax((a), (b), (cuComplex *)(c), (d), (e))
17: #define cublasXasum(a, b, c, d, e) cublasScasum((a), (b), (cuComplex *)(c), (d), (e))
18: #define cublasXgemv(a, b, c, d, e, f, g, h, i, j, k, l) cublasCgemv((a), (b), (c), (d), (cuComplex *)(e), (cuComplex *)(f), (g), (cuComplex *)(h), (i), (cuComplex *)(j), (cuComplex *)(k), (l))
19: #define cublasXgemm(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cublasCgemm((a), (b), (c), (d), (e), (f), (cuComplex *)(g), (cuComplex *)(h), (i), (cuComplex *)(j), (k), (cuComplex *)(l), (cuComplex *)(m), (n))
20: #define cublasXgeam(a, b, c, d, e, f, g, h, i, j, k, l, m) cublasCgeam((a), (b), (c), (d), (e), (cuComplex *)(f), (cuComplex *)(g), (h), (cuComplex *)(i), (cuComplex *)(j), (k), (cuComplex *)(l), (m))
21: #define cublasXgemvStridedBatched(a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p) \
22: cublasCgemvStridedBatched((a), (b), (c), (d), (const cuComplex *)(e), (const cuComplex *)(f), (g), (h), (const cuComplex *)(i), (j), (k), (const cuComplex *)(l), (cuComplex *)(m), (n), (o), (p))
23: #else /* complex double */
24: #define cublasXaxpy(a, b, c, d, e, f, g) cublasZaxpy((a), (b), (cuDoubleComplex *)(c), (cuDoubleComplex *)(d), (e), (cuDoubleComplex *)(f), (g))
25: #define cublasXscal(a, b, c, d, e) cublasZscal((a), (b), (cuDoubleComplex *)(c), (cuDoubleComplex *)(d), (e))
26: #define cublasXdotu(a, b, c, d, e, f, g) cublasZdotu((a), (b), (cuDoubleComplex *)(c), (d), (cuDoubleComplex *)(e), (f), (cuDoubleComplex *)(g))
27: #define cublasXdot(a, b, c, d, e, f, g) cublasZdotc((a), (b), (cuDoubleComplex *)(c), (d), (cuDoubleComplex *)(e), (f), (cuDoubleComplex *)(g))
28: #define cublasXswap(a, b, c, d, e, f) cublasZswap((a), (b), (cuDoubleComplex *)(c), (d), (cuDoubleComplex *)(e), (f))
29: #define cublasXnrm2(a, b, c, d, e) cublasDznrm2((a), (b), (cuDoubleComplex *)(c), (d), (e))
30: #define cublasIXamax(a, b, c, d, e) cublasIzamax((a), (b), (cuDoubleComplex *)(c), (d), (e))
31: #define cublasXasum(a, b, c, d, e) cublasDzasum((a), (b), (cuDoubleComplex *)(c), (d), (e))
32: #define cublasXgemv(a, b, c, d, e, f, g, h, i, j, k, l) cublasZgemv((a), (b), (c), (d), (cuDoubleComplex *)(e), (cuDoubleComplex *)(f), (g), (cuDoubleComplex *)(h), (i), (cuDoubleComplex *)(j), (cuDoubleComplex *)(k), (l))
33: #define cublasXgemm(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cublasZgemm((a), (b), (c), (d), (e), (f), (cuDoubleComplex *)(g), (cuDoubleComplex *)(h), (i), (cuDoubleComplex *)(j), (k), (cuDoubleComplex *)(l), (cuDoubleComplex *)(m), (n))
34: #define cublasXgeam(a, b, c, d, e, f, g, h, i, j, k, l, m) cublasZgeam((a), (b), (c), (d), (e), (cuDoubleComplex *)(f), (cuDoubleComplex *)(g), (h), (cuDoubleComplex *)(i), (cuDoubleComplex *)(j), (k), (cuDoubleComplex *)(l), (m))
35: #define cublasXgemvStridedBatched(a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p) \
36: cublasZgemvStridedBatched((a), (b), (c), (d), (const cuDoubleComplex *)(e), (const cuDoubleComplex *)(f), (g), (h), (const cuDoubleComplex *)(i), (j), (k), (const cuDoubleComplex *)(l), (cuDoubleComplex *)(m), (n), (o), (p))
38: #endif
39: #else /* real single */
40: #if defined(PETSC_USE_REAL_SINGLE)
41: #define cublasXaxpy cublasSaxpy
42: #define cublasXscal cublasSscal
43: #define cublasXdotu cublasSdot
44: #define cublasXdot cublasSdot
45: #define cublasXswap cublasSswap
46: #define cublasXnrm2 cublasSnrm2
47: #define cublasIXamax cublasIsamax
48: #define cublasXasum cublasSasum
49: #define cublasXgemv cublasSgemv
50: #define cublasXgemm cublasSgemm
51: #define cublasXgeam cublasSgeam
52: #define cublasXgemvStridedBatched cublasSgemvStridedBatched
53: #else /* real double */
54: #define cublasXaxpy cublasDaxpy
55: #define cublasXscal cublasDscal
56: #define cublasXdotu cublasDdot
57: #define cublasXdot cublasDdot
58: #define cublasXswap cublasDswap
59: #define cublasXnrm2 cublasDnrm2
60: #define cublasIXamax cublasIdamax
61: #define cublasXasum cublasDasum
62: #define cublasXgemv cublasDgemv
63: #define cublasXgemm cublasDgemm
64: #define cublasXgeam cublasDgeam
65: #define cublasXgemvStridedBatched cublasDgemvStridedBatched
66: #endif
67: #endif
69: #if defined(PETSC_USE_COMPLEX)
70: #if defined(PETSC_USE_REAL_SINGLE)
71: #define cusolverDnXpotrf(a, b, c, d, e, f, g, h) cusolverDnCpotrf((a), (b), (c), (cuComplex *)(d), (e), (cuComplex *)(f), (g), (h))
72: #define cusolverDnXpotrf_bufferSize(a, b, c, d, e, f) cusolverDnCpotrf_bufferSize((a), (b), (c), (cuComplex *)(d), (e), (f))
73: #define cusolverDnXpotrs(a, b, c, d, e, f, g, h, i) cusolverDnCpotrs((a), (b), (c), (d), (cuComplex *)(e), (f), (cuComplex *)(g), (h), (i))
74: #define cusolverDnXpotri(a, b, c, d, e, f, g, h) cusolverDnCpotri((a), (b), (c), (cuComplex *)(d), (e), (cuComplex *)(f), (g), (h))
75: #define cusolverDnXpotri_bufferSize(a, b, c, d, e, f) cusolverDnCpotri_bufferSize((a), (b), (c), (cuComplex *)(d), (e), (f))
76: #define cusolverDnXsytrf(a, b, c, d, e, f, g, h, i) cusolverDnCsytrf((a), (b), (c), (cuComplex *)(d), (e), (f), (cuComplex *)(g), (h), (i))
77: #define cusolverDnXsytrf_bufferSize(a, b, c, d, e) cusolverDnCsytrf_bufferSize((a), (b), (cuComplex *)(c), (d), (e))
78: #define cusolverDnXgetrf(a, b, c, d, e, f, g, h) cusolverDnCgetrf((a), (b), (c), (cuComplex *)(d), (e), (cuComplex *)(f), (g), (h))
79: #define cusolverDnXgetrf_bufferSize(a, b, c, d, e, f) cusolverDnCgetrf_bufferSize((a), (b), (c), (cuComplex *)(d), (e), (f))
80: #define cusolverDnXgetrs(a, b, c, d, e, f, g, h, i, j) cusolverDnCgetrs((a), (b), (c), (d), (cuComplex *)(e), (f), (g), (cuComplex *)(h), (i), (j))
81: #define cusolverDnXgeqrf_bufferSize(a, b, c, d, e, f) cusolverDnCgeqrf_bufferSize((a), (b), (c), (cuComplex *)(d), (e), (f))
82: #define cusolverDnXgeqrf(a, b, c, d, e, f, g, h, i) cusolverDnCgeqrf((a), (b), (c), (cuComplex *)(d), (e), (cuComplex *)(f), (cuComplex *)(g), (h), (i))
83: #define cusolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) cusolverDnCunmqr_bufferSize((a), (b), (c), (d), (e), (f), (cuComplex *)(g), (h), (cuComplex *)(i), (cuComplex *)(j), (k), (l))
84: #define cusolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cusolverDnCunmqr((a), (b), (c), (d), (e), (f), (cuComplex *)(g), (h), (cuComplex *)(i), (cuComplex *)(j), (k), (cuComplex *)(l), (m), (n))
85: #define cublasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) cublasCtrsm((a), (b), (c), (d), (e), (f), (g), (cuComplex *)(h), (cuComplex *)(i), (j), (cuComplex *)(k), (l))
86: #else /* complex double */
87: #define cusolverDnXpotrf(a, b, c, d, e, f, g, h) cusolverDnZpotrf((a), (b), (c), (cuDoubleComplex *)(d), (e), (cuDoubleComplex *)(f), (g), (h))
88: #define cusolverDnXpotrf_bufferSize(a, b, c, d, e, f) cusolverDnZpotrf_bufferSize((a), (b), (c), (cuDoubleComplex *)(d), (e), (f))
89: #define cusolverDnXpotrs(a, b, c, d, e, f, g, h, i) cusolverDnZpotrs((a), (b), (c), (d), (cuDoubleComplex *)(e), (f), (cuDoubleComplex *)(g), (h), (i))
90: #define cusolverDnXpotri(a, b, c, d, e, f, g, h) cusolverDnZpotri((a), (b), (c), (cuDoubleComplex *)(d), (e), (cuDoubleComplex *)(f), (g), (h))
91: #define cusolverDnXpotri_bufferSize(a, b, c, d, e, f) cusolverDnZpotri_bufferSize((a), (b), (c), (cuDoubleComplex *)(d), (e), (f))
92: #define cusolverDnXsytrf(a, b, c, d, e, f, g, h, i) cusolverDnZsytrf((a), (b), (c), (cuDoubleComplex *)(d), (e), (f), (cuDoubleComplex *)(g), (h), (i))
93: #define cusolverDnXsytrf_bufferSize(a, b, c, d, e) cusolverDnZsytrf_bufferSize((a), (b), (cuDoubleComplex *)(c), (d), (e))
94: #define cusolverDnXgetrf(a, b, c, d, e, f, g, h) cusolverDnZgetrf((a), (b), (c), (cuDoubleComplex *)(d), (e), (cuDoubleComplex *)(f), (g), (h))
95: #define cusolverDnXgetrf_bufferSize(a, b, c, d, e, f) cusolverDnZgetrf_bufferSize((a), (b), (c), (cuDoubleComplex *)(d), (e), (f))
96: #define cusolverDnXgetrs(a, b, c, d, e, f, g, h, i, j) cusolverDnZgetrs((a), (b), (c), (d), (cuDoubleComplex *)(e), (f), (g), (cuDoubleComplex *)(h), (i), (j))
97: #define cusolverDnXgeqrf_bufferSize(a, b, c, d, e, f) cusolverDnZgeqrf_bufferSize((a), (b), (c), (cuDoubleComplex *)(d), (e), (f))
98: #define cusolverDnXgeqrf(a, b, c, d, e, f, g, h, i) cusolverDnZgeqrf((a), (b), (c), (cuDoubleComplex *)(d), (e), (cuDoubleComplex *)(f), (cuDoubleComplex *)(g), (h), (i))
99: #define cusolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) cusolverDnZunmqr_bufferSize((a), (b), (c), (d), (e), (f), (cuDoubleComplex *)(g), (h), (cuDoubleComplex *)(i), (cuDoubleComplex *)(j), (k), (l))
100: #define cusolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cusolverDnZunmqr((a), (b), (c), (d), (e), (f), (cuDoubleComplex *)(g), (h), (cuDoubleComplex *)(i), (cuDoubleComplex *)(j), (k), (cuDoubleComplex *)(l), (m), (n))
101: #define cublasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) cublasZtrsm((a), (b), (c), (d), (e), (f), (g), (cuDoubleComplex *)(h), (cuDoubleComplex *)(i), (j), (cuDoubleComplex *)(k), (l))
102: #endif
103: #else /* real single */
104: #if defined(PETSC_USE_REAL_SINGLE)
105: #define cusolverDnXpotrf(a, b, c, d, e, f, g, h) cusolverDnSpotrf((a), (b), (c), (d), (e), (f), (g), (h))
106: #define cusolverDnXpotrf_bufferSize(a, b, c, d, e, f) cusolverDnSpotrf_bufferSize((a), (b), (c), (d), (e), (f))
107: #define cusolverDnXpotrs(a, b, c, d, e, f, g, h, i) cusolverDnSpotrs((a), (b), (c), (d), (e), (f), (g), (h), (i))
108: #define cusolverDnXpotri(a, b, c, d, e, f, g, h) cusolverDnSpotri((a), (b), (c), (d), (e), (f), (g), (h))
109: #define cusolverDnXpotri_bufferSize(a, b, c, d, e, f) cusolverDnSpotri_bufferSize((a), (b), (c), (d), (e), (f))
110: #define cusolverDnXsytrf(a, b, c, d, e, f, g, h, i) cusolverDnSsytrf((a), (b), (c), (d), (e), (f), (g), (h), (i))
111: #define cusolverDnXsytrf_bufferSize(a, b, c, d, e) cusolverDnSsytrf_bufferSize((a), (b), (c), (d), (e))
112: #define cusolverDnXgetrf(a, b, c, d, e, f, g, h) cusolverDnSgetrf((a), (b), (c), (d), (e), (f), (g), (h))
113: #define cusolverDnXgetrf_bufferSize(a, b, c, d, e, f) cusolverDnSgetrf_bufferSize((a), (b), (c), (d), (e), (f))
114: #define cusolverDnXgetrs(a, b, c, d, e, f, g, h, i, j) cusolverDnSgetrs((a), (b), (c), (d), (e), (f), (g), (h), (i), (j))
115: #define cusolverDnXgeqrf_bufferSize(a, b, c, d, e, f) cusolverDnSgeqrf_bufferSize((a), (b), (c), (float *)(d), (e), (f))
116: #define cusolverDnXgeqrf(a, b, c, d, e, f, g, h, i) cusolverDnSgeqrf((a), (b), (c), (float *)(d), (e), (float *)(f), (float *)(g), (h), (i))
117: #define cusolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) cusolverDnSormqr_bufferSize((a), (b), (c), (d), (e), (f), (float *)(g), (h), (float *)(i), (float *)(j), (k), (l))
118: #define cusolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cusolverDnSormqr((a), (b), (c), (d), (e), (f), (float *)(g), (h), (float *)(i), (float *)(j), (k), (float *)(l), (m), (n))
119: #define cublasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) cublasStrsm((a), (b), (c), (d), (e), (f), (g), (float *)(h), (float *)(i), (j), (float *)(k), (l))
120: #else /* real double */
121: #define cusolverDnXpotrf(a, b, c, d, e, f, g, h) cusolverDnDpotrf((a), (b), (c), (d), (e), (f), (g), (h))
122: #define cusolverDnXpotrf_bufferSize(a, b, c, d, e, f) cusolverDnDpotrf_bufferSize((a), (b), (c), (d), (e), (f))
123: #define cusolverDnXpotrs(a, b, c, d, e, f, g, h, i) cusolverDnDpotrs((a), (b), (c), (d), (e), (f), (g), (h), (i))
124: #define cusolverDnXpotri(a, b, c, d, e, f, g, h) cusolverDnDpotri((a), (b), (c), (d), (e), (f), (g), (h))
125: #define cusolverDnXpotri_bufferSize(a, b, c, d, e, f) cusolverDnDpotri_bufferSize((a), (b), (c), (d), (e), (f))
126: #define cusolverDnXsytrf(a, b, c, d, e, f, g, h, i) cusolverDnDsytrf((a), (b), (c), (d), (e), (f), (g), (h), (i))
127: #define cusolverDnXsytrf_bufferSize(a, b, c, d, e) cusolverDnDsytrf_bufferSize((a), (b), (c), (d), (e))
128: #define cusolverDnXgetrf(a, b, c, d, e, f, g, h) cusolverDnDgetrf((a), (b), (c), (d), (e), (f), (g), (h))
129: #define cusolverDnXgetrf_bufferSize(a, b, c, d, e, f) cusolverDnDgetrf_bufferSize((a), (b), (c), (d), (e), (f))
130: #define cusolverDnXgetrs(a, b, c, d, e, f, g, h, i, j) cusolverDnDgetrs((a), (b), (c), (d), (e), (f), (g), (h), (i), (j))
131: #define cusolverDnXgeqrf_bufferSize(a, b, c, d, e, f) cusolverDnDgeqrf_bufferSize((a), (b), (c), (double *)(d), (e), (f))
132: #define cusolverDnXgeqrf(a, b, c, d, e, f, g, h, i) cusolverDnDgeqrf((a), (b), (c), (double *)(d), (e), (double *)(f), (double *)(g), (h), (i))
133: #define cusolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) cusolverDnDormqr_bufferSize((a), (b), (c), (d), (e), (f), (double *)(g), (h), (double *)(i), (double *)(j), (k), (l))
134: #define cusolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) cusolverDnDormqr((a), (b), (c), (d), (e), (f), (double *)(g), (h), (double *)(i), (double *)(j), (k), (double *)(l), (m), (n))
135: #define cublasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) cublasDtrsm((a), (b), (c), (d), (e), (f), (g), (double *)(h), (double *)(i), (j), (double *)(k), (l))
136: #endif
137: #endif
138: #endif // PETSC_HAVE_CUDA
140: #if defined(PETSC_HAVE_HIP)
141: #if defined(PETSC_USE_COMPLEX)
142: #if defined(PETSC_USE_REAL_SINGLE)
143: #define hipsolverDnXpotrf(a, b, c, d, e, f, g, h) hipsolverCpotrf((a), (b), (c), (hipComplex *)(d), (e), (hipComplex *)(f), (g), (h))
144: #define hipsolverDnXpotrf_bufferSize(a, b, c, d, e, f) hipsolverCpotrf_bufferSize((a), (b), (c), (hipComplex *)(d), (e), (f))
145: #define hipsolverDnXpotrs(a, b, c, d, e, f, g, h, i, j, k) hipsolverCpotrs((a), (b), (c), (d), (hipComplex *)(e), (f), (hipComplex *)(g), (h), (hipComplex *)(i), (j), (k))
146: #define hipsolverDnXpotri(a, b, c, d, e, f, g, h) hipsolverCpotri((a), (b), (c), (hipComplex *)(d), (e), (hipComplex *)(f), (g), (h))
147: #define hipsolverDnXpotri_bufferSize(a, b, c, d, e, f) hipsolverCpotri_bufferSize((a), (b), (c), (hipComplex *)(d), (e), (f))
148: #define hipsolverDnXsytrf(a, b, c, d, e, f, g, h, i) hipsolverCsytrf((a), (b), (c), (hipComplex *)(d), (e), (f), (hipComplex *)(g), (h), (i))
149: #define hipsolverDnXsytrf_bufferSize(a, b, c, d, e) hipsolverCsytrf_bufferSize((a), (b), (hipComplex *)(c), (d), (e))
150: #define hipsolverDnXgetrf(a, b, c, d, e, f, g, h, i) hipsolverCgetrf((a), (b), (c), (hipComplex *)(d), (e), (hipComplex *)(f), (g), (h), (i))
151: #define hipsolverDnXgetrf_bufferSize(a, b, c, d, e, f) hipsolverCgetrf_bufferSize((a), (b), (c), (hipComplex *)(d), (e), (f))
152: #define hipsolverDnXgetrs(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverCgetrs((a), (b), (c), (d), (hipComplex *)(e), (f), (g), (hipComplex *)(h), (i), (hipComplex *)(j), (k), (l))
153: #define hipsolverDnXgeqrf_bufferSize(a, b, c, d, e, f) hipsolverCgeqrf_bufferSize((a), (b), (c), (hipComplex *)(d), (e), (f))
154: #define hipsolverDnXgeqrf(a, b, c, d, e, f, g, h, i) hipsolverCgeqrf((a), (b), (c), (hipComplex *)(d), (e), (hipComplex *)(f), (hipComplex *)(g), (h), (i))
155: #define hipsolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverCunmqr_bufferSize((a), (b), (c), (d), (e), (f), (hipComplex *)(g), (h), (hipComplex *)(i), (hipComplex *)(j), (k), (l))
156: #define hipsolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) hipsolverCunmqr((a), (b), (c), (d), (e), (f), (hipComplex *)(g), (h), (hipComplex *)(i), (hipComplex *)(j), (k), (hipComplex *)(l), (m), (n))
157: #define hipblasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) hipblasCtrsm((a), (b), (c), (d), (e), (f), (g), (hipComplex *)(h), (hipComplex *)(i), (j), (hipComplex *)(k), (l))
158: #else /* complex double */
159: #define hipsolverDnXpotrf(a, b, c, d, e, f, g, h) hipsolverZpotrf((a), (b), (c), (hipDoubleComplex *)(d), (e), (hipDoubleComplex *)(f), (g), (h))
160: #define hipsolverDnXpotrf_bufferSize(a, b, c, d, e, f) hipsolverZpotrf_bufferSize((a), (b), (c), (hipDoubleComplex *)(d), (e), (f))
161: #define hipsolverDnXpotrs(a, b, c, d, e, f, g, h, i, j, k) hipsolverZpotrs((a), (b), (c), (d), (hipDoubleComplex *)(e), (f), (hipDoubleComplex *)(g), (h), (hipDoubleComplex *)(i), (j), (k))
162: #define hipsolverDnXpotri(a, b, c, d, e, f, g, h) hipsolverZpotri((a), (b), (c), (hipDoubleComplex *)(d), (e), (hipDoubleComplex *)(f), (g), (h))
163: #define hipsolverDnXpotri_bufferSize(a, b, c, d, e, f) hipsolverZpotri_bufferSize((a), (b), (c), (hipDoubleComplex *)(d), (e), (f))
164: #define hipsolverDnXsytrf(a, b, c, d, e, f, g, h, i) hipsolverZsytrf((a), (b), (c), (hipDoubleComplex *)(d), (e), (f), (hipDoubleComplex *)(g), (h), (i))
165: #define hipsolverDnXsytrf_bufferSize(a, b, c, d, e) hipsolverZsytrf_bufferSize((a), (b), (hipDoubleComplex *)(c), (d), (e))
166: #define hipsolverDnXgetrf(a, b, c, d, e, f, g, h, i) hipsolverZgetrf((a), (b), (c), (hipDoubleComplex *)(d), (e), (hipDoubleComplex *)(f), (g), (h), (i))
167: #define hipsolverDnXgetrf_bufferSize(a, b, c, d, e, f) hipsolverZgetrf_bufferSize((a), (b), (c), (hipDoubleComplex *)(d), (e), (f))
168: #define hipsolverDnXgetrs(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverZgetrs((a), (b), (c), (d), (hipDoubleComplex *)(e), (f), (g), (hipDoubleComplex *)(h), (i), (hipDoubleComplex *)(j), (k), (l))
169: #define hipsolverDnXgeqrf_bufferSize(a, b, c, d, e, f) hipsolverZgeqrf_bufferSize((a), (b), (c), (hipDoubleComplex *)(d), (e), (f))
170: #define hipsolverDnXgeqrf(a, b, c, d, e, f, g, h, i) hipsolverZgeqrf((a), (b), (c), (hipDoubleComplex *)(d), (e), (hipDoubleComplex *)(f), (hipDoubleComplex *)(g), (h), (i))
171: #define hipsolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverZunmqr_bufferSize((a), (b), (c), (d), (e), (f), (hipDoubleComplex *)(g), (h), (hipDoubleComplex *)(i), (hipDoubleComplex *)(j), (k), (l))
172: #define hipsolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) hipsolverZunmqr((a), (b), (c), (d), (e), (f), (hipDoubleComplex *)(g), (h), (hipDoubleComplex *)(i), (hipDoubleComplex *)(j), (k), (hipDoubleComplex *)(l), (m), (n))
173: #define hipblasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) hipblasZtrsm((a), (b), (c), (d), (e), (f), (g), (hipDoubleComplex *)(h), (hipDoubleComplex *)(i), (j), (hipDoubleComplex *)(k), (l))
174: #endif
175: #else /* real single */
176: #if defined(PETSC_USE_REAL_SINGLE)
177: #define hipsolverDnXpotrf(a, b, c, d, e, f, g, h) hipsolverSpotrf((a), (b), (c), (float *)(d), (e), (float *)(f), (g), (h))
178: #define hipsolverDnXpotrf_bufferSize(a, b, c, d, e, f) hipsolverSpotrf_bufferSize((a), (b), (c), (float *)(d), (e), (f))
179: #define hipsolverDnXpotrs(a, b, c, d, e, f, g, h, i, j, k) hipsolverSpotrs((a), (b), (c), (d), (float *)(e), (f), (float *)(g), (h), (float *)(i), (j), (k))
180: #define hipsolverDnXpotri(a, b, c, d, e, f, g, h) hipsolverSpotri((a), (b), (c), (float *)(d), (e), (float *)(f), (g), (h))
181: #define hipsolverDnXpotri_bufferSize(a, b, c, d, e, f) hipsolverSpotri_bufferSize((a), (b), (c), (float *)(d), (e), (f))
182: #define hipsolverDnXsytrf(a, b, c, d, e, f, g, h, i) hipsolverSsytrf((a), (b), (c), (float *)(d), (e), (f), (float *)(g), (h), (i))
183: #define hipsolverDnXsytrf_bufferSize(a, b, c, d, e) hipsolverSsytrf_bufferSize((a), (b), (float *)(c), (d), (e))
184: #define hipsolverDnXgetrf(a, b, c, d, e, f, g, h, i) hipsolverSgetrf((a), (b), (c), (float *)(d), (e), (float *)(f), (g), (h), (i))
185: #define hipsolverDnXgetrf_bufferSize(a, b, c, d, e, f) hipsolverSgetrf_bufferSize((a), (b), (c), (float *)(d), (e), (f))
186: #define hipsolverDnXgetrs(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverSgetrs((a), (b), (c), (d), (float *)(e), (f), (g), (float *)(h), (i), (float *)(j), (k), (l))
187: #define hipsolverDnXgeqrf_bufferSize(a, b, c, d, e, f) hipsolverSgeqrf_bufferSize((a), (b), (c), (float *)(d), (e), (f))
188: #define hipsolverDnXgeqrf(a, b, c, d, e, f, g, h, i) hipsolverSgeqrf((a), (b), (c), (float *)(d), (e), (float *)(f), (float *)(g), (h), (i))
189: #define hipsolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverSormqr_bufferSize((a), (b), (c), (d), (e), (f), (float *)(g), (h), (float *)(i), (float *)(j), (k), (l))
190: #define hipsolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) hipsolverSormqr((a), (b), (c), (d), (e), (f), (float *)(g), (h), (float *)(i), (float *)(j), (k), (float *)(l), (m), (n))
191: #define hipblasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) hipblasStrsm((a), (b), (c), (d), (e), (f), (g), (float *)(h), (float *)(i), (j), (float *)(k), (l))
192: #else /* real double */
193: #define hipsolverDnXpotrf(a, b, c, d, e, f, g, h) hipsolverDpotrf((a), (b), (c), (double *)(d), (e), (double *)(f), (g), (h))
194: #define hipsolverDnXpotrf_bufferSize(a, b, c, d, e, f) hipsolverDpotrf_bufferSize((a), (b), (c), (double *)(d), (e), (f))
195: #define hipsolverDnXpotrs(a, b, c, d, e, f, g, h, i, j, k) hipsolverDpotrs((a), (b), (c), (d), (double *)(e), (f), (double *)(g), (h), (double *)(i), (j), (k))
196: #define hipsolverDnXpotri(a, b, c, d, e, f, g, h) hipsolverDpotri((a), (b), (c), (double *)(d), (e), (double *)(f), (g), (h))
197: #define hipsolverDnXpotri_bufferSize(a, b, c, d, e, f) hipsolverDpotri_bufferSize((a), (b), (c), (double *)(d), (e), (f))
198: #define hipsolverDnXsytrf(a, b, c, d, e, f, g, h, i) hipsolverDsytrf((a), (b), (c), (double *)(d), (e), (f), (double *)(g), (h), (i))
199: #define hipsolverDnXsytrf_bufferSize(a, b, c, d, e) hipsolverDsytrf_bufferSize((a), (b), (double *)(c), (d), (e))
200: #define hipsolverDnXgetrf(a, b, c, d, e, f, g, h, i) hipsolverDgetrf((a), (b), (c), (double *)(d), (e), (double *)(f), (g), (h), (i))
201: #define hipsolverDnXgetrf_bufferSize(a, b, c, d, e, f) hipsolverDgetrf_bufferSize((a), (b), (c), (double *)(d), (e), (f))
202: #define hipsolverDnXgetrs(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverDgetrs((a), (b), (c), (d), (double *)(e), (f), (g), (double *)(h), (i), (double *)(j), (k), (l))
203: #define hipsolverDnXgeqrf_bufferSize(a, b, c, d, e, f) hipsolverDgeqrf_bufferSize((a), (b), (c), (double *)(d), (e), (f))
204: #define hipsolverDnXgeqrf(a, b, c, d, e, f, g, h, i) hipsolverDgeqrf((a), (b), (c), (double *)(d), (e), (double *)(f), (double *)(g), (h), (i))
205: #define hipsolverDnXormqr_bufferSize(a, b, c, d, e, f, g, h, i, j, k, l) hipsolverDormqr_bufferSize((a), (b), (c), (d), (e), (f), (double *)(g), (h), (double *)(i), (double *)(j), (k), (l))
206: #define hipsolverDnXormqr(a, b, c, d, e, f, g, h, i, j, k, l, m, n) hipsolverDormqr((a), (b), (c), (d), (e), (f), (double *)(g), (h), (double *)(i), (double *)(j), (k), (double *)(l), (m), (n))
207: #define hipblasXtrsm(a, b, c, d, e, f, g, h, i, j, k, l) hipblasDtrsm((a), (b), (c), (d), (e), (f), (g), (double *)(h), (double *)(i), (j), (double *)(k), (l))
208: #endif
209: #endif
211: /* complex single */
212: #if defined(PETSC_USE_COMPLEX)
213: #if defined(PETSC_USE_REAL_SINGLE)
214: #define hipblasXaxpy(a, b, c, d, e, f, g) hipblasCaxpy((a), (b), (hipblasComplex *)(c), (hipblasComplex *)(d), (e), (hipblasComplex *)(f), (g))
215: #define hipblasXscal(a, b, c, d, e) hipblasCscal((a), (b), (hipblasComplex *)(c), (hipblasComplex *)(d), (e))
216: #define hipblasXdotu(a, b, c, d, e, f, g) hipblasCdotu((a), (b), (hipblasComplex *)(c), (d), (hipblasComplex *)(e), (f), (hipblasComplex *)(g))
217: #define hipblasXdot(a, b, c, d, e, f, g) hipblasCdotc((a), (b), (hipblasComplex *)(c), (d), (hipblasComplex *)(e), (f), (hipblasComplex *)(g))
218: #define hipblasXswap(a, b, c, d, e, f) hipblasCswap((a), (b), (hipblasComplex *)(c), (d), (hipblasComplex *)(e), (f))
219: #define hipblasXnrm2(a, b, c, d, e) hipblasScnrm2((a), (b), (hipblasComplex *)(c), (d), (e))
220: #define hipblasIXamax(a, b, c, d, e) hipblasIcamax((a), (b), (hipblasComplex *)(c), (d), (e))
221: #define hipblasXasum(a, b, c, d, e) hipblasScasum((a), (b), (hipblasComplex *)(c), (d), (e))
222: #define hipblasXgemv(a, b, c, d, e, f, g, h, i, j, k, l) hipblasCgemv((a), (b), (c), (d), (hipblasComplex *)(e), (hipblasComplex *)(f), (g), (hipblasComplex *)(h), (i), (hipblasComplex *)(j), (hipblasComplex *)(k), (l))
223: #define hipblasXgemm(a, b, c, d, e, f, g, h, i, j, k, l, m, n) hipblasCgemm((a), (b), (c), (d), (e), (f), (hipblasComplex *)(g), (hipblasComplex *)(h), (i), (hipblasComplex *)(j), (k), (hipblasComplex *)(l), (hipblasComplex *)(m), (n))
224: #define hipblasXgeam(a, b, c, d, e, f, g, h, i, j, k, l, m) hipblasCgeam((a), (b), (c), (d), (e), (hipblasComplex *)(f), (hipblasComplex *)(g), (h), (hipblasComplex *)(i), (hipblasComplex *)(j), (k), (hipblasComplex *)(l), (m))
225: #else /* complex double */
226: #define hipblasXaxpy(a, b, c, d, e, f, g) hipblasZaxpy((a), (b), (hipblasDoubleComplex *)(c), (hipblasDoubleComplex *)(d), (e), (hipblasDoubleComplex *)(f), (g))
227: #define hipblasXscal(a, b, c, d, e) hipblasZscal((a), (b), (hipblasDoubleComplex *)(c), (hipblasDoubleComplex *)(d), (e))
228: #define hipblasXdotu(a, b, c, d, e, f, g) hipblasZdotu((a), (b), (hipblasDoubleComplex *)(c), (d), (hipblasDoubleComplex *)(e), (f), (hipblasDoubleComplex *)(g))
229: #define hipblasXdot(a, b, c, d, e, f, g) hipblasZdotc((a), (b), (hipblasDoubleComplex *)(c), (d), (hipblasDoubleComplex *)(e), (f), (hipblasDoubleComplex *)(g))
230: #define hipblasXswap(a, b, c, d, e, f) hipblasZswap((a), (b), (hipblasDoubleComplex *)(c), (d), (hipblasDoubleComplex *)(e), (f))
231: #define hipblasXnrm2(a, b, c, d, e) hipblasDznrm2((a), (b), (hipblasDoubleComplex *)(c), (d), (e))
232: #define hipblasIXamax(a, b, c, d, e) hipblasIzamax((a), (b), (hipblasDoubleComplex *)(c), (d), (e))
233: #define hipblasXasum(a, b, c, d, e) hipblasDzasum((a), (b), (hipblasDoubleComplex *)(c), (d), (e))
234: #define hipblasXgemv(a, b, c, d, e, f, g, h, i, j, k, l) \
235: hipblasZgemv((a), (b), (c), (d), (hipblasDoubleComplex *)(e), (hipblasDoubleComplex *)(f), (g), (hipblasDoubleComplex *)(h), (i), (hipblasDoubleComplex *)(j), (hipblasDoubleComplex *)(k), (l))
236: #define hipblasXgemm(a, b, c, d, e, f, g, h, i, j, k, l, m, n) \
237: hipblasZgemm((a), (b), (c), (d), (e), (f), (hipblasDoubleComplex *)(g), (hipblasDoubleComplex *)(h), (i), (hipblasDoubleComplex *)(j), (k), (hipblasDoubleComplex *)(l), (hipblasDoubleComplex *)(m), (n))
238: #define hipblasXgeam(a, b, c, d, e, f, g, h, i, j, k, l, m) \
239: hipblasZgeam((a), (b), (c), (d), (e), (hipblasDoubleComplex *)(f), (hipblasDoubleComplex *)(g), (h), (hipblasDoubleComplex *)(i), (hipblasDoubleComplex *)(j), (k), (hipblasDoubleComplex *)(l), (m))
240: #endif
241: #else /* real single */
242: #if defined(PETSC_USE_REAL_SINGLE)
243: #define hipblasXaxpy hipblasSaxpy
244: #define hipblasXscal hipblasSscal
245: #define hipblasXdotu hipblasSdot
246: #define hipblasXdot hipblasSdot
247: #define hipblasXswap hipblasSswap
248: #define hipblasXnrm2 hipblasSnrm2
249: #define hipblasIXamax hipblasIsamax
250: #define hipblasXasum hipblasSasum
251: #define hipblasXgemv hipblasSgemv
252: #define hipblasXgemm hipblasSgemm
253: #define hipblasXgeam hipblasSgeam
254: #else /* real double */
255: #define hipblasXaxpy hipblasDaxpy
256: #define hipblasXscal hipblasDscal
257: #define hipblasXdotu hipblasDdot
258: #define hipblasXdot hipblasDdot
259: #define hipblasXswap hipblasDswap
260: #define hipblasXnrm2 hipblasDnrm2
261: #define hipblasIXamax hipblasIdamax
262: #define hipblasXasum hipblasDasum
263: #define hipblasXgemv hipblasDgemv
264: #define hipblasXgemm hipblasDgemm
265: #define hipblasXgeam hipblasDgeam
266: #endif
267: #endif
268: #endif // PETSC_HAVE_HIP
270: #endif // PETSCLEGACYCUPMBLAS_H