10 #ifndef PLASMA_DESCRIPTOR_H
11 #define PLASMA_DESCRIPTOR_H
13 #include "plasma_types.h"
14 #include "plasma_error.h"
79 static inline size_t plasma_element_size(
int type)
82 case PlasmaByte:
return 1;
83 case PlasmaInteger:
return sizeof(int);
84 case PlasmaRealFloat:
return sizeof(float);
85 case PlasmaRealDouble:
return sizeof(double);
86 case PlasmaComplexFloat:
return 2*
sizeof(float);
87 case PlasmaComplexDouble:
return 2*
sizeof(double);
93 static inline void *plasma_tile_addr_general(
plasma_desc_t A,
int m,
int n)
95 int mm = m + A.i/A.mb;
96 int nn = n + A.j/A.nb;
97 size_t eltsize = plasma_element_size(A.precision);
105 offset = A.mb*A.nb*(mm + (size_t)lm1 * nn);
107 offset = A.A12 + ((size_t)A.mb * (A.gn%A.nb) * mm);
110 offset = A.A21 + ((size_t)A.nb * (A.gm%A.mb) * nn);
114 return (
void*)((
char*)A.matrix + (offset*eltsize));
118 static inline void *plasma_tile_addr_triangle(
plasma_desc_t A,
int m,
int n)
120 int mm = m + A.i/A.mb;
121 int nn = n + A.j/A.nb;
122 size_t eltsize = plasma_element_size(A.precision);
130 if (A.type == PlasmaUpper) {
131 offset = A.mb*A.nb*(mm + (nn * (nn + 1))/2);
134 offset = A.mb*A.nb*((mm - nn) + (nn * (2*lm1 - (nn-1)))/2);
138 offset = A.A12 + ((size_t)A.mb * (A.gn%A.nb) * mm);
143 offset = A.A21 + ((size_t)A.nb * (A.gm%A.mb) * nn);
150 return (
void*)((
char*)A.matrix + (offset*eltsize));
154 static inline void *plasma_tile_addr_general_band(
plasma_desc_t A,
int m,
int n)
156 return plasma_tile_addr_general(A, (A.kut-1)+m-n, n);
160 static inline void *plasma_tile_addr(
plasma_desc_t A,
int m,
int n)
162 if (A.type == PlasmaGeneral) {
163 return plasma_tile_addr_general(A, m, n);
165 else if (A.type == PlasmaGeneralBand) {
166 return plasma_tile_addr_general_band(A, m, n);
168 else if (A.type == PlasmaUpper || A.type == PlasmaLower) {
169 return plasma_tile_addr_triangle(A, m, n);
172 plasma_fatal_error(
"invalid matrix type");
184 if (A.type == PlasmaGeneralBand) {
188 if (A.i/A.mb+k < A.gm/A.mb)
202 if (A.j/A.nb+k < A.gn/A.nb)
219 if ((A.i+A.m)%A.mb == 0)
222 return (A.i+A.m)%A.mb;
236 if ((A.j+A.n)%A.nb == 0)
239 return (A.j+A.n)%A.nb;
243 static inline int plasma_tile_mmain_band(
plasma_desc_t A,
int m,
int n)
245 return plasma_tile_mmain(A, (A.kut-1)+m-n);
249 int plasma_desc_general_create(plasma_enum_t dtyp,
int mb,
int nb,
250 int lm,
int ln,
int i,
int j,
int m,
int n,
253 int plasma_desc_general_band_create(plasma_enum_t dtyp, plasma_enum_t uplo,
254 int mb,
int nb,
int lm,
int ln,
255 int i,
int j,
int m,
int n,
int kl,
int ku,
258 int plasma_desc_triangular_create(plasma_enum_t dtyp, plasma_enum_t uplo,
int mb,
int nb,
259 int lm,
int ln,
int i,
int j,
int m,
int n,
264 int plasma_desc_general_init(plasma_enum_t precision,
void *matrix,
265 int mb,
int nb,
int lm,
int ln,
int i,
int j,
268 int plasma_desc_general_band_init(plasma_enum_t precision, plasma_enum_t uplo,
269 void *matrix,
int mb,
int nb,
int lm,
int ln,
270 int i,
int j,
int m,
int n,
int kl,
int ku,
273 int plasma_desc_triangular_init(plasma_enum_t precision, plasma_enum_t uplo,
void *matrix,
274 int mb,
int nb,
int lm,
int ln,
int i,
int j,
283 int plasma_descT_create(
plasma_desc_t A,
int ib, plasma_enum_t householder_mode,
290 #endif // PLASMA_DESCRIPTOR_H