急！急急需求n*m矩阵的协方差矩阵的特征值和特征向量的vc源代码！！

supercooling

各位大侠

supercooling

<>小弟在线等!</P><>谢谢!</P>

supercooling

<>只要给小弟一个求矩阵特征值和特征向量的的vc源代码就可以！</P>

knight

/*
  Name: QR Method
  Copyright: Knight Studio
  Author: Ran NI   
  Date: 10-07-04 11:46
  Description: QR Method for Solving Matrix's Eigenvalue
*/
#include &lt;math.h&gt;
#include &lt;vector&gt;
#include &lt;complex&gt;
#include &lt;iostream&gt;
#include &lt;fstream&gt;

using namespace std;

vector&lt;vector&lt;double&gt; &gt; Inverse(vector&lt;vector&lt;double&gt; &gt; a,bool &amp;judge)
{
         vector&lt;vector&lt;double&gt; &gt; b(a);
         if(b.size()!=b[0].size())
         {
                  judge=false;
                  return b;
         }
         else
         {
                  double temp=1;
                 for(int k=0;k&lt;(int)b.size();++k)
                  {
                          if(k&lt;((int)b.size()-1))
                          {
                                   int tk=k;
                                   double tem=fabs(b[k][k]);
                                   for(int i=k;i&lt;(int)b.size();++i)
                                   {
                                            if(fabs(b[k])&gt;tem)
                                            {
                                                      tk=i;
                                                      tem=fabs(b[k]);
                                            }
                                   }
                                   if(tk!=k)
                                   {
                                            swap(b[tk],b[k]);
                                   }
                                   if(b[k][k]==0)
                                   {
                                            temp=0;
                                   }
                                   else
                                   {
                                            for(int i=k+1;i&lt;(int)b.size();++i)
                                            {
                                                     tem=b[k];
                                                     for(int j=k;j&lt;(int)b.size();++j)
                                                     {
                                                              b[j]-=(b[k][j]*tem/b[k][k]);
                                                     }
                                            }
                                   }
                          }
                  }
                  if(temp!=0)
                  {
                          for(int i=0;i&lt;(int)b.size();++i)
                          {
                                  temp=temp*b;
                          }
                  }
                  if(temp==0)
                  {
                          judge=false;
                          return a;
                  }
                  else
                  {
                          vector&lt;vector&lt;double&gt; &gt; I;
                          for(int i=0;i&lt;(int)a.size();++i)
                          {
                                   vector&lt;double&gt; z(a.size(),0);
                                   z=1;
                                   I.push_back(z);
                          }
                          for(int k=0;k&lt;(int)a.size();++k)
                          {
                                   int tk=k;
                                   double tem=fabs(a[k][k]);
                                   for(int i=k;i&lt;(int)a.size();++i)
                                   {
                                            if(fabs(a[k])&gt;tem)
                                            {
                                                      tk=i;
                                                      tem=fabs(a[k]);
                                            }
                                   }
                                   if(tk!=k)
                                   {
                                            swap(a[tk],a[k]);
                                            swap(I[tk],I[k]);
                                   }
                                   tem=a[k][k];
                                   for(int i=0;i&lt;(int)a[k].size();++i)
                                   {
                                            a[k]/=tem;
                                            I[k]/=tem;
                                   }
                                   for(int i=0;i&lt;(int)a.size();++i)
                                   {
                                            if(i!=k)
                                            {
                                                     double temp=a[k];
                                                     for(int j=0;j&lt;(int)a.size();++j)
                                                     {
                                                              a[j]-=(a[k][j]*temp);
                                                              I[j]-=(I[k][j]*temp);
                                                     }
                                            }
                                   }
                          }
                          judge=true;
                          return I;
                  }
         }
}

vector&lt;vector&lt;double&gt; &gt; operator - (vector&lt;vector&lt;double&gt; &gt; a,vector&lt;vector&lt;double&gt; &gt; b)
{
         bool dd=true;
         if(a.size()!=b.size())
         {
                  dd=false;
         }
         else
         {
                 for(int i=0;i&lt;(int)a.size();++i)
                  {
                          if(a.size()!=b.size())
                          {
                                   dd=false;
                          }
                          else
                          {
                                   for(int j=0;j&lt;(int)a.size();++j)
                                   {
                                            a[j]-=b[j];
                                   }
                          }
                  }
         }
         vector&lt;vector&lt;double&gt; &gt; ss;
         if(dd==false)
         {
                  return ss;
         }
         else
         {
                  ss=a;
                  return ss;
         }
}

