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
|
/* This program demonstrates a bug fixed in librsb-1.2-rc7. */
#include <rsb.h>
#include <stdio.h>
int main(const int argc, char * const argv[])
{
rsb_err_t errval = RSB_ERR_NO_ERROR;
int i;
const int N = 2;
double complex x[N], y[N];
const double complex alpha = 1.0;
double complex VA[]={1,1,2};
int IA[]={0,0,1},JA[]={0,1,1},nnz=3;
struct rsb_mtx_t * mtxAp;
// matrix:
// (1,0) (1,0)
// (1,0) (1,0)
// represented as upper hermitian:
// (1,0) (1,0)
// (0,0) (1,0)
// multiplicand vector:
// (0,1)
// (0,0)
// result shall be:
// (0,1)
// (0,1)
errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS);
if(errval != RSB_ERR_NO_ERROR)
goto err;
mtxAp = rsb_mtx_alloc_from_coo_begin(2, RSB_NUMERICAL_TYPE_DOUBLE_COMPLEX, N, N, RSB_FLAG_DEFAULT_STORAGE_FLAGS|RSB_FLAG_UPPER_HERMITIAN , &errval);
errval = rsb_mtx_set_vals(mtxAp, VA, IA, JA, nnz, RSB_FLAG_NOFLAGS);
if(errval != RSB_ERR_NO_ERROR)
goto err;
errval = rsb_mtx_alloc_from_coo_end(&mtxAp);
if(errval != RSB_ERR_NO_ERROR)
goto err;
for (i = 0; i < N; i++) x[i] = y[i] = 0.0;
x[0] = I;
errval = rsb_spmv(RSB_TRANSPOSITION_N, &alpha, mtxAp, x, 1, NULL, y, 1);
if(errval != RSB_ERR_NO_ERROR)
goto err;
for (i = 0; i < N; i++)
printf("(%g,%g)\t", creal(y[i]), cimag(y[i]));
printf("\n");
if(cimag(y[1])!=1)
printf("Result seems incorrect -- are you using < librsb-1.2-rc7?\n");
else
printf("Result seems correct -- are you using >= librsb-1.2-rc7?\n");
rsb_mtx_free(mtxAp);
errval = rsb_lib_exit(RSB_NULL_EXIT_OPTIONS);
if(errval != RSB_ERR_NO_ERROR)
goto err;
return 0;
err:
rsb_perror(NULL,errval);
return -1;
}
|
