使用雅可比方法求矩阵A的全部特征值和特征向量C语言代码

以下是使用雅可比方法求矩阵A的全部特征值和特征向量的C语言代码：

#include <stdio.h>
#include <math.h>

#define N 3  //矩阵A的维数

int main()
{
    double A[N][N] = {{2,-1,0},{-1,2,-1},{0,-1,2}};  //矩阵A
    double eps = 1e-8;  //精度
    double sum;  //主对角线上元素的和
    int i,j,k,p,q;  //循环计数器
    double max;  //最大非主对角线元素的值
    double c,s,t;  //旋转角度sin、cos和tan

    double X[N][N];  //存放特征向量的矩阵
    for(i = 0; i < N; i++)  //初始化X为单位矩阵
        for(j = 0; j < N; j++)
            X[i][j] = (i == j ? 1 : 0);

    while(1)
    {
        //求出矩阵A的主对角线元素之和
        sum = 0;
        for(i = 0; i < N; i++)
            sum += A[i][i];

        //如果最大非主对角线元素的绝对值小于精度，则停止迭代
        max = 0;
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                if(i != j && fabs(A[i][j]) > max)
                {
                    max = fabs(A[i][j]);
                    p = i;
                    q = j;
                }
        if(max < eps)
            break;

        //计算旋转角度sin、cos和tan
        if(A[p][p] == A[q][q])
            t = 1;
        else
            t = (2*A[p][q]) / (A[p][p] - A[q][q]);
        double sign = (t >= 0 ? 1 : -1);
        c = 1 / sqrt(1 + t*t);
        s = sign * t * c;

        //构造Givens矩阵
        double G[N][N];
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                G[i][j] = (i == j ? 1 : 0);
        G[p][p] = c;
        G[p][q] = s;
        G[q][p] = -s;
        G[q][q] = c;

        //进行相似变换，更新矩阵A和特征向量矩阵X
        double temp[N][N];
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                temp[i][j] = 0;
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                for(k = 0; k < N; k++)
                    temp[i][j] += G[k][i] * A[k][j];
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                A[i][j] = 0;
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                for(k = 0; k < N; k++)
                    A[i][j] += temp[i][k] * G[k][j];
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                temp[i][j] = 0;
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                for(k = 0; k < N; k++)
                    temp[i][j] += G[k][i] * X[k][j];
        for(i = 0; i < N; i++)
            for(j = 0; j < N; j++)
                X[i][j] = temp[i][j];
    }

    //输出特征值和特征向量
    printf("Eigenvalues:\n");
    for(i = 0; i < N; i++)
        printf("%f\n", A[i][i]);
    printf("\nEigenvectors:\n");
    for(i = 0; i < N; i++)
    {
        for(j = 0; j < N; j++)
            printf("%f ", X[j][i]);
        printf("\n");
    }

    return 0;
}