EVSL  1.1.0
EigenValues Slicing Library
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros
evsl.h
Go to the documentation of this file.
1 
6 #ifndef EVSL_H
7 #define EVSL_H
8 
9 #include <complex.h>
10 #include <stdio.h>
11 #include "struct.h"
12 
13 /*- - - - - - - - - cheblanNr.c */
14 int ChebLanNr(double *intv, int maxit, double tol, double *vinit,
15  polparams *pol, int *nevOut, double **lamo, double **Wo,
16  double **reso, FILE *fstats);
17 
18 /*- - - - - - - - - cheblanTr.c */
19 int ChebLanTr(int lanm, int nev, double *intv, int maxit, double tol,
20  double *vinit, polparams *pol, int *nev2, double **vals,
21  double **W, double **resW, FILE *fstats);
22 
23 /*- - - - - - - - - chebpoly.c */
24 //
25 void set_pol_def(polparams *pol);
26 //
27 int find_pol(double *intv, polparams *pol);
28 //
29 void free_pol(polparams *pol);
30 
31 /*- - - - - - - - - chebsi.c */
32 int ChebSI(int nev, double *intv, int maxit, double tol, double *vinit,
33  polparams *pol, int *nevo, double **lamo, double **Yo, double **reso,
34  FILE *fstats);
35 
36 /*- - - - - - - - - - dos_utils.c */
37 void SetupPolRec(int n, int max_deg, double tol, double lmin, double lmax, BSolDataPol *data);
38 //
39 void SetupPolSqrt(int n, int max_deg, double tol, double lmin, double lmax, BSolDataPol *data);
40 //
41 void FreeBSolPolData(BSolDataPol *data);
42 //
43 void BSolPol(double *b, double *x, void *data);
44 //
45 void extrDiagCsr(csrMat *B, double *d);
46 //
47 void diagScalCsr(csrMat *A, double *d);
48 
49 /*- - - - - - - - - lanbounds.c */
50 int LanBounds(int msteps, double *v, double *lmin, double *lmax);
51 
52 /*- - - - - - - - - - landos.c */
53 // Computes the density of states (DOS, or spectral density)
54 int LanDos(const int nvec, int msteps, const int npts, double *xdos,
55  double *ydos, double *neig, const double *const intv);
56 
57 //*- - - - - - - - - landosG.c - Generalized lanDOS */
58 int LanDosG(const int nvec, int msteps, const int npts, double *xdos, double *ydos,
59  double *neig, const double *const intv);
60 
61 /*- - - - - - - - - lanTrbounds.c */
62 int LanTrbounds(int lanm, int maxit, double tol, double *vinit, int bndtype,
63  double *lammin, double *lammax, FILE *fstats);
64 
65 /*- - - -- - - - - - misc_la.c */
66 int scalEigVec(int n, int nev, double *Y, double* sqrtdiag);
67 
68 
69 /*- - - - - - - - - ratfilter.c */
70 //
71 void set_ratf_def(ratparams *rat);
72 //
73 int find_ratf(double *intv, ratparams *rat);
74 //
75 int set_ratf_solfunc(ratparams *rat, csrMat *A, csrMat *B, SolFuncC *funcs,
76  void **data);
77 //
78 void free_rat(ratparams *rat);
79 
80 /*- - - - - - - - - ratlanNr.c */
81 //
82 int RatLanNr(double *intv, int maxit, double tol, double *vinit, ratparams *rat,
83  int *nevOut, double **lamo, double **Wo, double **reso,
84  FILE *fstats);
85 
86 /*- - - - - - - - - ratlanTr.c */
87 //
88 int RatLanTr(int lanm, int nev, double *intv, int maxit, double tol,
89  double *vinit, ratparams *rat, int *nev2, double **vals,
90  double **W, double **resW, FILE *fstats);
91 
92 /*- - - - - - - - - spmat.c */
93 void csr_copy(csrMat *A, csrMat *B, int allocB);
94 // convert a COO matrix to a CSR matrix
95 int cooMat_to_csrMat(int cooidx, cooMat *coo, csrMat *csr);
96 // free a CSR
97 void free_csr(csrMat *csr);
98 // free a COO
99 void free_coo(cooMat *coo);
100 // matvec y = A*x
101 int matvec(char trans, csrMat *A, double *x, double *y);
102 
103 /*- - - - - - - - - evsl.c */
104 /* set an external matvec function */
105 int SetAMatvec(int n, MVFunc func, void *data);
106 int SetBMatvec(int n, MVFunc func, void *data);
107 /* unset an external matvec function */
108 int UnsetAMatvec();
109 int UnsetBMatvec();
110 int SetAMatrix(csrMat *A);
111 int SetBMatrix(csrMat *B);
112 int SetBSol(SolFuncR func, void *data);
113 int SetLTSol(SolFuncR func, void *data);
114 int SetASigmaBSol(ratparams *rat, SolFuncC *func, SolFuncC allf, void **data);
115 int SetStdEig();
116 int SetGenEig();
117 /* start EVSL */
118 int EVSLStart();
119 /* finalize EVSL */
120 int EVSLFinish();
121 void SetDiagScal(double *ds);
122 
123 /*- - - - - - - - - spslicer.c */
124 //
125 int spslicer(double *sli, double *mu, int Mdeg, double *intv, int n_int,
126  int npts);
127 //
128 int kpmdos(int Mdeg, int damping, int nvec, double *ab, double *mu,
129  double *ecnt);
130 
131 /*- - - - - - - - - - spslicer2.c */
132 void spslicer2(double *xi, double *yi, int n_int, int npts, double *sli);
133 
134 /*- - - - - - - - - timing.c */
135 //
136 double evsl_timer();
137 
138 /*- - - - - - - - - vect.c */
139 // generate a random vector
140 void rand_double(int n, double *v);
141 // generate a normally distributed random vector
142 void randn_double(int n, double *v);
143 // sort a vector
144 void sort_double(int n, double *v, int *ind);
145 //
146 void linspace(double a, double b, int num, double *arr);
147 
148 void StatsPrint(FILE *fstats);
149 
150 void StatsReset();
151 
152 #endif
void free_csr(csrMat *csr)
memory deallocation for csr matrix
Definition: spmat.c:91
void SetupPolRec(int n, int max_deg, double tol, double lmin, double lmax, BSolDataPol *data)
Definition: dos_utils.c:37
int EVSLFinish()
Finish EVSL.
Definition: evsl.c:48
void sort_double(int n, double *v, int *ind)
Definition: vect.c:92
int set_ratf_solfunc(ratparams *rat, csrMat *A, csrMat *B, SolFuncC *funcs, void **data)
int SetAMatrix(csrMat *A)
Set the matrix A.
Definition: evsl.c:68
void FreeBSolPolData(BSolDataPol *data)
Definition: dos_utils.c:53
void(* MVFunc)(double *x, double *y, void *data)
matvec function prototype
Definition: struct.h:96
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
void free_rat(ratparams *rat)
Definition: ratfilter.c:424
double evsl_timer()
evsl timer for mac
Definition: mactime.c:14
int LanDos(const int nvec, int msteps, const int npts, double *xdos, double *ydos, double *neig, const double *const intv)
Definition: landos.c:40
void SetupPolSqrt(int n, int max_deg, double tol, double lmin, double lmax, BSolDataPol *data)
Definition: dos_utils.c:45
void randn_double(int n, double *v)
Definition: vect.c:26
int ChebLanNr(double *intv, int maxit, double tol, double *vinit, polparams *pol, int *nevOut, double **lamo, double **Wo, double **reso, FILE *fstats)
Definition: cheblanNr.c:60
void extrDiagCsr(csrMat *B, double *d)
Definition: spmat.c:406
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
int ChebSI(int nev, double *intv, int maxit, double tol, double *vinit, polparams *pol, int *nevo, double **lamo, double **Yo, double **reso, FILE *fstats)
Chebyshev polynomial filtering Subspace Iteration.
Definition: chebsi.c:42
int UnsetBMatvec()
void diagScalCsr(csrMat *A, double *d)
Definition: spmat.c:390
int EVSLStart()
Initialize evslData.
Definition: evsl.c:27
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 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
int UnsetAMatvec()
int SetLTSol(SolFuncR func, void *data)
Set the solve routine for LT.
Definition: evsl.c:185
void csr_copy(csrMat *A, csrMat *B, int allocB)
copy a csr matrix A into B alloB: 0: will not allocate memory for B (have been alloced outside) 1: wi...
Definition: spmat.c:106
void rand_double(int n, double *v)
Definition: vect.c:11
void spslicer2(double *xi, double *yi, int n_int, int npts, double *sli)
Definition: spslice2.c:29
int LanTrbounds(int lanm, int maxit, double tol, double *vinit, int bndtype, double *lammin, double *lammax, FILE *fstats)
Lanczos process for eigenvalue bounds [Thick restart version].
Definition: lanTrbounds.c:38
void linspace(double a, double b, int num, double *arr)
Definition: vect.c:54
int matvec(char trans, csrMat *A, double *x, double *y)
structs used in evsl
void SetDiagScal(double *ds)
Set diagonal scaling matrix D.
Definition: evsl.c:197
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 find_ratf(double *intv, ratparams *rat)
Definition: ratfilter.c:375
void StatsReset()
Definition: stats.c:79
int SetStdEig()
Set the problem to standard eigenvalue problem.
Definition: evsl.c:148
int RatLanNr(double *intv, int maxit, double tol, double *vinit, ratparams *rat, int *nevOut, double **lamo, double **Wo, double **reso, FILE *fstats)
Rational filtering Lanczos process [NON-restarted version].
Definition: ratlanNr.c:54
int SetBMatrix(csrMat *B)
Set the B matrix.
Definition: evsl.c:83
int kpmdos(int Mdeg, int damping, int nvec, double *ab, double *mu, double *ecnt)
This function computes the coefficients of the density of states in the chebyshev basis...
Definition: spslice.c:32
int SetGenEig()
Set the problem to generalized eigenvalue problem.
Definition: evsl.c:158
void set_ratf_def(ratparams *rat)
Sets default values for ratparams struct.
Definition: ratfilter.c:345
int ChebLanTr(int lanm, int nev, double *intv, int maxit, double tol, double *vinit, polparams *pol, int *nev2, double **vals, double **W, double **resW, FILE *fstats)
Chebyshev polynomial filtering Lanczos process [Thick restart version].
Definition: cheblanTr.c:62
int cooMat_to_csrMat(int cooidx, cooMat *coo, csrMat *csr)
convert coo to csr
Definition: spmat.c:135
parameters for polynomial filter
Definition: struct.h:45
void set_pol_def(polparams *pol)
set default values for polparams struct.
Definition: chebpoly.c:18
void StatsPrint(FILE *fstats)
Definition: stats.c:14
int SetBSol(SolFuncR func, void *data)
Set the solve routine and the associated data for B.
Definition: evsl.c:132
int RatLanTr(int lanm, int nev, double *intv, int maxit, double tol, double *vinit, ratparams *rat, int *nev2, double **vals, double **W, double **resW, FILE *fstats)
RatLanTR filtering Lanczos process [Thick restart version].
Definition: ratlanTr.c:61
void free_coo(cooMat *coo)
memory deallocation for coo matrix
Definition: spmat.c:126
int scalEigVec(int n, int nev, double *Y, double *sqrtdiag)
void free_pol(polparams *pol)
Definition: chebpoly.c:488
int LanDosG(const int nvec, int msteps, const int npts, double *xdos, double *ydos, double *neig, const double *const intv)
Definition: landosG.c:43
int LanBounds(int msteps, double *v, double *lmin, double *lmax)
Definition: lanbounds.c:21
sparse matrix format: the coordinate (COO) format, 0-based
Definition: struct.h:16
int find_pol(double *intv, polparams *pol)
Sets the values in pol.
Definition: chebpoly.c:361
void BSolPol(double *b, double *x, void *data)
Definition: dos_utils.c:61
int spslicer(double *sli, double *mu, int Mdeg, double *intv, int n_int, int npts)
given the dos function defined by mu find a partitioning of sub-interval [a,b] of the spectrum so eac...
Definition: spslice.c:204
parameters for rational filter
Definition: struct.h:112