void QR(const Matrix& m, Matrix& Q,Matrix& R)
{
int n = 0;
int row = m.rows();
int col = m.cols();
if (row>col){ n =col; }
else{ n =row; }
Matrix mtemp = m;
Matrix q1(row, col);
Q.resize(row,row);
R.resize(row,col);
int i, j, k, r, l;
double temp, sum, dr, cr, hr;
double *ur = new double[n];
double *pr = new double[n];
double *wr = new double[n];
for (i = 0; i<n; i++)//定义单位矩阵
for (j = 0; j<n; j++)
{
if (i == j)Q[i][j] = 1;
else Q[i][j] = 0;
}
for (r = 0; r<n; r++)
{
temp = 0;
for (k = r; k<n; k++)
temp += fabs(mtemp[k][r]);
if (temp >= 0.0)
{
sum = 0;
for (k = r; k<n; k++)
sum += mtemp[k][r] * mtemp[k][r];
dr = sqrt(sum);
if (mtemp[r][r]>0.0) l = -1;
else l = 1;
cr = l*dr;
hr = cr*(cr - mtemp[r][r]);
for (i = 0; i<n; i++)//定义ur
{
if (i<r)ur[i] = 0;
if (i == r)ur[i] = mtemp[r][r] - cr;
if (i>r)ur[i] = mtemp[i][r];
}
for (i = 0; i<n; i++)//定义wr
{
sum = 0;
for (j = 0; j<n; j++)
sum += Q[i][j] * ur[j];
wr[i] = sum;
}
for (i = 0; i<row; i++)//定义qr
for (j = 0; j<col; j++)
{
q1[i][j] = Q[i][j] - wr[i] * ur[j] / hr;
}
for (i = 0; i<row; i++)//定义qr+1
for (j = 0; j<col; j++)
{
Q[i][j] = q1[i][j];
};
for (i = 0; i<col; i++)//定义pr
{
sum = 0;
for (j = 0; j<n; j++)
sum += mtemp[j][i] * ur[j];
pr[i] = sum / hr;
};
for (i = 0; i<row; i++)
for (j = 0; j<col; j++)
{
mtemp[i][j] = mtemp[i][j] - ur[i] * pr[j];
};
};
};
for (i = 0; i<row; i++)
for (j = 0; j<col; j++)
{
if (fabs(mtemp[i][j])<0.0) mtemp[i][j] = 0;
};
for (i = 0; i<row; i++)
for (j = 0; j<col; j++)
{
R[i][j] = mtemp[i][j];
}
}
什么问题呢?难道是分享资料不成!
哈,楼主还个代码写的还可以,就是不知道什么问题
有问题吗?还是只是来分享代码呢?