数值逼近课程设计报告样本.doc

资料内容仅供您学习参考，如有不当之处，请联系改正或者删除。

课程设计报告

课程名称 数 值 逼 近
专 业 信息与计算科学
班 级 计算092
姓 名 杜青
学 号
指导教师 秦新强、 胡钢
日 期 -07-01
理学院应用数学系

数值积分及其应用
报告1

一、 目的意义
(1)进一步熟悉掌握复化梯形公式。
(2)进一步掌握熟悉复化抛物线公式。
(3) 学会比较复化梯形公式和复化抛物线公式如何达到所要求的精度。
二、 内容要求
积分计算问题: 分别用复化梯形和复化Simpsom求积公式计算积分, 并比较计算量( 精度为10-8) 。
三、 问题解决的方法与算法
方法: 利用复化梯形和复化抛物线积分公式。
算法:
输入: 端点a、 b以及要计算的积分公式f(x);
输出: 积分f(x)在指定区间上的近似值以及当其达到所要求的精度时要做的等分数n的值。
Step1: 编写复化梯形公式程序。
Step2: 经过所要达到的精度作为条件, 算出要做的等分数以及对应的近视值。
Setp3: 编写复化抛物线积分公式程序。
Setp4: 经过所要达到的精度作为条件, 算出要做的等分数以及对应的近视值。
Setp5: 然后比较复化梯形和复化抛物线的所需等分数, 比较谁的精度比较高。
四、 计算程序
1.复化梯形
#include <stdio.h>
#include <math.h>
double f(double x)
{
double s;
s=13*(x-x*x)*exp(-1.5*x);
return s;
}
void main()
{
int n,i;
double h,m,y,a,b,t[3000];
printf("请输入端点的值a,b\n");
scanf("%lf",&a);
scanf("%lf",&b);
for(n=1;;n++)
{
h=(b-a)/n;
m=(f(a)+f(b))/2;
for(i=1;i<n;i++)
{
m+=f(a+i*h);
}
t[n]=m*h;
h=(b-a)/(n+1);
m=(f(a)+f(b))/2;
for(i=1;i<n+1;i++)
{
m+=f(a+i*h);
}
t[n+1]=m*h;
if(fabs(t[n+1]-t[n])<0.00000001) break;
}

printf("求得结果为n=%d",n);
printf("求得结果为:t[n+1]=%10.8f\n",t[n+1]);
}
2.复化抛物线
#include <stdio.h>
#include <math.h>
double f(double x)
{
double s;
s=13*(x-x*x)*exp(-1.5*x);
return s;
}
void main()
{
int i,n;
double h,m,p,q,x,s,a,b,t[1000];
printf("请输入端点的值a,b\n");
scanf("%lf",&a);
scanf("%lf",&b);
for(n=1;;n++)
{h=(b-a)/(2*n);
m=f(a)+f(b);
p=0;
q=0;
for(i=1;i<2*n;i++)
{ x=a+i*h;
if(i%2==0)
q=q+f(x);
else
p=p+f(x);
}
t[n]=h*(m+2*q+4*p)/3;
h=(b-a)/(2*(n+1));
m=f(a)+f(b);
p=0;
q=0;
for(i=1;i<2*(n+1);i++)
{ x=a+i*h;
if(i%2==0)
q=q+f(x);
else
p=p+f(x);
}
t[n+1]=h*(m+2*q+4*p)/3;
if(fabs(t[n+1]-t[n])<0.00000001) break;
}
printf("求得结果为:n=%d\n",n);
printf("求得结果为:%10.8f\n",t[n+1]);
}
五、 计算结果与分析
1.复化梯形的运行结果:

2.复化抛物线的运行结果:

分析与评价: 经过对复化梯形的运行结果和复化抛物线的运行结果的分析得到, 当其所要求的精度相同时, 复化抛物线的的等分数明显比复化梯形的等分数少, 因此能够得到复化抛物线的精度比复化梯形的精度高。
六、 参考文献
[1] 谭浩强. C语言程序设计[M]. 北京: 清华大学出版社, .
[2] 秦新强. 数值逼近, 西安:西安理工大学印刷厂, .

机械设计问题
报告2

一、 目的意义
(1)掌握分段线性插值多项式的算法实现。
(2)学会分析误差和精度, 熟练运用各类编程的语句。
(3)进一步提高用数值算法解决实际问题的能力。
二、 内容要求
机械设计问题: 万能拉拨机中有一个圆柱形凸轮(见图1), 其底圆半径R=30cm,凸轮的上端面不在同一平面上, 要根据从动杆位移变化的需要进行设计制造。

将底圆周长36等分为xi (i=0,1, …, 36), 每一圆弧段长为h=52.36mm,对应于每 一分点的柱高为yi (i=0,1, …, 36)。为方便, 将圆柱展开为平面, 柱面的的顶端成为图2所示的平面曲线, 并已知该曲线上的37个点的坐标( 表1) 。

