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 program is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License as published by
13  * the Free Software Foundation; either version 2 of the License, or
14  * (at your option) any later version.
15  *
16  * This program 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 General Public License for more details.
20  *
21  * You should have received a copy of the GNU General Public License along
22  * with this program; if not, write to the Free Software Foundation, Inc.,
23  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24  *
25  * @END LICENSE
26  */
27 
28 #ifndef _psi_src_lib_libmints_matrix_h_
29 #define _psi_src_lib_libmints_matrix_h_
30 
34 #include "psi4/libmints/typedefs.h"
36 
37 #include "psi4/pybind11.h"
38 
39 #include <cstdio>
40 #include <string>
41 #include <vector>
42 #include <memory>
43 
44 namespace psi {
45 
46 struct dpdfile2;
47 
48 class PSIO;
49 class Vector;
50 class SimpleVector;
51 class MatrixFactory;
52 class SimpleMatrix;
53 class Dimension;
54 class Molecule;
55 class Vector3;
56 
57 
60  ascending = 1,
63 };
64 
71 class Matrix : public std::enable_shared_from_this<Matrix> {
72 protected:
74  double ***matrix_;
76  int nirrep_;
82  std::string name_;
84  int symmetry_;
85 
87  void alloc();
89  void release();
90 
92  void copy_from(double ***);
93 
95  static double** matrix(int nrow, int ncol);
97  static void free(double** Block);
98 
99  void print_mat(const double *const *const a, int m, int n, std::string out) const;
100 
102  std::vector<int> numpy_shape_;
103 
104 public:
105 
107  Matrix();
113  Matrix(const std::string& name, int symmetry = 0);
115  Matrix(const Matrix& copy);
117  explicit Matrix(const SharedMatrix& copy);
119  explicit Matrix(const Matrix* copy);
127  Matrix(int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
136  Matrix(const std::string& name, int nirrep, const int *rowspi, const int *colspi, int symmetry = 0);
143  Matrix(int nirrep, int rows, const int *colspi);
144 
151  Matrix(int nirrep, const int* rowspi, int cols);
152 
161  Matrix(int rows, int cols);
171  Matrix(const std::string& name, int rows, int cols);
172 
178  Matrix(dpdfile2 *inFile);
179 
188  Matrix(const std::string& name, const Dimension& rows, const Dimension& cols, int symmetry = 0);
189 
197  Matrix(const Dimension& rows, const Dimension& cols, int symmetry = 0);
198 
200  virtual ~Matrix();
201 
211  void init(int nirrep, const int *rowspi, const int *colspi, const std::string& name = "", int symmetry = 0);
212 
213  void init(const Dimension& rowspi, const Dimension& colspi, const std::string& name = "", int symmetry = 0);
214 
216  SharedMatrix clone() const;
217 
221  static SharedMatrix create(const std::string& name,
222  const Dimension& rows,
223  const Dimension& cols);
224 
230  void copy(const SharedMatrix& cp);
231  void copy(const Matrix& cp);
232  void copy(const Matrix* cp);
239  static SharedMatrix horzcat(const std::vector<SharedMatrix >& mats);
240 
245  static SharedMatrix vertcat(const std::vector<SharedMatrix >& mats);
246 
257  SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn);
258 
260  void copy_to_row(int h, int row, double const * const data);
261 
262  enum SaveType {
266  };
267 
278  bool load(psi::PSIO* psio, unsigned int fileno, const std::string& tocentry, int nso);
279  bool load(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, const std::string& tocentry, int nso);
291  void load(psi::PSIO* const psio, unsigned int fileno, SaveType savetype=LowerTriangle);
292  void load(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, SaveType savetype=LowerTriangle);
300  void load(const std::string& filename);
301 
308  void load_mpqc(const std::string& filename);
309 
319  void save(const std::string& filename, bool append=true, bool saveLowerTriangle = true, bool saveSubBlocks=false);
330  void save(psi::PSIO* const psio, unsigned int fileno, SaveType savetype=LowerTriangle);
331  void save(std::shared_ptr<psi::PSIO>& psio, unsigned int fileno, SaveType savetype=LowerTriangle);
339  void set(double val);
340 
346  void set(const double * const tri);
347 
354  void set(const double * const * const sq);
364  void set(const double * const * const sq, int irrep);
373  void set(const SimpleMatrix * const sq);
374  void set(const std::shared_ptr<SimpleMatrix>& sq);
385  void set(int h, int m, int n, double val) { matrix_[h][m][n] = val; }
386 
394  void set(int m, int n, double val) { matrix_[0][m][n] = val; }
395 
402  void set_diagonal(const Vector * const vec);
403  void set_diagonal(const Vector& vec);
404  void set_diagonal(const std::shared_ptr<Vector>& vec);
415  double get(const int& h, const int& m, const int& n) const { return matrix_[h][m][n]; }
416 
424  double get(const int& m, const int& n) const { return matrix_[0][m][n]; }
425 
433  SharedVector get_row(int h, int m);
434 
442  SharedVector get_column(int h, int m);
443 
451  void set_row(int h, int m, SharedVector vec);
452 
460  void set_column(int h, int m, SharedVector vec);
461 
465  double pyget(const py::tuple& key);
469  void pyset(const py::tuple& key, double value);
470 
483  double** pointer(const int& h = 0) const { return matrix_[h]; }
484  const double** const_pointer(const int& h=0) const { return const_cast<const double**>(matrix_[h]); }
485 
498  double* get_pointer(const int& h = 0) const {
499  if(rowspi_[h]*(size_t)colspi_[h] > 0)
500  return &(matrix_[h][0][0]);
501  else
502  return 0;}
503  const double* get_const_pointer(const int& h=0) const {
504  if(rowspi_[h]*(size_t)colspi_[h] > 0)
505  return const_cast<const double*>(&(matrix_[h][0][0]));
506  else
507  return 0;}
508 
509  size_t size(const int &h=0) const { return colspi_[h] * (size_t)rowspi_[h]; }
510 
512  void apply_denominator(const Matrix* const);
514  void apply_denominator(const Matrix&);
516  void apply_denominator(const SharedMatrix&);
517 
523  double **to_block_matrix() const;
535  double *to_lower_triangle() const;
536 
542  SimpleMatrix *to_simple_matrix() const;
543 
549  void set_name(const std::string& name) { name_ = name; }
550 
554  const std::string& name() const { return name_; }
555 
557  void print_out() const { print("outfile"); }
558 
565  void print(std::string outfile = "outfile", const char *extra=NULL) const;
566 
568  void print_atom_vector(std::string OutFileRMR = "outfile");
569 
573  void print_to_mathematica();
574 
581  void eivprint(const Vector * const values, std::string out = "outfile");
583  void eivprint(const Vector& values, std::string out = "outfile");
585  void eivprint(const std::shared_ptr<Vector>& values, std::string out = "outfile");
586 
588  int rowdim(const int& h = 0) const { return rowspi_[h]; }
590  int coldim(const int& h = 0) const { return colspi_[h]; }
591 
593  const Dimension& rowspi() const {
594  return rowspi_;
595  }
597  int rowspi(const int& h) const {
598  return rowdim(h);
599  }
601  const Dimension& colspi() const {
602  return colspi_;
603  }
605  int colspi(const int& h) const {
606  return coldim(h);
607  }
609  const int& nirrep() const {
610  return nirrep_;
611  }
612 
614  int nrow() const {
615  int rows = 0;
616  for (int h=0; h<nirrep(); ++h)
617  rows += rowdim(h);
618  return rows;
619  }
620 
622  int ncol() const {
623  int cols = 0;
624  for (int h=0; h<nirrep(); ++h)
625  cols += coldim(h);
626  return cols;
627  }
628 
630  int max_nrow() const {
631  int row = 0;
632  for (int h=0; h<nirrep(); ++h)
633  if (row < rowdim(h))
634  row = rowdim(h);
635  return row;
636  }
637 
639  int max_ncol() const {
640  int col = 0;
641  for (int h=0; h<nirrep(); ++h)
642  if (col < coldim(h))
643  col = coldim(h);
644  return col;
645  }
646 
652  const int& symmetry() const {
653  return symmetry_;
654  }
655 
660  void symmetrize_gradient(std::shared_ptr<Molecule> mol);
661 
663  void identity();
665  void zero();
667  void zero_diagonal();
668 
669  // Math routines
671  double trace();
674 
676  void transpose_this();
677 
679  void add(const Matrix* const);
681  void add(const Matrix&);
683  void add(const SharedMatrix&);
684 
686  void subtract(const Matrix* const);
688  void subtract(const SharedMatrix&);
690  void accumulate_product(const Matrix* const, const Matrix* const);
691  void accumulate_product(const SharedMatrix&, const SharedMatrix&);
693  void scale(double);
695  double sum_of_squares();
697  double rms();
699  double absmax();
701  void add(int h, int m, int n, double val) {
702  #ifdef PSIDEBUG
703  if (m > rowspi_[h] || n > colspi_[h^symmetry_]) {
704  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
705  symmetry_, h, m, n);
706 
707  throw PSIEXCEPTION("What are you doing, Rob?");
708  }
709  #endif
710  matrix_[h][m][n] += val;
711  }
713  void add(int m, int n, double val) {
714  #ifdef PSIDEBUG
715  if (m > rowspi_[0] || n > colspi_[0^symmetry_]) {
716  outfile->Printf( "out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n",
717  symmetry_, 0, m, n);
718 
719  return;
720  }
721  #endif
722  matrix_[0][m][n] += val;
723  }
724 
726  for (int h=0; h<nirrep_; ++h) {
727  for (int i=0; i<rowspi_[h]; ++i) {
728  for (int j=0; j<i; ++j) {
729  matrix_[h][i][j] = matrix_[h][j][i] = (matrix_[h][i][j] + matrix_[h][j][i]);
730  }
731  }
732  }
733  }
734 
736  void scale_row(int h, int m, double a);
738  void scale_column(int h, int n, double a);
739 
746  void apply_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
747 
754  void remove_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
761  void transform(const SharedMatrix& L,
762  const SharedMatrix& F,
763  const SharedMatrix& R);
764 
767  void transform(const Matrix* const a, const Matrix* const transformer);
768  void transform(const SharedMatrix& a, const SharedMatrix& transformer);
770 
773  void transform(const Matrix* const transformer);
774  void transform(const SharedMatrix& transformer);
776 
779  void back_transform(const Matrix* const a, const Matrix* const transformer);
780  void back_transform(const SharedMatrix& a, const SharedMatrix& transformer);
782 
785  void back_transform(const Matrix* const transformer);
786  void back_transform(const SharedMatrix& transformer);
788 
790  double vector_dot(const Matrix* const rhs);
791  double vector_dot(const SharedMatrix& rhs);
792  double vector_dot(const Matrix& rhs);
793 
795 
803  void gemm(bool transa, bool transb, double alpha, const Matrix* const a, const Matrix* const b, double beta);
804  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const SharedMatrix& b, double beta);
805  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const Matrix& b, double beta);
806  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const SharedMatrix& b, double beta);
808 
810 
812  void gemm(const char& transa, const char& transb,
813  const std::vector<int>& m,
814  const std::vector<int>& n,
815  const std::vector<int>& k,
816  const double& alpha,
817  const SharedMatrix& a, const std::vector<int>& lda,
818  const SharedMatrix& b, const std::vector<int>& ldb,
819  const double& beta,
820  const std::vector<int>& ldc,
821  const std::vector<unsigned long>& offset_a = std::vector<unsigned long>(),
822  const std::vector<unsigned long>& offset_b = std::vector<unsigned long>(),
823  const std::vector<unsigned long>& offset_c = std::vector<unsigned long>());
824  void gemm(const char& transa, const char& transb,
825  const int& m,
826  const int& n,
827  const int& k,
828  const double& alpha,
829  const SharedMatrix& a, const int& lda,
830  const SharedMatrix& b, const int& ldb,
831  const double& beta,
832  const int& ldc,
833  const unsigned long& offset_a = 0,
834  const unsigned long& offset_b = 0,
835  const unsigned long& offset_c = 0);
837 
844  static SharedMatrix doublet(const SharedMatrix& A, const SharedMatrix& B, bool transA = false, bool transB = false);
845 
854  static SharedMatrix triplet(const SharedMatrix& A, const SharedMatrix& B, const SharedMatrix& C, bool transA = false, bool transB = false, bool transC = false);
855 
861  void axpy(double a, SharedMatrix X);
862 
868  SharedMatrix collapse(int dim = 0);
869 
872  void diagonalize(Matrix* eigvectors, Vector* eigvalues, diagonalize_order nMatz = ascending);
873  void diagonalize(SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
874  void diagonalize(SharedMatrix& eigvectors, Vector& eigvalues, diagonalize_order nMatz = ascending);
876 
879  void diagonalize(SharedMatrix& metric, SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
881 
884  void svd(SharedMatrix& U, SharedVector& S, SharedMatrix& V);
886 
893 
896  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_temps();
898 
901  std::tuple<SharedMatrix,SharedVector,SharedMatrix> svd_a_temps();
903 
906  SharedMatrix pseudoinverse(double condition = 0.0, bool* conditioned = NULL);
908 
915 
942  SharedMatrix partial_cholesky_factorize(double delta = 0.0, bool throw_if_negative = false);
943 
958  std::pair<SharedMatrix, SharedMatrix> partial_square_root(double delta = 0.0);
959 
965  void cholesky_factorize();
966 
971  void invert();
972 
977  void general_invert();
978 
1000  Dimension power(double alpha, double cutoff = 1.0E-12);
1001 
1009  void expm(int n = 2, bool scale = false);
1010 
1012  void swap_rows(int h, int i, int j);
1014  void swap_columns(int h, int i, int j);
1015 
1017  void hermitivitize();
1019  void copy_lower_to_upper();
1021  void copy_upper_to_lower();
1023  void zero_lower();
1025  void zero_upper();
1027  void zero_row(int h, int i);
1029  void zero_column(int h, int i);
1030 
1031  // Reference versions of the above functions
1033  void transform(const Matrix& a, const Matrix& transformer);
1035  void transform(const Matrix& transformer);
1037  void back_transform(const Matrix& a, const Matrix& transformer);
1039  void back_transform(const Matrix& transformer);
1040 
1046 
1062  bool schmidt_add_row(int h, int rows, Vector& v) throw();
1063  bool schmidt_add_row(int h, int rows, double* v) throw();
1065 
1067  void schmidt();
1068 
1071  void schmidt_orthog(SharedMatrix S, int n);
1072 
1079  Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double*res=0);
1080 
1088  void project_out(Matrix& v);
1089 
1091  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const Matrix& b, double beta);
1093  void diagonalize(Matrix& eigvectors, Vector& eigvalues, int nMatz = 1);
1094 
1097  double** operator[](int i) { return matrix_[i]; }
1098  double& operator()(int i, int j) { return matrix_[0][i][j]; }
1099  const double& operator()(int i, int j) const { return matrix_[0][i][j]; }
1100  double& operator()(int h, int i, int j) { return matrix_[h][i][j]; }
1101  const double& operator()(int h, int i, int j) const { return matrix_[h][i][j]; }
1103 
1104  // Serializable pure virtual functions:
1105  void send();
1106  void recv();
1107  void bcast(int broadcaster);
1111  void sum();
1112 
1114  void write_to_dpdfile2(dpdfile2 *outFile);
1115 
1120  bool equal(const Matrix& rhs);
1121  bool equal(const SharedMatrix& rhs);
1122  bool equal(const Matrix* rhs);
1124 
1129  bool equal_but_for_row_order(const Matrix& rhs, double TOL=1.0e-10);
1130  bool equal_but_for_row_order(const SharedMatrix& rhs, double TOL=1.0e-10);
1131  bool equal_but_for_row_order(const Matrix* rhs, double TOL=1.0e-10);
1133 
1138  void set_by_python_list(const py::list& data);
1139 
1143  void set_numpy_shape(std::vector<int> shape) { numpy_shape_ = shape; }
1144  std::vector<int> numpy_shape() { return numpy_shape_; }
1145  std::vector<py::buffer_info> array_interface();
1146 
1154  void rotate_columns(int h, int i, int j, double theta);
1155  friend class Vector;
1156 };
1157 
1158 }
1159 
1160 
1161 #endif // MATRIX_H
void send()
Definition: libmints/matrix.cc:3453
void diagonalize(Matrix *eigvectors, Vector *eigvalues, diagonalize_order nMatz=ascending)
Definition: libmints/matrix.cc:1853
std::vector< int > numpy_shape_
Numpy Shape.
Definition: libmints/matrix.h:102
Definition: libmints/matrix.h:264
const double & operator()(int h, int i, int j) const
Definition: libmints/matrix.h:1101
void set_diagonal(const Vector *const vec)
Definition: libmints/matrix.cc:648
double vector_dot(const Matrix *const rhs)
Returns the vector dot product of this by rhs.
Definition: libmints/matrix.cc:1826
bool add_and_orthogonalize_row(const SharedVector v)
Definition: libmints/matrix.cc:1682
std::vector< int > numpy_shape()
Definition: libmints/matrix.h:1144
int rowdim(const int &h=0) const
Returns the rows in irrep h.
Definition: libmints/matrix.h:588
int colspi(const int &h) const
Returns the columns per irrep array.
Definition: libmints/matrix.h:605
int max_nrow() const
Returns the row size of the largest block.
Definition: libmints/matrix.h:630
int * U
Definition: stringlist.cc:65
void zero()
Zeros this out.
Definition: libmints/matrix.cc:1023
#define TOL
Definition: detci/opdm.cc:52
Matrix()
Default constructor, zeros everything out.
Definition: libmints/matrix.cc:75
void alloc()
Allocates matrix_.
Definition: libmints/matrix.cc:519
bool load(psi::PSIO *psio, unsigned int fileno, const std::string &tocentry, int nso)
Definition: libmints/matrix.cc:3162
void set(int m, int n, double val)
Definition: libmints/matrix.h:394
double pyget(const py::tuple &key)
Definition: libmints/matrix.cc:3558
size_t size(const int &h=0) const
Definition: libmints/matrix.h:509
void schmidt_orthog(SharedMatrix S, int n)
static SharedMatrix horzcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:356
void add(int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:713
void set_by_python_list(const py::list &data)
Definition: libmints/matrix.cc:3573
virtual ~Matrix()
Destructor, frees memory.
Definition: libmints/matrix.cc:260
SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn)
Definition: libmints/matrix.cc:450
void back_transform(const Matrix *const a, const Matrix *const transformer)
Definition: libmints/matrix.cc:1322
void schmidt()
Definition: libmints/matrix.cc:1662
std::pair< SharedMatrix, SharedMatrix > partial_square_root(double delta=0.0)
Definition: libmints/matrix.cc:2305
void release()
Release matrix_.
Definition: libmints/matrix.cc:547
static double ** matrix(int nrow, int ncol)
allocate a block matrix – analogous to libciomr&#39;s block_matrix
Definition: libmints/matrix.cc:266
void pyset(const py::tuple &key, double value)
Definition: libmints/matrix.cc:3565
void set_name(const std::string &name)
Definition: libmints/matrix.h:549
bool schmidt_add_row(int h, int rows, Vector &v)
Definition: libmints/matrix.cc:1701
void eivprint(const Vector *const values, std::string out="outfile")
Definition: libmints/matrix.cc:933
void transpose_this()
In place transposition.
Definition: libmints/matrix.cc:1100
void transform(const SharedMatrix &L, const SharedMatrix &F, const SharedMatrix &R)
Definition: libmints/matrix.cc:1306
void load_mpqc(const std::string &filename)
Definition: libmints/matrix.cc:3264
void hermitivitize()
Definition: libmints/matrix.cc:2754
SharedMatrix to_block_sharedmatrix() const
Definition: libmints/matrix.cc:797
void zero_lower()
Definition: libmints/matrix.cc:2648
Definition: vector3.h:36
Definition: pointgrp.h:105
int ncol() const
Returns the total number of columns.
Definition: libmints/matrix.h:622
std::vector< py::buffer_info > array_interface()
Definition: libmints/matrix.cc:3604
void cholesky_factorize()
Definition: libmints/matrix.cc:2168
void write_to_dpdfile2(dpdfile2 *outFile)
Writes this to the dpdfile2 given.
Definition: libmints/matrix.cc:2943
static SharedMatrix create(const std::string &name, const Dimension &rows, const Dimension &cols)
Definition: libmints/matrix.cc:311
void swap_rows(int h, int i, int j)
Swap rows i and j.
Definition: libmints/matrix.cc:2158
const Dimension & colspi() const
Returns the columns per irrep array.
Definition: libmints/matrix.h:601
void svd(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:1948
void identity()
Set this to identity.
Definition: libmints/matrix.cc:1004
void copy_from(double ***)
Copies data from the passed matrix to this matrix_.
Definition: libmints/matrix.cc:560
double absmax()
Returns the absoluate maximum balue.
Definition: libmints/matrix.cc:1255
void symmetrize_gradient(std::shared_ptr< Molecule > mol)
Definition: libmints/matrix.cc:963
int nirrep_
Number of irreps.
Definition: libmints/matrix.h:76
void copy_upper_to_lower()
Definition: libmints/matrix.cc:2729
double & operator()(int i, int j)
Definition: libmints/matrix.h:1098
SharedMatrix collapse(int dim=0)
Definition: libmints/matrix.cc:1569
int symmetry_
Symmetry of this matrix (in most cases this will be 0 [totally symmetric])
Definition: libmints/matrix.h:84
double trace()
Returns the trace of this.
Definition: libmints/matrix.cc:1051
void zero_diagonal()
Zeros the diagonal.
Definition: libmints/matrix.cc:1037
SharedMatrix pseudoinverse(double condition=0.0, bool *conditioned=NULL)
Definition: libmints/matrix.cc:2061
void apply_symmetry(const SharedMatrix &a, const SharedMatrix &transformer)
Definition: libmints/matrix.cc:2786
void set(double val)
Definition: libmints/matrix.cc:570
void accumulate_product(const Matrix *const, const Matrix *const)
Multiplies the two arguments and adds their result to this.
Definition: libmints/matrix.cc:1191
double ** pointer(const int &h=0) const
Definition: libmints/matrix.h:483
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_a_temps()
Definition: libmints/matrix.cc:1934
void add(int h, int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:701
void gemm(bool transa, bool transb, double alpha, const Matrix *const a, const Matrix *const b, double beta)
Definition: libmints/matrix.cc:1427
void set_column(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:738
void svd_a(SharedMatrix &U, SharedVector &S, SharedMatrix &V)
Definition: libmints/matrix.cc:1995
SharedVector get_row(int h, int m)
Definition: libmints/matrix.cc:699
void zero_upper()
Definition: libmints/matrix.cc:2665
bool equal_but_for_row_order(const Matrix &rhs, double TOL=1.0e-10)
Definition: libmints/matrix.cc:3514
const double * get_const_pointer(const int &h=0) const
Definition: libmints/matrix.h:503
void subtract(const Matrix *const)
Subtracts a matrix from this.
Definition: libmints/matrix.cc:1150
Definition: libmints/matrix.h:62
std::tuple< SharedMatrix, SharedVector, SharedMatrix > svd_temps()
Definition: libmints/matrix.cc:1918
SharedVector get_column(int h, int m)
Definition: libmints/matrix.cc:713
void general_invert()
Definition: libmints/matrix.cc:2376
void print_atom_vector(std::string OutFileRMR="outfile")
Prints the matrix with atom and xyz styling.
Definition: libmints/matrix.cc:913
void element_add_mirror()
Definition: libmints/matrix.h:725
double * get_pointer(const int &h=0) const
Definition: libmints/matrix.h:498
diagonalize_order
Definition: libmints/matrix.h:58
void zero_row(int h, int i)
Definition: libmints/matrix.cc:2682
double rms()
Returns the rms of this.
Definition: libmints/matrix.cc:1238
double *** matrix_
Matrix data.
Definition: libmints/matrix.h:74
int rowspi(const int &h) const
Returns the rows per irrep array.
Definition: libmints/matrix.h:597
void expm(int n=2, bool scale=false)
Definition: libmints/matrix.cc:2484
SharedMatrix transpose()
Creates a new matrix which is the transpose of this.
Definition: libmints/matrix.cc:1068
void print_out() const
Python compatible printer.
Definition: libmints/matrix.h:557
bool equal(const Matrix &rhs)
Definition: libmints/matrix.cc:3477
const double ** const_pointer(const int &h=0) const
Definition: libmints/matrix.h:484
double * to_lower_triangle() const
Definition: libmints/matrix.cc:749
SharedMatrix clone() const
Creates an exact copy of the matrix and returns it.
Definition: libmints/matrix.cc:318
Definition: dimension.h:39
Definition: libmints/matrix.h:59
Makes using matrices just a little easier.
Definition: libmints/matrix.h:71
SaveType
Definition: libmints/matrix.h:262
void print_mat(const double *const *const a, int m, int n, std::string out) const
Definition: libmints/matrix.cc:812
static void free(double **Block)
free a (block) matrix – analogous to libciomr&#39;s free_block
Definition: libmints/matrix.cc:277
const double & operator()(int i, int j) const
Definition: libmints/matrix.h:1099
void recv()
Definition: libmints/matrix.cc:3457
void scale_row(int h, int m, double a)
Scale row m of irrep h by a.
Definition: libmints/matrix.cc:1213
const int & symmetry() const
Definition: libmints/matrix.h:652
Definition: libmints/matrix.h:263
Definition: libdpd/dpd.h:139
void save(const std::string &filename, bool append=true, bool saveLowerTriangle=true, bool saveSubBlocks=false)
Definition: libmints/matrix.cc:2985
void sum()
Definition: libmints/matrix.cc:3471
int nso
Definition: dx_write.cc:55
int coldim(const int &h=0) const
Returns the cols in irrep h.
Definition: libmints/matrix.h:590
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:501
void project_out(Matrix &v)
Definition: libmints/matrix.cc:1757
const Dimension & rowspi() const
Returns the rows per irrep array.
Definition: libmints/matrix.h:593
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:50
void rotate_columns(int h, int i, int j, double theta)
Definition: libmints/matrix.cc:3591
Definition: libmints/matrix.h:265
void zero_column(int h, int i)
Definition: libmints/matrix.cc:2693
void apply_denominator(const Matrix *const)
apply_denominators a matrix to this
Definition: libmints/matrix.cc:1165
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:2163
std::string name_
Name of the matrix.
Definition: libmints/matrix.h:82
void set(int h, int m, int n, double val)
Definition: libmints/matrix.h:385
void scale(double)
Scales this matrix.
Definition: libmints/matrix.cc:1202
void init(int nirrep, const int *rowspi, const int *colspi, const std::string &name="", int symmetry=0)
Definition: libmints/matrix.cc:283
SharedMatrix partial_cholesky_factorize(double delta=0.0, bool throw_if_negative=false)
Definition: libmints/matrix.cc:2194
static SharedMatrix doublet(const SharedMatrix &A, const SharedMatrix &B, bool transA=false, bool transB=false)
Definition: libmints/matrix.cc:1494
void set_numpy_shape(std::vector< int > shape)
Definition: libmints/matrix.h:1143
Dimension colspi_
Columns per irrep array.
Definition: libmints/matrix.h:80
std::shared_ptr< PsiOutStream > outfile
Definition: core.cc:99
SharedMatrix canonical_orthogonalization(double delta=0.0, SharedMatrix eigvec=SharedMatrix())
Definition: libmints/matrix.cc:2108
double ** to_block_matrix() const
Definition: libmints/matrix.cc:766
void copy_lower_to_upper()
Definition: libmints/matrix.cc:2704
Definition: pointgrp.h:105
void axpy(double a, SharedMatrix X)
Definition: libmints/matrix.cc:1550
Dimension rowspi_
Rows per irrep array.
Definition: libmints/matrix.h:78
void invert()
Definition: libmints/matrix.cc:2360
Definition: PsiFileImpl.h:38
void set_row(int h, int m, SharedVector vec)
Definition: libmints/matrix.cc:727
void add(const Matrix *const)
Adds a matrix to this.
Definition: libmints/matrix.cc:1130
Definition: libmints/matrix.h:61
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:1543
double & operator()(int h, int i, int j)
Definition: libmints/matrix.h:1100
const int & nirrep() const
Returns the number of irreps.
Definition: libmints/matrix.h:609
Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double *res=0)
Definition: libmints/matrix.cc:1670
double sum_of_squares()
Returns the sum of the squares of this.
Definition: libmints/matrix.cc:1223
const std::string & name() const
Definition: libmints/matrix.h:554
#define PSIEXCEPTION(message)
Definition: exception.h:47
void bcast(int broadcaster)
Definition: libmints/matrix.cc:3461
static SharedMatrix vertcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:403
Definition: vector.h:48
void print(std::string outfile="outfile", const char *extra=NULL) const
Definition: libmints/matrix.cc:861
void copy(const SharedMatrix &cp)
Definition: libmints/matrix.cc:514
void print_to_mathematica()
Definition: libmints/matrix.cc:883
SimpleMatrix * to_simple_matrix() const
std::shared_ptr< Vector > SharedVector
Definition: adc.h:52
void remove_symmetry(const SharedMatrix &a, const SharedMatrix &transformer)
Definition: libmints/matrix.cc:2842
int nrow() const
Returns the total number of rows.
Definition: libmints/matrix.h:614
double ** operator[](int i)
Definition: libmints/matrix.h:1097
Dimension power(double alpha, double cutoff=1.0E-12)
Definition: libmints/matrix.cc:2423
int max_ncol() const
Returns the column size of the largest block.
Definition: libmints/matrix.h:639
Definition: psio.hpp:196
Definition: libmints/matrix.h:60
void scale_column(int h, int n, double a)
Scale column n of irrep h by a.
Definition: libmints/matrix.cc:1218