void SolutionSlau::QRalgo(std::vector <std::vector <double>> &Q,
std::vector <std::vector <double>> &R, std::vector <double> &x)
{
int size = A.size();
R = A;
std::vector<std::vector<double> > E(size, std::vector<double>(size, 0));
for (int i = 0; i < size; ++i)
E[i][i] = 1;
for (int k = 0; k < size - 1; ++k)
{
double sum = 0;
for (int i = k; i < size; ++i)
{
sum += R[i][k] * R[i][k];
}
double betta = sign(-R[k][k])*sqrt(sum);
double mu = 1 / sqrt(2 * betta*(betta - R[k][k]));
std::vector<double> w;
for (int i = 0; i < size; ++i)
{
if (i < k)
w.push_back(0);
else if (i == k)
w.push_back(R[i][i] - betta);
else
w.push_back(R[i][k]);
}
std::vector<std::vector<double> > H(size, std::vector<double>(size, 0));
for (int i = 0; i < size; ++i)
{
for (int j = 0; j < size; ++j)
{
H[i][j] = E[i][j] - w[i] * w[j] * 2 * mu*mu;
}
}
std::vector<std::vector<double> > res(size, std::vector<double>(size, 0));
multiply(H, R, res);
R = res;
}
std::vector<std::vector<double> > temp(size, std::vector<double>(size, 0));
obrMatrix(A, E);
multiply(R, E, Q);
//transposition(Q);
std::vector<double> y(size);
for (int i = 0; i < size; ++i)
{
y[i] = 0;
for (int j = 0; j < size; ++j)
y[i] += Q[i][j] * b[j];
}
for (int i = size - 1; i >= 0; --i)
{
for (int j = i + 1; j < size; ++j)
{
y[i] -= R[i][j] * x[j];
}
x[i] = y[i] / R[i][i];
}
}