y C

A
B
yi yi

x0 O xi x17 x36 x

图1 凸轮模型 图2 展开曲线

表1: 测量数据表
xi
x0
x1
x2
x3
x4
x5
x6
x7
x8
yi
502.75
520.96
525
523.6
514.3
492
451
394.6
326.5
xi
X0
x10
x11
x12
x13
x14
x15
x16
x17- x36
yi
256.7
188.6
132.1
92.2
68.9
59.6
58.2
62.24
80.45-502.75
xi=jh, x0 =0, x36=1884.96mm, h=52.36mm。
是直线段, AB是曲线段, 为了数控加工, 需要计算出圆周上任一点处的柱高,试构造算法、 设计程序、 编程计算。
三．问题解决的方法和算法
方法: 利用分段线性插值和线性插值。
算法:
Step1: 输入题目所给的数据(x[i], y[i])。
Step2: 编写曲线段时要用的分段线性插值的程序。
Step3: 编写直线端时要用的线性插值的程序。
Step4: 任意输入一个插值点。
Step5: 经过判断条件计算插值点所对应的函数值。
四．程序设计
#include <stdio.h>
#include <math.h>
#define N 17
void main()
{ int i,n;
double X,Y;
double m1,m2,m3,m4;
double x[N+1];
doubley[N+1]={502.75,520.96,525,523.6,514.3,492,451,394.6,326.5,256.7,188.6,132.1,92.2,68.9,59.6,58.2,62.24,80.45};
double h=52.36;
for(i=0;i<N+1;i++)
{
x[i]=i*h;
printf("x[%d]=%lf\n",i,x[i]);
}
printf("请输入要插入的值X=");
scanf("%lf",&X);
if(X>=0&&X<=890.12)
{
for(i=0;i<N;i++)
{
if(X>=x[i]&&X<x[i+1])
{
m1=x[i];
m2=x[i+1];
m3=y[i];
m4=y[i+1];
Y=0.0;
Y=m3*(X-m2)/(m1-m2)+m4*(X-m1)/(m2-m1)+Y;
}
else if(X=x[i])
Y=y[i];
else
Y=y[i+1];
}
}
else if(X>=890.12&&X<=1884.96)
Y=(502.75-80.45)/(1884.96-890.12)*(X-890.12)+80.45;
else
printf("输入的插入点不在范围内,请重新输入\n");
printf("Y=%lf\n",Y);
}
五．对计算结果的分析和评价

此运行结果为第一种情况, 当插入点的值在曲线范围内所对应的x0~x17。

此运行结果为第二种情况, 当插入点的值在直线范围内所对应的x18~x36。

此运行结果为第三种情况, 当插入点的值不在题目所给的范围内时。
分析与评价:
运用数值算法解决实际问题时, 要注意各种误差。该题目运用的是分段线性插值和线性插值,
虽然存在误差, 可是对所要求的结果的影响不大, 因此本程序的结果是比较准确的。
六．参考文献
[1] 谭浩强. C语言程序设计[M]. 北京: 清华大学出版社, .
[2] 秦新强. 数值逼近, 西安:西安理工大学印刷厂, .

综合应用问题
报告3