vector&lt;vector&lt;double&gt; &gt; operator * (double a,vector&lt;vector&lt;double&gt; &gt; b)
{
         for(int i=0;i&lt;(int)b.size();++i)
         for(int j=0;j&lt;(int)b.size();++j)
         {
                  b[j]*=a;
         }
         return b;
}

vector&lt;double&gt; operator * (vector&lt;vector&lt;double&gt; &gt; a,vector&lt;double&gt; b)
{
         vector&lt;double&gt; c;
         for(int i=0;i&lt;(int)a.size();++i)
         {
                  double s=0;
                 for(int j=0;j&lt;(int)b.size();++j)
                  {
                          s+=a[j]*b[j];
                  }
                  c.push_back(s);
         }
         return c;
}

vector&lt;vector&lt;double&gt; &gt; operator * (vector&lt;vector&lt;double&gt; &gt; a,vector&lt;vector&lt;double&gt; &gt; b)
{
         if(a[0].size()!=b.size())
         {
                 vector&lt;vector&lt;double&gt; &gt; ss;
                  return ss;
         }
         else
         {
                 vector&lt;vector&lt;double&gt; &gt; ss(a.size());
                 for(int i=0;i&lt;(int)ss.size();++i)
                  {
                          for(int j=0;j&lt;(int)b[0].size();++j)
                          {
                                   double temp=0;
                                   for(int k=0;k&lt;(int)b.size();++k)
                                   {
                                             temp+=(a[k]*b[k][j]);
                                   }
                                   ss.push_back(temp);
                          }
                  }
                  return ss;
         }
}

vector&lt;vector&lt;double&gt; &gt; operator * (vector&lt;double&gt; a,vector&lt;double&gt; b)
{
         vector&lt;vector&lt;double&gt; &gt; result(a.size());
         for(int i=0;i&lt;(int)result.size();++i)
         {
                 for(int j=0;j&lt;(int)b.size();++j)
                  {
                          result.push_back(a*b[j]);
                  }
         }
         return result;
}

knight

int sgn(double x)
{
         if(x&gt;0)
         {
                  return 1;
         }else if(x==0)
         {
                  return 0;
         }else
         {
                  return -1;
         }
}

vector&lt;vector&lt;double&gt; &gt; Hessenberg(vector&lt;vector&lt;double&gt; &gt; A)
{
         for(int r=0;r&lt;(int)A.size()-2;++r)
         {
                 vector&lt;double&gt; ar(A.size(),0);
                 for(int i=0;i&lt;(int)A.size();++i)
                  {
                          ar=A[r];
                  }
                  double c=0;
                 for(int i=r+1;i&lt;(int)A.size();++i)
                  {
                          c+=pow(ar,2);
                  }
                  c=sqrt(c);
                  c=(-c*sgn(ar[r+1]));
                 double p=sqrt(2*c*(c-ar[r+1]));
                 vector&lt;double&gt; u(A.size(),0);
                 for(int i=r+1;i&lt;(int)A.size();++i)
                  {
                          if(i==r+1)
                          {
                                   u=(ar-c)/p;
                          }
                          else
                          {
                                   u=ar/p;
                          }
                  }
                 vector&lt;vector&lt;double&gt; &gt; I;
                 for(int i=0;i&lt;(int)A.size();++i)
                  {
                          vector&lt;double&gt; z(A.size(),0);
                          z=1;
                          I.push_back(z);
                  }
                 vector&lt;vector&lt;double&gt; &gt; H=I-2*(u*u);
                  bool s;
                  A=H*A*Inverse(H,s);
         }
         return A;
}

