EVSL  1.1.0
EigenValues Slicing Library
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros
internal_proto.h
Go to the documentation of this file.
1 
5 #ifndef INTERNAL_PROTO_H
6 #define INTERNAL_PROTO_H
7 
8 #include <stdio.h>
9 #include <complex.h>
10 #include "evsl.h"
11 #include "struct.h"
12 
13 /*- - - - - - - - - chebpoly.c */
14 //
15 int chebxCoefd(int m, double gam, int damping, double *mu);
16 //
17 int dampcf(int m, int damping, double *jac);
18 //
19 int chebxPltd(int m, double *mu, int n, double *xi, double *yi);
20 //
21 int ChebAv(polparams *pol, double *v, double *y, double *w);
22 //
23 void chext(polparams *pol, double aIn, double bIn);
24 
25 /*- - - - - - - - - dos_utils.c */
26 int apfun(const double c, const double h, const double* xi, double (*ffun)(double), const int npts, double* yi);
27 //
28 double rec(const double a);
29 //
30 double isqrt(const double a);
31 //
32 int pnav(double *mu, const int m, const double cc, const double dd, double *v,
33  double *y, double *w); // Really just ChebAv
34 //
35 int lsPol(double (*ffun)(double), BSolDataPol *pol);
36 
37 /*- - - - - - - - - dump.c */
38 //
39 void savemat(csrMat *A, const char *fn);
40 //
41 void savedensemat(double *A, int lda, int m, int n, const char *fn);
42 //
43 void save_vec(int n, const double *x, const char fn[]);
44 
45 /*- - - - - - - - - evsl.c */
46 //
47 
48 /*- - - - - - - - - misc_la.c */
49 //
50 int SymmTridEig(double *eigVal, double *eigVec, int n, const double *diag, const double *sdiag);
51 //
52 int SymmTridEigS(double *eigVal, double *eigVec, int n, double vl, double vu, int *nevO, const double *diag, const double *sdiag);
53 //
54 void SymEigenSolver(int n, double *A, int lda, double *Q, int ldq, double *lam);
55 //
56 void CGS_DGKS(int n, int k, int i_max, double *Q, double *v, double *nrmv, double *w);
57 //
58 void CGS_DGKS2(int n, int k, int i_max, double *Z, double *Q, double *v, double *w);
59 //
60 void orth(double *V, int n, int k, double *Vo, double *work);
61 
62 
63 /*- - - - - - - - - ratfilter.c */
64 //
65 void contQuad(int method, int n, complex double* zk);
66 //
67 void ratf2p2(int n, int *mulp, complex double *zk, complex double* alp, int m, double *z, double *x);
68 //
69 void pfe2(complex double s1, complex double s2, int k1, int k2, complex double* alp, complex double* bet);
70 //
71 complex double integg2(complex double s1, complex double s2, complex double* alp, int k1, complex double* bet, int k2, double a, double b);
72 //
73 void weights(int n, complex double* zk, int* pow, double lambda, complex double* omega);
74 //
75 int scaleweigthts(int n, double a, double b, complex double *zk, int* pow, complex double* omegaM);
76 
77 /*- - - - - - - - - ratlanNr.c */
78 void RatFiltApply(int n, ratparams *rat, double *b, double *x, double *w3);
79 
80 /*- - - - - - - - - - simpson.c */
81 void simpson(double *xi, double *yi, int npts);
82 
83 /*- - - - - - - - - spmat.c */
84 void matvec_csr(double *x, double *y, void *data);
85 
86 // memory allocation/reallocation for a CSR matrix
87 void csr_resize(int nrow, int ncol, int nnz, csrMat *csr);
88 //
89 void sortrow(csrMat *A);
90 int matadd(double alp, double bet, csrMat *A, csrMat *B, csrMat *C,
91  int *mapA, int *mapB);
92 /* sparse identity */
93 int speye(int n, csrMat *A);
94 /* extract upper triangular part of A */
95 void triuCsr(csrMat *A, csrMat *U);
96 
97 /*- - - - - - - - - timing.c */
98 int time_seeder();
99 
100 /*- - - - - - - - - vect.c */
101 void vecset(int n, double t, double *v);
102 void vec_perm(int n, int *p, double *x, double *y);
103 void vec_iperm(int n, int *p, double *x, double *y);
104 
105 /*------------------- inline functions */
110 static inline void matvec_A(double *x, double *y) {
111  CHKERR(!evsldata.Amv);
112  double tms = evsl_timer();
113  evsldata.Amv->func(x, y, evsldata.Amv->data);
114  double tme = evsl_timer();
115  evslstat.t_mvA += tme - tms;
116  evslstat.n_mvA ++;
117 }
118 
123 static inline void matvec_B(double *x, double *y) {
124  CHKERR(!evsldata.Bmv);
125  double tms = evsl_timer();
126  evsldata.Bmv->func(x, y, evsldata.Bmv->data);
127  double tme = evsl_timer();
128  evslstat.t_mvB += tme - tms;
129  evslstat.n_mvB ++;
130 }
131 
136 static inline void solve_B(double *x, double *y) {
137  CHKERR(!evsldata.Bsol);
138  double tms = evsl_timer();
139  evsldata.Bsol->func(x, y, evsldata.Bsol->data);
140  double tme = evsl_timer();
141  evslstat.t_svB += tme - tms;
142  evslstat.n_svB ++;
143 }
144 
149 static inline void solve_LT(double *x, double *y) {
151  double tms = evsl_timer();
153  double tme = evsl_timer();
154  evslstat.t_svLT += tme - tms;
155  evslstat.n_svLT ++;
156 }
157 
162 static inline void solve_ASigB(EVSLASIGMABSol *sol, int n,
163  double *br, double *bz,
164  double *xr, double *xz) {
165  double tms = evsl_timer();
166  (sol->func)(n, br, bz, xr, xz, sol->data);
167  double tme = evsl_timer();
168  evslstat.t_svASigB+= tme - tms;
169  evslstat.n_svASigB ++;
170 }
171 
172 /*- - - - - - - - - - check if an interval is valid */
173 static inline int check_intv(double *intv, FILE *fstats) {
174  /* intv[4]: ( intv[0], intv[1] ) is the inteval of interest
175  * ( intv[2], intv[3] ) is the spectrum bounds
176  * return 0: ok
177  * < 0: interval is invalid
178  */
179  double a=intv[0], b=intv[1], lmin=intv[2], lmax=intv[3];
180  if (a >= b) {
181  if (fstats) {
182  fprintf(fstats, " error: invalid interval (%e, %e)\n", a, b);
183  }
184  return -1;
185  }
186 
187  if (a >= lmax || b <= lmin) {
188  if (fstats) {
189  fprintf(fstats, " error: interval (%e, %e) is outside (%e %e) \n",
190  a, b, lmin, lmax);
191  }
192  return -2;
193  }
194 
195  return 0;
196 }
197 
198 #endif
void csr_resize(int nrow, int ncol, int nnz, csrMat *csr)
memory allocation for csr matrix
Definition: spmat.c:79
size_t n_svASigB
Definition: struct.h:218
EVSLMatvec * Amv
Definition: struct.h:167
void CGS_DGKS(int n, int k, int i_max, double *Q, double *v, double *nrmv, double *w)
Classical GS reortho with Daniel, Gragg, Kaufman, Stewart test.
Definition: misc_la.c:195
void save_vec(int n, const double *x, const char fn[])
Definition: dumps.c:38
#define CHKERR(ierr)
Definition: def.h:19
int SymmTridEig(double *eigVal, double *eigVec, int n, const double *diag, const double *sdiag)
compute all eigenvalues and eigenvectors of a symmetric tridiagonal matrix
Definition: misc_la.c:36
int speye(int n, csrMat *A)
return an identity matrix of dimension n
Definition: spmat.c:357
void vecset(int n, double t, double *v)
Definition: vect.c:48
double t_mvB
Definition: struct.h:203
EVSLBSol * Bsol
Definition: struct.h:169
int chebxPltd(int m, double *mu, int n, double *xi, double *yi)
function yi = chebxPltd computes yi = p_mu (xi),
Definition: chebpoly.c:106
void * data
Definition: struct.h:104
void * data
Definition: struct.h:138
void savemat(csrMat *A, const char *fn)
Definition: dumps.c:33
void vec_iperm(int n, int *p, double *x, double *y)
Definition: vect.c:127
void simpson(double *xi, double *yi, int npts)
Definition: simpson.c:27
size_t n_svLT
Definition: struct.h:217
void * data
Definition: struct.h:156
EVSLMatvec * Bmv
Definition: struct.h:168
real *8 function ffun(x, y, z)
Definition: functns.f90:40
complex double integg2(complex double s1, complex double s2, complex double *alp, int k1, complex double *bet, int k2, double a, double b)
Integration of 1/[(z-s1)^k1 (z-s2)^k2] from a to b.
Definition: ratfilter.c:129
int pnav(double *mu, const int m, const double cc, const double dd, double *v, double *y, double *w)
Definition: dos_utils.c:109
void contQuad(int method, int n, complex double *zk)
Compute the locations of the poles.
Definition: ratfilter.c:23
double isqrt(const double a)
Definition: dos_utils.c:17
void SymEigenSolver(int n, double *A, int lda, double *Q, int ldq, double *lam)
interface to LAPACK SYMMETRIC EIGEN-SOLVER
Definition: misc_la.c:156
double t_mvA
Definition: struct.h:202
evslStat evslstat
global statistics of EVSL
Definition: evsl.c:21
void * data
Definition: struct.h:147
void ratf2p2(int n, int *mulp, complex double *zk, complex double *alp, int m, double *z, double *x)
Compute the function value of the multiple pole rational filter at real locations.
Definition: ratfilter.c:69
double rec(const double a)
Definition: dos_utils.c:15
int SymmTridEigS(double *eigVal, double *eigVec, int n, double vl, double vu, int *nevO, const double *diag, const double *sdiag)
compute eigenvalues and eigenvectors of a symmetric tridiagonal matrix in a slice ...
Definition: misc_la.c:93
SolFuncR func
Definition: struct.h:155
void pfe2(complex double s1, complex double s2, int k1, int k2, complex double *alp, complex double *bet)
Get the fraction expansion of 1/[(z-s1)^k1 (z-s2)^k2].
Definition: ratfilter.c:95
structs used in evsl
int matadd(double alp, double bet, csrMat *A, csrMat *B, csrMat *C, int *mapA, int *mapB)
matrix addition C = alp * A + bet * B
Definition: spmat.c:307
user-provided function and data for solving (A - SIGMA*B) x = b
Definition: struct.h:102
sparse matrix format: the compressed sparse row (CSR) format, 0-based
Definition: struct.h:31
void chext(polparams *pol, double aIn, double bIn)
Determines polynomial for end interval cases.
Definition: chebpoly.c:171
int dampcf(int m, int damping, double *jac)
Computes damping coefficient for cheb. expansions.
Definition: chebpoly.c:40
This file contains function prototypes and constant definitions for EVSL.
double t_svASigB
Definition: struct.h:206
void CGS_DGKS2(int n, int k, int i_max, double *Z, double *Q, double *v, double *w)
Classical GS reortho. No test. just do i_max times used in generalized ev problems.
Definition: misc_la.c:235
int apfun(const double c, const double h, const double *xi, double(*ffun)(double), const int npts, double *yi)
Definition: dos_utils.c:83
void orth(double *V, int n, int k, double *Vo, double *work)
Orthogonalize columns of n-by-k matrix V.
Definition: misc_la.c:269
EVSLLTSol * LTsol
Definition: struct.h:170
MVFunc func
Definition: struct.h:137
void triuCsr(csrMat *A, csrMat *U)
Definition: spmat.c:424
void sortrow(csrMat *A)
Sort each row of a csr by increasing column order By double transposition.
Definition: spmat.c:56
size_t n_mvB
Definition: struct.h:215
int lsPol(double(*ffun)(double), BSolDataPol *pol)
Definition: dos_utils.c:168
void vec_perm(int n, int *p, double *x, double *y)
Definition: vect.c:114
int chebxCoefd(int m, double gam, int damping, double *mu)
evslData evsldata
global variable of EVSL
Definition: evsl.c:15
void RatFiltApply(int n, ratparams *rat, double *b, double *x, double *w3)
Apply rational filter R to a vetor b.
Definition: ratfilter.c:448
parameters for polynomial filter
Definition: struct.h:45
double t_svLT
Definition: struct.h:205
int ChebAv(polparams *pol, double *v, double *y, double *w)
Computes y=P(A) v, where pn is a Cheb. polynomial expansion
Definition: chebpoly.c:510
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 time_seeder()
Uses the timer to generate a seed to be used for srand.
Definition: mactime.c:30
size_t n_mvA
Definition: struct.h:214
double t_svB
Definition: struct.h:204
SolFuncR func
Definition: struct.h:146
void weights(int n, complex double *zk, int *pow, double lambda, complex double *omega)
Compute the LS weight for each multiple pole.
Definition: ratfilter.c:169
double evsl_timer()
evsl timer for mac
Definition: mactime.c:14
SolFuncC func
Definition: struct.h:103
void savedensemat(double *A, int lda, int m, int n, const char *fn)
Definition: dumps.c:49
int scaleweigthts(int n, double a, double b, complex double *zk, int *pow, complex double *omegaM)
Compute the weights and pole locations on [a, b].
Definition: ratfilter.c:313
parameters for rational filter
Definition: struct.h:112
size_t n_svB
Definition: struct.h:216