一、 目的意义
(1)掌握最小二乘法的基本算法。
(2)经过最小二乘法能够拟合多次多项式的算法。
(3)熟练运用编程所要用的各种语句, 该算法主要学会矩阵的乘法以及如何解方程组的程序的实现。
二、 内容要求
6.多项式拟合问题:考虑将[0, 1] 30 等分节点, 用多项式y = 1+ x +L+ x 5生成数据, 再分别用多项式拟合求其3 次、 5 次、 10 次、 15 次插值多项式, 并分析误差产生原因。
三．算法和方法
方法: 最小二乘法。
算法:
Step1: 计算和输入题目所给的数据对(x[i],y[i])。
Step2: 输入系数矩阵和系数矩阵的转置, 两个矩阵相乘得到新矩阵。
Step3: 系数矩阵的转置与函数值矩阵相乘, 得到新矩阵。
Step4: 得到所要求的方程组解方程组。
Step5: 得到多项式的各项系数。即能够写出所要求的拟合多项式。
四．程序设计。
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <process.h>
#define N 7//N个点
#define T 15 //T次拟合
#define W 1//权函数
#define PRECISION 0.00001
float pow_n(float a,int n)
{
int i;
if(n==0)
return(1);
float res=a;
for(i=1;i<n;i++)
{
res*=a;
}
return(res);
}
void mutiple(float a[][N],float b[][T+1],float c[][T+1])
{
float res=0;
int i,j,k;
for(i=0;i<T+1;i++)
for(j=0;j<T+1;j++)
{
res=0;
for(k=0;k<N;k++)
{
res+=a[i][k]*b[k][j];
c[i][j]=res;
}
}
}
void matrix_trans(float a[][T+1],float b[][N])
{
int i,j;
for(i=0;i<N;i++)
{
for(j=0;j<T+1;j++)
{
b[j][i]=a[i][j];
}
}
}
void init(float x_y[][2],int n)
{
int i;
printf("请输入%d个已知点: \n",N);
for(i=0;i<n;i++)
{
printf("(x%d y%d):",i,i);
scanf("%f %f",&x_y[i][0],&x_y[i][1]);
}
}
void get_A(float matrix_A[][T+1],float x_y[][2],int n)
{
int i,j;
for(i=0;i<N;i++)
{
for(j=0;j<T+1;j++)
{
matrix_A[i][j]=W*pow_n(x_y[i][0],j);
}
}
}
void print_array(float array[][T+1],int n)
{
int i,j;
for(i=0;i<n;i++)
{
for(j=0;j<T+1;j++)
{
printf("%-g",array[i][j]);
}
printf("\n");
}
}
void convert(float argu[][T+2],int n)
{
int i,j,k,p,t;
float rate,temp;
for(i=1;i<n;i++)
{
for(j=i;j<n;j++)
{
if(argu[i-1][i-1]==0)
{
for(p=i;p<n;p++)
{
if(argu[p][i-1]!=0)
break;
}
if(p==n)
{
printf("方程组无解!\n");
exit(0);
}
for(t=0;t<n+1;t++)
{
temp=argu[i-1][t];
argu[i-1][t]=argu[p][t];
argu[p][t]=temp;
}
}
rate=argu[j][i-1]/argu[i-1][i-1];
for(k=i-1;k<n+1;k++)
{
argu[j][k]-=argu[i-1][k]*rate;
if(fabs(argu[j][k])<=PRECISION)
argu[j][k]=0;
}
}
}
}
void compute(float argu[][T+2],int n,float root[])
{
int i,j;
float temp;
for(i=n-1;i>=0;i--)
{
temp=argu[i][n];
for(j=n-1;j>i;j--)
{
temp-=argu[i][j]*root[j];
}
root[i]=temp/argu[i][i];
}
}
void get_y(float trans_A[][N],float x_y[][2],float y[],int n)
{
int i,j;
float temp;
for(i=0;i<n;i++)
{
temp=0;
for(j=0;j<N;j++)
{
temp+=trans_A[i][j]*x_y[j][1];
}
y[i]=temp;
}
}
void cons_formula(float coef_A[][T+1],float y[],float coef_form[][T+2])
{
int i,j;
for(i=0;i<T+1;i++)
{
for(j=0;j<T+2;j++)
{
if(j==T+1)
coef_form[i][j]=y[i];
else
coef_form[i][j]=coef_A[i][j];
}
}
}
void print_root(float a[],int n)
{
int i,j;
printf("%d个点的%d次拟合的多项式系数为:\n",N,T);
for(i=0;i<n;i++)
{
printf("a[%d]=%g,",i+1,a[i]);
}
printf("\n");
printf("拟合曲线方程为:\ny(x)=%g",a[0]);
for(i=1;i<n;i++)
{
printf(" + %g",a[i]);
for(j=0;j<i;j++)
{
printf("*X");
}
}
printf("\n");
}
void process()
{
float x_y[N][2],matrix_A[N][T+1],trans_A[T+1][N],coef_A[T+1][T+1],coef_formu[T+1][T+2],y[T+1],a[T+1];
init(x_y,N);
get_A(matrix_A,x_y,N);
printf("矩阵A为: \n");
print_array(matrix_A,N);
matrix_trans(matrix_A,trans_A);
mutiple(trans_A,matrix_A,coef_A);
printf("法矩阵为: \n");
print_array(coef_A,T+1);
get_y(trans_A,x_y,y,T+1);
cons_formula(coef_A,y,coef_formu);
convert(coef_formu,T+1);
compute(coef_formu,T+1,a);
print_root(a,T+1);
}
void main()
{
process();
}
五．计算结果与分析
三次多项式时:

五次多项式时:

十次多项式和十五次多项式时:

分析与评价:
根据运行结果可得, 已知的数据对数只能够拟合成三次多项式和五次多项式, 不能拟合成十次多项式和十五次多项式, 因为其多项式系数作为未知量的方程组无解。
其误差产生原因是在计算过程中数据的四舍五入以及再输入数据时本身就输入的是近似值。
六．参考文献
[1] 谭浩强. C语言程序设计[M]. 北京: 清华大学出版社, .
[2] 秦新强. 数值逼近, 西安:西安理工大学印刷厂, .