-
Notifications
You must be signed in to change notification settings - Fork 0
/
matvec.h
92 lines (71 loc) · 3.44 KB
/
matvec.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
#ifndef MATVEC_H
#define MATVEC_H
#include <Eigen/Dense>
void sparseAvec(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
// a1.mtx 9*9 diag
void a1vec(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
void b1vec(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& bvecs);
// toys
// template <typename Derived>
void avec(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
void bvec(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& bvecs);
// void precnd(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs);
// preconditioner for avec
void mprec(int n, int m, const Eigen::MatrixXd& x, Eigen::MatrixXd& px, double shift);
void tridiagA_precnd_a(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs, double shift);
void avec_Si2(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
void precnd_Si2(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs, double shift);
void avec_Na5(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
void precnd_Na5(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs, double shift);
void avec_Si5H12(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs);
// void precnd_Si5H12(int n, int m, const Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs);
// template <typename Derived>
// void avec(int n, int m, Eigen::DenseBase<Derived>& vecs, Eigen::DenseBase<Derived>& avecs){
// // void avec(int n, int m, Eigen::MatrixXd& vecs, Eigen::MatrixXd& avecs){
// // check vecs and avecs are of size(n,m)
// if(vecs.rows() != n || vecs.cols() != m){
// std::cerr << "vecs must be of size (n,m)"; return;
// }
// if(avecs.rows() != n || avecs.cols() != m){
// std::cerr << "bvecs must be of size (n,m)"; return;
// }
// // Assume 'a' is a global variable or class member matrix that is already defined and initialized.
// static Eigen::MatrixXd a(n, n);
// for (int i = 0; i < n; ++i) {
// a(i, i) = static_cast<double>(i + 1) + 1.0;
// for (int j = 0; j < i; ++j) {
// a(j,i) = 1.0 / static_cast<double>(i + j);
// a(i,j) = a(j,i);
// }
// }
// // Perform matrix-vector multiplication for each column of x
// for (int icol = 0; icol < m; ++icol) {
// avecs.col(icol) = a * vecs.col(icol);
// }
// }
// template <typename Derived>
// void bvec(int n, int m, Eigen::DenseBase<Derived>& vecs, Eigen::DenseBase<Derived>& bvecs){
// // void bvec(int n, int m, Eigen::MatrixXd& vecs, Eigen::MatrixXd& bvecs){
// // check vecs and bvecs are of size(n,m)
// if(vecs.rows() != n || vecs.cols() != m){
// std::cerr << "vecs must be of size (n,m)"; return;
// }
// if(bvecs.rows() != n || bvecs.cols() != m){
// std::cerr << "bvecs must be of size (n,m)"; return;
// }
// bvecs = vecs;
// }
// // identity preconditioner
// template <typename Derived>
// void precnd(int n, int m, Eigen::DenseBase<Derived>& vecs, Eigen::DenseBase<Derived>& tvecs){
// // void precnd(int n, int m, Eigen::MatrixXd& vecs, Eigen::MatrixXd& tvecs){
// // check vecs and tvecs are of size(n,m)
// if(vecs.rows() != n || vecs.cols() != m){
// std::cerr << "vecs must be of size (n,m)"; return;
// }
// if(tvecs.rows() != n || tvecs.cols() != m){
// std::cerr << "tvecs must be of size (n,m)"; return;
// }
// tvecs = vecs;
// }
#endif // MATVEC_H