1、数值分析第一次大作业数值分析计算实习报告第一题院 系: 机械工程及自动化学院_ 学 号: _姓 名: _ _2017年11月7日一、算法设计方案1、求,和的值1)利用幂法计算出矩阵按模最大的特征值,设其为。2)令矩阵(I为单位矩阵),同样利用幂法计算出矩阵按模最大的特征值。3)令。由计算过程可知和分别为矩阵所有特征值按大小排序后,序列两端的值。即,。4) 利用反幂法计算。其中,反幂法每迭代一次都要求解线性方程组,由于矩阵A为带状矩阵,故可用三角分解法解带状线性方程组的方法求解得到。2、求的与数最接近的特征值1) 令矩阵,利用反幂法计算出矩阵按模最小的特征值,则。3、求的(谱范数)条件数cond
2、()2和行列式det1) ,前文已算出和,直接带入即可。2) 反幂法计算时,已经对矩阵A进行过Doolittle分解,得到A=LU。而L为对角线上元素全为1的下三角矩阵,U为上三角矩阵,可知,即有。最后,为节省存储量,需对矩阵进行压缩,将中带内元素存储为数组。二、源程序代码#include#include#include#includeusing namespace std;#define N 501#define K 39#define r 2#define s 2#define EPSI 1.0e-12/求两个整数中的最大值int int_max2(int a, int b) return
3、(ab ? a : b);/求两个整数中的最小值int int_min2(int a, int b) return(ab) t = a; else t = b; if (tc) t = c; return(t);/定义向量内积double dianji(double x, double y) double sum = 0; for (int i = 0; iN; i+) sum = sum + xi * yi; return(sum);/计算两个数之间的相对误差double erro(double lamd0, double lamd1) double e, d, l; e = fabs(la
4、md1 - lamd0); d = fabs(lamd1); l = e / d; return(l);/矩阵A的压缩存储初始化成Cvoid init_c(double cN) int i, j; for (i = 0; ir + s + 1; i+) for (j = 0; jN; j+) if (i = 0 | i = 4) cij = -0.064; else if (i = 1 | i = 3) cij = 0.16; else cij = (1.64 - 0.024*(j + 1)*sin(0.2*(j + 1) - 0.64*exp(0.1 / (j + 1);/矩阵复制void
5、fuzhi_c(double c_constN, double cN) int i, j; for (i = 0; ir + s + 1; i+) for (j = 0; jN; j+) cij = c_constij;/LU三角分解void LUDet_c(double c_constN, double c_LUN) double sum; int k, i, j; fuzhi_c(c_const, c_LU); for (k = 1; k = N; k+) for (j = k; j = int_min2(k + s, N); j+) sum = 0; for (i = int_max3(
6、1, k - r, j - s); i = k - 1; i+) sum += c_LUk - i + si - 1 * c_LUi - j + sj - 1; c_LUk - j + sj - 1 -= sum; for (j = k + 1; j = int_min2(k + r, N); j+) sum = 0; for (i = int_max3(1, j - r, k - s); i = k - 1; i+) sum += c_LUj - i + si - 1 * c_LUi - k + sk - 1; c_LUj - k + sk - 1 = (c_LUj - k + sk - 1
7、 - sum) / c_LUsk - 1; /三角分解法解带状线性方程组void jiefc(double c_constN, double b_const, double x) int i, j; double bN, c_LUr + s + 1N, sum; for (i = 0; iN; i+) bi = b_consti; LUDet_c(c_const, c_LU); for (i = 2; i = N; i+) sum = 0; for (j = int_max2(i - 2, 1); j = 1; i-) sum = 0; for (j = i + 1; j = int_min2
8、(i + 2, N); j+) sum += c_LUi - j + 2j - 1 * xj - 1; xi - 1 = (bi - 1 - sum) / c_LU2i - 1; /幂法求按模最大特征值double mifa_c(double c_constN) double uN, yN; double sum, length_u, beta0, beta1; int i, j; for (i = 0; iN; i+) /迭代初始向量 ui = 0.5; length_u = sqrt(dianji(u, u); for (i = 0; iN; i+) yi = ui / length_u;
9、 for (i = 1; i = N; i+) sum = 0; for (j = int_max2(i - 2, 1); j = int_min2(i + 2, N); j+) sum = sum + c_consti - j + 2j - 1 * yj - 1; ui - 1 = sum; beta1 = dianji(u, y); do beta0 = beta1; length_u = sqrt(dianji(u, u); for (i = 0; iN; i+) yi = ui / length_u; for (i = 1; i = N; i+) sum = 0; for (j = i
10、nt_max2(i - 2, 1); j = EPSI); return(beta1);/反幂法求按模最小特征值double fmifa_c(double c_constN) double uN, yN; double length_u, beta0, beta1; int i; for (i = 0; iN; i+) /迭代初始向量 ui = 0.5; length_u = sqrt(dianji(u, u); for (i = 0; iN; i+) yi = ui / length_u; jiefc(c_const, y, u); beta1 = dianji(y, u); do beta
11、0 = beta1; length_u = sqrt(dianji(u, u); for (i = 0; i= EPSI); beta1 = 1 / beta1; return(beta1);/计算lamd_1、lamd_501、lamd_svoid calculate1(double c_constN, double &lamd_1, double &lamd_501, double &lamd_s) int i; double lamd_mifa0, lamd_mifa1, cr + s + 1N; lamd_mifa0 = mifa_c(c_const); fuzhi_c(c_const
12、, c); for (i = 0; iN; i+) c2i = c2i - lamd_mifa0; lamd_mifa1 = mifa_c(c) + lamd_mifa0; if (lamd_mifa0lamd_mifa1) lamd_1 = lamd_mifa0; lamd_501 = lamd_mifa1; else lamd_501 = lamd_mifa0; lamd_1 = lamd_mifa1; lamd_s = fmifa_c(c_const);/平移+反幂法求最接近u_k的特征值 void calculate2(double c_constN, double lamd_1, double lamd_501, double lamd_k) i
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1