Loading [MathJax]/extensions/tex2jax.js
SPI
SLEPc PETSc Interface is used for rapid development and intuitive matrix operations similar to MatLab or NumPy style of operations
•All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
SPIMat.hpp
Go to the documentation of this file.
1 #ifndef SPIMAT_H
2 #define SPIMAT_H
3 #include <iostream>
4 #include <vector>
5 #include <petscksp.h>
6 #include <petscmath.h>
7 #include <slepceps.h>
8 #include <slepcpep.h>
9 #include <string>
10 #include <tuple>
11 #include "SPIVec.hpp"
12 
13 namespace SPI{
14  template <class T>
15  using Block2D = std::vector<std::vector<T>>;
16 
17  struct SPIMat{
18  PetscInt rows;
19  PetscInt cols;
20 
21  // Constructors
22  SPIMat(std::string _name="SPIMat"); // constructor with no arguments (no initialization)
23  SPIMat(const std::vector<SPIVec> &A, std::string _name="SPIMat"); // constructor using a vector of column vectors
24  SPIMat(const SPIMat &A, std::string _name="SPIMat"); // constructor using another SPIMat
25  SPIMat(PetscInt rowscols, std::string _name="SPIMat"); // constructor with one arguement to make square matrix
26  SPIMat(PetscInt rowsm, PetscInt colsn, std::string _name="SPIMat"); // constructor of rectangular matrix
27 
28  Mat mat;
29  PetscErrorCode ierr;
30 
31  // flags
32  PetscBool flag_init=PETSC_FALSE;
33  std::string name;
34 
35  PetscInt Init(PetscInt m,PetscInt n, std::string name="SPIMat"); // initialize the matrix of size m by n
36  SPIMat& set(PetscInt m, PetscInt n,const PetscScalar v); // set a scalar value at position row m and column n
37  SPIMat& set_col(const PetscInt col,const SPIVec &v); // set a column into a matrix
38  SPIMat& add(PetscInt m, PetscInt n,const PetscScalar v); // add a scalar value at position row m and column n
39  // () operators
40  PetscScalar operator()(PetscInt m, PetscInt n, PetscBool global=PETSC_FALSE); // get local value at row m, column n
41  PetscScalar operator()(PetscInt m, PetscInt n, PetscBool global=PETSC_FALSE) const; // get local value at row m, column n
42  SPIMat& operator()(PetscInt m, PetscInt n,const PetscScalar v); // set operator the same as set function
43  SPIMat& operator()(PetscInt m, PetscInt n,const double v); // set operator the same as set function
44  SPIMat& operator()(PetscInt m, PetscInt n,const int v); // set operator the same as set function
45  SPIMat& operator()(PetscInt m, PetscInt n,const SPIMat &Asub, InsertMode addv=INSERT_VALUES); // set submatrix into matrix at row m, col n
46  SPIMat& operator()(); // assemble the matrix
47  // +- operators
48  SPIMat& operator+=(const SPIMat &X); // MatAXPY, Y = 1.*X + Y operation
49  SPIMat& axpy(const PetscScalar a, const SPIMat &X); // MatAXPY function call to add a*X to the current mat
50  SPIMat operator+(const SPIMat &X); // Y + X operation
51  SPIMat& operator-=(const SPIMat &X); // Y = -1.*X + Y operation
52  SPIMat operator-(const SPIMat &X); // Y - X operation
53  SPIMat operator-() const; // -X operation
54  // * operators
55  SPIMat operator*(const PetscScalar a); // Y*a operation
56  SPIMat operator*(const double a); // Y*a operation
57  SPIVec operator*(const SPIVec &x); // A*x operation to return a vector
58  SPIMat& operator*=(const PetscScalar a); // Y = Y*a operation
59  SPIMat& operator*=(const double a); // Y = Y*a operation
60  SPIMat& operator/=(const PetscScalar a); // Y = Y/a operation
61  SPIMat operator/(const PetscScalar a); // Z = Y/a operation
62  SPIMat operator/(const SPIMat &A); // Z = Y/A elementwise operation
63  SPIMat operator*(const SPIMat &A); // Y*A operation
64  // = operator
65  SPIMat& operator=(const SPIMat &A); // Y=X with initialization of Y using MatConvert
66  // overload % for element wise multiplication
67  //SPIMat operator%(SPIMat A);
68  // Transpose functions
69  PetscInt T(SPIMat &A); // A = Transpose(*this.mat) operation with initialization of A
70  SPIMat& T(); // Transpose the current mat
71  // conjugate
72  PetscInt H(SPIMat &A); // A = Hermitian Transpose(*this.mat) operation with initialization of A (tranpose and complex conjugate)
73  SPIMat& H(); // Hermitian Transpose the current mat
74  SPIVec H(const SPIVec &q); // Hermitian Transpose and multiply by vector
75  SPIMat& conj(); // elemenwise conjugate current matrix
76  SPIMat& real(); // take the real part of the matrix (alters current matrix)
77  SPIVec diag(); // get diagonal of matrix
78  SPIMat& zero_row(const PetscInt row); // zero a row
79  SPIMat& eye_row(const PetscInt row); // zero a row and put 1 in the diagonal entry
80  SPIMat& zero_row_full(const PetscInt row); // zero a row
81  SPIMat& zero_rows(std::vector<PetscInt> rows); // zero every row
82  SPIMat& eye_rows(std::vector<PetscInt> rows); // zero every row
83  SPIVec col(const PetscInt i); // get column vector using MatGetColumnVector(mat,vec,i);
84  PetscInt print(); // print mat to screen using PETSC_VIEWER_STDOUT_WORLD
85 
86  ~SPIMat(); // destructor to delete memory
87 
88  };
89  SPIMat operator*(const PetscScalar a, const SPIMat A); // a*A operation to be equivalent to A*a
90  SPIMat operator*(const SPIMat A, const PetscScalar a); // A*a operation to be equivalent to A*a
91  SPIVec operator/(const SPIVec &b, const SPIMat &A); // Solve linear system, Ax=b using b/A notation
92  SPIMat operator^(const PetscScalar a, const SPIMat &A); // Y = a^A operation
93  //SPIVec solve(const SPIVec &b, const SPIMat &A); // Solve linear system, Ax=b using solve(A,b) notation
94  SPIVec solve(const SPIMat &A, const SPIVec &b); // Solve linear system, Ax=b using solve(A,b) notation
95  SPIMat eye(const PetscInt n); // create, form, and return identity matrix of size n
96  SPIMat inv(const SPIMat &A); // get inverse of matrix by solving A*Ainv = B using MatMatSolve
97  SPIMat zeros(const PetscInt m,const PetscInt n); // create, form, and return zero matrix of size mxn
98  SPIMat diag(const SPIVec &diag,const PetscInt k=0); // set diagonal of matrix
99  SPIMat kron(const SPIMat &A, const SPIMat &B); // set kronecker inner product of two matrices
100  std::tuple<std::vector<PetscReal>,std::vector<SPIVec>,std::vector<SPIVec>> svd(const SPIMat &A); // solve general SVD problem of A = U*E*V^H
101  SPIVec lstsq(const SPIMat &A, SPIVec &y); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector)
102  SPIVec lstsq(const std::vector<SPIVec> &A, SPIVec &y); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector)
103  std::tuple<PetscScalar,SPIVec,SPIVec> eig(const SPIMat &A, const SPIMat &B, const PetscScalar target,const PetscReal tol=-1,const PetscInt max_iter=-1); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector)
104  std::tuple<PetscScalar,SPIVec> eig_right(const SPIMat &A, const SPIMat &B, const PetscScalar target,const PetscReal tol=-1,const PetscInt max_iter=-1); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector)
105  std::tuple<PetscScalar,SPIVec, SPIVec> eig_init(const SPIMat &A, const SPIMat &B, const PetscScalar target,const SPIVec &ql, const SPIVec &qr, PetscReal tol=-1,const PetscInt max_iter=-1); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector) using initial subspace from q
106  std::tuple<PetscScalar, SPIVec> eig_init_right(const SPIMat &A, const SPIMat &B, const PetscScalar target, const SPIVec &qr, PetscReal tol=-1,const PetscInt max_iter=-1); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector) using initial subspace from q
107  std::tuple<std::vector<PetscScalar>, std::vector<SPIVec>> eig_init_rights(const SPIMat &A, const SPIMat &B, const std::vector<PetscScalar> targets, const std::vector<SPIVec> &qrs, PetscReal tol=-1,const PetscInt max_iter=-1); // solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector) using initial subspace from q
108  std::tuple<PetscScalar,SPIVec> polyeig(const std::vector<SPIMat> &As, const PetscScalar target,const PetscReal tol=-1.,const PetscInt max_iter=-1); // solve general polynomial eigenvalue problem of (A0 + A1x + A2x^2 ...) = 0 and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector)
109  std::tuple<PetscScalar,SPIVec> polyeig_init(const std::vector<SPIMat> &As, const PetscScalar target, const SPIVec &qr, const PetscReal tol=-1.,const PetscInt max_iter=-1); // solve general polynomial eigenvalue problem of (A0 + A1x + A2x^2 ...) = 0 and return a tuple of tie(PetscScalar alpha, SPIVec eig_vector) using initial subspace from q
110  //SPIMat block(const SPIMat Blocks[], const PetscInt rows,const PetscInt cols); // set block matrices using an input array of size rows*cols. Fills rows first
111  //SPIMat block(const std::vector<std::vector<SPIMat>> Blocks); // set block matrices using an input array of size rows*cols. Fills rows first
112  SPIMat block(const Block2D<SPIMat> Blocks); // set block matrices using an input array of size rows*cols. Fills rows first
113  std::tuple<SPIMat,SPIMat> meshgrid(SPIVec &x, SPIVec &y); // create meshgrid from two grids using ij indexing. i.e. X(i,j) = x(i) and Y(i,j) = y(j)
114  PetscInt save(const SPIMat &A, const std::string filename); // save matrix to filename to binary format
115  PetscInt save(const std::vector<SPIMat> &A, const std::string filename); // save matrix to filename to binary format
116  PetscInt load(SPIMat &A, const std::string filename); // load matrix to filename from binary format
117  PetscInt load(std::vector<SPIMat> &A, const std::string filename); // load matrix to filename from binary format
118  PetscInt draw(const SPIMat &A); // draw nonzero structure and wait at command line input
119  template <class T> SPIMat _Function_on_each_element( T (*f)(T const&), const SPIMat &A); // function handle template for operations
120  SPIMat sin(const SPIMat &A); // sin of matrix
121  SPIMat cos(const SPIMat &A); // cos of matrix
122  SPIMat acos(const SPIMat &A); // cos of matrix
123  SPIMat tan(const SPIMat &A); // tan of matrix
124  SPIMat abs(const SPIMat &A); // abs of matrix
125  SPIMat operator%(const SPIMat &A,const SPIMat &B); // A*B pointwise operation
126  SPIMat orthogonalize(const std::vector<SPIVec> &x); // create orthonormal basis from array of vectors
127 }
128 
129 
130 
131 #endif // SPIMAT_H
SPI::SPIMat::flag_init
PetscBool flag_init
flag if it has been initialized
Definition: SPIMat.hpp:32
SPI::operator^
SPIMat operator^(const PetscScalar a, const SPIMat &A)
Y=a^A operation.
Definition: SPIMat.cpp:578
SPI::sin
SPIMat sin(const SPIMat &A)
take the sin of each element in a matrix
Definition: SPIMat.cpp:2025
SPI::polyeig
std::tuple< PetscScalar, SPIVec > polyeig(const std::vector< SPIMat > &As, const PetscScalar target, const PetscReal tol=-1., const PetscInt max_iter=-1)
solve general polynomial eigenvalue problem of (A0 + A1*alpha + A2*alpha^2 + ...)*x = 0 and return a ...
Definition: SPIMat.cpp:1686
SPI::polyeig_init
std::tuple< PetscScalar, SPIVec > polyeig_init(const std::vector< SPIMat > &As, const PetscScalar target, const SPIVec &qr, const PetscReal tol=-1., const PetscInt max_iter=-1)
solve general polynomial eigenvalue problem of (A0 + A1*alpha + A2*alpha^2 + ...)*x = 0 and return a ...
Definition: SPIMat.cpp:1754
SPI::SPIMat::print
PetscInt print()
print mat to screen using PETSC_VIEWER_STDOUT_WORLD
Definition: SPIMat.cpp:498
SPI::SPIMat::rows
PetscInt rows
number of rows in mat
Definition: SPIMat.hpp:18
SPI::SPIMat::zero_rows
SPIMat & zero_rows(std::vector< PetscInt > rows)
set rows to zero
Definition: SPIMat.cpp:475
SPI::SPIMat::Init
PetscInt Init(PetscInt m, PetscInt n, std::string name="SPIMat")
Definition: SPIMat.cpp:55
SPI::SPIMat::operator+=
SPIMat & operator+=(const SPIMat &X)
MatAXPY, Y = 1.*X + Y operation.
Definition: SPIMat.cpp:225
SPI::SPIMat::mat
Mat mat
petsc Mat data
Definition: SPIMat.hpp:28
SPI::eig_right
std::tuple< PetscScalar, SPIVec > eig_right(const SPIMat &A, const SPIMat &B, const PetscScalar target, const PetscReal tol=-1, const PetscInt max_iter=-1)
solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha,...
Definition: SPIMat.cpp:1072
SPI::lstsq
SPIVec lstsq(const SPIMat &A, SPIVec &y)
solve least squares problem of A*x=y for skinny A matrix using SVD
Definition: SPIMat.cpp:905
SPI::Block2D
std::vector< std::vector< T > > Block2D
Definition: SPIMat.hpp:15
SPI::SPIMat::name
std::string name
Matrix name.
Definition: SPIMat.hpp:33
SPI::SPIMat::operator/=
SPIMat & operator/=(const PetscScalar a)
Y = Y/a operation.
Definition: SPIMat.cpp:327
SPI::SPIVec
Definition: SPIVec.hpp:13
SPI::SPIMat::operator-=
SPIMat & operator-=(const SPIMat &X)
Y = -1.*X + Y operation.
Definition: SPIMat.cpp:256
SPI::SPIMat::operator()
SPIMat & operator()()
assmelbe the matrix
Definition: SPIMat.cpp:218
SPI
Definition: SPIbaseflow.hpp:16
SPI::save
PetscInt save(const SPIMat &A, const std::string filename)
save matrix to filename in binary format (see Petsc documentation for format ) Format is (from Petsc ...
Definition: SPIMat.cpp:1924
SPI::SPIMat::ierr
PetscErrorCode ierr
ierr for various routines and operators
Definition: SPIMat.hpp:29
SPI::SPIMat::add
SPIMat & add(PetscInt m, PetscInt n, const PetscScalar v)
add a scalar value at position row m and column n
Definition: SPIMat.cpp:97
SPI::kron
SPIMat kron(const SPIMat &A, const SPIMat &B)
set kronecker inner product of two matrices
Definition: SPIMat.cpp:780
SPI::SPIMat::SPIMat
SPIMat(std::string _name="SPIMat")
constructor with no arguments (no initialization)
Definition: SPIMat.cpp:13
SPI::SPIMat::axpy
SPIMat & axpy(const PetscScalar a, const SPIMat &X)
MatAXPY function call to add a*X to the current mat.
Definition: SPIMat.cpp:232
SPI::SPIMat::operator-
SPIMat operator-() const
-X operation
Definition: SPIMat.cpp:278
SPI::SPIMat::zero_row_full
SPIMat & zero_row_full(const PetscInt row)
set a row to zero using dense format
Definition: SPIMat.cpp:465
SPI::SPIMat::diag
SPIVec diag()
get diagonal of matrix
Definition: SPIMat.cpp:445
SPI::_Function_on_each_element
SPIMat _Function_on_each_element(T(*f)(T const &), const SPIMat &A)
take the function of each element in a matrix, e.g. (*f)(A(i,j)) for each i,j
Definition: SPIMat.cpp:2011
SPI::SPIMat::operator*=
SPIMat & operator*=(const PetscScalar a)
Y = Y*a operation.
Definition: SPIMat.cpp:320
SPI::eye
SPIMat eye(const PetscInt n)
create, form, and return identity matrix of size n
Definition: SPIMat.cpp:697
SPI::inv
SPIMat inv(const SPIMat &A)
Solve linear system A*Ainv=B using MatMatSolve.
Definition: SPIMat.cpp:603
SPI::orthogonalize
std::vector< SPIVec > orthogonalize(std::vector< SPIVec > &x, SPIgrid1D &grid)
Definition: SPIgrid.cpp:1323
SPI::draw
PetscInt draw(const SPIMat &A)
draw nonzero structure of matrix
Definition: SPIMat.cpp:1994
SPI::SPIMat::~SPIMat
~SPIMat()
destructor to delete memory
Definition: SPIMat.cpp:508
SPI::tan
SPIMat tan(const SPIMat &A)
take the tan of each element in a matrix
Definition: SPIMat.cpp:2031
SPI::diag
SPIMat diag(const SPIVec &diag, const PetscInt k=0)
set diagonal of matrix
Definition: SPIMat.cpp:728
SPI::SPIMat::T
SPIMat & T()
Transpose the current mat.
Definition: SPIMat.cpp:405
SPI::SPIMat::operator+
SPIMat operator+(const SPIMat &X)
Y + X operation.
Definition: SPIMat.cpp:241
SPI::abs
SPIMat abs(const SPIMat &A)
take the abs of each element in a matrix
Definition: SPIMat.cpp:2045
SPI::SPIMat::eye_row
SPIMat & eye_row(const PetscInt row)
set a row to zero and set 1 in diagonal entry
Definition: SPIMat.cpp:458
SPIVec.hpp
SPI::eig_init_rights
std::tuple< std::vector< PetscScalar >, std::vector< SPIVec > > eig_init_rights(const SPIMat &A, const SPIMat &B, const std::vector< PetscScalar > targets, const std::vector< SPIVec > &qrs, PetscReal tol=-1, const PetscInt max_iter=-1)
solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha,...
Definition: SPIMat.cpp:1521
SPI::solve
SPIVec solve(const SPIMat &A, const SPIVec &b)
Solve linear system, Ax=b using solve(A,b) notation.
Definition: SPIMat.cpp:596
SPI::SPIMat::operator*
SPIMat operator*(const PetscScalar a)
Y*a operation.
Definition: SPIMat.cpp:285
SPI::eig_init
std::tuple< PetscScalar, SPIVec, SPIVec > eig_init(const SPIMat &A, const SPIMat &B, const PetscScalar target, const SPIVec &ql, const SPIVec &qr, PetscReal tol=-1, const PetscInt max_iter=-1)
solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha,...
Definition: SPIMat.cpp:1217
SPI::SPIMat::real
SPIMat & real()
take the real part of the vector
Definition: SPIMat.cpp:440
SPI::SPIMat::H
SPIMat & H()
Hermitian Transpose the current mat.
Definition: SPIMat.cpp:422
SPI::load
PetscInt load(SPIMat &A, const std::string filename)
load matrix from filename from binary format (works with save(SPIMat,std::string) function
Definition: SPIMat.cpp:1966
SPI::SPIMat
Definition: SPIMat.hpp:17
SPI::operator%
SPIMat operator%(const SPIMat &A, const SPIMat &B)
take the elementwise multiplication of each element in a matrix
Definition: SPIMat.cpp:2033
SPI::block
SPIMat block(const Block2D< SPIMat > Blocks)
set block matrices using an input array of size rows*cols. Fills rows first
Definition: SPIMat.cpp:1857
SPI::SPIMat::conj
SPIMat & conj()
elemenwise conjugate current matrix
Definition: SPIMat.cpp:435
SPI::SPIMat::operator/
SPIMat operator/(const PetscScalar a)
Z = Y/a operation.
Definition: SPIMat.cpp:334
SPI::acos
SPIMat acos(const SPIMat &A)
take the arccos of each element in a matrix
Definition: SPIMat.cpp:2029
SPI::SPIMat::set
SPIMat & set(PetscInt m, PetscInt n, const PetscScalar v)
Definition: SPIMat.cpp:77
SPI::operator*
SPIMat operator*(const PetscScalar a, const SPIMat A)
a*A operation to be equivalent to A*a
Definition: SPIMat.cpp:517
SPI::eig_init_right
std::tuple< PetscScalar, SPIVec > eig_init_right(const SPIMat &A, const SPIMat &B, const PetscScalar target, const SPIVec &qr, PetscReal tol=-1, const PetscInt max_iter=-1)
solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha,...
Definition: SPIMat.cpp:1368
SPI::SPIMat::operator=
SPIMat & operator=(const SPIMat &A)
Y=X with initialization of Y using MatConvert.
Definition: SPIMat.cpp:367
SPI::meshgrid
std::tuple< SPIMat, SPIMat > meshgrid(SPIVec &x, SPIVec &y)
Definition: SPIMat.cpp:1894
SPI::SPIMat::col
SPIVec col(const PetscInt i)
get column vector
Definition: SPIMat.cpp:489
SPI::eig
std::tuple< PetscScalar, SPIVec, SPIVec > eig(const SPIMat &A, const SPIMat &B, const PetscScalar target, const PetscReal tol=-1, const PetscInt max_iter=-1)
solve general eigenvalue problem of Ax = kBx and return a tuple of tie(PetscScalar alpha,...
Definition: SPIMat.cpp:928
SPI::cos
SPIMat cos(const SPIMat &A)
take the cos of each element in a matrix
Definition: SPIMat.cpp:2027
SPI::SPIMat::zero_row
SPIMat & zero_row(const PetscInt row)
set a row to zero
Definition: SPIMat.cpp:451
SPI::operator/
SPIVec operator/(const SPIVec &b, const SPIMat &A)
Solve linear system, Ax=b using b/A notation.
Definition: SPIMat.cpp:535
SPI::SPIMat::set_col
SPIMat & set_col(const PetscInt col, const SPIVec &v)
Definition: SPIMat.cpp:86
SPI::svd
std::tuple< std::vector< PetscReal >, std::vector< SPIVec >, std::vector< SPIVec > > svd(const SPIMat &A)
solve SVD problem of A=U*E*V^H for skinny A matrix
Definition: SPIMat.cpp:851
SPI::zeros
SPIMat zeros(const PetscInt m, const PetscInt n)
create, form, and return zeros matrix of size mxn
Definition: SPIMat.cpp:708
SPI::SPIMat::cols
PetscInt cols
number of columns in mat
Definition: SPIMat.hpp:19
SPI::SPIMat::eye_rows
SPIMat & eye_rows(std::vector< PetscInt > rows)
set rows to zero and set main diagonal to 1
Definition: SPIMat.cpp:482
Generated by   doxygen 1.8.17