28#ifndef _math_scmat_matrix_h
29#define _math_scmat_matrix_h
36#include <math/scmat/abstract.h>
42class SymmSCMatrixdouble;
43class DiagSCMatrixdouble;
45class SCMatrixBlockIter;
46class SCMatrixRectBlock;
47class SCMatrixLTriBlock;
48class SCMatrixDiagBlock;
49class SCVectorSimpleBlock;
94 void set_element(
int i,
double val)
const;
95 void accumulate_element(
int i,
double val)
const;
96 double get_element(
int)
const;
102 double maxabs()
const;
105 void normalize()
const;
106 void randomize()
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;
121 void print(std::ostream&out)
const;
122 void print(
const char*title=0,
123 std::ostream&out=
ExEnv::out0(),
int precision=10)
const;
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);
195 void accumulate_row(
const RefSCVector&,
int)
const;
196 void accumulate_column(
const RefSCVector&,
int)
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;
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,
228 double trace()
const;
299 void set_element(
int,
int,
double)
const;
300 void accumulate_element(
int,
int,
double)
const;
301 double get_element(
int,
int)
const;
303 RefSCMatrix get_subblock(
int br,
int er,
int bc,
int ec);
305 void assign_subblock(
const RefSCMatrix&,
int br,
int er,
int bc,
int ec);
307 void accumulate_subblock(
const RefSCMatrix&,
int,
int,
int,
int);
313 void accumulate_symmetric_outer_product(
const RefSCVector&)
const;
315 void accumulate_symmetric_product(
const RefSCMatrix&)
const;
316 void accumulate_symmetric_sum(
const RefSCMatrix&)
const;
319 SCMatrix::Transform = SCMatrix::NormalTransform)
const;
321 SCMatrix::Transform = SCMatrix::NormalTransform)
const;
325 void randomize()
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;
340 double trace()
const;
344 void print(std::ostream&)
const;
345 void print(
const char*title=0,
415 void set_element(
int,
double)
const;
416 void accumulate_element(
int,
double)
const;
417 double get_element(
int)
const;
418 void randomize()
const;
420 void scale(
double)
const;
421 void assign(
double)
const;
422 void assign(
const double*)
const;
423 void convert(
double*)
const;
434 double trace()
const;
435 void print(std::ostream&)
const;
436 void print(
const char*title=0,
467 double operator=(
double a);
484 double operator=(
double a);
501 double operator=(
double a);
518 double operator=(
double a);
526#ifdef INLINE_FUNCTIONS
527#include <math/scmat/matrix_i.h>
The SymmSCMatrix class is the abstract base class for diagonal double valued matrices.
Definition abstract.h:503
static std::ostream & out0()
Return an ostream that writes from node 0.
The RefDiagSCMatrix class is a smart pointer to an DiagSCMatrix specialization.
Definition matrix.h:380
RefDiagSCMatrix()
Initializes the matrix pointer to 0.
RefDiagSCMatrix(const RefDiagSCMatrix &m)
Make this and m refer to the same SCMatrix.
RefDiagSCMatrix operator+(const RefDiagSCMatrix &) const
Matrix addition and subtraction.
RefSCMatrix operator*(const RefSCMatrix &) const
Multiply this by a matrix and return a matrix.
RefDiagSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0.
RefDiagSCMatrix i() const
Return the inverse of this.
RefDiagSCMatrix(const RefSCDimension &, const Ref< SCMatrixKit > &)
Create a diagonal matrix with dimension d by d.
RefDiagSCMatrix & operator=(DiagSCMatrix *m)
Make this refer to m.
int nblock() const
If this matrix is blocked return the number of blocks.
RefDiagSCMatrix & operator=(const RefDiagSCMatrix &m)
Make this and m refer to the same matrix.
RefDiagSCMatrix gi() const
Return the generalized inverse of this.
double determ() const
Returns the determinant of the referenced matrix.
RefDiagSCMatrix block(int i) const
If this matrix is blocked return block i.
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
RefDiagSCMatrix(DiagSCMatrix *m)
Make this refer to m.
DiagSCMatrixdouble operator()(int i) const
Assign and examine matrix elements.
The RefSCDimension class is a smart pointer to an SCDimension specialization.
Definition dim.h:156
The RefSCMatrix class is a smart pointer to an SCMatrix specialization.
Definition matrix.h:135
RefSCMatrix()
Initializes the matrix pointer to 0.
void svd(const RefSCMatrix &U, const RefDiagSCMatrix &sigma, const RefSCMatrix &V)
Compute the singular value decomposition, this = U sigma V.t().
RefSCMatrix(const RefSCMatrix &m)
Make this and m refer to the same SCMatrix.
double determ() const
Returns the determinant of the referenced matrix.
RefSCMatrix operator*(const RefSCMatrix &) const
Multiply this by a matrix and return a matrix.
RefSCMatrix operator-(const RefSCMatrix &) const
Matrix subtraction.
double solve_lin(const RefSCVector &v) const
Solves this x = v.
RefSCMatrix gi() const
Return the generalized inverse of this.
RefSCMatrix operator+(const RefSCMatrix &) const
Matrix addition.
RefSCMatrix i() const
Return the inverse of this.
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
SCMatrixdouble operator()(int i, int j) const
Assign and examine matrix elements.
RefSCMatrix t() const
Return the transpose of this.
RefSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0.
RefSCVector operator*(const RefSCVector &) const
Multiply this by a vector and return a vector.
RefSCMatrix operator*(double) const
Multiply this by a scalar and return the result.
RefSCMatrix block(int i) const
If this matrix is blocked return block i.
RefSCMatrix(SCMatrix *m)
Make this refer to m.
RefSCMatrix & operator=(const RefSCMatrix &m)
Make this and m refer to the same matrix.
int nblock() const
If this matrix is blocked return the number of blocks.
RefSCMatrix(const RefSCDimension &d1, const RefSCDimension &d2, const Ref< SCMatrixKit > &)
Create a vector with dimension d1 by d2.
RefSCMatrix & operator=(SCMatrix *m)
Make this refer to m.
The RefSCVector class is a smart pointer to an SCVector specialization.
Definition matrix.h:55
RefSCVector(const RefSCDimension &dim, const Ref< SCMatrixKit > &)
Create a vector with dimension dim.
SCVectordouble operator[](int) const
Return an l-value that can be used to assign or retrieve an element.
RefSCVector(const RefSCVector &v)
Make this and v refer to the same SCVector.
RefSCVector operator-(const RefSCVector &a) const
Subtract two vectors.
RefSCVector & operator=(const RefSCVector &v)
Make this and v refer to the same SCVector.
SCVectordouble operator()(int) const
Return an l-value that can be used to assign or retrieve an element.
RefSymmSCMatrix symmetric_outer_product() const
The outer product of this with itself is a symmetric matrix.
RefSCMatrix outer_product(const RefSCVector &v) const
Return the outer product between this and v.
RefSCVector(SCVector *v)
Make this refer to v.
void restore(StateIn &)
Restores the matrix from StateIn object. The vector must have been initialized already.
RefSCVector & operator=(SCVector *v)
Make this refer to v.
RefSCVector operator+(const RefSCVector &a) const
Add two vectors.
RefSCVector operator*(double) const
Scale a vector.
RefSCVector()
Initializes the vector pointer to 0.
The RefSymmSCMatrix class is a smart pointer to an SCSymmSCMatrix specialization.
Definition matrix.h:261
RefSymmSCMatrix i() const
Return the inverse of this.
RefSymmSCMatrix & operator=(const RefSymmSCMatrix &m)
Make this and m refer to the same matrix.
RefSymmSCMatrix operator+(const RefSymmSCMatrix &) const
Matrix addition and subtraction.
double solve_lin(const RefSCVector &) const
Solves this x = v.
RefSCMatrix eigvecs() const
Returns the eigenvectors of the reference matrix.
RefSymmSCMatrix(SymmSCMatrix *m)
Make this refer to m.
RefSymmSCMatrix(const RefSCDimension &d, const Ref< SCMatrixKit > &)
Create a vector with dimension d by d.
RefSCVector operator*(const RefSCVector &a) const
Multiply this by a vector and return a vector.
int nblock() const
If this matrix is blocked return the number of blocks.
SymmSCMatrixdouble operator()(int i, int j) const
Assign and examine matrix elements.
RefSymmSCMatrix & operator=(SymmSCMatrix *m)
Make this refer to m.
RefSymmSCMatrix()
Initializes the matrix pointer to 0.
RefSymmSCMatrix(const RefSymmSCMatrix &m)
Make this and m refer to the same SCMatrix.
RefSCMatrix operator*(const RefSCMatrix &) const
Multiply this by a matrix and return a matrix.
RefDiagSCMatrix eigvals() const
Returns the eigenvalues of the reference matrix.
void accumulate_transform(const RefSCMatrix &a, const RefSymmSCMatrix &b, SCMatrix::Transform=SCMatrix::NormalTransform) const
Add a * b * a.t() to this.
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...
RefSymmSCMatrix block(int i) const
If this matrix is blocked return block i.
void restore(StateIn &)
Restores the matrix from StateIn object. The matrix must have been initialized already.
RefSymmSCMatrix clone() const
These call the SCMatrix members of the same name after checking for references to 0.
double determ() const
Returns the determinant of the referenced matrix.
RefSymmSCMatrix gi() const
Return the generalized inverse of this.
A template class that maintains references counts.
Definition ref.h:332
SCMatrix & operator*() const
Returns a C++ reference to the reference counted object.
Definition ref.h:390
The SCMatrix class is the abstract base class for general double valued n by m matrices.
Definition abstract.h:195
The SCVector class is the abstract base class for double valued vectors.
Definition abstract.h:97
Restores objects that derive from SavableState.
Definition statein.h:70
Serializes objects that derive from SavableState.
Definition stateout.h:61
The SymmSCMatrix class is the abstract base class for symmetric double valued matrices.
Definition abstract.h:364