-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathserialSchroedinger.cpp
44 lines (36 loc) · 1.49 KB
/
serialSchroedinger.cpp
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
#include "io.hpp"
#include <chrono>
int main() {
unsigned length = 1000;
double a,b,h;
Mat mat = Mat(length-2);
a = 0;
b = 1;
h = (b-a)/(length-1);
for(unsigned i = 0; i < length-2; i++) {
if(i == 0) {
mat.insert(i,i)= 2/pow(h,2);
mat.insert(i,i+1) = -1/pow(h,2);
}
else if(i == length -3) {
mat.insert(i,i-1) = -1/pow(h,2);
mat.insert(i,i)= 2/pow(h,2);
}
else {
mat.insert(i,i-1) = -1/pow(h,2);
mat.insert(i,i)= 2/pow(h,2);
mat.insert(i,i+1) = -1/pow(h,2);
}
}
EigenValue solver;
solver.setEpsilon(20);
std::chrono::steady_clock::time_point begin = std::chrono::steady_clock::now();
solver.solveJacobi(20,&mat);
std::chrono::steady_clock::time_point end = std::chrono::steady_clock::now();
// std::cout << "EigenValues: error relativ" << std::endl;
std::vector<double> buff;
for(unsigned i = 0; i < length-2; i++) buff.push_back(solver.getEigenValue(i));
std::sort(buff.begin(), buff.end());
for(unsigned i = 0; i < length-2; i++) std::cout << 100*(buff[i]-pow((i+1)*M_PI,2))/pow((i+1)*M_PI,2)<< " %" << std::endl;
std::cout << "Time difference = " << std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count() << " Milliseconds" << std::endl;
}