C语言课程设计 (3)

数值分析上机

实验报告

东北大学软件

2011年

姓名闫三木学号20093514 班级软件工程0908 指导教师邵新慧实验名称曲线拟合最小二乘法线性方程组迭代法

开设学期2010 – 2011第二学期

评定成绩

评定人签字

评定日期

课题三 解线性方程组的迭代法

一． 实验内容

1、设线性方程组

?????????????????????

???

???

?????--------------------------1368

2438100

412029137264221234179111016103524312053621775868323376162449113151

20130123122400

1

0563568000

121324??????

???

???

????????????????????10987654321x x x x x x x x x x =??

??

?

??

???

?

??

???

????????????????-2119381346323125

x *= ( 1, -1, 0, 1, 2, 0, 3, 1, -1, 2 )T 2、设对称正定阵系数阵线方程组

??????????????????????????----------------------19243360021411035204111443343104221812334161206538114140231212200420424??????????????????????????87654321x x x x x x x x = ??

??

?

???

????

??????????????---4515229232060

x * = ( 1, -1, 0, 2, 1, -1, 0, 2 )T 3、三对角形线性方程组

????????????????????????????????------------------4100000000141000000001410000000014100000000141000000001410000000014100000000141000000001410000000014???????????????

?????????????????10987654321x x x x x x x x x x = ?????

??

???

?

?????

????????????????----5541412621357

x *= ( 2, 1, -3, 0, 1, -2, 3, 0, 1, -1 )T

试分别选用Jacobi 迭代法，Gauss-Seidol 迭代法和SOR 方法计算其解。

二． 实验要求

1、体会迭代法求解线性方程组，并能与消去法做以比较；

2、分别对不同精度要求，由迭代次数体会该迭代法的收敛快慢；

3、对方程组2，3使用SOR 方法时，选取松弛因子ω=0.8，0.9，1，1.1，1.2等，试看对算法收敛性的影响，并能找出你所选用的松弛因子的最佳者；

4、给出各种算法的设计程序和计算结果。

三． 目的和意义

1、通过上机计算体会迭代法求解线性方程组的特点，并能和消去法比较；

2、运用所学的迭代法算法，解决各类线性方程组，编出算法程序；

3、体会上机计算时，终止步骤∞

