C#，码海拾贝（35）——求“实对称矩阵““特征值与特征向量“的“雅可比法“之C#源代码

在C#中，我们可以使用雅可比法（Jacobi method）来求解实对称矩阵的特征值与特征向量。下面是一个简单的示例代码：

csharpusing System;

class JacobiMethod{
 static void Main()
 {
 // 定义一个实对称矩阵 double[,] matrix = {
 {4, -2,2},
 {-2,3, -1},
 {2, -1,3}
 };

 // 调用Jacobi方法求解特征值与特征向量 Jacobi(matrix);
 }

 static void Jacobi(double[,] matrix)
 {
 int n = matrix.GetLength(0);
 double[,] eigenVectors = new double[n, n];
 double[] eigenValues = new double[n];

 // 初始化特征向量矩阵为单位矩阵 for (int i =0; i < n; i++)
 {
 eigenVectors[i, i] =1;
 }

 // 迭代求解 for (int iter =0; iter <100; iter++)
 {
 // 寻找最大非对角元素 int p =0, q =1;
 double max = Math.Abs(matrix[p, q]);
 for (int i =0; i < n; i++)
 {
 for (int j = i +1; j < n; j++)
 {
 if (Math.Abs(matrix[i, j]) > max)
 {
 max = Math.Abs(matrix[i, j]);
 p = i;
 q = j;
 }
 }
 }

 // 计算旋转角度 double theta =0.5 * Math.Atan2(2 * matrix[p, q], matrix[q, q] - matrix[p, p]);

 // 构造旋转矩阵 double[,] rotation = new double[n, n];
 for (int i =0; i < n; i++)
 {
 rotation[i, i] =1;
 }
 rotation[p, p] = Math.Cos(theta);
 rotation[p, q] = -Math.Sin(theta);
 rotation[q, p] = Math.Sin(theta);
 rotation[q, q] = Math.Cos(theta);

 // 更新矩阵和特征向量 matrix = MatrixMultiply(MatrixTranspose(rotation), MatrixMultiply(matrix, rotation));
 eigenVectors = MatrixMultiply(eigenVectors, rotation);
 }

 // 提取特征值和特征向量 for (int i =0; i < n; i++)
 {
 eigenValues[i] = matrix[i, i];
 }

 // 输出结果 Console.WriteLine("Eigenvalues:");
 foreach (var value in eigenValues)
 {
 Console.WriteLine(value);
 }

 Console.WriteLine("Eigenvectors:");
 for (int i =0; i < n; i++)
 {
 for (int j =0; j < n; j++)
 {
 Console.Write(eigenVectors[i, j] + " ");
 }
 Console.WriteLine();
 }
 }

 static double[,] MatrixMultiply(double[,] a, double[,] b)
 {
 int n = a.GetLength(0);
 double[,] result = new double[n, n];
 for (int i =0; i < n; i++)
 {
 for (int j =0; j < n; j++)
 {
 for (int k =0; k < n; k++)
 {
 result[i, j] += a[i, k] * b[k, j];
 }
 }
 }
 return result;
 }

 static double[,] MatrixTranspose(double[,] matrix)
 {
 int n = matrix.GetLength(0);
 double[,] result = new double[n, n];
 for (int i =0; i < n; i++)
 {
 for (int j =0; j < n; j++)
 {
 result[i, j] = matrix[j, i];
 }
 }
 return result;
 }
}