vector&lt;vector&lt;double&gt; &gt; QR(vector&lt;vector&lt;double&gt; &gt; A,vector&lt;vector&lt;double&gt; &gt; &amp;Q)
{
         vector&lt;vector&lt;double&gt; &gt; I;
         for(int i=0;i&lt;(int)A.size();++i)
         {
                 vector&lt;double&gt; z(A.size(),0);
                  z=1;
                  I.push_back(z);
         }         
         for(int i=0;i&lt;(int)A.size()-1;++i)
         {
                 double theta=atan(A[i+1]/A);
                 vector&lt;vector&lt;double&gt; &gt; P;
                 for(int r=0;r&lt;(int)A.size();++r)
                  {
                         vector&lt;double&gt; z(A.size(),0);
                         z[r]=1;
                         P.push_back(z);
                  }
                  P[i+1]=sin(theta);
                  P=cos(theta);
                  P[i+1][i+1]=cos(theta);
                  P[i+1]=(-sin(theta));
                  I=P*I;
                 vector&lt;double&gt; aa(A);
                 vector&lt;double&gt; aa1(A[i+1]);
                 for(int j=i;j&lt;(int)A.size();++j)
                  {
                          A[j]=aa[j]*cos(theta)+aa1[j]*sin(theta);
                          A[i+1][j]=(-aa[j]*sin(theta)+aa1[j]*cos(theta));
                  }
         }
         bool s;
         Q=Inverse(I,s);
         return A;
}

void Print(vector&lt;vector&lt;double&gt; &gt; a,string s)
{
         ofstream out(s.c_str());
         for(int i=0;i&lt;(int)a.size();++i)
         {
                 for(int j=0;j&lt;(int)a.size();++j)
                  {
                          out&lt;&lt;a[j]&lt;&lt;"\t";
                  }
                  out&lt;&lt;endl;
         }
}

double Delta(vector&lt;vector&lt;double&gt; &gt; a)
{
        double ss=0;
        for(int i=0;i&lt;(int)a.size();++i)
        {
               for(int j=0;j&lt;(int)a.size();++j)
                {
                       if(i!=j)
                       {
                                if(fabs(a[j])&gt;ss)
                                {
                                       ss=fabs(a[j]);
                                }
                       }
                }
        }
        return ss;
}

double Delta1(vector&lt;vector&lt;double&gt; &gt; A)
{
        double d=0;
        for(int i=0;i&lt;(int)A.size()-1;++i)
        {
                if(fabs(A[i+1])&gt;d)
                {
                        d=fabs(A[i+1]);
                }
                if(fabs(A[i+1])&gt;d)
                {
                        d=fabs(A[i+1]);
                }
        }
        return d;
}

vector&lt;complex&lt;double&gt; &gt; Namta(vector&lt;vector&lt;double&gt; &gt; A,double delta)
{
        double delta1=0;
        while((Delta(A)&gt;=delta)&amp;&amp;(fabs(Delta1(A)-delta1)&gt;=delta))
        {
                 delta1=Delta1(A);
                 A=Hessenberg(A);
                vector&lt;vector&lt;double&gt; &gt; Q;
                vector&lt;vector&lt;double&gt; &gt; R;
                 R=QR(A,Q);
                 A=R*Q;
        }
        string s("An.txt");
        Print(A,s);
        vector&lt;complex&lt;double&gt; &gt; namta;
        if(Delta(A)&lt;delta)
        {
                for(int i=0;i&lt;(int)A.size();++i)
                 {
                         complex&lt;double&gt; dd(A,0);
                         namta.push_back(dd);
                 }
        }
        else
        {
                 int r=0;
                 while(r&lt;A.size()-1)
                 {
                        if(fabs(A[r+1][r])&gt;=delta)
                        {
                                double b=-(A[r][r]+A[r+1][r+1]);
                                double c=(A[r][r]*A[r+1][r+1]-A[r][r+1]*A[r+1][r]);
                                if(pow(b,2)-4*c&lt;0)
                                {
                                        complex&lt;double&gt; d1(-b/2,sqrt(4*c-pow(b,2))/2);
                                        complex&lt;double&gt; d2(-b/2,-sqrt(4*c-pow(b,2))/2);
                                        namta.push_back(d1);
                                        namta.push_back(d2);
                                }
                                else
                                {
                                        complex&lt;double&gt; d1(-b/2+sqrt(pow(b,2)-4*c)/2,0);
                                        complex&lt;double&gt; d2(-b/2-sqrt(pow(b,2)-4*c)/2,0);
                                        namta.push_back(d1);
                                        namta.push_back(d2);
                                }
                                r+=2;
                        }
                        else
                        {
                                complex&lt;double&gt; d(A[r][r],0);
                                namta.push_back(d);
                                ++r;
                        }
                 }
                 if(r==A.size()-1)
                 {
                        complex&lt;double&gt; d(A[r][r],0);
                        namta.push_back(d);
                 }
        }
        return namta;
}                 

