MPQC  2.3.1
matrix.h
1 //
2 // matrix.h
3 //
4 // Copyright (C) 1996 Limit Point Systems, Inc.
5 //
6 // Author: Curtis Janssen <cljanss@limitpt.com>
7 // Maintainer: LPS
8 //
9 // This file is part of the SC Toolkit.
10 //
11 // The SC Toolkit is free software; you can redistribute it and/or modify
12 // it under the terms of the GNU Library General Public License as published by
13 // the Free Software Foundation; either version 2, or (at your option)
14 // any later version.
15 //
16 // The SC Toolkit is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU Library General Public License for more details.
20 //
21 // You should have received a copy of the GNU Library General Public License
22 // along with the SC Toolkit; see the file COPYING.LIB. If not, write to
23 // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
24 //
25 // The U.S. Government is granted a limited license as per AL 91-7.
26 //
27 
28 #ifndef _math_scmat_matrix_h
29 #define _math_scmat_matrix_h
30 #ifdef __GNUC__
31 #pragma interface
32 #endif
33 
34 #include <iostream>
35 
36 #include <math/scmat/abstract.h>
37 
38 namespace sc {
39 
40 class SCVectordouble;
41 class SCMatrixdouble;
42 class SymmSCMatrixdouble;
43 class DiagSCMatrixdouble;
44 
45 class SCMatrixBlockIter;
46 class SCMatrixRectBlock;
47 class SCMatrixLTriBlock;
48 class SCMatrixDiagBlock;
49 class SCVectorSimpleBlock;
50 
51 class RefSCMatrix;
52 class RefSymmSCMatrix;
55 class RefSCVector: public Ref<SCVector> {
56  // standard overrides
57  public:
60  RefSCVector();
62  RefSCVector(const RefSCVector& v);
65  // don't allow automatic conversion from any reference to a
66  // described class
67  ~RefSCVector();
72 
73  // vector specific members
74  public:
77  RefSCVector(const RefSCDimension& dim,const Ref<SCMatrixKit>&);
78 
80  SCVectordouble operator()(int) const;
82  SCVectordouble operator[](int) const;
84  RefSCVector operator+(const RefSCVector&a) const;
86  RefSCVector operator-(const RefSCVector&a) const;
88  RefSCVector operator*(double) const;
90  RefSCMatrix outer_product(const RefSCVector& v) const;
93 
94  void set_element(int i,double val) const;
95  void accumulate_element(int i,double val) const;
96  double get_element(int) const;
97  int n() const;
98  RefSCDimension dim() const;
99  Ref<SCMatrixKit> kit() const;
100  RefSCVector clone() const;
101  RefSCVector copy() const;
102  double maxabs() const;
103  double scalar_product(const RefSCVector&) const;
104  double dot(const RefSCVector&) const;
105  void normalize() const;
106  void randomize() const;
107  void assign(const RefSCVector& v) const;
108  void assign(double val) const;
109  void assign(const double* v) const;
110  void convert(double*) const;
111  void scale(double val) const;
112  void accumulate(const RefSCVector& v) const;
113  void accumulate_product(const RefSymmSCMatrix&, const RefSCVector&);
114  void accumulate_product(const RefSCMatrix&, const RefSCVector&);
115  void element_op(const Ref<SCElementOp>& op) const;
116  void element_op(const Ref<SCElementOp2>&,
117  const RefSCVector&) const;
118  void element_op(const Ref<SCElementOp3>&,
119  const RefSCVector&,
120  const RefSCVector&) const;
121  void print(std::ostream&out) const;
122  void print(const char*title=0,
123  std::ostream&out=ExEnv::out0(), int precision=10) const;
124  void save(StateOut&);
126  void restore(StateIn&);
127 };
128 RefSCVector operator*(double,const RefSCVector&);
129 
130 class RefSymmSCMatrix;
131 class RefDiagSCMatrix;
135 class RefSCMatrix: public Ref<SCMatrix> {
136  // standard overrides
137  public:
140  RefSCMatrix();
142  RefSCMatrix(const RefSCMatrix& m);
144  RefSCMatrix(SCMatrix* m);
145  ~RefSCMatrix();
149  RefSCMatrix& operator=(const RefSCMatrix& m);
150 
151  // matrix specific members
152  public:
155  RefSCMatrix(const RefSCDimension& d1,const RefSCDimension& d2,
156  const Ref<SCMatrixKit>&);
157 
159  RefSCVector operator*(const RefSCVector&) const;
160 
162  RefSCMatrix operator*(const RefSCMatrix&) const;
163  RefSCMatrix operator*(const RefSymmSCMatrix&) const;
164  RefSCMatrix operator*(const RefDiagSCMatrix&) const;
165 
167  RefSCMatrix operator*(double) const;
168 
170  RefSCMatrix operator+(const RefSCMatrix&) const;
172  RefSCMatrix operator-(const RefSCMatrix&) const;
173 
175  RefSCMatrix t() const;
177  RefSCMatrix i() const;
179  RefSCMatrix gi() const;
180 
183  RefSCMatrix clone() const;
184  RefSCMatrix copy() const;
185 
186  RefSCMatrix get_subblock(int br, int er, int bc, int ec);
187  void assign_subblock(const RefSCMatrix&, int br, int er, int bc, int ec,
188  int source_br = 0, int source_bc = 0);
189  void accumulate_subblock(const RefSCMatrix&, int, int, int, int,
190  int source_br = 0, int source_bc = 0);
191  RefSCVector get_row(int) const;
192  RefSCVector get_column(int) const;
193  void assign_row(const RefSCVector&, int) const;
194  void assign_column(const RefSCVector&, int) const;
195  void accumulate_row(const RefSCVector&, int) const;
196  void accumulate_column(const RefSCVector&, int) const;
197 
198  void accumulate_outer_product(const RefSCVector&,const RefSCVector&) const;
199  void accumulate_product(const RefSCMatrix&,const RefSCMatrix&) const;
200  void assign(const RefSCMatrix&) const;
201  void scale(double) const;
202  void randomize() const;
203  void assign(double) const;
204  void assign(const double*) const;
205  void assign(const double**) const;
206  void convert(double*) const;
207  void convert(double**) const;
208  void accumulate(const RefSCMatrix&) const;
209  void accumulate(const RefSymmSCMatrix&) const;
210  void accumulate(const RefDiagSCMatrix&) const;
211  void element_op(const Ref<SCElementOp>&) const;
212  void element_op(const Ref<SCElementOp2>&,
213  const RefSCMatrix&) const;
214  void element_op(const Ref<SCElementOp3>&,
215  const RefSCMatrix&,
216  const RefSCMatrix&) const;
217  int nrow() const;
218  int ncol() const;
219  RefSCDimension rowdim() const;
220  RefSCDimension coldim() const;
221  Ref<SCMatrixKit> kit() const;
222  void set_element(int,int,double) const;
223  void accumulate_element(int,int,double) const;
224  double get_element(int,int) const;
225  void print(std::ostream&) const;
226  void print(const char*title=0,
227  std::ostream&out=ExEnv::out0(), int =10) const;
228  double trace() const;
229  void save(StateOut&);
231  void restore(StateIn&);
232 
237  void svd(const RefSCMatrix &U,
238  const RefDiagSCMatrix &sigma,
239  const RefSCMatrix &V);
241  double solve_lin(const RefSCVector& v) const;
243  double determ() const;
245  SCMatrixdouble operator()(int i,int j) const;
246 
250  int nblock() const;
254  RefSCMatrix block(int i) const;
255 };
257 RefSCMatrix operator*(double,const RefSCMatrix&);
258 
261 class RefSymmSCMatrix: public Ref<SymmSCMatrix> {
262  // standard overrides
263  public:
266  RefSymmSCMatrix();
271  ~RefSymmSCMatrix();
276 
277  // matrix specific members
278  public:
283  RefSCMatrix operator*(const RefSCMatrix&) const;
284  RefSCMatrix operator*(const RefSymmSCMatrix&) const;
286  RefSCVector operator*(const RefSCVector&a) const;
287  RefSymmSCMatrix operator*(double) const;
290  RefSymmSCMatrix operator-(const RefSymmSCMatrix&) const;
292  RefSymmSCMatrix i() const;
294  RefSymmSCMatrix gi() const;
297  RefSymmSCMatrix clone() const;
298  RefSymmSCMatrix copy() const;
299  void set_element(int,int,double) const;
300  void accumulate_element(int,int,double) const;
301  double get_element(int,int) const;
302 
303  RefSCMatrix get_subblock(int br, int er, int bc, int ec);
304  RefSymmSCMatrix get_subblock(int br, int er);
305  void assign_subblock(const RefSCMatrix&, int br, int er, int bc, int ec);
306  void assign_subblock(const RefSymmSCMatrix&, int br, int er);
307  void accumulate_subblock(const RefSCMatrix&, int, int, int, int);
308  void accumulate_subblock(const RefSymmSCMatrix&, int, int);
309  RefSCVector get_row(int);
310  void assign_row(const RefSCVector&, int);
311  void accumulate_row(const RefSCVector&, int);
312 
313  void accumulate_symmetric_outer_product(const RefSCVector&) const;
314  double scalar_product(const RefSCVector&) const;
315  void accumulate_symmetric_product(const RefSCMatrix&) const;
316  void accumulate_symmetric_sum(const RefSCMatrix&) const;
318  void accumulate_transform(const RefSCMatrix&a,const RefSymmSCMatrix&b,
319  SCMatrix::Transform = SCMatrix::NormalTransform) const;
320  void accumulate_transform(const RefSCMatrix&a,const RefDiagSCMatrix&b,
321  SCMatrix::Transform = SCMatrix::NormalTransform) const;
323  const RefSymmSCMatrix&b) const;
324 
325  void randomize() const;
326  void assign(const RefSymmSCMatrix&) const;
327  void scale(double) const;
328  void assign(double) const;
329  void assign(const double*) const;
330  void assign(const double**) const;
331  void convert(double*) const;
332  void convert(double**) const;
333  void accumulate(const RefSymmSCMatrix&) const;
334  void element_op(const Ref<SCElementOp>&) const;
335  void element_op(const Ref<SCElementOp2>&,
336  const RefSymmSCMatrix&) const;
337  void element_op(const Ref<SCElementOp3>&,
338  const RefSymmSCMatrix&,
339  const RefSymmSCMatrix&) const;
340  double trace() const;
341  int n() const;
342  RefSCDimension dim() const;
343  Ref<SCMatrixKit> kit() const;
344  void print(std::ostream&) const;
345  void print(const char*title=0,
346  std::ostream&out=ExEnv::out0(), int =10) const;
347  void save(StateOut&);
349  void restore(StateIn&);
350 
352  double solve_lin(const RefSCVector&) const;
354  double determ() const;
356  RefDiagSCMatrix eigvals() const;
358  RefSCMatrix eigvecs() const;
362  void diagonalize(const RefDiagSCMatrix& eigvals,
363  const RefSCMatrix& eigvecs) const;
365  SymmSCMatrixdouble operator()(int i,int j) const;
369  int nblock() const;
373  RefSymmSCMatrix block(int i) const;
374 };
376 RefSymmSCMatrix operator*(double,const RefSymmSCMatrix&);
377 
380 class RefDiagSCMatrix: public Ref<DiagSCMatrix> {
381  // standard overrides
382  public:
385  RefDiagSCMatrix();
390  ~RefDiagSCMatrix();
395 
396  // matrix specific members
397  public:
402  RefSCMatrix operator*(const RefSCMatrix&) const;
403  RefDiagSCMatrix operator*(double) const;
406  RefDiagSCMatrix operator-(const RefDiagSCMatrix&) const;
408  RefDiagSCMatrix i() const;
410  RefDiagSCMatrix gi() const;
413  RefDiagSCMatrix clone() const;
414  RefDiagSCMatrix copy() const;
415  void set_element(int,double) const;
416  void accumulate_element(int,double) const;
417  double get_element(int) const;
418  void randomize() const;
419  void assign(const RefDiagSCMatrix&) const;
420  void scale(double) const;
421  void assign(double) const;
422  void assign(const double*) const;
423  void convert(double*) const;
424  void accumulate(const RefDiagSCMatrix&) const;
425  void element_op(const Ref<SCElementOp>&) const;
426  void element_op(const Ref<SCElementOp2>&,
427  const RefDiagSCMatrix&) const;
428  void element_op(const Ref<SCElementOp3>&,
429  const RefDiagSCMatrix&,
430  const RefDiagSCMatrix&) const;
431  int n() const;
432  RefSCDimension dim() const;
433  Ref<SCMatrixKit> kit() const;
434  double trace() const;
435  void print(std::ostream&) const;
436  void print(const char*title=0,
437  std::ostream&out=ExEnv::out0(), int =10) const;
438  void save(StateOut&);
440  void restore(StateIn&);
442  double determ() const;
444  DiagSCMatrixdouble operator()(int i) const;
448  int nblock() const;
452  RefDiagSCMatrix block(int i) const;
453 };
455 RefDiagSCMatrix operator*(double,const RefDiagSCMatrix&);
456 
458  friend class RefSCVector;
459  private:
460  RefSCVector vector;
461  int i;
462 
463  SCVectordouble(SCVector*,int);
464  public:
466  ~SCVectordouble();
467  double operator=(double a);
468  double operator=(const SCVectordouble&);
469  operator double();
470  double val() const;
471 };
472 
474  friend class RefSCMatrix;
475  private:
476  RefSCMatrix matrix;
477  int i;
478  int j;
479 
480  SCMatrixdouble(SCMatrix*,int,int);
481  public:
483  ~SCMatrixdouble();
484  double operator=(double a);
485  double operator=(const SCMatrixdouble&);
486  operator double();
487  double val() const;
488 };
489 
491  friend class RefSymmSCMatrix;
492  private:
493  RefSymmSCMatrix matrix;
494  int i;
495  int j;
496 
498  public:
501  double operator=(double a);
502  double operator=(const SymmSCMatrixdouble&);
503  operator double();
504  double val() const;
505 };
506 
508  friend class RefDiagSCMatrix;
509  private:
510  RefDiagSCMatrix matrix;
511  int i;
512  int j;
513 
515  public:
518  double operator=(double a);
519  double operator=(const DiagSCMatrixdouble&);
520  operator double();
521  double val() const;
522 };
523 
524 }
525 
526 #ifdef INLINE_FUNCTIONS
527 #include <math/scmat/matrix_i.h>
528 #endif
529 
530 #endif
531 
532 // Local Variables:
533 // mode: c++
534 // c-file-style: "CLJ"
535 // End:
SCVector & operator*() const
Returns a C++ reference to the reference counted object.
Definition: ref.h:390
void restore(StateIn &)
Restores the matrix from StateIn object. The vector must have been initialized already.
RefSCMatrix eigvecs() const
Returns the eigenvectors of the reference matrix.
RefDiagSCMatrix & operator=(DiagSCMatrix *m)
Make this refer to m.
RefDiagSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0...
RefSymmSCMatrix block(int i) const
If this matrix is blocked return block i.
RefSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0...
RefSymmSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0...
Definition: matrix.h:473
RefDiagSCMatrix()
Initializes the matrix pointer to 0.
RefSCVector operator-(const RefSCVector &a) const
Subtract two vectors.
The RefSCVector class is a smart pointer to an SCVector specialization.
Definition: matrix.h:55
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
Definition: matrix.h:507
Serializes objects that derive from SavableState.
Definition: stateout.h:61
int nblock() const
If this matrix is blocked return the number of blocks.
RefDiagSCMatrix i() const
Return the inverse of this.
RefSymmSCMatrix & operator=(SymmSCMatrix *m)
Make this refer to m.
RefSymmSCMatrix gi() const
Return the generalized inverse of this.
The SCMatrix class is the abstract base class for general double valued n by m matrices.
Definition: abstract.h:195
RefSymmSCMatrix()
Initializes the matrix pointer to 0.
void accumulate_transform(const RefSCMatrix &a, const RefSymmSCMatrix &b, SCMatrix::Transform=SCMatrix::NormalTransform) const
Add a * b * a.t() to this.
The RefDiagSCMatrix class is a smart pointer to an DiagSCMatrix specialization.
Definition: matrix.h:380
double determ() const
Returns the determinant of the referenced matrix.
double solve_lin(const RefSCVector &v) const
Solves this x = v.
SymmSCMatrixdouble operator()(int i, int j) const
Assign and examine matrix elements.
A template class that maintains references counts.
Definition: ref.h:332
RefSCMatrix operator-(const RefSCMatrix &) const
Matrix subtraction.
int nblock() const
If this matrix is blocked return the number of blocks.
SCVectordouble operator()(int) const
Return an l-value that can be used to assign or retrieve an element.
The RefSymmSCMatrix class is a smart pointer to an SCSymmSCMatrix specialization. ...
Definition: matrix.h:261
double solve_lin(const RefSCVector &) const
Solves this x = v.
void diagonalize(const RefDiagSCMatrix &eigvals, const RefSCMatrix &eigvecs) const
Sets eigvals to the eigenvalues and eigvecs to the eigenvalues and eigenvectors of the referenced mat...
The SymmSCMatrix class is the abstract base class for diagonal double valued matrices.
Definition: abstract.h:503
RefSymmSCMatrix operator+(const RefSymmSCMatrix &) const
Matrix addition and subtraction.
RefSCMatrix & operator=(SCMatrix *m)
Make this refer to m.
RefDiagSCMatrix block(int i) const
If this matrix is blocked return block i.
double determ() const
Returns the determinant of the referenced matrix.
SCMatrixdouble operator()(int i, int j) const
Assign and examine matrix elements.
Definition: mpqcin.h:13
RefSCVector()
Initializes the vector pointer to 0.
The RefSCDimension class is a smart pointer to an SCDimension specialization.
Definition: dim.h:156
Restores objects that derive from SavableState.
Definition: statein.h:70
The SCVector class is the abstract base class for double valued vectors.
Definition: abstract.h:97
RefSCVector & operator=(SCVector *v)
Make this refer to v.
DiagSCMatrixdouble operator()(int i) const
Assign and examine matrix elements.
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
RefSCMatrix t() const
Return the transpose of this.
RefSCMatrix()
Initializes the matrix pointer to 0.
static std::ostream & out0()
Return an ostream that writes from node 0.
The RefSCMatrix class is a smart pointer to an SCMatrix specialization.
Definition: matrix.h:135
The SymmSCMatrix class is the abstract base class for symmetric double valued matrices.
Definition: abstract.h:364
RefSymmSCMatrix symmetric_outer_product() const
The outer product of this with itself is a symmetric matrix.
void svd(const RefSCMatrix &U, const RefDiagSCMatrix &sigma, const RefSCMatrix &V)
Compute the singular value decomposition, this = U sigma V.t().
RefSCMatrix block(int i) const
If this matrix is blocked return block i.
RefDiagSCMatrix eigvals() const
Returns the eigenvalues of the reference matrix.
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
SCVectordouble operator[](int) const
Return an l-value that can be used to assign or retrieve an element.
RefSCMatrix i() const
Return the inverse of this.
Definition: matrix.h:457
RefSCMatrix operator+(const RefSCMatrix &) const
Matrix addition.
Definition: matrix.h:490
RefSymmSCMatrix i() const
Return the inverse of this.
RefSCVector operator+(const RefSCVector &a) const
Add two vectors.
RefDiagSCMatrix operator+(const RefDiagSCMatrix &) const
Matrix addition and subtraction.
double determ() const
Returns the determinant of the referenced matrix.
RefSCMatrix gi() const
Return the generalized inverse of this.
RefSCMatrix outer_product(const RefSCVector &v) const
Return the outer product between this and v.
int nblock() const
If this matrix is blocked return the number of blocks.
RefDiagSCMatrix gi() const
Return the generalized inverse of this.

Generated at Fri Feb 16 2018 01:48:56 for MPQC 2.3.1 using the documentation package Doxygen 1.8.14.