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 
33 #include "psi4/libmints/typedefs.h"
35 
36 #include <cstdio>
37 #include <string>
38 #include <vector>
39 #include <memory>
40 
41 namespace psi {
42 
43 struct dpdfile2;
44 
45 class PSIO;
46 class Vector;
47 class SimpleMatrix;
48 class Dimension;
49 class Molecule;
50 class Vector3;
51 
52 
55  ascending = 1,
58 };
59 
66 class Matrix : public std::enable_shared_from_this<Matrix> {
67 protected:
69  double ***matrix_;
71  int nirrep_;
77  std::string name_;
79  int symmetry_;
80 
82  void alloc();
84  void release();
85 
87  void copy_from(double ***);
88 
90  static double** matrix(int nrow, int ncol);
92  static void free(double** Block);
93 
94  void print_mat(const double *const *const a, int m, int n, std::string out) const;
95 
97  std::vector<int> numpy_shape_;
98 
99 public:
100 
102  Matrix();
108  Matrix(const std::string& name, int symmetry = 0);
110  Matrix(const Matrix& copy);
112  explicit Matrix(const SharedMatrix& copy);
114  explicit Matrix(const Matrix* copy);
122  Matrix(int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
131  Matrix(const std::string& name, int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
138  Matrix(int nirrep, int rows, const int *colspi);
139 
146  Matrix(int nirrep, const int* rowspi, int cols);
147 
156  Matrix(int rows, int cols);
166  Matrix(const std::string& name, int rows, int cols);
167 
173  Matrix(dpdfile2 *inFile);
174 
183  Matrix(const std::string& name, const Dimension& rows, const Dimension& cols, int symmetry = 0);
184 
192  Matrix(const Dimension& rows, const Dimension& cols, int symmetry = 0);
193 
195  virtual ~Matrix();
196 
206  void init(int nirrep, const int *rowspi, const int *colspi, const std::string& name = "", int symmetry = 0);
207 
208  void init(const Dimension& rowspi, const Dimension& colspi, const std::string& name = "", int symmetry = 0);
209 
211  SharedMatrix clone() const;
212 
216  static SharedMatrix create(const std::string& name,
217  const Dimension& rows,
218  const Dimension& cols);
219 
225  void copy(const SharedMatrix& cp);
226  void copy(const Matrix& cp);
227  void copy(const Matrix* cp);
234  static SharedMatrix horzcat(const std::vector<SharedMatrix >& mats);
235 
240  static SharedMatrix vertcat(const std::vector<SharedMatrix >& mats);
241 
252  SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn);
253 
255  void copy_to_row(int h, int row, double const * const data);
256 
257  enum SaveType {
261  };
262 
273  bool load(psi::PSIO* psio, size_t fileno, const std::string& tocentry, int nso);
274  bool load(std::shared_ptr<psi::PSIO>& psio, size_t fileno, const std::string& tocentry, int nso);
286  void load(psi::PSIO* const psio, size_t fileno, SaveType savetype=LowerTriangle);
287  void load(std::shared_ptr<psi::PSIO>& psio, size_t fileno, SaveType savetype=LowerTriangle);
295  void load(const std::string& filename);
296 
303  void load_mpqc(const std::string& filename);
304 
314  void save(const std::string& filename, bool append=true, bool saveLowerTriangle = true, bool saveSubBlocks=false);
325  void save(psi::PSIO* const psio, size_t fileno, SaveType savetype=LowerTriangle);
326  void save(std::shared_ptr<psi::PSIO>& psio, size_t fileno, SaveType savetype=LowerTriangle);
334  void set(double val);
335 
341  void set(const double * const tri);
342 
349  void set(const double * const * const sq);
359  void set(const double * const * const sq, int irrep);
368  void set(const SimpleMatrix * const sq);
369  void set(const std::shared_ptr<SimpleMatrix>& sq);
380  void set(int h, int m, int n, double val) { matrix_[h][m][n] = val; }
381 
389  void set(int m, int n, double val) { matrix_[0][m][n] = val; }
390 
397  void set_diagonal(const Vector * const vec);
398  void set_diagonal(const Vector& vec);
399  void set_diagonal(const std::shared_ptr<Vector>& vec);
410  double get(const int& h, const int& m, const int& n) const { return matrix_[h][m][n]; }
411 
419  double get(const int& m, const int& n) const { return matrix_[0][m][n]; }
420 
428  SharedVector get_row(int h, int m);
429 
437  SharedVector get_column(int h, int m);
438 
446  void set_row(int h, int m, SharedVector vec);
447 
455  void set_column(int h, int m, SharedVector vec);
456 
464  SharedMatrix get_block(const Slice& rows,const Slice& cols);
465 
473  void set_block(const Slice& rows,const Slice& cols,SharedMatrix block);
474 
487  double** pointer(const int& h = 0) const { return matrix_[h]; }
488  const double** const_pointer(const int& h=0) const { return const_cast<const double**>(matrix_[h]); }
489 
502  double* get_pointer(const int& h = 0) const {
503  if(rowspi_[h]*(size_t)colspi_[h] > 0)
504  return &(matrix_[h][0][0]);
505  else
506  return 0;}
507  const double* get_const_pointer(const int& h=0) const {
508  if(rowspi_[h]*(size_t)colspi_[h] > 0)
509  return const_cast<const double*>(&(matrix_[h][0][0]));
510  else
511  return 0;}
512 
513  size_t size(const int &h=0) const { return colspi_[h] * (size_t)rowspi_[h]; }
514 
516  void apply_denominator(const Matrix* const);
518  void apply_denominator(const Matrix&);
520  void apply_denominator(const SharedMatrix&);
521 
527  double **to_block_matrix() const;
539  double *to_lower_triangle() const;
540 
546  SimpleMatrix *to_simple_matrix() const;
547 
553  void set_name(const std::string& name) { name_ = name; }
554 
558  const std::string& name() const { return name_; }
559 
561  void print_out() const { print("outfile"); }
562 
569  void print(std::string outfile = "outfile", const char *extra=NULL) const;
570 
572  void print_atom_vector(std::string out_fname = "outfile");
573 
577  void print_to_mathematica();
578 
585  void eivprint(const Vector * const values, std::string out = "outfile");
587  void eivprint(const Vector& values, std::string out = "outfile");
589  void eivprint(const std::shared_ptr<Vector>& values, std::string out = "outfile");
590 
592  int rowdim(const int& h = 0) const { return rowspi_[h]; }
594  int coldim(const int& h = 0) const { return colspi_[h]; }
595 
597  const Dimension& rowspi() const {
598  return rowspi_;
599  }
601  int rowspi(const int& h) const {
602  return rowdim(h);
603  }
605  const Dimension& colspi() const {
606  return colspi_;
607  }
609  int colspi(const int& h) const {
610  return coldim(h);
611  }
613  const int& nirrep() const {
614  return nirrep_;
615  }
616 
618  int nrow() const {
619  int rows = 0;
620  for (int h=0; h<nirrep(); ++h)
621  rows += rowdim(h);
622  return rows;
623  }
624 
626  int ncol() const {
627  int cols = 0;
628  for (int h=0; h<nirrep(); ++h)
629  cols += coldim(h);
630  return cols;
631  }
632 
634  int max_nrow() const {
635  int row = 0;
636  for (int h=0; h<nirrep(); ++h)
637  if (row < rowdim(h))
638  row = rowdim(h);
639  return row;
640  }
641 
643  int max_ncol() const {
644  int col = 0;
645  for (int h=0; h<nirrep(); ++h)
646  if (col < coldim(h))
647  col = coldim(h);
648  return col;
649  }
650 
656  const int& symmetry() const {
657  return symmetry_;
658  }
659 
664  void symmetrize_gradient(std::shared_ptr<Molecule> mol);
665 
670  void symmetrize_hessian(std::shared_ptr<Molecule> mol);
671 
673  void identity();
675  void zero();
677  void zero_diagonal();
678 
679  // Math routines
681  double trace();
684 
686  void transpose_this();
687 
689  void add(const Matrix* const);
691  void add(const Matrix&);
693  void add(const SharedMatrix&);
694 
696  void subtract(const Matrix* const);
698  void subtract(const SharedMatrix&);
700  void accumulate_product(const Matrix* const, const Matrix* const);
701  void accumulate_product(const SharedMatrix&, const SharedMatrix&);
703  void scale(double);
705  double sum_of_squares();
707  double rms();
709  double absmax();
711  void add(int h, int m, int n, double val) {
712  #ifdef PSIDEBUG
713  if (m > rowspi_[h] || n > colspi_[h^symmetry_]) {
714  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
715  symmetry_, h, m, n);
716 
717  throw PSIEXCEPTION("What are you doing, Rob?");
718  }
719  #endif
720  matrix_[h][m][n] += val;
721  }
723  void add(int m, int n, double val) {
724  #ifdef PSIDEBUG
725  if (m > rowspi_[0] || n > colspi_[0^symmetry_]) {
726  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
727  symmetry_, 0, m, n);
728 
729  return;
730  }
731  #endif
732  matrix_[0][m][n] += val;
733  }
734 
736  for (int h=0; h<nirrep_; ++h) {
737  for (int i=0; i<rowspi_[h]; ++i) {
738  for (int j=0; j<i; ++j) {
739  matrix_[h][i][j] = matrix_[h][j][i] = (matrix_[h][i][j] + matrix_[h][j][i]);
740  }
741  }
742  }
743  }
744 
746  void scale_row(int h, int m, double a);
748  void scale_column(int h, int n, double a);
749 
756  void apply_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
757 
764  void remove_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
771  void transform(const SharedMatrix& L,
772  const SharedMatrix& F,
773  const SharedMatrix& R);
774 
777  void transform(const Matrix* const a, const Matrix* const transformer);
778  void transform(const SharedMatrix& a, const SharedMatrix& transformer);
780 
783  void transform(const Matrix* const transformer);
784  void transform(const SharedMatrix& transformer);
786 
789  void back_transform(const Matrix* const a, const Matrix* const transformer);
790  void back_transform(const SharedMatrix& a, const SharedMatrix& transformer);
792 
795  void back_transform(const Matrix* const transformer);
796  void back_transform(const SharedMatrix& transformer);
798 
800  double vector_dot(const Matrix* const rhs);
801  double vector_dot(const SharedMatrix& rhs);
802  double vector_dot(const Matrix& rhs);
803 
805 
813  void gemm(bool transa, bool transb, double alpha, const Matrix* const a, const Matrix* const b, double beta);
814  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const SharedMatrix& b, double beta);
815  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const Matrix& b, double beta);
816  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const SharedMatrix& b, double beta);
818 
820 
822  void gemm(const char& transa, const char& transb,
823  const std::vector<int>& m,
824  const std::vector<int>& n,
825  const std::vector<int>& k,
826  const double& alpha,
827  const SharedMatrix& a, const std::vector<int>& lda,
828  const SharedMatrix& b, const std::vector<int>& ldb,
829  const double& beta,
830  const std::vector<int>& ldc,
831  const std::vector<unsigned long>& offset_a = std::vector<unsigned long>(),
832  const std::vector<unsigned long>& offset_b = std::vector<unsigned long>(),
833  const std::vector<unsigned long>& offset_c = std::vector<unsigned long>());
834  void gemm(const char& transa, const char& transb,
835  const int& m,
836  const int& n,
837  const int& k,
838  const double& alpha,
839  const SharedMatrix& a, const int& lda,
840  const SharedMatrix& b, const int& ldb,
841  const double& beta,
842  const int& ldc,
843  const unsigned long& offset_a = 0,
844  const unsigned long& offset_b = 0,
845  const unsigned long& offset_c = 0);
847 
854  static SharedMatrix doublet(const SharedMatrix& A, const SharedMatrix& B, bool transA = false, bool transB = false);
855 
864  static SharedMatrix triplet(const SharedMatrix& A, const SharedMatrix& B, const SharedMatrix& C, bool transA = false, bool transB = false, bool transC = false);
865 
871  void axpy(double a, SharedMatrix X);
872 
878  SharedMatrix collapse(int dim = 0);
879 
882  void diagonalize(Matrix* eigvectors, Vector* eigvalues, diagonalize_order nMatz = ascending);
883  void diagonalize(SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
884  void diagonalize(SharedMatrix& eigvectors, Vector& eigvalues, diagonalize_order nMatz = ascending);
886 
889  void diagonalize(SharedMatrix& metric, SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
891 
894  void svd(SharedMatrix& U, SharedVector& S, SharedMatrix& V);
896 
903 
906  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_temps();
908 
911  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_a_temps();
913 
916  SharedMatrix pseudoinverse(double condition, int &nremoved);
918 
925 
952  SharedMatrix partial_cholesky_factorize(double delta = 0.0, bool throw_if_negative = false);
953 
968  std::pair<SharedMatrix, SharedMatrix> partial_square_root(double delta = 0.0);
969 
975  void cholesky_factorize();
976 
981  void invert();
982 
987  void general_invert();
988 
1010  Dimension power(double alpha, double cutoff = 1.0E-12);
1011 
1019  void expm(int n = 2, bool scale = false);
1020 
1022  void swap_rows(int h, int i, int j);
1024  void swap_columns(int h, int i, int j);
1025 
1027  void hermitivitize();
1029  void copy_lower_to_upper();
1031  void copy_upper_to_lower();
1033  void zero_lower();
1035  void zero_upper();
1037  void zero_row(int h, int i);
1039  void zero_column(int h, int i);
1040 
1041  // Reference versions of the above functions
1043  void transform(const Matrix& a, const Matrix& transformer);
1045  void transform(const Matrix& transformer);
1047  void back_transform(const Matrix& a, const Matrix& transformer);
1049  void back_transform(const Matrix& transformer);
1050 
1056 
1072  bool schmidt_add_row(int h, int rows, Vector& v) throw();
1073  bool schmidt_add_row(int h, int rows, double* v) throw();
1075 
1077  void schmidt();
1078 
1081  void schmidt_orthog(SharedMatrix S, int n);
1082 
1089  Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double*res=0);
1090 
1098  void project_out(Matrix& v);
1099 
1101  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const Matrix& b, double beta);
1103  void diagonalize(Matrix& eigvectors, Vector& eigvalues, int nMatz = 1);
1104 
1107  double** operator[](int i) { return matrix_[i]; }
1108  double& operator()(int i, int j) { return matrix_[0][i][j]; }
1109  const double& operator()(int i, int j) const { return matrix_[0][i][j]; }
1110  double& operator()(int h, int i, int j) { return matrix_[h][i][j]; }
1111  const double& operator()(int h, int i, int j) const { return matrix_[h][i][j]; }
1113 
1115  void write_to_dpdfile2(dpdfile2 *outFile);
1116 
1121  bool equal(const Matrix& rhs);
1122  bool equal(const SharedMatrix& rhs);
1123  bool equal(const Matrix* rhs);
1125 
1130  bool equal_but_for_row_order(const Matrix& rhs, double TOL=1.0e-10);
1131  bool equal_but_for_row_order(const SharedMatrix& rhs, double TOL=1.0e-10);
1132  bool equal_but_for_row_order(const Matrix* rhs, double TOL=1.0e-10);
1134 
1138  void set_numpy_shape(std::vector<int> shape) { numpy_shape_ = shape; }
1139  std::vector<int> numpy_shape() { return numpy_shape_; }
1140 
1148  void rotate_columns(int h, int i, int j, double theta);
1149  friend class Vector;
1150 };
1151 
1152 }
1153 
1154 
1155 #endif // MATRIX_H
void diagonalize(Matrix *eigvectors, Vector *eigvalues, diagonalize_order nMatz=ascending)
Definition: libmints/matrix.cc:1991
std::vector< int > numpy_shape_
Numpy Shape.
Definition: libmints/matrix.h:97
Definition: libmints/matrix.h:259
const double & operator()(int h, int i, int j) const
Definition: libmints/matrix.h:1111
void set_diagonal(const Vector *const vec)
Definition: libmints/matrix.cc:651
double vector_dot(const Matrix *const rhs)
Returns the vector dot product of this by rhs.
Definition: libmints/matrix.cc:1964
bool add_and_orthogonalize_row(const SharedVector v)
Definition: libmints/matrix.cc:1820
std::vector< int > numpy_shape()
Definition: libmints/matrix.h:1139
int rowdim(const int &h=0) const
Returns the rows in irrep h.
Definition: libmints/matrix.h:592
int colspi(const int &h) const
Returns the columns per irrep array.
Definition: libmints/matrix.h:609
int max_nrow() const
Returns the row size of the largest block.
Definition: libmints/matrix.h:634
int * U
Definition: stringlist.cc:66
void zero()
Zeros this out.
Definition: libmints/matrix.cc:1161
#define TOL
Definition: detci/opdm.cc:57
Matrix()
Default constructor, zeros everything out.
Definition: libmints/matrix.cc:75
void alloc()
Allocates matrix_.
Definition: libmints/matrix.cc:522
void set(int m, int n, double val)
Definition: libmints/matrix.h:389
size_t size(const int &h=0) const
Definition: libmints/matrix.h:513
void schmidt_orthog(SharedMatrix S, int n)
static SharedMatrix horzcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:359
void add(int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:723
virtual ~Matrix()
Destructor, frees memory.
Definition: libmints/matrix.cc:263
SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn)
Definition: libmints/matrix.cc:453
void back_transform(const Matrix *const a, const Matrix *const transformer)
Definition: libmints/matrix.cc:1460
void schmidt()
Definition: libmints/matrix.cc:1800
std::pair< SharedMatrix, SharedMatrix > partial_square_root(double delta=0.0)
Definition: libmints/matrix.cc:2438
void release()
Release matrix_.
Definition: libmints/matrix.cc:550
static double ** matrix(int nrow, int ncol)
allocate a block matrix – analogous to libciomr&#39;s block_matrix
Definition: libmints/matrix.cc:269
void set_name(const std::string &name)
Definition: libmints/matrix.h:553
bool schmidt_add_row(int h, int rows, Vector &v)
Definition: libmints/matrix.cc:1839
void eivprint(const Vector *const values, std::string out="outfile")
Definition: libmints/matrix.cc:1009
void transpose_this()
In place transposition.
Definition: libmints/matrix.cc:1238
void transform(const SharedMatrix &L, const SharedMatrix &F, const SharedMatrix &R)
Definition: libmints/matrix.cc:1444
void load_mpqc(const std::string &filename)
Definition: libmints/matrix.cc:3375
void hermitivitize()
Definition: libmints/matrix.cc:2887
SharedMatrix to_block_sharedmatrix() const
Definition: libmints/matrix.cc:873
void zero_lower()
Definition: libmints/matrix.cc:2781
Definition: vector3.h:37
Definition: pointgrp.h:106
int ncol() const
Returns the total number of columns.
Definition: libmints/matrix.h:626
void cholesky_factorize()
Definition: libmints/matrix.cc:2301
void write_to_dpdfile2(dpdfile2 *outFile)
Writes this to the dpdfile2 given.
Definition: libmints/matrix.cc:3076
static SharedMatrix create(const std::string &name, const Dimension &rows, const Dimension &cols)
Definition: libmints/matrix.cc:314
void swap_rows(int h, int i, int j)
Swap rows i and j.
Definition: libmints/matrix.cc:2291
const Dimension & colspi() const
Returns the columns per irrep array.
Definition: libmints/matrix.h:605
void print_atom_vector(std::string out_fname="outfile")
Prints the matrix with atom and xyz styling.
Definition: libmints/matrix.cc:989
void svd(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:2086
void identity()
Set this to identity.
Definition: libmints/matrix.cc:1142
void copy_from(double ***)
Copies data from the passed matrix to this matrix_.
Definition: libmints/matrix.cc:563
double absmax()
Returns the absoluate maximum balue.
Definition: libmints/matrix.cc:1393
void symmetrize_gradient(std::shared_ptr< Molecule > mol)
Definition: libmints/matrix.cc:1039
int nirrep_
Number of irreps.
Definition: libmints/matrix.h:71
void copy_upper_to_lower()
Definition: libmints/matrix.cc:2862
double & operator()(int i, int j)
Definition: libmints/matrix.h:1108
SharedMatrix collapse(int dim=0)
Definition: libmints/matrix.cc:1707
int symmetry_
Symmetry of this matrix (in most cases this will be 0 [totally symmetric])
Definition: libmints/matrix.h:79
double trace()
Returns the trace of this.
Definition: libmints/matrix.cc:1189
void zero_diagonal()
Zeros the diagonal.
Definition: libmints/matrix.cc:1175
void apply_symmetry(const SharedMatrix &a, const SharedMatrix &transformer)
Definition: libmints/matrix.cc:2919
void set(double val)
Definition: libmints/matrix.cc:573
void accumulate_product(const Matrix *const, const Matrix *const)
Multiplies the two arguments and adds their result to this.
Definition: libmints/matrix.cc:1329
double ** pointer(const int &h=0) const
Definition: libmints/matrix.h:487
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_a_temps()
Definition: libmints/matrix.cc:2072
void add(int h, int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:711
void gemm(bool transa, bool transb, double alpha, const Matrix *const a, const Matrix *const b, double beta)
Definition: libmints/matrix.cc:1565
Slicing for Matrices and Vectors objects.
Definition: dimension.h:120
void set_column(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:741
void svd_a(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:2132
SharedVector get_row(int h, int m)
Definition: libmints/matrix.cc:702
void zero_upper()
Definition: libmints/matrix.cc:2798
bool equal_but_for_row_order(const Matrix &rhs, double TOL=1.0e-10)
Definition: libmints/matrix.cc:3601
const double * get_const_pointer(const int &h=0) const
Definition: libmints/matrix.h:507
void subtract(const Matrix *const)
Subtracts a matrix from this.
Definition: libmints/matrix.cc:1288
Definition: libmints/matrix.h:57
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_temps()
Definition: libmints/matrix.cc:2056
SharedVector get_column(int h, int m)
Definition: libmints/matrix.cc:716
void general_invert()
Definition: libmints/matrix.cc:2509
void element_add_mirror()
Definition: libmints/matrix.h:735
double * get_pointer(const int &h=0) const
Definition: libmints/matrix.h:502
diagonalize_order
Definition: libmints/matrix.h:53
void zero_row(int h, int i)
Definition: libmints/matrix.cc:2815
double rms()
Returns the rms of this.
Definition: libmints/matrix.cc:1376
double *** matrix_
Matrix data.
Definition: libmints/matrix.h:69
int rowspi(const int &h) const
Returns the rows per irrep array.
Definition: libmints/matrix.h:601
long int size_t
Definition: sortintegrals.cc:40
void symmetrize_hessian(std::shared_ptr< Molecule > mol)
Definition: libmints/matrix.cc:1079
void expm(int n=2, bool scale=false)
Definition: libmints/matrix.cc:2617
bool load(psi::PSIO *psio, size_t fileno, const std::string &tocentry, int nso)
Definition: libmints/matrix.cc:3295
SharedMatrix transpose()
Creates a new matrix which is the transpose of this.
Definition: libmints/matrix.cc:1206
void print_out() const
Python compatible printer.
Definition: libmints/matrix.h:561
bool equal(const Matrix &rhs)
Definition: libmints/matrix.cc:3564
const double ** const_pointer(const int &h=0) const
Definition: libmints/matrix.h:488
double * to_lower_triangle() const
Definition: libmints/matrix.cc:825
SharedMatrix clone() const
Creates an exact copy of the matrix and returns it.
Definition: libmints/matrix.cc:321
Definition: dimension.h:38
Definition: libmints/matrix.h:54
Makes using matrices just a little easier.
Definition: libmints/matrix.h:66
SaveType
Definition: libmints/matrix.h:257
void print_mat(const double *const *const a, int m, int n, std::string out) const
Definition: libmints/matrix.cc:888
static void free(double **Block)
free a (block) matrix – analogous to libciomr&#39;s free_block
Definition: libmints/matrix.cc:280
const double & operator()(int i, int j) const
Definition: libmints/matrix.h:1109
void scale_row(int h, int m, double a)
Scale row m of irrep h by a.
Definition: libmints/matrix.cc:1351
const int & symmetry() const
Definition: libmints/matrix.h:656
Definition: libmints/matrix.h:258
Definition: libdpd/dpd.h:140
SharedMatrix get_block(const Slice &rows, const Slice &cols)
Definition: libmints/matrix.cc:752
void save(const std::string &filename, bool append=true, bool saveLowerTriangle=true, bool saveSubBlocks=false)
Definition: libmints/matrix.cc:3118
int nso
Definition: dx_write.cc:56
int coldim(const int &h=0) const
Returns the cols in irrep h.
Definition: libmints/matrix.h:594
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:504
void project_out(Matrix &v)
Definition: libmints/matrix.cc:1895
const Dimension & rowspi() const
Returns the rows per irrep array.
Definition: libmints/matrix.h:597
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:49
void rotate_columns(int h, int i, int j, double theta)
Definition: libmints/matrix.cc:3645
Definition: libmints/matrix.h:260
void zero_column(int h, int i)
Definition: libmints/matrix.cc:2826
void apply_denominator(const Matrix *const)
apply_denominators a matrix to this
Definition: libmints/matrix.cc:1303
int delta(const int i, const int j)
Definition: bend.cc:175
void swap_columns(int h, int i, int j)
Swap cols i and j.
Definition: libmints/matrix.cc:2296
std::string name_
Name of the matrix.
Definition: libmints/matrix.h:77
void set(int h, int m, int n, double val)
Definition: libmints/matrix.h:380
void scale(double)
Scales this matrix.
Definition: libmints/matrix.cc:1340
void init(int nirrep, const int *rowspi, const int *colspi, const std::string &name="", int symmetry=0)
Definition: libmints/matrix.cc:286
SharedMatrix partial_cholesky_factorize(double delta=0.0, bool throw_if_negative=false)
Definition: libmints/matrix.cc:2327
static SharedMatrix doublet(const SharedMatrix &A, const SharedMatrix &B, bool transA=false, bool transB=false)
Definition: libmints/matrix.cc:1632
void set_numpy_shape(std::vector< int > shape)
Definition: libmints/matrix.h:1138
Dimension colspi_
Columns per irrep array.
Definition: libmints/matrix.h:75
std::shared_ptr< PsiOutStream > outfile
Definition: core.cc:102
SharedMatrix canonical_orthogonalization(double delta=0.0, SharedMatrix eigvec=SharedMatrix())
Definition: libmints/matrix.cc:2241
double ** to_block_matrix() const
Definition: libmints/matrix.cc:842
void copy_lower_to_upper()
Definition: libmints/matrix.cc:2837
Definition: pointgrp.h:106
SharedMatrix pseudoinverse(double condition, int &nremoved)
Definition: libmints/matrix.cc:2198
void axpy(double a, SharedMatrix X)
Definition: libmints/matrix.cc:1688
Dimension rowspi_
Rows per irrep array.
Definition: libmints/matrix.h:73
void invert()
Definition: libmints/matrix.cc:2493
void set_row(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:730
void add(const Matrix *const)
Adds a matrix to this.
Definition: libmints/matrix.cc:1268
Definition: libmints/matrix.h:56
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:1681
double & operator()(int h, int i, int j)
Definition: libmints/matrix.h:1110
const int & nirrep() const
Returns the number of irreps.
Definition: libmints/matrix.h:613
Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double *res=0)
Definition: libmints/matrix.cc:1808
double sum_of_squares()
Returns the sum of the squares of this.
Definition: libmints/matrix.cc:1361
const std::string & name() const
Definition: libmints/matrix.h:558
#define PSIEXCEPTION(message)
Definition: exception.h:48
static SharedMatrix vertcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:406
Definition: vector.h:48
void print(std::string outfile="outfile", const char *extra=NULL) const
Definition: libmints/matrix.cc:937
void copy(const SharedMatrix &cp)
Definition: libmints/matrix.cc:517
void print_to_mathematica()
Definition: libmints/matrix.cc:959
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:2975
int nrow() const
Returns the total number of rows.
Definition: libmints/matrix.h:618
double ** operator[](int i)
Definition: libmints/matrix.h:1107
void set_block(const Slice &rows, const Slice &cols, SharedMatrix block)
Definition: libmints/matrix.cc:785
Dimension power(double alpha, double cutoff=1.0E-12)
Definition: libmints/matrix.cc:2556
int max_ncol() const
Returns the column size of the largest block.
Definition: libmints/matrix.h:643
Definition: psio.hpp:197
Definition: libmints/matrix.h:55
void scale_column(int h, int n, double a)
Scale column n of irrep h by a.
Definition: libmints/matrix.cc:1356