1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
|
#include BLAS_HEADER
#include <algorithm>
#include <interfaces/blas/C/xerbla.h>
#include <ulmblas/level1/copy.h>
#include <ulmblas/level1extensions/gecopy.h>
#include <ulmblas/level2/sylmv.h>
//#define SCATTER
#ifdef SCATTER
#define SCATTER_INCROWA 2
#define SCATTER_INCCOLA 3
#define SCATTER_INCX 4
#define SCATTER_INCY 5
#endif
extern "C" {
void
ULMBLAS(dsymv)(enum UpLo upLo,
int n,
double alpha,
const double *A,
int ldA,
const double *x,
int incX,
double beta,
double *y,
int incY)
{
//
// Test the input parameters
//
int info = 0;
if (upLo!=Upper && upLo!=Lower) {
info = 1;
} else if (n<0) {
info = 2;
} else if (ldA<std::max(1,n)) {
info = 5;
} else if (incX==0) {
info = 7;
} else if (incY==0) {
info = 10;
}
if (info!=0) {
ULMBLAS(xerbla)("DSYMV ", &info);
}
#ifndef SCATTER
//
// Start the operations.
//
if (incX<0) {
x -= incX*(n-1);
}
if (incY<0) {
y -= incY*(n-1);
}
if (upLo==Lower) {
ulmBLAS::sylmv(n, alpha, A, 1, ldA, x, incX, beta, y, incY);
} else {
ulmBLAS::sylmv(n, alpha, A, ldA, 1, x, incX, beta, y, incY);
}
#else
//
// Scatter operands
//
double *A_ = new double[ldA*n*SCATTER_INCROWA*SCATTER_INCCOLA];
double *x_ = new double[n*incX*SCATTER_INCX];
double *y_ = new double[n*incY*SCATTER_INCY];
ulmBLAS::gecopy(n, n,
A, 1, ldA,
A_, SCATTER_INCROWA, ldA*SCATTER_INCCOLA);
ulmBLAS::copy(n, x, incX, x_, incX*SCATTER_INCX);
ulmBLAS::copy(n, y, incY, y_, incY*SCATTER_INCY);
//
// Start the operations.
//
if (upLo==Lower) {
ulmBLAS::sylmv(n, alpha,
A_, SCATTER_INCROWA, ldA*SCATTER_INCCOLA,
x_, incX*SCATTER_INCX,
beta,
y_, incY*SCATTER_INCY);
} else {
ulmBLAS::sylmv(n, alpha,
A_, ldA*SCATTER_INCCOLA, SCATTER_INCROWA,
x_, incX*SCATTER_INCX,
beta,
y_, incY*SCATTER_INCY);
}
ulmBLAS::copy(n, y_, incY*SCATTER_INCY, y, incY);
//
// Gather result
//
delete [] A_;
delete [] x_;
delete [] y_;
#endif
}
} // extern "C"
|