EVSL  1.1.0
EigenValues Slicing Library
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros
evsl.c
Go to the documentation of this file.
1 #include <stdlib.h>
2 #include "def.h"
3 #include "struct.h"
4 #include "internal_proto.h"
5 
16 
22 
27 int EVSLStart() {
28  evsldata.n = 0;
29  evsldata.ifGenEv = 0;
30  evsldata.Amv = NULL;
31  evsldata.Bmv = NULL;
32  evsldata.Bsol = NULL;
33  evsldata.LTsol = NULL;
34  evsldata.ds = NULL;
35 
36  StatsReset();
37 
38  /* do a dummy timer call to improve accuracy of timing */
39  evsl_timer();
40 
41  return 0;
42 }
43 
48 int EVSLFinish() {
49  if (evsldata.Amv) {
50  free(evsldata.Amv);
51  }
52  if (evsldata.Bmv) {
53  free(evsldata.Bmv);
54  }
55  if (evsldata.Bsol) {
56  free(evsldata.Bsol);
57  }
58  if (evsldata.LTsol) {
59  free(evsldata.LTsol);
60  }
61  return 0;
62 }
63 
68 int SetAMatrix(csrMat *A) {
69  evsldata.n = A->ncols;
70  if (!evsldata.Amv) {
71  Calloc(evsldata.Amv, 1, EVSLMatvec);
72  }
73  evsldata.Amv->func = matvec_csr;
74  evsldata.Amv->data = (void *) A;
75 
76  return 0;
77 }
78 
83 int SetBMatrix(csrMat *B) {
84  evsldata.n = B->ncols;
85  if (!evsldata.Bmv) {
86  Calloc(evsldata.Bmv, 1, EVSLMatvec);
87  }
88  evsldata.Bmv->func = matvec_csr;
89  evsldata.Bmv->data = (void *) B;
90 
91  return 0;
92 }
93 
100 int SetAMatvec(int n, MVFunc func, void *data) {
101  evsldata.n = n;
102  if (!evsldata.Amv) {
103  Calloc(evsldata.Amv, 1, EVSLMatvec);
104  }
105  evsldata.Amv->func = func;
106  evsldata.Amv->data = data;
107 
108  return 0;
109 }
110 
117 int SetBMatvec(int n, MVFunc func, void *data) {
118  evsldata.n = n;
119  if (!evsldata.Bmv) {
120  Calloc(evsldata.Bmv, 1, EVSLMatvec);
121  }
122  evsldata.Bmv->func = func;
123  evsldata.Bmv->data = data;
124 
125  return 0;
126 }
127 
128 
132 int SetBSol(SolFuncR func, void *data) {
133  if (!evsldata.Bsol) {
134  Calloc(evsldata.Bsol, 1, EVSLBSol);
135  }
136 
137  evsldata.Bsol->func = func;
138  evsldata.Bsol->data = data;
139 
140  return 0;
141 }
142 
143 
148 int SetStdEig() {
149  evsldata.ifGenEv = 0;
150 
151  return 0;
152 }
153 
158 int SetGenEig() {
159  evsldata.ifGenEv = 1;
160 
161  return 0;
162 }
163 
168 int SetASigmaBSol(ratparams *rat, SolFuncC *func, SolFuncC allf, void **data) {
169  int i,num;
170  num = rat->num;
171  Calloc(rat->ASIGBsol, num, EVSLASIGMABSol);
172 
173  for (i=0; i<num; i++) {
174  rat->ASIGBsol[i].func = func ? func[i] : allf;
175  rat->ASIGBsol[i].data = data[i];
176  }
177 
178  return 0;
179 }
180 
181 
185 int SetLTSol(SolFuncR func, void *data) {
186  if (!evsldata.LTsol) {
187  Calloc(evsldata.LTsol, 1, EVSLLTSol);
188  }
189  evsldata.LTsol->func = func;
190  evsldata.LTsol->data = data;
191  return 0;
192 }
193 
197 void SetDiagScal(double *ds) {
198  evsldata.ds = ds;
199 }
200 
EVSLMatvec * Amv
Definition: struct.h:167
int num
Definition: struct.h:114
void(* MVFunc)(double *x, double *y, void *data)
matvec function prototype
Definition: struct.h:96
defs in EVSL
int SetGenEig()
Set the problem to generalized eigenvalue problem.
Definition: evsl.c:158
int ifGenEv
Definition: struct.h:166
void(* SolFuncC)(int n, double *br, double *bz, double *xr, double *xz, void *data)
linear solver function prototype: [complex version] which is used for solving system with A-SIGMA B n...
Definition: struct.h:86
int ncols
Definition: struct.h:32
int SetAMatrix(csrMat *A)
Set the matrix A.
Definition: evsl.c:68
EVSLBSol * Bsol
Definition: struct.h:169
int SetBMatvec(int n, MVFunc func, void *data)
Set the user-input matvec routine and the associated data for B. Save them in evsldata.
Definition: evsl.c:117
int SetBMatrix(csrMat *B)
Set the B matrix.
Definition: evsl.c:83
void * data
Definition: struct.h:104
void * data
Definition: struct.h:138
user-provided function for solving L^{T} x = b
Definition: struct.h:154
user-provided function and data for solving B x = b
Definition: struct.h:145
int EVSLStart()
Initialize evslData.
Definition: evsl.c:27
int n
Definition: struct.h:165
void * data
Definition: struct.h:156
#define Calloc(base, nmem, type)
Definition: def.h:32
This file contains function prototypes and constant definitions internally used in EVSL...
EVSLMatvec * Bmv
Definition: struct.h:168
EVSLASIGMABSol * ASIGBsol
Definition: struct.h:127
evslStat evslstat
global statistics of EVSL
Definition: evsl.c:21
void * data
Definition: struct.h:147
SolFuncR func
Definition: struct.h:155
structs used in evsl
user-provided function and data for solving (A - SIGMA*B) x = b
Definition: struct.h:102
int SetLTSol(SolFuncR func, void *data)
Set the solve routine for LT.
Definition: evsl.c:185
int SetStdEig()
Set the problem to standard eigenvalue problem.
Definition: evsl.c:148
void(* SolFuncR)(double *b, double *x, void *data)
function prototype for applying the solve B x = b
Definition: struct.h:91
sparse matrix format: the compressed sparse row (CSR) format, 0-based
Definition: struct.h:31
int EVSLFinish()
Finish EVSL.
Definition: evsl.c:48
int SetBSol(SolFuncR func, void *data)
Set the solve routine and the associated data for B.
Definition: evsl.c:132
EVSLLTSol * LTsol
Definition: struct.h:170
MVFunc func
Definition: struct.h:137
void SetDiagScal(double *ds)
Set diagonal scaling matrix D.
Definition: evsl.c:197
double * ds
Definition: struct.h:171
wrapper of all global variables in EVSL
Definition: struct.h:164
evslData evsldata
global variable of EVSL
Definition: evsl.c:15
void StatsReset()
Definition: stats.c:79
void matvec_csr(double *x, double *y, void *data)
matvec for a CSR matrix, y = A*x. void *data points to csrMat, compatible form with EVSLMatvec (see s...
Definition: spmat.c:247
int SetASigmaBSol(ratparams *rat, SolFuncC *func, SolFuncC allf, void **data)
Set the solve routine and the associated data for A-SIGMA*B if func == NULL, set all functions to be ...
Definition: evsl.c:168
timing and memory statistics of EVSL
Definition: struct.h:196
SolFuncR func
Definition: struct.h:146
int SetAMatvec(int n, MVFunc func, void *data)
Set the user-input matvec routine and the associated data for A. Save them in evsldata.
Definition: evsl.c:100
user-provided Mat-Vec function and data for y = A * x or y = B * x
Definition: struct.h:136
double evsl_timer()
evsl timer for mac
Definition: mactime.c:14
SolFuncC func
Definition: struct.h:103
parameters for rational filter
Definition: struct.h:112