+-)()1(k k x x < ε 或k >（予给的迭代次

数），对迭代法敛散性的意义；

4、体会初始解 x )0(，松弛因子的选取，对计算结果的影响。

四． 流程图设计 1.主要部分流程

2.Jacobi 算法流程

3.GS 算法流程图

4.SOR算法流程

五．结构程序设计

#include

using namespace std;

#define N 40

const int n=10;

int jacobi(fl oat *p,float b[],fl oat X[],fl oat x[],int n); int GS(float *p,fl oat b[],float X[],float x[],int n); int SOR(float *p,fl oat b[],float X[],float x[],int n); void print(fl oat *a,int r);

void main()

{

fl oat A[10][10]={4,2,-3,-1,2,1,0,0,0,0,

8,6,-5,-3,6,5,0,1,0,0,

4,2,-2,-1,3,2,-1,1,9,4,

0,-2,1,5,-1,3,-1,1,9,4,

-4,2,6,-1,6,7,-3,3,2,3,

0,2,-1,3,-4,2,5,3,0,1,

16,10,-11,-9,17,34,2,-1,2,2,

4,6,2,-7,13,9,2,0,12,4,

0,0,-1,8,-3,-24,-8,6,3,-1};

fl oat a[10]={7,5,-13,2,6,-12,14,-4,5,-5};

fl oat X1[10]={0,0,0,0,0,0,0,0,0,0};

fl oat x1[10];

fl oat B[8][8]={4,2,-4,0,2,4,0,0,

2,2,-1,-2,1,3,2,0,

-4,-1,14,1,-8,-3,5,6,

0,-2,1,6,-1,-4,-3,3,

2,1,-8,-1,22,4,-10,-3,

4,3,-3,-4,4,11,1,-4,

0,2,5,-3,-10,1,14,2,

0,0,6,3,-3,-4,2,19};

fl oat b[8]={0,-6,6,23,11,-22,-15,45};

fl oat X2[8]={0,0,0,0,0,0,0,0};

fl oat x2[8];

fl oat C[10][10]={4,-1,0,0,0,0,0,0,0,0,

-1,4,-1,0,0,0,0,0,0,0,

0,-1,4,-1,0,0,0,0,0,0,

0,0,-1,4,-1,0,0,0,0,0,

0,0,0,-1,4,-1,0,0,0,0,

0,0,0,0,-1,4,-1,0,0,0,

0,0,0,0,0,-1,4,-1,0,0,

0,0,0,0,0,0,-1,4,-1,0,

0,0,0,0,0,0,0,-1,4,-1,

0,0,0,0,0,0,0,0,-1,4};

fl oat c[10]={7,5,-13,2,6,-12,14,-4,5,-5};

fl oat x3[10];

fl oat X3[10]={0,0,0,0,0,0,0,0,0,0};

fl oat *p[3];

p[0]=&A[0][0];

p[1]=&B[0][0];

p[2]=&C[0][0];

cout<<"Jacobi迭代法解第1个方程:"<

cout<<"Jacobi迭代法解第2个方程:"<

cout<<"Jacobi迭代法解第3个方程:"<

cout<<"Gauss-Seidel迭代法解第1个方程:"<

cout<<"Gauss-Seidel迭代法解第2个方程:"<

cout<<"Gauss-Seidel迭代法解第3个方程:"<

cout<<"SOR迭代法解第1个方程:"<

SOR(p[0],a,X1,x1,10);

cout<<"SOR迭代法解第2个方程:"<

SOR(p[1],b,X2,x2,8);

cout<<"SOR迭代法解第3个方程:"<

SOR(p[2],c,X3,x3,10);

}

int jacobi(fl oat *p,float b[],fl oat X[],fl oat x[],int n) {

int k,i,j;

fl oat m,R,r,e;

cout<<"请输入精度e:";

cin>>e;

for(k=0;k

{

R=0;

for(i=0;i

{

m=0;

for(j=0;j

{

m=m+(*(p+i*n+j))*X[j];

}

x[i]=X[i]+(b[i]-m)/(*(p+i*n+i));

r=x[i]-X[i];

if(r<0)

r=X[i]-x[i];

if(r>R)

R=r;

}

print(x,n);

cout<<"迭代次数为:"<

return k;

}

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

X[j]=x[j];

}

print(x,n);

cout<<"迭代次数为:"<

cout<<"方程解发散,无法用Jacobi方法解此方程!"<

}

int GS(float *p,fl oat b[],float X[],float x[],int n)

{

int i,j,k;

fl oat t,R,r,e;

cout<<"请输入精度e:";

cin>>e;

for(k=0;k

{

for(i=0;i

{

t=0;

for(j=0;j

{

if(j

t+=(*(p+i*n+j))*x[j];

if(j>i)

t+=(*(p+i*n+j))*X[j];

}

x[i]=(b[i]-t)/(*(p+i*n+i));

}

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

{

r=x[i]-X[i];

if(r<0)

r=X[i]-x[i];

}

if(R<=e&&R>0)

{

print(x,n);

cout<<"迭代次数为:"<

return k;

}

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

X[j]=x[j];

}

print(x,n);

cout<<"迭代次数为:"<

cout<<"方程解发散,无法用Gauss-Seidel方法解此方程!"<

}

int SOR(float *p,fl oat b[],float X[],float x[],int n)

{

int i,j,k;

fl oat t,R,r,e,w;

cout<<"请输入松弛因子w(0

cin>>w;

cout<<"请输入精度e:";

cin>>e;

for(i=0;i

x[i]=X[i];

for(k=0;k

{

R=0;

for(i=0;i

{

t=0;

for(j=0;j

t+=(*(p+i*n+j))*x[j];

r=w*(b[i]-t)/(*(p+i*n+i));

x[i]+=r;

if(r<0)

r=-r;

}

if(R<=e&&R>0)

{

print(x,n);

cout<<"迭代次数:"<

return k;

}

}

print(x,n);

cout<<"迭代次数:"<

cout<<"方程解发散,无法使用SOR方法解此方程!!!!"<

}

void print(fl oat *a,int n)

{

int j;

fl oat *t=a;

cout<<"x=( ";

for(j=0;j

cout<<*(t+j)<<",";

cout<<*(t+j)<<")"<

}

六．结果讨论和分析

程序运行截图：

七．小结和体会

通过这次实验，我对线性方程组的迭代算法有了更深的理解。迭代算法给计算机解决问题提供了有效的途径。使得一些人为难以计算出来的复杂线性方程组通过迭代算法交由计算机来进行计算。

同时我也发现了，J，GS和SOR迭代方法各自的特点，对课堂学习的知识有了更好的理解和认识。

课题八曲线拟合的最小二乘法

一．实验内容

利用数据拟合的最小二乘法从一组数据中找出其规律性，并给出其数学模型的近似表达式。

问题提出：

在某冶炼过程中，通过实验检测得到含碳量与时间关系的数据

如下，求含碳量y与时间t内关系的拟合曲线。

用最小二乘法进行三次多项式的曲线拟合，并计算均方误差二．实验要求

1.用最小二乘法进行三次多项式的曲线拟合；

2.计算yj与y(tj)的误差；j=1,2,3….

3.另外选取一个拟合函数，进行拟合效果的比较；

4.绘制出曲线拟合的图形。

三．目的和意义

1.掌握最小二乘法应用方法

2.理解最小二乘法进行曲线拟合的意义

四．结构程序设计

#include

#include

#include

using namespace std;

#define deviation 0.0000000001

#define N 20

class fit

{

protected:

int m, n;

double *a, *x, *y;

double dt;

public:

//构造函数

fit(int m, int n)

{

this->m = m;

this->n = n;

dt = 0;

x = new double[m];//一共有m个数据点，数据点的x值

y = new double[m];

a = new double[n+1];

void input();

//进行最小二乘曲线拟合

void fitting();

double value(double, int);

double* equation(double **, double *);

void deviat();

void output();

};

//从文件读入数据

void fit::input()

{

int i;

char fileName[20];//保存文件名

//提示用户输入用户名n

cout<<"*********请输入数据文件名称:*********\n"; cin>>fileName;

//从文件读取数据点

ifstream fin(fileName);

if(!fin)

{

cout<<"不能正常打开文件"<

exit(1);

}

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

{

fin>>x[i];

fin>>y[i];

}

fin.close();//关闭文件流

}

//进行最小二乘曲线拟合

void fit::fitting()

{

double *vy = new double[m+1];

double *f = new double[n+1];

int i, j, k;

double **vx = new double*[n+1];

double **coe = new double*[n+1];

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

{

vx[i] = new double[m+1];

coe[i] = new double[n+1];

}

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

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

{

vx[i][j] = value(x[j], i);

}

}

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

{

vy[j] = y[j];

}

//计算方程组的系数矩阵

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

{

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

{

coe[i][j] = 0;

for(k = 0; k < m; k++)

{

coe[i][j] += vx[i][k] * vx[j][k];

}

}

}

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

{

f[i] = 0;

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

{

f[i] += vx[i][j] * vy[j];

}

}

//解正则方程组，并将结果存储到数组a中

a = equation(coe, f);

//计算均方误差，并将结果存储到统计变量dt中deviat();

}

//计算基函数的函数值

double fit::value(double x, int m)

{

double result = 1;

for(int i = 0; i < m; i++)

{

result *= x;

}

return result;

}

{

double *tx = new double[n+1];

long double r, r1;

int w = 1.5;

int i, j, k;

//初始化初始解

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

{

tx[i] = 0;

}

for(k = 0; k < N; k++)

{

r = 0;

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

{

//计算误差

r1 = b[i];

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

{

r1 -= coe[i][j] * tx[j];

}

r1 /= coe[i][i];

r1 *= w;

//取最大误差

if(fabs(r1) > fabs(r))

{

r = r1;

}

tx[i] += r1;

}

//如果达到精度要求，则结束迭代

if(fabs(r) <= deviation)

{

break;

}

}

//返回计算结果

return tx;

}

//计算误差

void fit::deviat()

{

double *funcValue = new double[m]; double sum = 0;

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

{

funcValue[j] = 0;

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

{

funcValue[j] += a[i] * value(x[j], i);

}

sum += (y[j] - funcValue[j]) * (y[j] - funcValue[j]);

}

this->dt = sqrt(sum);

}

//输出最终拟合结果，并保存到指定文件

void fit::output()

{

int i;

//输出各项系数

cout<<"*********拟合函数的各项系数为(按各项次数从小到大排列):*********\n"<

{

cout<

}

cout<<"*********均方误差为：*********\n"<

}

//main函数

int main()

{

int m, n;

cout<<"*********请输入数据点的个数：*********\n";

cin>>m;

cout<<"*********请输入需要进行拟合的次数：*********\n";

cin>>n;

fit express(m, n);

express.input();

express.fitting();

express.output();

system("Pause");

return 1;

}

五．结果讨论和分析

程序运行结果：

运用MATLAB绘制的拟合图像

六．小结和体会

通过这次实验，我对最小二乘法进行的曲线拟合有了更多的认识。通过几种不同形式的计算，由运行结果可知，三次多项式的均方误差小于二次多项式的均方误差，二次多项式的均方误差小于四次多项式的均方误差，因此，三次多项式的拟合曲线相对较好。

通过这个实验我也理解到了，数值分析是一个工具学科，它

教给了我们分析和解决数值计算问题的方法。使我从中得到很多关于算法的思想，从中受益匪浅。