ConicBundle
LanczMaxEig.hxx
Go to the documentation of this file.
1 
2 
3 #ifndef CH_MATRIX_CLASSES__LANCZMAXEIG_HXX
4 #define CH_MATRIX_CLASSES__LANCZMAXEIG_HXX
5 
15 #ifndef CH_MATRIX_CLASSES__LANCZOS_HXX
16 #include "lanczos.hxx"
17 #endif
18 #include "clock.hxx"
19 
20 namespace CH_Matrix_Classes {
21 
22 
26 
27 
35 class LanczMaxEig:public Lanczos,protected Memarrayuser
36 {
37 private:
38  int ierr;
49 
50  //--- global variables to allow reastarting
55 
62 
72 
74 
75  int stop_above;
77 
80 
82 
88 
90  std::ostream* myout;
91 
93  void error(const char *s){
94  if (myout)
95  (*myout)<<"Lanczos Error: "<<s<<std::endl;
96  }
97 
99  int guess_extremes(Integer nproposed);
100 
102  int dhua(Integer nreig,Integer nproposed, Integer maxmult);
103 
105  int bklanc(Integer neigfound,Integer blocksz,Integer s,
106  const Matrix &d,Matrix &C,Matrix &X,
107  Matrix &e,Matrix &u,Matrix &v);
108 
110  int bkqrlanc(Integer neigfound,Integer blocksz,Integer s,
111  const Matrix &d,Matrix &C,Matrix &X,
112  Matrix &e,Matrix &u,Matrix &v);
113 
115  int bklanccheb(Integer neigfound,Integer blocksz,Integer s,
116  Matrix &C,Matrix &X,Matrix &e,Matrix &v);
117 
119  int bkqrlanccheb(Integer neigfound,Integer blocksz,Integer s,
120  Matrix &C,Matrix &X,Matrix &e,Matrix &v);
121 
123  int err(Integer neigfound,Integer blocksz,const Matrix& X,Matrix &e);
124 
126  int cnvtst(Integer neigfound,Integer blocksz,Real& errc,Real eps,
127  const Matrix& d,const Matrix &e,Indexmatrix& newconv);
128 
130  int eigen(Integer neigfound,Integer blocksz,Integer sbsz,Matrix& C,
131  Matrix& d,Matrix& u,Matrix& v,Real& af);
132 
134  int sectn(Matrix& X,Integer neigfound,Integer blocksz,
135  Matrix& C,Matrix& d,Matrix& u,Matrix& v,
136  Real& af);
137 
139  int rotate_extremes(Integer neigfound,Integer sbs,Matrix& d,
140  const Matrix& C,Matrix &X,Matrix& v);
141 
143  int rotate(Integer neigfound,Integer sbs,Integer l,
144  const Matrix& C,Matrix &X,Matrix& v);
145 
147  int resort(Integer neigfound,const Indexmatrix& newconv,
148  Matrix &X,Matrix& d,Matrix& e,Matrix& v);
149 
151  int random(Matrix& X,Integer j);
152 
154  Integer orthog(Integer offset,Integer blocksz,Matrix &X,Matrix& B);
155 
157  int blockcheby(Integer col_offset,const Matrix& X,Matrix& v);
158 
160  Real scalarcheby(Real xval);
161 
163  int tred2(Integer n,const Matrix& C,Matrix& u,Matrix& v,Matrix& Z);
164 
166  int tql2(Integer n,Matrix &u,Matrix &v, Matrix& Z);
167 
168 public:
170  LanczMaxEig();
172  ~LanczMaxEig();
173 
174 
181 
183  void set_mineig(Real ie){minval=ie;}
185  void set_maxmult(Integer mop){maxop=mop;}
187  void set_maxiter(Integer mi){maxiter=mi;}
189  void set_relprec(Real relprec){eps=relprec;}
191  void set_nchebit(Integer cheb){choicencheb=cheb;}
193  void set_nblockmult(Integer nb){choicenbmult=nb;}
194 
196  void enable_stop_above(Real ub){stop_above=1; upper_bound=ub;}
198  void disable_stop_above(){stop_above=0;}
199 
201  void set_retlanvecs(Integer nl){retlanvecs=nl;}
202 
204  int get_lanczosvecs(Matrix& val,Matrix& vecs) const;
205 
207  Real get_relprec(void){return eps;} //relative precision
208 
210  int get_err() const{return ierr;}
212  Integer get_iter() const{return iter;}
214  Integer get_nmult() const{return nmult;}
215 
217 
219  int compute(const Lanczosmatrix* bigmat,
220  Matrix& eigval,
221  Matrix& eigvec,
222  Integer nreig,
223  Integer in_blocksz=0,
224  Integer maxcol=0
225  );
226 
227 
228 
233 
235  void set_out(std::ostream* o=0,int pril=1)
236  {myout=o; print_level=pril;}
237 
238 
240  std::ostream& save(std::ostream& out) const;
241 
243  std::istream& restore(std::istream& in);
245 };
246 
248 
249 }
250 
251 #endif
252 
int Integer
all integer numbers in calculations and indexing are of this type
Definition: matop.hxx:40
Real maxval
current approximation of largest eigenvalue
Definition: LanczMaxEig.hxx:59
Integer orthog(Integer offset, Integer blocksz, Matrix &X, Matrix &B)
orthonormalize columns X(:,offset:offset+l_blocksz-1) to all previous columns
std::istream & restore(std::istream &in)
restore the data from in where it was stored by save()
Integer guessmult
number of matrix vector multiplications to guess interval
Definition: LanczMaxEig.hxx:42
Integer ncalls
counts number of calls to LanczMaxEig
Definition: LanczMaxEig.hxx:78
void set_maxmult(Integer mop)
set an upper bound on the number of matrix vector multiplications
Definition: LanczMaxEig.hxx:185
double Real
all real numbers in calculations are of this type
Definition: matop.hxx:50
Integer mymaxj
maximum amount of storage (columns) provided in X and C
Definition: LanczMaxEig.hxx:79
int compute(const Lanczosmatrix *bigmat, Matrix &eigval, Matrix &eigvec, Integer nreig, Integer in_blocksz=0, Integer maxcol=0)
the main routine: compute the nreig maximum eigenvalues of the matrix specified by bigmat ...
Integer nmult
number of single vector multiplications with matrix
Definition: LanczMaxEig.hxx:54
std::ostream & save(std::ostream &out) const
save all data in out so that the current state can be recovered completely by restore() ...
Matrix v
temporary matrix
Definition: LanczMaxEig.hxx:69
Integer get_iter() const
returns the number of restarts of the last call
Definition: LanczMaxEig.hxx:212
void set_out(std::ostream *o=0, int pril=1)
set output stream and level of detail of log output (for debugging)
Definition: LanczMaxEig.hxx:235
int rotate_extremes(Integer neigfound, Integer sbs, Matrix &d, const Matrix &C, Matrix &X, Matrix &v)
rotate the eigenvectors of the largest and smallest eigenvalues of the tridiagonal matrix...
Matrix class for integral values of type Integer
Definition: indexmat.hxx:195
int cnvtst(Integer neigfound, Integer blocksz, Real &errc, Real eps, const Matrix &d, const Matrix &e, Indexmatrix &newconv)
check convergence of maximum Ritz value/vector
int stop_above
1 if algorithm is to stop after upper bound is exceeded
Definition: LanczMaxEig.hxx:75
Integer neigfound
number of eigenvalues known, these are in the first neigfound columns of X
Definition: LanczMaxEig.hxx:51
Integer nlanczvecs
number of columns of storage matrix X carrying "meaningful" Ritz vectors
Definition: LanczMaxEig.hxx:47
A Lanczos method allowing spectral transformation by Chebycheff polynomials and premature termination...
Definition: LanczMaxEig.hxx:35
int ierr
error return code
Definition: LanczMaxEig.hxx:38
void set_mineig(Real ie)
set a guess on the value of the smallest
Definition: LanczMaxEig.hxx:183
int dhua(Integer nreig, Integer nproposed, Integer maxmult)
main routine performing the Lancozs iterations
int random(Matrix &X, Integer j)
assign a random vector to column j of X
Abstract base class for supplying the input matrix for Lanzcosmethods.
Definition: lanczos.hxx:33
int blockcheby(Integer col_offset, const Matrix &X, Matrix &v)
apply spectral transformation by using a Chebycheff polynomial on the block multiplications ...
CH_Tools::GB_rand randgen
local random number generator
Definition: LanczMaxEig.hxx:81
void set_maxiter(Integer mi)
set an upper bpound on the number of restarts
Definition: LanczMaxEig.hxx:187
CH_Tools::Microseconds time_iter
time spent in last lanczos iteration
Definition: LanczMaxEig.hxx:86
Integer choicenbmult
user&#39;s choice for number of block multiplications (<0 -> automatic determination, min 6) ...
Definition: LanczMaxEig.hxx:45
const Lanczosmatrix * bigmatrix
pointer giving the (virtual) input matrix
Definition: LanczMaxEig.hxx:73
Real mcheps
machine precision (computed in constructor)
Definition: LanczMaxEig.hxx:58
int sectn(Matrix &X, Integer neigfound, Integer blocksz, Matrix &C, Matrix &d, Matrix &u, Matrix &v, Real &af)
compute matrix for eigenvalue computation into C
CH_Tools::Microseconds time_sum
time spent in this call to compute()
Definition: LanczMaxEig.hxx:87
Matrix e
error bounds
Definition: LanczMaxEig.hxx:66
allows measuring time difference to its initialization time in Microseconds
Definition: clock.hxx:282
int bklanc(Integer neigfound, Integer blocksz, Integer s, const Matrix &d, Matrix &C, Matrix &X, Matrix &e, Matrix &u, Matrix &v)
block lanzos multiplication without Chebycheff spectral transformation
void enable_stop_above(Real ub)
allow the algorithm to stop as soon as the maximum Ritz value exceeds the value ub ...
Definition: LanczMaxEig.hxx:196
~LanczMaxEig()
destructor, nothing particular
extra long integer number for expressing and computing time measurements in microseconds.
Definition: clock.hxx:46
device independent random number generator based on long int with seed
Definition: gb_rand.hxx:28
Integer iter
iteration counter for Lanczos restarts (interval guess+computation)
Definition: LanczMaxEig.hxx:53
void set_retlanvecs(Integer nl)
set an upper bound on the number of vectors returned in get_lanczosvecs()
Definition: LanczMaxEig.hxx:201
Matrix minvec
temporary matrix for guessing minimal eigenvalue
Definition: LanczMaxEig.hxx:71
Abstract interface to Lanzcos methods for computing a few extremal eigenvalues given via a Lanczosmat...
Definition: lanczos.hxx:53
int eigen(Integer neigfound, Integer blocksz, Integer sbsz, Matrix &C, Matrix &d, Matrix &u, Matrix &v, Real &af)
compute eigenvalues of current (block) tridiagonalization
Real get_relprec(void)
returns current relative precision requirement
Definition: LanczMaxEig.hxx:207
void set_relprec(Real relprec)
set relative precision requirement for termination
Definition: LanczMaxEig.hxx:189
int guess_extremes(Integer nproposed)
compute a guess for minimum and maximum eigenvalue
int bkqrlanccheb(Integer neigfound, Integer blocksz, Integer s, Matrix &C, Matrix &X, Matrix &e, Matrix &v)
block lanzos multiplication with complete QR orthogonalization and Chebycheff spectral transformation...
Integer choicencheb
user&#39;s choice for number of block Chebychev iterations (<0 -> automatic determination) ...
Definition: LanczMaxEig.hxx:43
Integer nchebit
number of block Chebychev iterations within one iteration
Definition: LanczMaxEig.hxx:44
int get_lanczosvecs(Matrix &val, Matrix &vecs) const
returns the Lanczos vectors of the last call with their Ritz values
LanczMaxEig()
intialize all to default values
Matrix Classes and Linear Algebra. See Matrix Classes (namespace CH_Matrix_Classes) for a quick intro...
Definition: PSCOracle.hxx:20
Real scalarcheby(Real xval)
compute the same polynomial as in blockcheby but for the scalar value xval
CH_Tools::Microseconds time_mult
for each restart, the time spent in lanczosmult
Definition: LanczMaxEig.hxx:84
int print_level
level of iteration information that should be displayed
Definition: LanczMaxEig.hxx:89
Real minval
current approximation of smallest eigenvalue
Definition: LanczMaxEig.hxx:60
Real polval
the Chebychef polynomial will have this value at maxval
Definition: LanczMaxEig.hxx:61
int bklanccheb(Integer neigfound, Integer blocksz, Integer s, Matrix &C, Matrix &X, Matrix &e, Matrix &v)
block lanzos multiplication with Chebycheff spectral transformation
Integer retlanvecs
user&#39;s choice for number of returend Ritz vectors (<0 -> nlanczvecs)
Definition: LanczMaxEig.hxx:48
std::ostream * myout
everything is output to *myout, default: myout=&cout (may be 0 for no output)
Definition: LanczMaxEig.hxx:90
Integer maxguessiter
upper bound on number of restarts to guess spectral interval
Definition: LanczMaxEig.hxx:41
Matrix u
temporary matrix
Definition: LanczMaxEig.hxx:68
Header declaring the abstract classes CH_Matrix_Classes::Lanczosmatrix and CH_Matrix_Classes::Lanczos...
Integer maxop
upper bound on matrix vector multiplications
Definition: LanczMaxEig.hxx:39
Matrix class for real values of type Real
Definition: matrix.hxx:74
int get_err() const
returns the error code of the last call
Definition: LanczMaxEig.hxx:210
Integer nblockmult
number of blockmultiplications in one restart
Definition: LanczMaxEig.hxx:46
void error(const char *s)
output of error messages
Definition: LanczMaxEig.hxx:93
int bkqrlanc(Integer neigfound, Integer blocksz, Integer s, const Matrix &d, Matrix &C, Matrix &X, Matrix &e, Matrix &u, Matrix &v)
block lanzos multiplication, with complete QR orthogonalization, without Chebycheff spectral transfor...
int err(Integer neigfound, Integer blocksz, const Matrix &X, Matrix &e)
compute norms of deviations of the Ritz vectors from being eigenvectors
int rotate(Integer neigfound, Integer sbs, Integer l, const Matrix &C, Matrix &X, Matrix &v)
rotate the lanzos vectors to Ritz vectors
Matrix C
X^tAX, intermediate eigenvalue computations, orthogonalizations, etc.
Definition: LanczMaxEig.hxx:64
Matrix X
provides storage for the Lanczos vectors during computation
Definition: LanczMaxEig.hxx:63
All derived classes share a common Memarray memory manager, which is generated with the first user an...
Definition: memarray.hxx:117
void disable_stop_above()
do not allow premature termination as in enable_stop_above()
Definition: LanczMaxEig.hxx:198
Real errc
error accumulation
Definition: LanczMaxEig.hxx:56
Integer get_nmult() const
returns the number of matrix-vector multiplications of the last call
Definition: LanczMaxEig.hxx:214
CH_Tools::Clock myclock
for time measurements
Definition: LanczMaxEig.hxx:83
Real eps
relative precision
Definition: LanczMaxEig.hxx:57
Matrix d
diagonals (Ritz values)
Definition: LanczMaxEig.hxx:65
int resort(Integer neigfound, const Indexmatrix &newconv, Matrix &X, Matrix &d, Matrix &e, Matrix &v)
rotate the lanzos vectors to Ritz vectors
Header declaring and (inline) implementing the classes CH_Tools::Microseconds and CH_Tools::Clock as ...
CH_Tools::Microseconds time_mult_sum
sum over time_mult for all restarts
Definition: LanczMaxEig.hxx:85
Matrix Xqr
for complete orthogonalization with Householder QR
Definition: LanczMaxEig.hxx:67
Integer maxiter
upper bound on number of restarts
Definition: LanczMaxEig.hxx:40
void set_nblockmult(Integer nb)
set maximum number of block multiplications within one restart
Definition: LanczMaxEig.hxx:193
void set_nchebit(Integer cheb)
set the degree of the Chebycheff polynomial for the spectral transformation
Definition: LanczMaxEig.hxx:191
Real upper_bound
stop if current maximum Ritz value exceeds this value
Definition: LanczMaxEig.hxx:76
int tred2(Integer n, const Matrix &C, Matrix &u, Matrix &v, Matrix &Z)
tridiagonalize a symmetric blockdiagonal matrix
Matrix w
temporary matrix
Definition: LanczMaxEig.hxx:70
Integer blocksz
working block of lanczos vectors X(:,neigfound:neigfound+blocksz-1)
Definition: LanczMaxEig.hxx:52
int tql2(Integer n, Matrix &u, Matrix &v, Matrix &Z)
compute the eigenvalues of a tridiagonal matrix