VC++6新手,这个算法在TC2.0下经验证没有问题,拿到MFC中应用就会出错,怎么回事啊??
一个矩阵QR分解的算法,已经在TC2.0下验证了没问题,为了应用改到MFC中就出现异常.
具体表现在设置矩阵的有效列数时,当小于8时分解一切正常,等于8时一按button1按钮启动运算程序就自行退出,当等于9时运行出错,显示 "0x00000000 "指令引用的 "0x00000000 "内存.该内存不能为 "read ".
我不是学计算机专业的,这个问题已经困扰我好久,项目没法进行,各位高手救命啊!!!
程序如下:
void CTestDlg::OnButton1()
{
double B[5][10]={{1,2,1,-1,5,9,4,-1,2,3},{1,1,-1,2,9,-6,7,4,5,8},{-1,0,0,1,2,5,4,7,6,0},{1,-2,3,4,-5,6,-7,8,9,5},{2,3,4,5,-6,7,-8,9,10,-1}};
int n=9;// *****************就是这个值，小于等于7是完全正常，大于7时就会出错！！！！！！*************************
double Bet[10]={0};
double V[5][10]={0};
double H[10][10]={0};
double Q[10][10]={0};
double Q_temp[10][10]={0};
double R[10][5]={0};
double temp=0.0;
for(int i=0;i <n;i++) //初始化Q_temp为单位阵
Q_temp[i][i]=1;
for(int k=0;k <n-1;k++) //QR分解
{
Householder(B[k]+k,n-k,Bet+k,V[k]+k); //调用函数计算H阵的各参数
for(int i=0;i <n;i++) //初始化H为单位阵
{
for(int j=0;j <n;j++)
{
if(i==j) H[i][j]=1;
else H[i][j]=0;
}
}
for(i=0;i <n;i++) //计算H阵
{
for(int j=0;j <n;j++)
{
if(i+k <n && j+k <n) H[i+k][j+k]-=Bet[k]*V[k][k+i]*V[k][k+j];
else break;
}
}
for(i=0;i <n;i++) //计算Q阵(Q_temp*H)
{
for(int j=0;j <n;j++)
{
temp=0;
for(int p=0;p <n;p++)
temp+=Q_temp[i][p]*H[p][j];
Q[i][j]=temp;
}
}
if(k <n-1)
{
for(i=0;i <n;i++) //更新Q_temp
{
for(int j=0;j <n;j++)
Q_temp[i][j]=Q[i][j];
}
}
for(i=0;i <n;i++) //计算H*B直至得到R阵
{
for(int j=0;j <5;j++) //5
{
temp=0;
for(int p=0;p <n;p++)
temp+=H[i][p]*B[j][p];
R[i][j]=temp;
}
}
if(k <n-1)
{
for(i=0;i <n;i++) //更新B
{
for(int j=0;j <5;j++) //5
B[j][i]=R[i][j];
}
}
}
void Householder(double *x,int n,double *bet,double *v) //Householder变换函数
{
int i;
double alph,norm,norm_temp;
norm=fabs(x[0]);
for(i=1;i <n;i++)
if(norm <fabs(x[i])) norm=fabs(x[i]); //x的无穷范数
if(norm==0) //全零报错
{
//cout < < "Fail, because x is zero. " < <endl;
//exit(0);
}
for(i=0;i <n;i++) //无穷范数规范化
{
x[i]/=norm;
}
norm_temp=norm; //x发生变化，暂存以恢复x
norm=0.0;
for(i=1;i <n;i++) //求除首元素外的平方和
norm+=x[i]*x[i];
v[0]=1.0; //构造v，首元素为零，其余同规范化的x
for(i=1;i <n;i++)
v[i]=x[i];
if(norm==0) *bet=0.0;
else
{
alph=sqrt(x[0]*x[0]+norm); //计算x的2范数
if(x[0] <=0) v[0]=x[0]-alph;
else v[0]=-norm/(x[0]+alph);
*bet=2.0*v[0]*v[0]/(norm+v[0]*v[0]);
alph=v[0];
for(i=0;i <n;i++) //以2范数规范化v
v[i]=v[i]/alph;
}
for(i=0;i <n;i++) //恢复x
x[i]*=norm_temp;
}
[解决办法]
double B[5][10]={{...
...
for(int k=0;k <n-1;k++) //QR分解
{
Householder(B[k]+k,n-k,Bet+k,V[k]+k); //调用函数计算H阵的各参数
当 n > 6 的时候，
Householder(B[k]+k,n-k,Bet+k,V[k]+k); 中的 B[k] 存在越界。
