Psi4
libmints/matrix.h
Go to the documentation of this file.
1 /*
2  * @BEGIN LICENSE
3  *
4  * Psi4: an open-source quantum chemistry software package
5  *
6  * Copyright (c) 2007-2017 The Psi4 Developers.
7  *
8  * The copyrights for code used from other parties are included in
9  * the corresponding files.
10  *
11  * This file is part of Psi4.
12  *
13  * Psi4 is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU Lesser General Public License as published by
15  * the Free Software Foundation, version 3.
16  *
17  * Psi4 is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU Lesser General Public License for more details.
21  *
22  * You should have received a copy of the GNU Lesser General Public License along
23  * with Psi4; if not, write to the Free Software Foundation, Inc.,
24  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
25  *
26  * @END LICENSE
27  */
28 
29 #ifndef _psi_src_lib_libmints_matrix_h_
30 #define _psi_src_lib_libmints_matrix_h_
31 
35 #include "psi4/libmints/typedefs.h"
37 
38 #include "psi4/pybind11.h"
39 
40 #include <cstdio>
41 #include <string>
42 #include <vector>
43 #include <memory>
44 
45 namespace psi {
46 
47 struct dpdfile2;
48 
49 class PSIO;
50 class Vector;
51 class SimpleMatrix;
52 class Dimension;
53 class Molecule;
54 class Vector3;
55 
56 
59  ascending = 1,
62 };
63 
70 class Matrix : public std::enable_shared_from_this<Matrix> {
71 protected:
73  double ***matrix_;
75  int nirrep_;
81  std::string name_;
83  int symmetry_;
84 
86  void alloc();
88  void release();
89 
91  void copy_from(double ***);
92 
94  static double** matrix(int nrow, int ncol);
96  static void free(double** Block);
97 
98  void print_mat(const double *const *const a, int m, int n, std::string out) const;
99 
101  std::vector<int> numpy_shape_;
102 
103 public:
104 
106  Matrix();
112  Matrix(const std::string& name, int symmetry = 0);
114  Matrix(const Matrix& copy);
116  explicit Matrix(const SharedMatrix& copy);
118  explicit Matrix(const Matrix* copy);
126  Matrix(int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
135  Matrix(const std::string& name, int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
142  Matrix(int nirrep, int rows, const int *colspi);
143 
150  Matrix(int nirrep, const int* rowspi, int cols);
151 
160  Matrix(int rows, int cols);
170  Matrix(const std::string& name, int rows, int cols);
171 
177  Matrix(dpdfile2 *inFile);
178 
187  Matrix(const std::string& name, const Dimension& rows, const Dimension& cols, int symmetry = 0);
188 
196  Matrix(const Dimension& rows, const Dimension& cols, int symmetry = 0);
197 
199  virtual ~Matrix();
200 
210  void init(int nirrep, const int *rowspi, const int *colspi, const std::string& name = "", int symmetry = 0);
211 
212  void init(const Dimension& rowspi, const Dimension& colspi, const std::string& name = "", int symmetry = 0);
213 
215  SharedMatrix clone() const;
216 
220  static SharedMatrix create(const std::string& name,
221  const Dimension& rows,
222  const Dimension& cols);
223 
229  void copy(const SharedMatrix& cp);
230  void copy(const Matrix& cp);
231  void copy(const Matrix* cp);
238  static SharedMatrix horzcat(const std::vector<SharedMatrix >& mats);
239 
244  static SharedMatrix vertcat(const std::vector<SharedMatrix >& mats);
245 
256  SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn);
257 
259  void copy_to_row(int h, int row, double const * const data);
260 
261  enum SaveType {
265  };
266 
277  bool load(psi::PSIO* psio, unsigned int fileno, const std::string& tocentry, int nso);
278  bool load(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, const std::string& tocentry, int nso);
290  void load(psi::PSIO* const psio, unsigned int fileno, SaveType savetype=LowerTriangle);
291  void load(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, SaveType savetype=LowerTriangle);
299  void load(const std::string& filename);
300 
307  void load_mpqc(const std::string& filename);
308 
318  void save(const std::string& filename, bool append=true, bool saveLowerTriangle = true, bool saveSubBlocks=false);
329  void save(psi::PSIO* const psio, unsigned int fileno, SaveType savetype=LowerTriangle);
330  void save(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, SaveType savetype=LowerTriangle);
338  void set(double val);
339 
345  void set(const double * const tri);
346 
353  void set(const double * const * const sq);
363  void set(const double * const * const sq, int irrep);
372  void set(const SimpleMatrix * const sq);
373  void set(const std::shared_ptr<SimpleMatrix>& sq);
384  void set(int h, int m, int n, double val) { matrix_[h][m][n] = val; }
385 
393  void set(int m, int n, double val) { matrix_[0][m][n] = val; }
394 
401  void set_diagonal(const Vector * const vec);
402  void set_diagonal(const Vector& vec);
403  void set_diagonal(const std::shared_ptr<Vector>& vec);
414  double get(const int& h, const int& m, const int& n) const { return matrix_[h][m][n]; }
415 
423  double get(const int& m, const int& n) const { return matrix_[0][m][n]; }
424 
432  SharedVector get_row(int h, int m);
433 
441  SharedVector get_column(int h, int m);
442 
450  void set_row(int h, int m, SharedVector vec);
451 
459  void set_column(int h, int m, SharedVector vec);
460 
468  SharedMatrix get_block(const Slice& rows,const Slice& cols);
469 
477  void set_block(const Slice& rows,const Slice& cols,SharedMatrix block);
478 
482  double pyget(const py::tuple& key);
486  void pyset(const py::tuple& key, double value);
487 
500  double** pointer(const int& h = 0) const { return matrix_[h]; }
501  const double** const_pointer(const int& h=0) const { return const_cast<const double**>(matrix_[h]); }
502 
515  double* get_pointer(const int& h = 0) const {
516  if(rowspi_[h]*(size_t)colspi_[h] > 0)
517  return &(matrix_[h][0][0]);
518  else
519  return 0;}
520  const double* get_const_pointer(const int& h=0) const {
521  if(rowspi_[h]*(size_t)colspi_[h] > 0)
522  return const_cast<const double*>(&(matrix_[h][0][0]));
523  else
524  return 0;}
525 
526  size_t size(const int &h=0) const { return colspi_[h] * (size_t)rowspi_[h]; }
527 
529  void apply_denominator(const Matrix* const);
531  void apply_denominator(const Matrix&);
533  void apply_denominator(const SharedMatrix&);
534 
540  double **to_block_matrix() const;
552  double *to_lower_triangle() const;
553 
559  SimpleMatrix *to_simple_matrix() const;
560 
566  void set_name(const std::string& name) { name_ = name; }
567 
571  const std::string& name() const { return name_; }
572 
574  void print_out() const { print("outfile"); }
575 
582  void print(std::string outfile = "outfile", const char *extra=NULL) const;
583 
585  void print_atom_vector(std::string OutFileRMR = "outfile");
586 
590  void print_to_mathematica();
591 
598  void eivprint(const Vector * const values, std::string out = "outfile");
600  void eivprint(const Vector& values, std::string out = "outfile");
602  void eivprint(const std::shared_ptr<Vector>& values, std::string out = "outfile");
603 
605  int rowdim(const int& h = 0) const { return rowspi_[h]; }
607  int coldim(const int& h = 0) const { return colspi_[h]; }
608 
610  const Dimension& rowspi() const {
611  return rowspi_;
612  }
614  int rowspi(const int& h) const {
615  return rowdim(h);
616  }
618  const Dimension& colspi() const {
619  return colspi_;
620  }
622  int colspi(const int& h) const {
623  return coldim(h);
624  }
626  const int& nirrep() const {
627  return nirrep_;
628  }
629 
631  int nrow() const {
632  int rows = 0;
633  for (int h=0; h<nirrep(); ++h)
634  rows += rowdim(h);
635  return rows;
636  }
637 
639  int ncol() const {
640  int cols = 0;
641  for (int h=0; h<nirrep(); ++h)
642  cols += coldim(h);
643  return cols;
644  }
645 
647  int max_nrow() const {
648  int row = 0;
649  for (int h=0; h<nirrep(); ++h)
650  if (row < rowdim(h))
651  row = rowdim(h);
652  return row;
653  }
654 
656  int max_ncol() const {
657  int col = 0;
658  for (int h=0; h<nirrep(); ++h)
659  if (col < coldim(h))
660  col = coldim(h);
661  return col;
662  }
663 
669  const int& symmetry() const {
670  return symmetry_;
671  }
672 
677  void symmetrize_gradient(std::shared_ptr<Molecule> mol);
678 
683  void symmetrize_hessian(std::shared_ptr<Molecule> mol);
684 
686  void identity();
688  void zero();
690  void zero_diagonal();
691 
692  // Math routines
694  double trace();
697 
699  void transpose_this();
700 
702  void add(const Matrix* const);
704  void add(const Matrix&);
706  void add(const SharedMatrix&);
707 
709  void subtract(const Matrix* const);
711  void subtract(const SharedMatrix&);
713  void accumulate_product(const Matrix* const, const Matrix* const);
714  void accumulate_product(const SharedMatrix&, const SharedMatrix&);
716  void scale(double);
718  double sum_of_squares();
720  double rms();
722  double absmax();
724  void add(int h, int m, int n, double val) {
725  #ifdef PSIDEBUG
726  if (m > rowspi_[h] || n > colspi_[h^symmetry_]) {
727  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
728  symmetry_, h, m, n);
729 
730  throw PSIEXCEPTION("What are you doing, Rob?");
731  }
732  #endif
733  matrix_[h][m][n] += val;
734  }
736  void add(int m, int n, double val) {
737  #ifdef PSIDEBUG
738  if (m > rowspi_[0] || n > colspi_[0^symmetry_]) {
739  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
740  symmetry_, 0, m, n);
741 
742  return;
743  }
744  #endif
745  matrix_[0][m][n] += val;
746  }
747 
749  for (int h=0; h<nirrep_; ++h) {
750  for (int i=0; i<rowspi_[h]; ++i) {
751  for (int j=0; j<i; ++j) {
752  matrix_[h][i][j] = matrix_[h][j][i] = (matrix_[h][i][j] + matrix_[h][j][i]);
753  }
754  }
755  }
756  }
757 
759  void scale_row(int h, int m, double a);
761  void scale_column(int h, int n, double a);
762 
769  void apply_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
770 
777  void remove_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
784  void transform(const SharedMatrix& L,
785  const SharedMatrix& F,
786  const SharedMatrix& R);
787 
790  void transform(const Matrix* const a, const Matrix* const transformer);
791  void transform(const SharedMatrix& a, const SharedMatrix& transformer);
793 
796  void transform(const Matrix* const transformer);
797  void transform(const SharedMatrix& transformer);
799 
802  void back_transform(const Matrix* const a, const Matrix* const transformer);
803  void back_transform(const SharedMatrix& a, const SharedMatrix& transformer);
805 
808  void back_transform(const Matrix* const transformer);
809  void back_transform(const SharedMatrix& transformer);
811 
813  double vector_dot(const Matrix* const rhs);
814  double vector_dot(const SharedMatrix& rhs);
815  double vector_dot(const Matrix& rhs);
816 
818 
826  void gemm(bool transa, bool transb, double alpha, const Matrix* const a, const Matrix* const b, double beta);
827  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const SharedMatrix& b, double beta);
828  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const Matrix& b, double beta);
829  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const SharedMatrix& b, double beta);
831 
833 
835  void gemm(const char& transa, const char& transb,
836  const std::vector<int>& m,
837  const std::vector<int>& n,
838  const std::vector<int>& k,
839  const double& alpha,
840  const SharedMatrix& a, const std::vector<int>& lda,
841  const SharedMatrix& b, const std::vector<int>& ldb,
842  const double& beta,
843  const std::vector<int>& ldc,
844  const std::vector<unsigned long>& offset_a = std::vector<unsigned long>(),
845  const std::vector<unsigned long>& offset_b = std::vector<unsigned long>(),
846  const std::vector<unsigned long>& offset_c = std::vector<unsigned long>());
847  void gemm(const char& transa, const char& transb,
848  const int& m,
849  const int& n,
850  const int& k,
851  const double& alpha,
852  const SharedMatrix& a, const int& lda,
853  const SharedMatrix& b, const int& ldb,
854  const double& beta,
855  const int& ldc,
856  const unsigned long& offset_a = 0,
857  const unsigned long& offset_b = 0,
858  const unsigned long& offset_c = 0);
860 
867  static SharedMatrix doublet(const SharedMatrix& A, const SharedMatrix& B, bool transA = false, bool transB = false);
868 
877  static SharedMatrix triplet(const SharedMatrix& A, const SharedMatrix& B, const SharedMatrix& C, bool transA = false, bool transB = false, bool transC = false);
878 
884  void axpy(double a, SharedMatrix X);
885 
891  SharedMatrix collapse(int dim = 0);
892 
895  void diagonalize(Matrix* eigvectors, Vector* eigvalues, diagonalize_order nMatz = ascending);
896  void diagonalize(SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
897  void diagonalize(SharedMatrix& eigvectors, Vector& eigvalues, diagonalize_order nMatz = ascending);
899 
902  void diagonalize(SharedMatrix& metric, SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
904 
907  void svd(SharedMatrix& U, SharedVector& S, SharedMatrix& V);
909 
916 
919  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_temps();
921 
924  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_a_temps();
926 
929  SharedMatrix pseudoinverse(double condition, int &nremoved);
931 
938 
965  SharedMatrix partial_cholesky_factorize(double delta = 0.0, bool throw_if_negative = false);
966 
981  std::pair<SharedMatrix, SharedMatrix> partial_square_root(double delta = 0.0);
982 
988  void cholesky_factorize();
989 
994  void invert();
995 
1000  void general_invert();
1001 
1023  Dimension power(double alpha, double cutoff = 1.0E-12);
1024 
1032  void expm(int n = 2, bool scale = false);
1033 
1035  void swap_rows(int h, int i, int j);
1037  void swap_columns(int h, int i, int j);
1038 
1040  void hermitivitize();
1042  void copy_lower_to_upper();
1044  void copy_upper_to_lower();
1046  void zero_lower();
1048  void zero_upper();
1050  void zero_row(int h, int i);
1052  void zero_column(int h, int i);
1053 
1054  // Reference versions of the above functions
1056  void transform(const Matrix& a, const Matrix& transformer);
1058  void transform(const Matrix& transformer);
1060  void back_transform(const Matrix& a, const Matrix& transformer);
1062  void back_transform(const Matrix& transformer);
1063 
1069 
1085  bool schmidt_add_row(int h, int rows, Vector& v) throw();
1086  bool schmidt_add_row(int h, int rows, double* v) throw();
1088 
1090  void schmidt();
1091 
1094  void schmidt_orthog(SharedMatrix S, int n);
1095 
1102  Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double*res=0);
1103 
1111  void project_out(Matrix& v);
1112 
1114  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const Matrix& b, double beta);
1116  void diagonalize(Matrix& eigvectors, Vector& eigvalues, int nMatz = 1);
1117 
1120  double** operator[](int i) { return matrix_[i]; }
1121  double& operator()(int i, int j) { return matrix_[0][i][j]; }
1122  const double& operator()(int i, int j) const { return matrix_[0][i][j]; }
1123  double& operator()(int h, int i, int j) { return matrix_[h][i][j]; }
1124  const double& operator()(int h, int i, int j) const { return matrix_[h][i][j]; }
1126 
1127  // Serializable pure virtual functions:
1128  void send();
1129  void recv();
1130  void bcast(int broadcaster);
1134  void sum();
1135 
1137  void write_to_dpdfile2(dpdfile2 *outFile);
1138 
1143  bool equal(const Matrix& rhs);
1144  bool equal(const SharedMatrix& rhs);
1145  bool equal(const Matrix* rhs);
1147 
1152  bool equal_but_for_row_order(const Matrix& rhs, double TOL=1.0e-10);
1153  bool equal_but_for_row_order(const SharedMatrix& rhs, double TOL=1.0e-10);
1154  bool equal_but_for_row_order(const Matrix* rhs, double TOL=1.0e-10);
1156 
1161  void set_by_python_list(const py::list& data);
1162 
1166  void set_numpy_shape(std::vector<int> shape) { numpy_shape_ = shape; }
1167  std::vector<int> numpy_shape() { return numpy_shape_; }
1168  std::vector<py::buffer_info> array_interface();
1169 
1177  void rotate_columns(int h, int i, int j, double theta);
1178  friend class Vector;
1179 };
1180 
1181 }
1182 
1183 
1184 #endif // MATRIX_H
void send()
Definition: libmints/matrix.cc:3584
void diagonalize(Matrix *eigvectors, Vector *eigvalues, diagonalize_order nMatz=ascending)
Definition: libmints/matrix.cc:1989
std::vector< int > numpy_shape_
Numpy Shape.
Definition: libmints/matrix.h:101
Definition: libmints/matrix.h:263
const double & operator()(int h, int i, int j) const
Definition: libmints/matrix.h:1124
void set_diagonal(const Vector *const vec)
Definition: libmints/matrix.cc:649
double vector_dot(const Matrix *const rhs)
Returns the vector dot product of this by rhs.
Definition: libmints/matrix.cc:1962
bool add_and_orthogonalize_row(const SharedVector v)
Definition: libmints/matrix.cc:1818
std::vector< int > numpy_shape()
Definition: libmints/matrix.h:1167
int rowdim(const int &h=0) const
Returns the rows in irrep h.
Definition: libmints/matrix.h:605
int colspi(const int &h) const
Returns the columns per irrep array.
Definition: libmints/matrix.h:622
int max_nrow() const
Returns the row size of the largest block.
Definition: libmints/matrix.h:647
int * U
Definition: stringlist.cc:66
void zero()
Zeros this out.
Definition: libmints/matrix.cc:1159
#define TOL
Definition: detci/opdm.cc:53
Matrix()
Default constructor, zeros everything out.
Definition: libmints/matrix.cc:76
void alloc()
Allocates matrix_.
Definition: libmints/matrix.cc:520
bool load(psi::PSIO *psio, unsigned int fileno, const std::string &tocentry, int nso)
Definition: libmints/matrix.cc:3293
void set(int m, int n, double val)
Definition: libmints/matrix.h:393
double pyget(const py::tuple &key)
Definition: libmints/matrix.cc:3689
size_t size(const int &h=0) const
Definition: libmints/matrix.h:526
void schmidt_orthog(SharedMatrix S, int n)
static SharedMatrix horzcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:357
void add(int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:736
void set_by_python_list(const py::list &data)
Definition: libmints/matrix.cc:3704
virtual ~Matrix()
Destructor, frees memory.
Definition: libmints/matrix.cc:261
SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn)
Definition: libmints/matrix.cc:451
void back_transform(const Matrix *const a, const Matrix *const transformer)
Definition: libmints/matrix.cc:1458
void schmidt()
Definition: libmints/matrix.cc:1798
std::pair< SharedMatrix, SharedMatrix > partial_square_root(double delta=0.0)
Definition: libmints/matrix.cc:2436
void release()
Release matrix_.
Definition: libmints/matrix.cc:548
static double ** matrix(int nrow, int ncol)
allocate a block matrix – analogous to libciomr&#39;s block_matrix
Definition: libmints/matrix.cc:267
void pyset(const py::tuple &key, double value)
Definition: libmints/matrix.cc:3696
void set_name(const std::string &name)
Definition: libmints/matrix.h:566
bool schmidt_add_row(int h, int rows, Vector &v)
Definition: libmints/matrix.cc:1837
void eivprint(const Vector *const values, std::string out="outfile")
Definition: libmints/matrix.cc:1007
void transpose_this()
In place transposition.
Definition: libmints/matrix.cc:1236
void transform(const SharedMatrix &L, const SharedMatrix &F, const SharedMatrix &R)
Definition: libmints/matrix.cc:1442
void load_mpqc(const std::string &filename)
Definition: libmints/matrix.cc:3395
void hermitivitize()
Definition: libmints/matrix.cc:2885
SharedMatrix to_block_sharedmatrix() const
Definition: libmints/matrix.cc:871
void zero_lower()
Definition: libmints/matrix.cc:2779
Definition: vector3.h:37
Definition: pointgrp.h:106
int ncol() const
Returns the total number of columns.
Definition: libmints/matrix.h:639
std::vector< py::buffer_info > array_interface()
Definition: libmints/matrix.cc:3735
void cholesky_factorize()
Definition: libmints/matrix.cc:2299
void write_to_dpdfile2(dpdfile2 *outFile)
Writes this to the dpdfile2 given.
Definition: libmints/matrix.cc:3074
static SharedMatrix create(const std::string &name, const Dimension &rows, const Dimension &cols)
Definition: libmints/matrix.cc:312
void swap_rows(int h, int i, int j)
Swap rows i and j.
Definition: libmints/matrix.cc:2289
const Dimension & colspi() const
Returns the columns per irrep array.
Definition: libmints/matrix.h:618
void svd(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:2084
void identity()
Set this to identity.
Definition: libmints/matrix.cc:1140
void copy_from(double ***)
Copies data from the passed matrix to this matrix_.
Definition: libmints/matrix.cc:561
double absmax()
Returns the absoluate maximum balue.
Definition: libmints/matrix.cc:1391
void symmetrize_gradient(std::shared_ptr< Molecule > mol)
Definition: libmints/matrix.cc:1037
int nirrep_
Number of irreps.
Definition: libmints/matrix.h:75
void copy_upper_to_lower()
Definition: libmints/matrix.cc:2860
double & operator()(int i, int j)
Definition: libmints/matrix.h:1121
SharedMatrix collapse(int dim=0)
Definition: libmints/matrix.cc:1705
int symmetry_
Symmetry of this matrix (in most cases this will be 0 [totally symmetric])
Definition: libmints/matrix.h:83
double trace()
Returns the trace of this.
Definition: libmints/matrix.cc:1187
void zero_diagonal()
Zeros the diagonal.
Definition: libmints/matrix.cc:1173
void apply_symmetry(const SharedMatrix &a, const SharedMatrix &transformer)
Definition: libmints/matrix.cc:2917
void set(double val)
Definition: libmints/matrix.cc:571
void accumulate_product(const Matrix *const, const Matrix *const)
Multiplies the two arguments and adds their result to this.
Definition: libmints/matrix.cc:1327
double ** pointer(const int &h=0) const
Definition: libmints/matrix.h:500
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_a_temps()
Definition: libmints/matrix.cc:2070
void add(int h, int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:724
void gemm(bool transa, bool transb, double alpha, const Matrix *const a, const Matrix *const b, double beta)
Definition: libmints/matrix.cc:1563
Slicing for Matrices and Vectors objects.
Definition: dimension.h:120
void set_column(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:739
void svd_a(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:2130
SharedVector get_row(int h, int m)
Definition: libmints/matrix.cc:700
void zero_upper()
Definition: libmints/matrix.cc:2796
bool equal_but_for_row_order(const Matrix &rhs, double TOL=1.0e-10)
Definition: libmints/matrix.cc:3645
const double * get_const_pointer(const int &h=0) const
Definition: libmints/matrix.h:520
void subtract(const Matrix *const)
Subtracts a matrix from this.
Definition: libmints/matrix.cc:1286
Definition: libmints/matrix.h:61
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_temps()
Definition: libmints/matrix.cc:2054
SharedVector get_column(int h, int m)
Definition: libmints/matrix.cc:714
void general_invert()
Definition: libmints/matrix.cc:2507
void print_atom_vector(std::string OutFileRMR="outfile")
Prints the matrix with atom and xyz styling.
Definition: libmints/matrix.cc:987
void element_add_mirror()
Definition: libmints/matrix.h:748
double * get_pointer(const int &h=0) const
Definition: libmints/matrix.h:515
diagonalize_order
Definition: libmints/matrix.h:57
void zero_row(int h, int i)
Definition: libmints/matrix.cc:2813
double rms()
Returns the rms of this.
Definition: libmints/matrix.cc:1374
double *** matrix_
Matrix data.
Definition: libmints/matrix.h:73
int rowspi(const int &h) const
Returns the rows per irrep array.
Definition: libmints/matrix.h:614
void symmetrize_hessian(std::shared_ptr< Molecule > mol)
Definition: libmints/matrix.cc:1077
void expm(int n=2, bool scale=false)
Definition: libmints/matrix.cc:2615
SharedMatrix transpose()
Creates a new matrix which is the transpose of this.
Definition: libmints/matrix.cc:1204
void print_out() const
Python compatible printer.
Definition: libmints/matrix.h:574
bool equal(const Matrix &rhs)
Definition: libmints/matrix.cc:3608
const double ** const_pointer(const int &h=0) const
Definition: libmints/matrix.h:501
double * to_lower_triangle() const
Definition: libmints/matrix.cc:823
SharedMatrix clone() const
Creates an exact copy of the matrix and returns it.
Definition: libmints/matrix.cc:319
Definition: dimension.h:38
Definition: libmints/matrix.h:58
Makes using matrices just a little easier.
Definition: libmints/matrix.h:70
SaveType
Definition: libmints/matrix.h:261
void print_mat(const double *const *const a, int m, int n, std::string out) const
Definition: libmints/matrix.cc:886
static void free(double **Block)
free a (block) matrix – analogous to libciomr&#39;s free_block
Definition: libmints/matrix.cc:278
const double & operator()(int i, int j) const
Definition: libmints/matrix.h:1122
void recv()
Definition: libmints/matrix.cc:3588
void scale_row(int h, int m, double a)
Scale row m of irrep h by a.
Definition: libmints/matrix.cc:1349
const int & symmetry() const
Definition: libmints/matrix.h:669
Definition: libmints/matrix.h:262
Definition: libdpd/dpd.h:140
SharedMatrix get_block(const Slice &rows, const Slice &cols)
Definition: libmints/matrix.cc:750
void save(const std::string &filename, bool append=true, bool saveLowerTriangle=true, bool saveSubBlocks=false)
Definition: libmints/matrix.cc:3116
void sum()
Definition: libmints/matrix.cc:3602
int nso
Definition: dx_write.cc:56
int coldim(const int &h=0) const
Returns the cols in irrep h.
Definition: libmints/matrix.h:607
void copy_to_row(int h, int row, double const *const data)
Copies data to the row specified. Assumes data is of correct length.
Definition: libmints/matrix.cc:502
void project_out(Matrix &v)
Definition: libmints/matrix.cc:1893
const Dimension & rowspi() const
Returns the rows per irrep array.
Definition: libmints/matrix.h:610
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:49
void rotate_columns(int h, int i, int j, double theta)
Definition: libmints/matrix.cc:3722
Definition: libmints/matrix.h:264
void zero_column(int h, int i)
Definition: libmints/matrix.cc:2824
void apply_denominator(const Matrix *const)
apply_denominators a matrix to this
Definition: libmints/matrix.cc:1301
int delta(const int i, const int j)
Definition: bend.cc:176
void swap_columns(int h, int i, int j)
Swap cols i and j.
Definition: libmints/matrix.cc:2294
std::string name_
Name of the matrix.
Definition: libmints/matrix.h:81
void set(int h, int m, int n, double val)
Definition: libmints/matrix.h:384
void scale(double)
Scales this matrix.
Definition: libmints/matrix.cc:1338
void init(int nirrep, const int *rowspi, const int *colspi, const std::string &name="", int symmetry=0)
Definition: libmints/matrix.cc:284
SharedMatrix partial_cholesky_factorize(double delta=0.0, bool throw_if_negative=false)
Definition: libmints/matrix.cc:2325
static SharedMatrix doublet(const SharedMatrix &A, const SharedMatrix &B, bool transA=false, bool transB=false)
Definition: libmints/matrix.cc:1630
void set_numpy_shape(std::vector< int > shape)
Definition: libmints/matrix.h:1166
Dimension colspi_
Columns per irrep array.
Definition: libmints/matrix.h:79
std::shared_ptr< PsiOutStream > outfile
Definition: core.cc:101
SharedMatrix canonical_orthogonalization(double delta=0.0, SharedMatrix eigvec=SharedMatrix())
Definition: libmints/matrix.cc:2239
double ** to_block_matrix() const
Definition: libmints/matrix.cc:840
void copy_lower_to_upper()
Definition: libmints/matrix.cc:2835
Definition: pointgrp.h:106
SharedMatrix pseudoinverse(double condition, int &nremoved)
Definition: libmints/matrix.cc:2196
void axpy(double a, SharedMatrix X)
Definition: libmints/matrix.cc:1686
Dimension rowspi_
Rows per irrep array.
Definition: libmints/matrix.h:77
void invert()
Definition: libmints/matrix.cc:2491
Definition: PsiFileImpl.h:39
void set_row(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:728
void add(const Matrix *const)
Adds a matrix to this.
Definition: libmints/matrix.cc:1266
Definition: libmints/matrix.h:60
static SharedMatrix triplet(const SharedMatrix &A, const SharedMatrix &B, const SharedMatrix &C, bool transA=false, bool transB=false, bool transC=false)
Definition: libmints/matrix.cc:1679
double & operator()(int h, int i, int j)
Definition: libmints/matrix.h:1123
const int & nirrep() const
Returns the number of irreps.
Definition: libmints/matrix.h:626
Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double *res=0)
Definition: libmints/matrix.cc:1806
double sum_of_squares()
Returns the sum of the squares of this.
Definition: libmints/matrix.cc:1359
const std::string & name() const
Definition: libmints/matrix.h:571
#define PSIEXCEPTION(message)
Definition: exception.h:48
void bcast(int broadcaster)
Definition: libmints/matrix.cc:3592
static SharedMatrix vertcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:404
Definition: vector.h:50
void print(std::string outfile="outfile", const char *extra=NULL) const
Definition: libmints/matrix.cc:935
void copy(const SharedMatrix &cp)
Definition: libmints/matrix.cc:515
void print_to_mathematica()
Definition: libmints/matrix.cc:957
SimpleMatrix * to_simple_matrix() const
std::shared_ptr< Vector > SharedVector
Definition: adc.h:51
void remove_symmetry(const SharedMatrix &a, const SharedMatrix &transformer)
Definition: libmints/matrix.cc:2973
int nrow() const
Returns the total number of rows.
Definition: libmints/matrix.h:631
double ** operator[](int i)
Definition: libmints/matrix.h:1120
void set_block(const Slice &rows, const Slice &cols, SharedMatrix block)
Definition: libmints/matrix.cc:783
Dimension power(double alpha, double cutoff=1.0E-12)
Definition: libmints/matrix.cc:2554
int max_ncol() const
Returns the column size of the largest block.
Definition: libmints/matrix.h:656
Definition: psio.hpp:197
Definition: libmints/matrix.h:59
void scale_column(int h, int n, double a)
Scale column n of irrep h by a.
Definition: libmints/matrix.cc:1354