现将小生研一数值分析作业摘录如下,上机的代老师要求我们分为四个部分:解题理论依据、计算程序、计算结果截图和问题讨论,因公式输入不便所以第一部分的公式略去,需要的可去第一篇博客“研究生数值分析课后题(上机编程)-1 ”下载附件。因为本人编程水平有限,故使用的是最基本的函数,但保证原创和结果正确,所有程序在VC++6.0上验证通过,还望各位兄台不吝赐教。
第四题:牛顿法求方程近似根
1.解题理论依据或方法应用条件:
依据牛顿迭代法进行计算
2.计算程序(使用软件:VC):
#include<math.h>
void main()
{double x0,x1=1.9,f0,f1;
for(x0=1.0;fabs(x0-x1)>=0.00001;)
{x0=x1;
f0=pow(x0,7)-28*pow(x0,4)+14;
f1=7*pow(x0,6)-112*pow(x0,3);
if(f1==0)break;
else x1=x0-f0/f1;}
if(f1==0)printf("出错\n");
else printf("方程的解为:%.16lf\n",x1);
}
3.计算结果
4.问题讨论(误差分析、上机出现情况等)
(1)在做除法运算的时候,应对除数是否等于零进行判断。
(2)编程的过程出现问题如下:①在编写程序过程中,不小心将输入换成了全角输入造成错误②忘记头文件中增加#include<math.h>,导致pow函数无法使用。
第五题:用Romberg算法求积分近似值 1.解题理论依据或方法应用条件: 依据Romberg法列出相应的方程,先求出第一列的相应元素,然后利用Romberg中三角法求出后面列的元素,将得到的第一行元素取差的绝对值,当小于误差要求的时候可以停止运算 2.计算程序(使用软件:VC): 方法一: #include<stdio.h> #include<math.h> #define f(x) pow(3,(x))*pow((x),1.4)*(5*(x)+7)*sin((x)*(x)) void main() {double t[9][9]; double sum=0; int j=1,i=1,m=0,k=1; double a=1,b=3; t[0][1]=(b-a)/2*(f(a)+f(b)); sum=f(a+(2*1-1)*(b-a)/pow(2,1)); t[1][1]=0.5*(t[0][1]+sum*(b-a)/pow(2,0)); t[0][2]=(pow(4,1)*t[1][1]-t[0][1])/(pow(4,1)-1); /*先求出t[0][1],t[0][2]的值*/ for(j=1;fabs(t[0][j+1]-t[0][j])>0.00001;) /*判断误差是否符合条件*/ {j++; sum=0; /*累加和清零*/ for(i=1;i<=pow(2,j-1);i++) sum+=f(a+(2*i-1)*(b-a)/pow(2,j)); /*求出式子中的累加清零*/ t[j][1]=0.5*(t[j-1][1]+sum*(b-a)/pow(2,j-1)); /*求出第一列元素的值*/ for(k=j;k>=1;k--) {m++;t[k-1][m+1]=(pow(4,m)*t[k][m]-t[k-1][m])/(pow(4,m)-1);} /*依次用三角形法求出其他元素的值*/ m=0; } printf("积分的值为:%.16lf\n",t[0][j+1]); } 方法二:(逐列求值,较浪费电脑资源、时间) #include<stdio.h> #include<math.h> #define f(x) pow(3,(x))*pow((x),1.4)*(5*(x)+7)*sin((x)*(x)) void main() {double t[9][9]; double sum; int j=1,i=1,m=1,k=1,l=1; double a=1,b=3; t[0][1]=(b-a)/2*(f(a)+f(b)); for(l=1;l<=7;l++) {sum=0; for(i=1;i<=pow(2,l-1);i++) sum+=f(a+(2*i-1)*(b-a)/pow(2,l)); t[l][1]=0.5*(t[l-1][1]+sum*(b-a)/pow(2,l-1)); } for(m=1;m<=7;m++) for(k=1;k<=(7-m+1);k++) t[k-1][m+1]=(pow(4,m)*t[k][m]-t[k-1][m])/(pow(4,m)-1); for(j=1;fabs(t[0][j+1]-t[0][j])>0.00001;j++); printf("积分的值为:%.16lf\n",t[0][j+1]); } 3.计算结果 4.问题讨论(误差分析、上机出现情况等) (1)这道题目是前三道题目中耗时最长、但也是学到最多东西的题目,一开始使用的是第二个方法进行编程,VC报错,一直未找到哪里出错,最后发现原来是将 int j=1,i=1,m=1,k=1,l=1; 放在了 t[0][1]=(b-a)/2*(f(a)+f(b)); 之后,VC规定必须将所有的定义放在程序开头,这是自己以前所不知道的。 (2)改掉上面的错误后,系统不再报错,但是数值不对,逐行看程序发现原来是sum未清零,导致后面的数据出错。 (3)改掉上面的错误后,得到的数据依然不对,纠缠很久以及翻阅C语言书籍发现原来是#define f(x) pow(3,(x))*pow((x),1.4)*(5*(x)+7)*sin((x)*(x)) 语句中,X必须加上括号,这样才能保证运算顺序正确。 (4)这样程序的结果正确了,但是第二种方法较浪费电脑资源、时间,故对原程序更改有了第一种利用循环、用多少数据求多少数据的方法。 第六题:定步长四阶Runge-Kutta法求方程组 1.解题理论依据或方法应用条件: 依据定步长四阶Runge-Kutta法进行计算 2.计算程序(使用软件:VC): #include<stdio.h> main() {float y[4][205]={{0},{0},{0},{0}}; /*定义第一列值全为零*/ int i,j; float k1,k2,k3,k4; float h=0.0005; for(j=1;(j*0.0005)<=0.1;j++) for(i=1;i<=3;i++) /*本解法未使用0行,故从1开始取值*/ {if(i=1) {k1=1;k2=1;k3=1;k4=1; y[1][j]=y[1][j-1]+h/6*(k1+2*k2+2*k3+k4);} /*计算第一行数值*/ if(i=2) {k1=y[3][j-1]; k2=y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1]); k3=y[3][j-1]+0.5*h*(1000-1000*(y[2][j-1]+0.5*h*y[3][j-1])-100*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1]))); k4=y[3][j-1]+h*(1000-1000*(y[2][j-1]+0.5*h*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1])))-100*(y[3][j-1]+0.5*h*(1000-1000*(y[2][j-1]+0.5*h*y[3][j-1])-100*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1]))))); y[2][j]=y[2][j-1]+h/6*(k1+2*k2+2*k3+k4);} /*计算第二行数值*/ if(i=3) {k1=1000-1000*y[2][j-1]-100*y[3][j-1]; k2=1000-1000*(y[2][j-1]+0.5*h*y[3][j-1])-100*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1])); k3=1000-1000*(y[2][j-1]+0.5*h*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1])))-100*(y[3][j-1]+0.5*h*(1000-1000*(y[2][j-1]+0.5*h*y[3][j-1])-100*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1])))); k4=1000-1000*(y[2][j-1]+h*(y[3][j-1]+0.5*h*(1000-1000*(y[2][j-1]+0.5*h*y[3][j-1])-100*(y[3][j-1]+0.5*h*(1000-1000*y[2][j-1]-100*y[3][j-1])))))-100*(y[3][j-1]+h*k3); y[3][j]=y[3][j-1]+h/6*(k1+2*k2+2*k3+k4);} /*计算第三行数值*/ } printf("y1(0.025),y2(0.025),y3(0.025)值分别为:%f,%f,%f\n",y[1][50],y[2][50],y[3][50]); printf("y1(0.045),y2(0.045),y3(0.045)值分别为:%f,%f,%f\n",y[1][90],y[2][90],y[3][90]); printf("y1(0.085),y2(0.085),y3(0.085)值分别为:%f,%f,%f\n",y[1][170],y[2][170],y[3][170]); printf("y1(0.100),y2(0.100),y3(0.100)值分别为:%f,%f,%f\n",y[1][200],y[2][200],y[3][200]); } 3.计算结果 4.问题讨论(误差分析、上机出现情况等) (1)本题的编写过程相对前一题较为顺利,过程中只出现一个较大问题,原程序中我在for循环体内使用的是if(){} else if(){} else (){} 的结构,无错,但运行不出结果,改用if (){} if (){} if (){} 后出现结果,至今仍未找清楚原因。 (2)另一个小问题是,在计算第三行数据时,由于要使用前一列二行的数据,刚开始直接将复制进来,稍后发现这时的K值是二行中的数据,放在第三行中将出错,由于发现及时,未对后面程序编写造成影响。(故后续将K值都由具体的算式替代,避免出错)