实验9 Jacobi旋转法.docVIP

《数值分析》实验9 一．实验名称：Jacobi旋转法 二、实验目的： （1） 掌握求解实对称矩阵的全部特征值及其特征向量的Jacobi旋转法； （2） 了解Jacobi过关法。 三、实验要求 按照题目要求完成实验内容 写出相应的C语言程序 给出实验结果，结果分析 写出相应的实验报告 四、实验步骤 1、Jacobi旋转法的算法 1）在矩阵A的非对角元素中选取绝对值最大的元素，记为； 2）对确定的，用确定出； 3）产生平面旋转阵，得到 （p75：3.3.1~3.3.5）； 4）如果近似为对角阵，迭代终止；否则返回1）。 5）对角线上的元素近似为A的特征值；所有平面旋转矩阵的乘积矩阵的列向量为对应的特征向量。 2、实验题目： 用求以下矩阵的特征值及其特征向量迭代终止误差为1e-5： ， 1程序： #includestdio.h #includemath.h void main(){ int i, j, k, n = 3, ii, jj, kk; double aa, tan2x, x, pi = 3.1415926, err = 1e-3, sum_abs_aij = 1, t; double c[3][3], p[3][3], p1[3][3], p2[3][3], pp[3][3]; double a[3][3]={2,-1,0,-1,2,-1,0,-1,2}; for (i = 0; i n; i++) for (j = 0; j n; j++){ if (j == i) p[i][j] = 1; else p[i][j] = 0; } p1 = p; for (kk = 1; kk1e4sum_abs_aijerr; kk++){ c = a; aa = a[0][1]; ii = 0; jj = 1; sum_abs_aij = 0; for (i = 0; i n; i++){ for (j = i + 1; j n; j++){ if (fabs(aa) fabs(a[i][j])) { ii = i; jj = j; aa = a[i][j]; } sum_abs_aij += pow(a[i][j], 2); } } if ((a[ii][ii] - a[jj][jj]) == 0) x = pi / 2; else { tan2x = -2 * a[ii][jj] / (a[ii][ii] - a[jj][jj]); x = atan(tan2x); } x /= 2; p2 = p; p2[ii][ii] = p2[jj][jj] = cos(x); p2[ii][jj] = sin(x); p2[jj][ii] = -sin(x); for (i = 0; i n; i++) { for (j = 0; j n; j++) { t = 0; for (k = 0; k n; k++) t += p1[i][k] * p2[k][j]; pp[i][j] = t; } } p1 = pp; c[ii][ii] = a[ii][ii] * pow(cos(x), 2) + a[jj][jj] * pow(sin(x), 2) - a[ii][jj] * sin(2 * x); c[jj][jj] = a[ii][ii] * pow(sin(x), 2) + a[jj][jj] * pow(cos(x), 2) + a[ii][jj] * sin(2 * x); c[ii][jj] = c[jj][ii] = (a[ii][ii] - a[jj][jj]) / 2 * sin(2 * x) + a[ii][jj] * cos(2 * x); for (k = 0; k n; k++) if (k != jjk != ii) c[ii][k] = c[k][ii] = a[ii][k] * cos(x) - a[jj][k] * sin(x); for (k = 0; k n; k++) if (k != jjk != ii) c[jj][k] = c[k][jj] = a[jj][k] * cos(x) + a[ii][k] * sin(x); a = c; } kk--; printf(迭代次数kk=%d,\nA_%d=\n, kk, kk,n,n);//A_