void Print(vector&lt;complex&lt;double&gt; &gt; a)
{
        for(int i=0;i&lt;(int)a.size();++i)
        {
                cout&lt;&lt;a&lt;&lt;endl;
        }
}

void Print(vector&lt;double&gt; a)
{
        for(int i=0;i&lt;(int)a.size();++i)
        {
                cout&lt;&lt;a&lt;&lt;"\t";
        }
        cout&lt;&lt;endl;
}

double Delta(vector&lt;double&gt; a,vector&lt;double&gt; b)
{
        double x=0;
        for(int i=0;i&lt;(int)a.size();++i)
        {
                if(fabs(a-b)&gt;x)
                {
                       x=fabs(a-b);
                }
        }
        return x;
}

double Max(vector&lt;double&gt; a)
{
        double x=0;
        int n=0;
        for(int i=0;i&lt;(int)a.size();++i)
        {
                if(fabs(a)&gt;x)
                {
                       x=fabs(a);
                       n=i;
                }
        }
        return a[n];
}

vector&lt;double&gt; operator / (vector&lt;double&gt; a,double b)
{
        for(int i=0;i&lt;(int)a.size();++i)
        {
                a/=b;
        }
        return a;
}

vector&lt;double&gt; ComputeVector(vector&lt;vector&lt;double&gt; &gt; A,complex&lt;double&gt; namta,double delta)
{
        vector&lt;vector&lt;double&gt; &gt; I;
        for(int i=0;i&lt;(int)A.size();++i)
        {
                vector&lt;double&gt; z(A.size(),0);
                z=1;
                I.push_back(z);
        }
        A=A-namta.real()*I;
        string s("A-1.txt");
        bool d;
        A=Inverse(A,d);
        Print(A,s);
        vector&lt;double&gt; z(A.size(),1);
        vector&lt;double&gt; z1;
        do
        {
                z1=z;
                vector&lt;double&gt; y=A*z;
                double m=Max(y);
                z=y/m;
        }while(Delta(z,z1)&gt;=delta);
        return z;
}

int main()
{
        ifstream ina("A.txt");
        cout&lt;&lt;"输入矩阵的维数："&lt;&lt;endl;
        int n;
        cin&gt;&gt;n;
        vector&lt;vector&lt;double&gt; &gt; A(n);
        for(int i=0;i&lt;n;++i)
        {
                 for(int j=0;j&lt;n;++j)
                 {
                          double x;
                          if(ina&gt;&gt;x)
                          {
                                   A.push_back(x);
                          }
                 }
        }
        cout&lt;&lt;"输入误差:"&lt;&lt;endl;
        double delta;
        cin&gt;&gt;delta;
        string s("namta.txt");
        vector&lt;complex&lt;double&gt; &gt; namta=Namta(A,delta);
        cout&lt;&lt;"特征值："&lt;&lt;endl;
        Print(namta);
        cout&lt;&lt;namta[1]&lt;&lt;"所对应的特征向量："&lt;&lt;endl;
        Print(ComputeVector(A,namta[1],0.00001));
        system("AUSE");
}

knight

以上是我一年前写的QR方法的程序，写得比较简单，不过算法实现都有了，你只要看看，改改就可以用了，用了一些新规范,也许VC6会编译不通过,用VC2003或者更新的编译器应该没有问题.

supercooling

<>各位大侠，小弟现在急需求n*m矩阵的协方差矩阵的特征值和特征向量的vc源代码，请高手慷慨解囊，</P><><a href="mailtcjhjj@sohu.com" target="_blank" >cjhjj@sohu.com</A></P><>感激涕零！！</P>

supercooling

<>小弟在线等!</P><>谢谢!</P>

supercooling

<>谢谢！</P><>感激涕零！！</P>