计算方法编程题实验报告Word格式文档下载.docx
- 文档编号:21307312
- 上传时间:2023-01-29
- 格式:DOCX
- 页数:18
- 大小:69.79KB
计算方法编程题实验报告Word格式文档下载.docx
《计算方法编程题实验报告Word格式文档下载.docx》由会员分享,可在线阅读,更多相关《计算方法编程题实验报告Word格式文档下载.docx(18页珍藏版)》请在冰豆网上搜索。
);
exit
(1);
}
a[0][0]=12;
a[0][1]=-3;
a[0][2]=3;
a[1][0]=-18;
a[1][1]=3;
a[1][2]=-1;
a[2][0]=1;
a[2][1]=1;
a[2][2]=1;
b[0]=15;
b[1]=-15;
b[2]=6;
if(!
GS(n,a,b,ep))
不可以用高斯消去法求解\n"
exit(0);
printf("
该方程组的解为:
\n"
for(i=0;
i<
3;
i++)
x%d=%.2f\n"
i,b[i]);
TwoArrayFree(a);
free(b);
}
intGS(n,a,b,ep)
intn;
double**a;
double*b;
doubleep;
inti,j,k,l;
doublet;
for(k=1;
k<
=n;
k++)
for(l=k;
l<
l++)
if(fabs(a[l-1][k-1])>
ep)
break;
elseif(l==n)
return(0);
if(l!
=k)
{
for(j=k;
j<
j++)
{
t=a[k-1][j-1];
a[k-1][j-1]=a[l-1][j-1];
a[l-1][j-1]=t;
}
t=b[k-1];
b[k-1]=b[l-1];
b[l-1]=t;
}
t=1/a[k-1][k-1];
for(j=k+1;
a[k-1][j-1]=t*a[k-1][j-1];
b[k-1]*=t;
for(i=k+1;
for(j=k+1;
a[i-1][j-1]-=a[i-1][k-1]*a[k-1][j-1];
b[i-1]-=a[i-1][k-1]*b[k-1];
for(i=n-1;
i>
=1;
i--)
for(j=i+1;
b[i-1]-=a[i-1][j-1]*b[j-1];
return
(1);
double**TwoArrayAlloc(intr,intc)
double*x,**y;
intn;
x=(double*)calloc(r*c,sizeof(double));
y=(double**)calloc(r,sizeof(double*));
for(n=0;
n<
=r-1;
++n)
y[n]=&
x[c*n];
return(y);
voidTwoArrayFree(double**x)
free(x[0]);
free(x);
【实验结果】
P40/3.2.3
线性方程组为:
12x1-3x2+3x3=15
-18x1+3x2-x3=-15
x1+x2+x3=6
正确结果为:
x1=1.000
x2=2.000
x3=3.000
实验二最小二乘法
1.曲线拟合的最小二乘法的基本思路和拟合步骤
2.能根据给定的函数值表构造出次数不相同的拟合多项式
conio.h>
process.h>
#defineN5//N个点
#defineT3//T次拟合
#defineW1//权函数
#definePRECISION0.00001
floatpow_n(floata,intn)
inti;
floatres;
if(n==0)
res=a;
for(i=1;
n;
res*=a;
return(res);
voidmutiple(floata[][N],floatb[][T+1],floatc[][T+1])
floatres=0;
inti,j,k;
for(i=0;
T+1;
for(j=0;
res=0;
for(k=0;
N;
res+=a[i][k]*b[k][j];
c[i][j]=res;
voidmatrix_trans(floata[][T+1],floatb[][N])
inti,j;
b[j][i]=a[i][j];
voidinit(floatx_y[][2],intn)
printf("
请输入%d个已知点:
N);
(x%dy%d):
"
i,i);
scanf("
%f%f"
&
x_y[i][0],&
x_y[i][1]);
voidget_A(floatmatrix_A[][T+1],floatx_y[][2],intn)
matrix_A[i][j]=W*pow_n(x_y[i][0],j);
voidprint_array(floatarray[][T+1],intn)
%-g"
array[i][j]);
voidconvert(floatargu[][T+2],intn)
inti,j,k,p,t;
floatrate,temp;
for(j=i;
if(argu[i-1][i-1]==0)
for(p=i;
p<
p++)
if(argu[p][i-1]!
=0)
break;
if(p==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;
argu[j][k]-=argu[i-1][k]*rate;
if(fabs(argu[j][k])<
=PRECISION)
argu[j][k]=0;
voidcompute(floatargu[][T+2],intn,floatroot[])
floattemp;
for(i=n-1;
=0;
temp=argu[i][n];
for(j=n-1;
j>
i;
j--)
temp-=argu[i][j]*root[j];
root[i]=temp/argu[i][i];
voidget_y(floattrans_A[][N],floatx_y[][2],floaty[],intn)
temp=0;
temp+=trans_A[i][j]*x_y[j][1];
y[i]=temp;
voidcons_formula(floatcoef_A[][T+1],floaty[],floatcoef_form[][T+2])
T+2;
if(j==T+1)
coef_form[i][j]=y[i];
else
coef_form[i][j]=coef_A[i][j];
voidprint_root(floata[],intn)
%d个点的%d次拟合的多项式系数为:
N,T);
a[%d]=%g,"
i+1,a[i]);
拟合曲线方程为:
\ny(x)=%g"
a[0]);
+%g"
a[i]);
*X"
voidprocess()
floatx_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);
矩阵A为:
print_array(matrix_A,N);
matrix_trans(matrix_A,trans_A);
mutiple(trans_A,matrix_A,coef_A);
法矩阵为:
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);
process();
实验三龙贝格方法
1.理解龙贝格方法的基本思路
2.用龙贝格方法设计算法,编程求解一个数值积分的问题。
(程序参考:
#include"
stdio.h"
conio.h"
math.h"
#defineN101
#defineCtrl0.000000001
#definefunction(a,b,d,x)a/(b+pow(x,d))
intmain()
doublea,b,d,up,down,t[N]={0},s[N]={0},c[N]={0},r[N]={0},h=0,sum=0,fa=0,fb=0;
inti=0,j=0,goon=0;
label:
本程序只支持函数是a/(b+x^d)请输入a,b,d:
\n"
scanf("
%lf%lf%lf"
a,&
b,&
d);
请输入积分上下限(a,b):
%lf%lf"
up,&
down);
Theresultsareasfollows!
fa=a/(b+pow(down,d));
fb=a/(b+pow(up,d));
t[0]=0.5*(fa+fb);
t[0]=%lf\n"
t[0]);
for(i=1;
20;
sum=0;
h=(up-down)/pow(2,i);
//h认为在每次计算中都是变化的//
h=%lf\n"
h);
for(j=1;
pow(2,i);
j+=2)
sum+=a/(b+pow(down+j*h,d));
//printf("
sum=%lf\n"
sum);
t[i]=0.5*t[i-1]+h*sum;
s[i-1]=(4*t[i]-t[i-1])/3.0;
if(i>
=2)
c[i-2]=(16*s[i-1]-s[i-2])/15.0;
=3)
r[i-3]=(64*c[i-2]-c[i-3])/63.0;
t[%lf]=%lf!
pow(2,i),t[i]);
Risasfollows!
i<
15;
r[%lf]=%lf,\n"
pow(2,i),r[i]);
Cisasfollows!
c[%lf]=%lf,\n"
pow(2,i),c[i]);
Sisasfollows!
s[%lf]=%lf,\n"
pow(2,i),s[i]);
Continuetouse,press1!
%d"
goon);
if(goon==1)
gotolabel;
else
return0;
getch();
实验四龙格-库塔法
1.用标准四阶龙格-库塔方法设计算法
2.编程解微分方程初值问题;
#include<
/*n表示几等分,n+1表示他输出的个数*/
intRungeKutta(doubley0,doublea,doubleb,intn,double*x,double*y,intstyle,double(*function)(double,double))
doubleh=(b-a)/n,k1,k2,k3,k4;
inti;
//x=(double*)malloc((n+1)*sizeof(double));
//y=(double*)malloc((n+1)*sizeof(double));
x[0]=a;
y[0]=y0;
switch(style)
case2:
for(i=0;
x[i+1]=x[i]+h;
k1=function(x[i],y[i]);
k2=function(x[i]+h/2,y[i]+h*k1/2);
y[i+1]=y[i]+h*k2;
break;
case3:
i++){
k3=function(x[i]+h,y[i]-h*k1+2*h*k2);
y[i+1]=y[i]+h*(k1+4*k2+k3)/6;
case4:
k3=function(x[i]+h/2,y[i]+h*k2/2);
k4=function(x[i]+h,y[i]+h*k3);
y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6;
default:
return0;
return1;
doublefunction(doublex,doubley)
returny-2*x/y;
//例子求y'
=y-2*x/y(0<
x<
1);
y0=1;
doublex[6],y[6];
用二阶龙格-库塔方法\n"
RungeKutta(1,0,1,5,x,y,2,function);
for(inti=0;
6;
x[%d]=%f,y[%d]=%f\n"
i,x[i],i,y[i]);
用三阶龙格-库塔方法\n"
RungeKutta(1,0,1,5,x,y,3,function);
用四阶龙格-库塔方法\n"
RungeKutta(1,0,1,5,x,y,4,function);
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 编程 实验 报告