雅克比迭代法(2022年)

发布时间:2022-07-01 12:35:04

下面是小编为大家整理的雅克比迭代法(2022年),供大家参考。

雅克比迭代法(2022年)

 

 实验项目名称

 运用插值法 实验成绩

  实验者

 江骏 专业班级

 软件 0803 组别

 同组者

 第一部分:实验分析与设计

 一、实验内容描述 实验日期

 年 月 日

 ⑴ 研究用

 Jacobi

 迭代法与

 Gauss- - Seidel

 迭代法解下列方程组

 Ax=b

 的收敛性,通过上机计算,验证分析是否正确,并观察右端项对迭代收敛是否有影响,比较两法的收敛速度

 ⑵ 松弛因子对超松弛因子迭代法收敛速度的影响,要求对不同的阶数进行迭代

 ⑶ 观察欧拉显式方法的收敛性

 ⑷ 观察欧拉隐式方法的收敛性

 ⑸ 写出实验报告

 二、实验基本原理与设计

 ①Jacobi

 迭代法:

 ∑ n n

 a a

 x x

 =b

 (i=1,2,…,n)

 j=1

 ij

 j j

 i i

  ②Gauss- - Seidel

 迭代法 :x x

 (K+1)

 =[b

 - - ∑ i i- -1 1

 a x

 (k+1)

 - - ∑ n n

 a x (k)

 ]/a

 i i

 i i

 j=1

 ij

  j j

 j=i+1

 ij

  j j

 ij

  ③ 超松弛因子迭代法:

 Dx (k+1)

 =Dx (k)

 +Ex

 (k+1)

 +Fx (k)

 +b- - Dx

 (k)

  X X (k+1) =x (k) +D - -1 1

 R R

  ④ 三角分解法:利用三角矩阵将方程组化解为两个方程组,从而简化运算。

 ⑤ 欧拉显式方法:在任意节 点 t t

 =t

 +(n+1)h

 处, u(t

 ) ) 的近似值由

 Euler

 公式给出:

 n+1

 0 0

 n+1

  u u

 =u

 +hf(t

 ,u

 ) )

 n+1

 n n

 n n

 n n

  ⑥ 欧拉隐式方法:改进 的 Euler

 公式

 y =y

 +hf(x

 ,y

 ) )

 n+1

 n n

 n n

 n n

  h y y

 =y

 + +

 [f(x

 ,y

 )+f(x

 , ,

 y )]

 n+1

 n n

 2 n n

 n n

 n+1

  n+1

  三、主要仪器设备及耗材

 WindowsXp VC++

 第二部分:实验调试与结果分析

 一、调试过程

  ①Jacobi 迭代法#include <> #define N 10

 float ABS(float,float); int

 main (void)

 { {

 int i, j, n;

 float a[N][N], b[N];

 float x[N], y[N]; float e,total;

 printf ("Please input the dimension:"); scanf ("%d", &n);

 printf ("Please input the Coefficient Matrix: "); for (i = 0; i < n; i++)

 for (j = 0; j < n; j++) scanf ("%f", &a[i][j]);

 printf ("Please input the Vector: "); for (i = 0; i < n; i++)

 scanf ("%f", &b[i]);

 printf ("Please input the initial vector: "); for (i = 0; i < n; i++)

 scanf ("%f", &x[i]);

 do

 { {

 for (i = 0; i < n; i++)

 { {

 total = ;

 for (j = 0; j < n; j++)

 { {

 if (i != j)

 total += a[i][j] * x[j];

 } }

 y[i] = (b[i] - -

 total) / a[i][i];

 } }

 e=;

 for (j = 0; j < n; j++)

 e = e + ABS(x[j], y[j]);

 printf("%f\ \ n",e);

 for (i = 0; i < n; i++)

 x[i] = y[i];

 } }

 while (e > ;

 for (i = 0; i < n; i++)

 printf ("x%d=%f\ \ t", i, x[i]); printf ("\ \ n");

 return 0;

 } }

 float ABS(float x,float y)

 { {

 int total;

 if(x<y) total=y- - x; else total=x- - y; return total;

 } }

 结 果 :

 流程图:

  开始

 阶数 n,系数矩阵,初始向量

 for(i=0;i<n;i+

 Y

 i=j

  N

 total

 +=

 e e

 = =

 e e

 + +

 结束

  ②Gauss- - Seidel

 迭代法#include<>

 int

 main()

 { {

 double x[3] = {0, 0, 0};

 double a[3][3] = {6,2,- - 1,1,4,- - 2,- - 3,1,4};

 double y[3] = {- - 3,2,4};

 double d[3][3],g[3]; int round = 5, i,j; for (i=0; i<3; ++i)

 { {

 g[i] = y[i] / a[i][i]; for (j=0; j<3; ++j)

 { {

 d[i][j] = i==j ? 0 : - - a[i][j]/a[i][i];

 } }

 } }

 while (round -- ) { for (i=0; i<3;

 ++i)

 { {

 x[i] = g[i];

 for (j=0; j<3; ++j)

 { {

 x[i] += d[i][j] * x[j];

 } }

 printf("%lf", x[i]);

 } }

 printf("\ \ n");

 } }

 } }

 结 果 :

 流程图:

 开始

 系数矩阵

 for

 (i=0;

 i<3;

 ++i) d[i][j]

 = =

 i==j

 ? ?

 0 0

 : :

 for

 (j=0;

 j<3;

 ++j)

 结束

 ③ 超松弛因子迭代法#include <iostream> #include <cmath> using namespace std;

 float *one_array_malloc(int n);

 float **two_array_malloc(int m,int n); float matrix_category(float* x,int n); int main(){ const int MAX=100;

 int n,i,j,k; float** a; float* x_0; float* x_k;

 float

 precision;

 float w; cout<<" 输

 入

 精

 度

 e:"; cin>>precision;

 cout<<endl<<" 输入系数矩阵的 阶数,N N :

 "; cin>>n; a=two_array_malloc(n,n+1); cout<<endl<<" 输入增广矩阵的各值:\ \ n"; for(i=0;i<n;i++)

 { {

 for(j=0;j<n+1;j++)

 { {

 cin>>a[i][j];

 } }

 } }

 x_0=one_array_malloc(n); cout<<endl<<" 输入初始向量: :\ \ n"; for(i=0;i<n;i++)

 { {

 cin>>x_0[i]; } x_k=one_array_malloc(n);

 cout<<" 输入松弛因子w w

 (1<w<2) :\ \ n"; cin>>w;

 float temp; for(k=0;k<MAX;k++)

 { {

 for(i=0;i<n;i++)

 { {

 temp=0; for(j=0;j<i;j++)

 { {

 temp=temp+a[i][j]*x_k[j];

 } }

 x_k[i]=a[i][n]- - temp; temp=0; for(j=i+1;j<n;j++)

 { {

 temp=temp+a[i][j]*x_0[j];

 } }

 x_k[i]=(x_k[i]- - temp)/a[i][i];

 x_k[i]=(1- - w)*x_0[i]+w*x_k[i];

 } }

 for(i=0;i<n;i++)

 { {

 x_0[i]=x_k[i]- - x_0[i];

 } if(matrix_category(x_0,n)<precision)

 { {

 break;

 } }

 else

 { {

 for(i=0;i<n;i++)

 { {

 x_0[i]=x_k[i];

 } }

 } }

 } }

 if(MAX==k)

 { {

 cout<<" 迭代不收敛\ \ n";

 } cout<<" 迭代次数为:

 "<<k<<endl; cout<<" 解

 向

 量

 为

 :\ \ n"; for(i=0;i<n;i++)

 { {

 cout<<"x"<<i<<":

 "<<x_k[i]<<endl;

 } return

 0;

 } }

 float *one_array_malloc(int n)

 { {

 float

 *a;

 a=(float *)malloc(sizeof(float)*n); return a;

 } }

 float **two_array_malloc(int m,int n)

 { {

 float **a; int

 i;

 a=(float **)malloc(m*sizeof(float *)); for (i=0;i<m;i++)

 { {

 a[i]=(float

 *)malloc(n*sizeof(float));

 } }

 return a;

 } }

 float matrix_category(float* x,int

 n)

 { int i; float temp=0;

 for(i=0;i<n;i++)

 { {

 temp=temp+fabs(x[i]);

 } }

 return temp;

 } }

 结果:

 流程图:

 开始

 系数矩阵阶, 精度,增广阵,初始

  for(k=0;k<MAX;

 Y

 matrix_category(x

 N

 for(i=0;i<n;i++)

 a=(float

 *)malloc(sizeof

 结束

 ④ 三角分解法#include<iostream> using namespace std; int main()

 { {

 const int N=100;

 static

 double

 a[N][N],b[N];

 int i,j,k,num,p; double m,t,q;

 cout<<" 请输入矩阵阶数:

 "; cin>>num;

 ."<<num<<" :"<<endl; for(i=1;i<=num;i++)

 { {

 for(j=1;j<=num;j++)

 { {

 cout<<"a["<<i<<"]["<<j<<"]=";

 cin>>a[i][j];

 } }

 } }

 cout<<"now input the matrix b[i],i=1..."<<num<<" :"<<endl; for(i=1;i<=num;i++)

 { {

 cout<<"b["<<i<<"]="; cin>>b[i];

 } }

 t=0;

 for(i=1;i<=num;i++)

 { {

 m=0;

 for(j=i;j<=num+1;j++)

 { {

 for(k=1;k<=i- - 1;k++)

 m=m+a[i][k]*a[k][j];

 a[i][j]=

 a[i][j]- - m;

 } }

 for(j=i+1;j<=num;j++)

 { {

 for(k=1;k<=i- - 1;k++) t=t+a[j][k]*a[k][i];

 a[j][i]=(a[j][i]- - t)/a[i][i];

 } }

 } }

 a[num][num+1]=a[num][num+1]/a[num][num]; q=0;

 for(k=num- - 1;k>=1;k --) )

 { {

 for(j=k+1;j<=num;j++)

 q=q+a[k][j]*a[j][num+1];

 a[k][num+1]=(a[k][num+1]- - q)/a[k][k];

 } }

 cout<<" 方

 程

 组

 的

 解

 为

 :

 "; for(i=1;i<=num;i++) cout<<"x["<<i<<"]="<<a[i][num+1]<<endl; return 0;

 } }

 结 果 :

 流程图:

 开始

  矩阵阶数

 for(k=1;k<=i- - 1;k++) m=m+a[i][k]*a[k][j];

 a[k][num+1]=(a[k][

 结束

  ⑤ 欧拉显式方法#include <> #include <> #include <>

 double f(double x,double y)

 { {

 return x*pow(y,3);

 } }

 int main()

 { {

 int i;

 double x,y,y0=1,dx=; double xx[11];

 double euler[11],euler_2[11];

 double temp;

 double f(double x,double y); for (i=0;i<11;i++) xx[i]=i*dx;

 euler[0]=y0;

 for (i=1,x=0;i<11;i++,x+=dx) euler[i]=euler[i- - 1]+dx*f(x,euler[i- - 1]); euler_2[0]=y0;

 for (i=1,x=0;i<11;i++,x+=dx)

 { {

 temp=euler_2[i- - 1]+dx*f(x,euler_2[i- - 1]);

 euler_2[i]=euler_2[i- - 1]+dx*(f(x,euler_2[i- - 1])+f(x+dx,temp))/2;

 } }

 for (i=0,x=0;i<11;i++,x+=dx) printf("x=%lf\ \ teluer=%lf\ \ teuler_2=%lf\ \ taccu=%lf\ \ n",x,euler[i],euler_2[i],pow(1+x*x,3)); getch();

 } }

 结 果 :

 流程图:

 ⑥ 欧拉隐式方法

 #include<> #include<>

 #define f(x,y)(x+y ) int

 main()

 {int m; int i;

 double a,b,y0;

 double xn,yn,xn1,yn1,yn1b; double h;

 printf("\ \ nInput the begin and end of x:"); scanf("%lf%lf",&a,&b);

 printf("Input the y value at %f:",a); scanf("%lf",&y0);

 printf("Input m value[divide(%f,%f)]:",a,b); scanf("%d",&m);

 if(m<=0)

 {printf("Please input a number larger than 1.\ \ n"); return(1);

 } }

 h=(b- - a)/m; xn=a;yn=y0; for(i=1;i<=m;i++)

 {xn1=xn+h; yn1b=yn+h*f(xn,yn);

 yn1=yn+h/2*(f(xn,yn)+f(xn1,yn1b)); printf("x%d=% f,y%d=%f\ \ n",i,xn1,i,yn1); xn=xn1;

 yn=yn1;

 }return(0);

 } }

 结 果 :

 流程图:

 二、实验结果及分析 雅克比迭代法的全部分量都是间接利用,由于新值比旧值好,Gauss-Seidel 迭代法直接对算出

 的分量加以使用,这样使得迭代的收敛情况有所改善,而超松弛法是在 Gauss-Seidel 迭代法的 基础上,应用加速收敛思想而构造的一种方法;欧拉公式和梯形公式在计算上有明显的区别, 欧拉公式的特点是可以有 U 直接计算 U ,也就是显式的,而梯形公式右端也有 U ,必须通过解 N n+1 n+1

 方程才能求得,这类是隐式的。从算法的角度看显式远比隐式方便,隐式通常用迭代法求解。

 三、实验小节、建议及体会 任何一种方法都有它的优点和局限性,只有找到一种方法能够解决我们的问题,又尽量减少出错,减少 计算量就是符合这个问题的好方法。

 良好的编程习惯有助于检查错误,也可以增强程序的可阅读性。

推荐访问:雅克比迭代法 迭代法 雅克

版权所有:众一秘书网 2005-2024 未经授权禁止复制或建立镜像[众一秘书网]所有资源完全免费共享

Powered by 众一秘书网 © All Rights Reserved.。备案号: 辽ICP备05005627号-1