-
Notifications
You must be signed in to change notification settings - Fork 0
/
Gradpo.cpp
59 lines (50 loc) · 1.02 KB
/
Gradpo.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
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#include "Gradpo.h"
#include <iostream>
#include <math.h>
#include <iostream>
#include <Eigen>
using Eigen::MatrixXd;
using namespace Eigen;
using namespace std;
//constructeur
Gradpo::Gradpo(MatrixXd A,VectorXd b) : A_(A),b_(b)
{
}
//gradient à pas optimal
VectorXd Gradpo::Solve()const
{
MatrixXd A=A_;
int n = A.rows();
VectorXd x(n);
for (int i = 0; i < n; i++)
{
x[i]=0.;
}
VectorXd r(n),b(b_);
//cout<<b<<endl;
r = b -A*x;
double alpha;
VectorXd rSuivant(n);
VectorXd xSuivant(n);
VectorXd z(n);
int j = 0;
double beta=r.norm();
int nb_iterat_=0;
while (beta > pow(10,-10))
{
z=A*r;
alpha= (r.dot(r) ) / (z.dot(r)) ;
xSuivant= x + alpha*r;
rSuivant=r-alpha*z;
x=xSuivant ;
//cout<<x<<endl;
//cout<<"----------------------------------------"<<endl;
r=rSuivant;
beta=r.norm();
nb_iterat_=nb_iterat_ +1;
j++;
}
cout<<"le nombre d'itération : "<<nb_iterat_<<endl;
cout<<"----------------Gradient à pas optimal------------------------"<<endl;
return x;
}