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-2019 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 #pragma once
30 
31 #include <string>
32 #include <vector>
33 #include <memory>
34 
36 
37 #include "dimension.h"
38 
39 namespace psi {
40 
41 struct dpdfile2;
42 struct dpdbuf4;
43 
44 class PSIO;
45 class Vector;
46 using SharedVector = std::shared_ptr<Vector>;
47 class Dimension;
48 class Molecule;
49 class Vector3;
50 class Matrix;
51 using SharedMatrix = std::shared_ptr<Matrix>;
52 
54 
55 namespace linalg {
60 PSI_API
61 SharedMatrix horzcat(const std::vector<SharedMatrix>& mats);
62 
67 PSI_API
68 SharedMatrix vertcat(const std::vector<SharedMatrix>& mats);
69 
76 PSI_API
77 SharedMatrix doublet(const SharedMatrix& A, const SharedMatrix& B, bool transA = false, bool transB = false);
78 
87 PSI_API
88 SharedMatrix triplet(const SharedMatrix& A, const SharedMatrix& B, const SharedMatrix& C, bool transA = false,
89  bool transB = false, bool transC = false);
90 
91 namespace detail {
95 PSI_API
96 double** matrix(int nrow, int ncol);
97 
101 PSI_API
102 void free(double** Block);
103 } // namespace detail
104 } // namespace linalg
105 
112 class PSI_API Matrix : public std::enable_shared_from_this<Matrix> {
113  protected:
115  double*** matrix_;
117  int nirrep_;
123  std::string name_;
126 
128  void alloc();
130  void release();
131 
133  void copy_from(double***);
134 
135  void print_mat(const double* const* const a, int m, int n, std::string out) const;
136 
138  std::vector<int> numpy_shape_;
139 
140  public:
142  Matrix();
148  Matrix(const std::string& name, int symmetry = 0);
150  Matrix(const Matrix& copy);
152  explicit Matrix(const SharedMatrix& copy);
154  explicit Matrix(const Matrix* copy);
162  Matrix(int nirrep, const int* rowspi, const int* colspi, int symmetry = 0);
171  Matrix(const std::string& name, int nirrep, const int* rowspi, const int* colspi, int symmetry = 0);
178  Matrix(int nirrep, int rows, const int* colspi);
179 
186  Matrix(int nirrep, const int* rowspi, int cols);
187 
195  Matrix(int rows, int cols);
204  Matrix(const std::string& name, int rows, int cols);
205 
211  Matrix(dpdfile2* inFile);
212 
218  Matrix(dpdbuf4 *inBuf);
219 
220 
229  Matrix(const std::string& name, const Dimension& rows, const Dimension& cols, int symmetry = 0);
230 
238  Matrix(const Dimension& rows, const Dimension& cols, int symmetry = 0);
239 
241  virtual ~Matrix();
242 
252  void init(int nirrep, const int* rowspi, const int* colspi, const std::string& name = "", int symmetry = 0);
253 
254  void init(const Dimension& rowspi, const Dimension& colspi, const std::string& name = "", int symmetry = 0);
255 
257  SharedMatrix clone() const;
258 
264  void copy(const SharedMatrix& cp);
265  void copy(const Matrix& cp);
266  void copy(const Matrix* cp);
279  SharedMatrix matrix_3d_rotation(Vector3 axis, double phi, bool Sn);
280 
282  void copy_to_row(int h, int row, double const* const data);
283 
284  enum SaveType { Full, SubBlocks, LowerTriangle };
285 
296  bool load(psi::PSIO* psio, size_t fileno, const std::string& tocentry, int nso);
297  bool load(std::shared_ptr<psi::PSIO>& psio, size_t fileno, const std::string& tocentry, int nso);
309  void load(psi::PSIO* const psio, size_t fileno, SaveType savetype = LowerTriangle);
310  void load(std::shared_ptr<psi::PSIO>& psio, size_t fileno, SaveType savetype = LowerTriangle);
318  void load(const std::string& filename);
319 
326  void load_mpqc(const std::string& filename);
327 
338  void save(const std::string& filename, bool append = true, bool saveLowerTriangle = true,
339  bool saveSubBlocks = false);
350  void save(psi::PSIO* const psio, size_t fileno, SaveType savetype = LowerTriangle);
351  void save(std::shared_ptr<psi::PSIO>& psio, size_t fileno, SaveType savetype = LowerTriangle);
359  void set(double val);
360 
366  void set(const double* const tri);
367 
374  void set(const double* const* const sq);
384  void set(const double* const* const sq, int irrep);
395  void set(int h, int m, int n, double val) { matrix_[h][m][n] = val; }
396 
404  void set(int m, int n, double val) { matrix_[0][m][n] = val; }
405 
412  void set_diagonal(const Vector* const vec);
413  void set_diagonal(const Vector& vec);
414  void set_diagonal(const std::shared_ptr<Vector>& vec);
425  double get(const int& h, const int& m, const int& n) const { return matrix_[h][m][n]; }
426 
434  double get(const int& m, const int& n) const { return matrix_[0][m][n]; }
435 
443  SharedVector get_row(int h, int m);
444 
452  SharedVector get_column(int h, int m);
453 
461  void set_row(int h, int m, SharedVector vec);
462 
470  void set_column(int h, int m, SharedVector vec);
471 
479  SharedMatrix get_block(const Slice& rows, const Slice& cols);
480 
488  void set_block(const Slice& rows, const Slice& cols, SharedMatrix block);
489 
502  double** pointer(const int& h = 0) const { return matrix_[h]; }
503  const double** const_pointer(const int& h = 0) const { return const_cast<const double**>(matrix_[h]); }
504 
517  double* get_pointer(const int& h = 0) const {
518  if (rowspi_[h] * (size_t)colspi_[h] > 0)
519  return &(matrix_[h][0][0]);
520  else
521  return nullptr;
522  }
523  const double* get_const_pointer(const int& h = 0) const {
524  if (rowspi_[h] * (size_t)colspi_[h] > 0)
525  return const_cast<const double*>(&(matrix_[h][0][0]));
526  else
527  return nullptr;
528  }
529 
530  size_t size(const int& h = 0) const { return colspi_[h] * (size_t)rowspi_[h]; }
531 
533  void apply_denominator(const Matrix* const);
535  void apply_denominator(const Matrix&);
537  void apply_denominator(const SharedMatrix&);
538 
544  double** to_block_matrix() const;
550  SharedMatrix to_block_sharedmatrix() const;
556  double* to_lower_triangle() const;
557 
563  void set_name(const std::string& name) { name_ = name; }
564 
568  const std::string& name() const { return name_; }
569 
571  void print_out() const { print("outfile"); }
572 
579  void print(std::string outfile = "outfile", const char* extra = nullptr) const;
580 
582  void print_atom_vector(std::string out_fname = "outfile");
583 
587  void print_to_mathematica();
588 
595  void eivprint(const Vector* const values, std::string out = "outfile");
597  void eivprint(const Vector& values, std::string out = "outfile");
599  void eivprint(const std::shared_ptr<Vector>& values, std::string out = "outfile");
600 
602  int rowdim(const int& h = 0) const { return rowspi_[h]; }
604  int coldim(const int& h = 0) const { return colspi_[h]; }
605 
607  const Dimension& rowspi() const { return rowspi_; }
609  int rowspi(const int& h) const { return rowdim(h); }
611  const Dimension& colspi() const { return colspi_; }
613  int colspi(const int& h) const { return coldim(h); }
615  const int& nirrep() const { return nirrep_; }
616 
618  int nrow() const {
619  int rows = 0;
620  for (int h = 0; h < nirrep(); ++h) rows += rowdim(h);
621  return rows;
622  }
623 
625  int ncol() const {
626  int cols = 0;
627  for (int h = 0; h < nirrep(); ++h) cols += coldim(h);
628  return cols;
629  }
630 
632  int max_nrow() const {
633  int row = 0;
634  for (int h = 0; h < nirrep(); ++h)
635  if (row < rowdim(h)) row = rowdim(h);
636  return row;
637  }
638 
640  int max_ncol() const {
641  int col = 0;
642  for (int h = 0; h < nirrep(); ++h)
643  if (col < coldim(h)) col = coldim(h);
644  return col;
645  }
646 
652  const int& symmetry() const { return symmetry_; }
653 
658  void symmetrize_gradient(std::shared_ptr<Molecule> mol);
659 
664  void symmetrize_hessian(std::shared_ptr<Molecule> mol);
665 
667  void identity();
669  void zero();
671  void zero_diagonal();
672 
673  // Math routines
675  double trace();
677  SharedMatrix transpose();
678 
680  void transpose_this();
681 
683  void add(const Matrix* const);
685  void add(const Matrix&);
687  void add(const SharedMatrix&);
688 
690  void subtract(const Matrix* const);
692  void subtract(const SharedMatrix&);
694  void accumulate_product(const Matrix* const, const Matrix* const);
695  void accumulate_product(const SharedMatrix&, const SharedMatrix&);
697  void scale(double);
699  double sum_of_squares();
701  double rms();
703  double absmax();
705  void add(int h, int m, int n, double val) {
706 #ifdef PSIDEBUG
707  if (m > rowspi_[h] || n > colspi_[h ^ symmetry_]) {
708  outfile->Printf("out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n", symmetry_, h, m, n);
709 
710  throw PSIEXCEPTION("What are you doing, Rob?");
711  }
712 #endif
713  matrix_[h][m][n] += val;
714  }
716  void add(int m, int n, double val) {
717 #ifdef PSIDEBUG
718  if (m > rowspi_[0] || n > colspi_[0 ^ symmetry_]) {
719  outfile->Printf("out of bounds: symmetry_ = %d, h = %d, m = %d, n = %d\n", symmetry_, 0, m, n);
720 
721  return;
722  }
723 #endif
724  matrix_[0][m][n] += val;
725  }
726 
728  for (int h = 0; h < nirrep_; ++h) {
729  for (int i = 0; i < rowspi_[h]; ++i) {
730  for (int j = 0; j < i; ++j) {
731  matrix_[h][i][j] = matrix_[h][j][i] = (matrix_[h][i][j] + matrix_[h][j][i]);
732  }
733  }
734  }
735  }
736 
738  void scale_row(int h, int m, double a);
740  void scale_column(int h, int n, double a);
741 
747  void apply_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
748 
754  void remove_symmetry(const SharedMatrix& a, const SharedMatrix& transformer);
761  void transform(const SharedMatrix& L, const SharedMatrix& F, const SharedMatrix& R);
762 
765  void transform(const Matrix* const a, const Matrix* const transformer);
766  void transform(const SharedMatrix& a, const SharedMatrix& transformer);
768 
771  void transform(const Matrix* const transformer);
772  void transform(const SharedMatrix& transformer);
774 
777  void back_transform(const Matrix* const a, const Matrix* const transformer);
778  void back_transform(const SharedMatrix& a, const SharedMatrix& transformer);
780 
783  void back_transform(const Matrix* const transformer);
784  void back_transform(const SharedMatrix& transformer);
786 
788  double vector_dot(const Matrix* const rhs);
789  double vector_dot(const SharedMatrix& rhs);
790  double vector_dot(const Matrix& rhs);
791 
793 
801  void gemm(bool transa, bool transb, double alpha, const Matrix* const a, const Matrix* const b, double beta);
802  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const SharedMatrix& b, double beta);
803  void gemm(bool transa, bool transb, double alpha, const SharedMatrix& a, const Matrix& b, double beta);
804  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const SharedMatrix& b, double beta);
806 
808 
810  void gemm(const char& transa, const char& transb, const std::vector<int>& m, const std::vector<int>& n,
811  const std::vector<int>& k, const double& alpha, const SharedMatrix& a, const std::vector<int>& lda,
812  const SharedMatrix& b, const std::vector<int>& ldb, const double& beta, const std::vector<int>& ldc,
813  const std::vector<unsigned long>& offset_a = std::vector<unsigned long>(),
814  const std::vector<unsigned long>& offset_b = std::vector<unsigned long>(),
815  const std::vector<unsigned long>& offset_c = std::vector<unsigned long>());
816  void gemm(const char& transa, const char& transb, const int& m, const int& n, const int& k, const double& alpha,
817  const SharedMatrix& a, const int& lda, const SharedMatrix& b, const int& ldb, const double& beta,
818  const int& ldc, const unsigned long& offset_a = 0, const unsigned long& offset_b = 0,
819  const unsigned long& offset_c = 0);
821 
827  void axpy(double a, SharedMatrix X);
828 
834  SharedMatrix collapse(int dim = 0);
835 
838  void diagonalize(Matrix* eigvectors, Vector* eigvalues, diagonalize_order nMatz = ascending);
839  void diagonalize(SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues, diagonalize_order nMatz = ascending);
840  void diagonalize(SharedMatrix& eigvectors, Vector& eigvalues, diagonalize_order nMatz = ascending);
842 
846  void diagonalize(SharedMatrix& metric, SharedMatrix& eigvectors, std::shared_ptr<Vector>& eigvalues,
847  diagonalize_order nMatz = ascending);
849 
852  void svd(SharedMatrix& U, SharedVector& S, SharedMatrix& V);
854 
859  void svd_a(SharedMatrix& U, SharedVector& S, SharedMatrix& V);
861 
864  std::tuple<SharedMatrix, SharedVector, SharedMatrix> svd_temps();
866 
869  std::tuple<SharedMatrix, SharedVector, SharedMatrix> svd_a_temps();
871 
874  SharedMatrix pseudoinverse(double condition, int& nremoved);
876 
882  SharedMatrix canonical_orthogonalization(double delta = 0.0, SharedMatrix eigvec = SharedMatrix());
883 
910  SharedMatrix partial_cholesky_factorize(double delta = 0.0, bool throw_if_negative = false);
911 
926  std::pair<SharedMatrix, SharedMatrix> partial_square_root(double delta = 0.0);
927 
933  void cholesky_factorize();
934 
939  void invert();
940 
945  void general_invert();
946 
968  Dimension power(double alpha, double cutoff = 1.0E-12);
969 
977  void expm(int n = 2, bool scale = false);
978 
980  void swap_rows(int h, int i, int j);
982  void swap_columns(int h, int i, int j);
983 
985  void hermitivitize();
987  void copy_lower_to_upper();
989  void copy_upper_to_lower();
991  void zero_lower();
993  void zero_upper();
995  void zero_row(int h, int i);
997  void zero_column(int h, int i);
998 
999  // Reference versions of the above functions
1001  void transform(const Matrix& a, const Matrix& transformer);
1003  void transform(const Matrix& transformer);
1005  void back_transform(const Matrix& a, const Matrix& transformer);
1007  void back_transform(const Matrix& transformer);
1008 
1013  bool add_and_orthogonalize_row(const SharedVector v);
1014 
1030  bool schmidt_add_row(int h, int rows, Vector& v);
1031  bool schmidt_add_row(int h, int rows, double* v) noexcept;
1033 
1035  void schmidt();
1036 
1039  void schmidt_orthog(SharedMatrix S, int n);
1040 
1047  Dimension schmidt_orthog_columns(SharedMatrix S, double tol, double* res = nullptr);
1048 
1056  void project_out(Matrix& v);
1057 
1059  void gemm(bool transa, bool transb, double alpha, const Matrix& a, const Matrix& b, double beta);
1061  void diagonalize(Matrix& eigvectors, Vector& eigvalues, int nMatz = 1);
1062 
1065  double** operator[](int i) { return matrix_[i]; }
1066  double& operator()(int i, int j) { return matrix_[0][i][j]; }
1067  const double& operator()(int i, int j) const { return matrix_[0][i][j]; }
1068  double& operator()(int h, int i, int j) { return matrix_[h][i][j]; }
1069  const double& operator()(int h, int i, int j) const { return matrix_[h][i][j]; }
1071 
1073  void write_to_dpdfile2(dpdfile2* outFile);
1074 
1076  void write_to_dpdbuf4(dpdbuf4 *outBuf);
1077 
1082  bool equal(const Matrix& rhs, double TOL = 1.0e-10);
1083  bool equal(const SharedMatrix& rhs, double TOL = 1.0e-10);
1084  bool equal(const Matrix* rhs, double TOL = 1.0e-10);
1086 
1091  bool equal_but_for_row_order(const Matrix& rhs, double TOL = 1.0e-10);
1092  bool equal_but_for_row_order(const SharedMatrix& rhs, double TOL = 1.0e-10);
1093  bool equal_but_for_row_order(const Matrix* rhs, double TOL = 1.0e-10);
1095 
1099  void set_numpy_shape(std::vector<int> shape) { numpy_shape_ = shape; }
1100  std::vector<int> numpy_shape() { return numpy_shape_; }
1101 
1109  void rotate_columns(int h, int i, int j, double theta);
1110 
1112  "Using `Matrix::matrix` instead of `linalg::detail::matrix` is deprecated, and in 1.4 it will "
1113  "stop working")
1114  static double** matrix(int nrow, int ncol) { return linalg::detail::matrix(nrow, ncol); }
1115 
1117  "Using `Matrix::free` instead of `linalg::detail::free` is deprecated, and in 1.4 it will "
1118  "stop working")
1119  static void free(double** Block) { linalg::detail::free(Block); }
1120 
1122  "Using `Matrix::create` instead of `auto my_mat = std::make_shared<Matrix>(name, rows, cols);` "
1123  "is deprecated, and in 1.4 it will "
1124  "stop working")
1125  static SharedMatrix create(const std::string& name, const Dimension& rows, const Dimension& cols) {
1126  return std::make_shared<Matrix>(name, rows, cols);
1127  }
1128 
1130  "Using `Matrix::horzcat` instead of `linalg::horzcat` is deprecated, and in 1.4 it will "
1131  "stop working")
1132  static SharedMatrix horzcat(const std::vector<SharedMatrix>& mats) { return linalg::horzcat(mats); }
1133 
1135  "Using `Matrix::vertcat` instead of `linalg::vertcat` is deprecated, and in 1.4 it will "
1136  "stop working")
1137  static SharedMatrix vertcat(const std::vector<SharedMatrix>& mats) { return linalg::vertcat(mats); }
1138 
1140  "Using `Matrix::doublet` instead of `doublet` is deprecated, and in 1.4 it will "
1141  "stop working")
1142  static SharedMatrix doublet(const SharedMatrix& A, const SharedMatrix& B, bool transA = false,
1143  bool transB = false) {
1144  return linalg::doublet(A, B, transA, transB);
1145  }
1146 
1148  "Using `Matrix::triplet` instead of `triplet` is deprecated, and in 1.4 it will "
1149  "stop working")
1150  static SharedMatrix triplet(const SharedMatrix& A, const SharedMatrix& B, const SharedMatrix& C,
1151  bool transA = false, bool transB = false, bool transC = false) {
1152  return linalg::triplet(A, B, C, transA, transB, transC);
1153  }
1154 };
1155 } // namespace psi
SharedMatrix triplet(const SharedMatrix &A, const SharedMatrix &B, const SharedMatrix &C, bool transA, bool transB, bool transC)
Definition: libmints/matrix.cc:3306
std::vector< int > numpy_shape_
Numpy Shape.
Definition: libmints/matrix.h:138
Definition: libmints/matrix.h:284
const double & operator()(int h, int i, int j) const
Definition: libmints/matrix.h:1069
std::vector< int > numpy_shape()
Definition: libmints/matrix.h:1100
int rowdim(const int &h=0) const
Returns the rows in irrep h.
Definition: libmints/matrix.h:602
int colspi(const int &h) const
Returns the columns per irrep array.
Definition: libmints/matrix.h:613
double ** matrix(int nrow, int ncol)
allocate a block matrix – analogous to libciomr&#39;s block_matrix
Definition: libmints/matrix.cc:3315
int max_nrow() const
Returns the row size of the largest block.
Definition: libmints/matrix.h:632
int * U
Definition: stringlist.cc:66
#define TOL
Definition: detci/opdm.cc:58
void set(int m, int n, double val)
Definition: libmints/matrix.h:404
void set_block(Slice rows, Slice cols, SharedMatrix block)
size_t size(const int &h=0) const
Definition: libmints/matrix.h:530
void schmidt(double **A, int rows, int cols, std::string out_fname)
Definition: schmidt.cc:60
PSI_DEPRECATED("Using `Matrix::free` instead of `linalg::detail::free` is deprecated, and in 1.4 it will ""stop working") static void free(double **Block)
Definition: libmints/matrix.h:1116
PSI_DEPRECATED("Using `Matrix::vertcat` instead of `linalg::vertcat` is deprecated, and in 1.4 it will ""stop working") static SharedMatrix vertcat(const std
Definition: libmints/matrix.h:1134
void add(int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:716
SharedMatrix doublet(const SharedMatrix &A, const SharedMatrix &B, bool transA, bool transB)
Definition: libmints/matrix.cc:3296
PSI_DEPRECATED("Using `Matrix::horzcat` instead of `linalg::horzcat` is deprecated, and in 1.4 it will ""stop working") static SharedMatrix horzcat(const std
Definition: libmints/matrix.h:1129
void set_name(const std::string &name)
Definition: libmints/matrix.h:563
Definition: matrixtmp.h:42
Definition: vector3.h:40
Definition: pointgrp.h:104
int ncol() const
Returns the total number of columns.
Definition: libmints/matrix.h:625
const Dimension & colspi() const
Returns the columns per irrep array.
Definition: libmints/matrix.h:611
int nirrep_
Number of irreps.
Definition: libmints/matrix.h:117
double & operator()(int i, int j)
Definition: libmints/matrix.h:1066
int symmetry_
Symmetry of this matrix (in most cases this will be 0 [totally symmetric])
Definition: libmints/matrix.h:125
PSI_API void print_mat(double **a, int rows, int cols, std::string out)
Definition: print_mat.cc:52
std::shared_ptr< PsiOutStream > outfile
Definition: core.cc:106
double ** pointer(const int &h=0) const
Definition: libmints/matrix.h:502
void add(int h, int m, int n, double val)
Add val to an element of this.
Definition: libmints/matrix.h:705
Slicing for Matrices and Vectors objects.
Definition: dimension.h:125
const double * get_const_pointer(const int &h=0) const
Definition: libmints/matrix.h:523
Definition: libmints/matrix.h:53
void element_add_mirror()
Definition: libmints/matrix.h:727
double * get_pointer(const int &h=0) const
Definition: libmints/matrix.h:517
diagonalize_order
Definition: libmints/matrix.h:53
double *** matrix_
Matrix data.
Definition: libmints/matrix.h:115
int rowspi(const int &h) const
Returns the rows per irrep array.
Definition: libmints/matrix.h:609
long int size_t
Definition: sortintegrals.cc:41
SharedMatrix vertcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:3250
void print_out() const
Python compatible printer.
Definition: libmints/matrix.h:571
const double ** const_pointer(const int &h=0) const
Definition: libmints/matrix.h:503
Definition: dimension.h:40
Definition: libmints/matrix.h:53
Makes using matrices just a little easier.
Definition: libmints/matrix.h:112
SaveType
Definition: libmints/matrix.h:284
PSI_DEPRECATED("Using `Matrix::create` instead of `auto my_mat = std::make_shared<Matrix>(name, rows, cols);` ""is deprecated, and in 1.4 it will ""stop working") static SharedMatrix create(const std
Definition: libmints/matrix.h:1121
const double & operator()(int i, int j) const
Definition: libmints/matrix.h:1067
const int & symmetry() const
Definition: libmints/matrix.h:652
Definition: libdpd/dpd.h:140
int nso
Definition: dx_write.cc:56
int coldim(const int &h=0) const
Returns the cols in irrep h.
Definition: libmints/matrix.h:604
const Dimension & rowspi() const
Returns the rows per irrep array.
Definition: libmints/matrix.h:607
std::shared_ptr< Matrix > SharedMatrix
Definition: adc.h:49
int delta(const int i, const int j)
Definition: bend.cc:175
std::string name_
Name of the matrix.
Definition: libmints/matrix.h:123
#define PSI_API
Definition: pragma.h:155
void set(int h, int m, int n, double val)
Definition: libmints/matrix.h:395
Definition: libdpd/dpd.h:105
void set_numpy_shape(std::vector< int > shape)
Definition: libmints/matrix.h:1099
Dimension colspi_
Columns per irrep array.
Definition: libmints/matrix.h:121
Definition: pointgrp.h:104
Dimension rowspi_
Rows per irrep array.
Definition: libmints/matrix.h:119
SharedMatrix horzcat(const std::vector< SharedMatrix > &mats)
Definition: libmints/matrix.cc:3204
Definition: libmints/matrix.h:53
double & operator()(int h, int i, int j)
Definition: libmints/matrix.h:1068
const int & nirrep() const
Returns the number of irreps.
Definition: libmints/matrix.h:615
const std::string & name() const
Definition: libmints/matrix.h:568
#define PSIEXCEPTION(message)
Definition: exception.h:48
#define PSI_DEPRECATED(msg)
Definition: pragma.h:164
Definition: vector.h:44
void free(double **Block)
free a (block) matrix – analogous to libciomr&#39;s free_block
Definition: libmints/matrix.cc:3325
std::shared_ptr< Vector > SharedVector
Definition: adc.h:51
int nrow() const
Returns the total number of rows.
Definition: libmints/matrix.h:618
double ** operator[](int i)
Definition: libmints/matrix.h:1065
int max_ncol() const
Returns the column size of the largest block.
Definition: libmints/matrix.h:640
Definition: psio.hpp:195
Definition: libmints/matrix.h:53