计算方法与实习上机题问题详解

实习题1

1用两种不容的顺序计算

644834.11000

1

2

≈∑=-n n

，分析误差的变化

（1）顺序计算 源代码：

#include #include void main() {

double sum=0; int n=1; while(1) {

sum=sum+(1/pow(n,2));

if(n%1000==0)printf("sun[%d]=%-30f",n,sum); if(n>=10000)break; n++; }

printf("sum[%d]=%f\n",n,sum); }

结果：

（2）逆序计算 源代码：

#include #include void main() {

double sum=0; int n=10000; while(1) {

sum=sum+(1/pow(n,2));

if(n%1000==0)

printf("sum[%d]=%-30f",n,sum); if(n<=1)break; n--; }

printf("sum[%d]=%f\n",n,sum); }

结果：

2已知连分数

))//(.../(322101

n n b a a b a b a b f ++++

=

利用下面的方法计算f:

1

1

)0,...,2,1(,d f n n i d a b d b d i i i i n n =--=+==++

写一个程序，读入n,n n b a ,,计算并打印f

源代码：

#include #include void main() {

int i=0,n;

float a[1024],b[1024],d[1024]; printf("please input n,n="); scanf("%d",&n);

printf("\nplease input a[1] to a[n]:\n"); for(i=1;i<=n;i++) {

printf("a[%d]=",i);

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

}

printf("\nplease input b[0] to b[n]:\n"); for(i=0;i<=n;i++) {

printf("b[%d]=",i); scanf("%f",&b[i]); }

d[n]=b[n];

for(i=n-1;i>=0;i--)

d[i]=b[i]+a[i+1]/d[i+1]; printf("\nf=%f\n",d[0]); }

结果：

3给出一个有效的算法和一个无效的算法计算积分

?=+=10

)10,...1,0(14n dx x x y n n

源代码：

#include #include main() {

double y_0=(1/4.0)*log(5),y_1;

double y_2=(1.0/55.0+1.0/11.0)/2,y_3; int n=1,m=10;

printf("有效算法输出结果：\n"); printf("y[0]=%-20f",y_0);

while(1) {

y_1=1.0/(4*n)+y_0/(-4.0); printf("y[%d]=%-20f",n,y_1); if(n>=10) break; y_0=y_1; n++;

if(n%3==0)

printf("\n"); }

printf("\n 无效算法的输出结果：\n"); printf("y[10]=%-20f",y_2); while(1) {

y_3=1.0/n-4.0*y_2;

printf("y[%d]=%-20f",m-1,y_3); if(m<=1)

break; y_2=y_3; m--;

if(m%2==0) printf("\n"); } }

结果：

4设∑

=-=

N

j N j S 2

211

，已知其精确值为

)1

1123(21+--N N （1）编制按从小到大顺序计算N S 的程序 （2）编制按从小达到的顺序计算N S 的程序

（3）按两种顺序分别计算30000100001000,,S S S ,并指出有效位数

源代码：

#include main() {

int N;

double SN[30000];

SN[30000]=(3.0/2.0-1.0/30000.0-1/30001.0)/2.0; for(N=30000;N>=2;N--)

SN[N-1]=SN[N]-1.0/(N*N-1); printf("从大到小顺序计算： \nSN[1000]=%f\nSN[10000]=%f\nSN[30000]=%f\n",SN[1000],SN[10000],SN[30000]); SN[2]=(3.0/2-1.0/2.0-1/3.0)/2.0; for(N=3;N<=30000;N++)

SN[N]=SN[N-1]+1.0/(N*N-1); printf("从小到大顺序计算： \nSN[1000]=%f\nSN[10000]=%f\nSN[30000]=%f\n",SN[1000],SN[10000],SN[30000]); }

结果：

实习题2 1.用牛顿法求下列方程的根

2=-x e x

01=-x xe 02lg =-+x x

源代码：

#include #include

typedef float (*p)(float ); float ff1(float x) {

return x*x-exp(x); }

float ff2(float x) {

return x*exp(x)-1; }

float ff3(float x) {

return log(x)+x-2; }

float answer(float(*p)(float)) {

int k=2;

float m=1,n=-1,x2,a,b,c; if (p==ff3)n=2;

printf("x[0] = %.4f, x[1] = %.4f, ",m,n); while (1) {

if (fabs(m-n)<1e-4) break; a=p(n)*(n-m); b=p(n)-p(m); c=a/b;

x2=n-c; m = n; n = x2;

printf("x[%d] = %.4f, ",k,x2); k++;

if (k%3==0) printf("\n"); }

if (k%3!=0) printf("\n");

printf("iteration times: %d, roots: %.4f\n ",k-2,n);

return 0;

}

main()

{

printf("x*x-exp(x),\n");

answer(ff1);

printf("x*exp(x)-1,\n");

answer(ff2);

printf("lg(x)+x-2,\n");

answer(ff3);

return 0;

}

结果：

2.编写一个割线法的程序，求解上述各方程

源代码：

#include

#include

float gexian(float,float);

float f(float);

main()

{

int i,j;

float x1=2.2;

float x2=2,x3;

scanf("%d",&i);

if(i==1)

printf("%f",x1); else if(i==2)

printf("%f",x2); else {

for(j=3;j<=i;j++) {

x3=gexian(x1,x2); x1=x2; x2=x3; }

printf("%f",gexian(x1,x2)); } }

float f(float x) {

return (x*x-exp(x)); }

float gexian(float x1,float x2) {

return (x2-(f(x2)/(f(x2)-f(x1)))*(x2-x1)); }

结果：

实习题3 1用列主元消去法解下列方程组；

??????

?=++=-++--=+--=--+4

3443233312)1(42143214

3214321x x x x x x x x x x x x x x x ??????

?=++--=++-=-+--=-+-4

341220332282)2(4321

3

214

3214321x x x x x x x x x x x x x x x 源程序：

#include

#include

void ColPivot(float*,int,float[]);

void ColPivot(float*c,int n,float x[])

{

int i,j,t,k;

float p;

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

{

k=i;

for(j=i+1;j<=n-1;j++)

if(fabs(*(c+j*(n+1)+i))>(fabs(*(c+k*(n+1)+i))))

k=j;

if(k!=i)

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

{

p=*(c+i*(n+1)+j);

*(c+i*(n+1)+j)=*(c+k*(n+1)+j);

*(c+k*(n+1)+j)=p;

}

for(j=i+1;j<=n-1;j++)

{

p=(*(c+j*(n+1)+i))/(*(c+i*(n+1)+i));

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

*(c+j*(n+1)+t)-=p*(*(c+i*(n+1)+t));

}

}

for(i=n-1;i>=0;i--)

{

for(j=n-1;j>=i+1;j--)

(*(c+i*(n+1)+n))-=x[j]*(*(c+i*(n+1)+j));

x[i]=*(c+i*(n+1)+n)/(*(c+i*(n+1)+i));

}

}

void main()

{

int i;

float x[4];

float c[4][5]={1,1,0,3,4,2,1,-1,1,1,3,-1,-1,3,-3,-1,2,3,-1,4};

ColPivot(c[0],4,x);

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

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

}

结果：

第（1）题

第（2）题

2、

源代码：

#include

void main()

{

float x[4];

int i;

float a[4][5]={48,-24,0,-12,4,-24,24,12,12,4,0,6,20,2,-2,-6,6,2,16,-2};

void DirectLU(float*,int,float[]);

DirectLU(a[0],4,x);

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

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

}

void DirectLU(float*u,int n,float x[])

{

int i,r,k;

for(r=0;r<=n-1;r++)

{

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

for(k=0;k<=r-1;k++)

*(u+r*(n+1)+i)-=*(u+r*(n+1)+k)*(*(u+k*(n+1)+i));

for(i=r+1;i<=n-1;i++)

{

for(k=0;k<=r-1;k++)

*(u+i*(n+1)+r)-=*(u+i*(n+1)+k)*(*(u+k*(n+1)+r));

*(u+i*(n+1)+r)/=*(u+r*(n+1)+r);

}

}

for(i=n-1;i>=0;i--)

{

for(r=n-1;r>=i+1;r--)

*(u+i*(n+1)+n)-=*(u+i*(n+1)+r)*x[r];

x[i]=*(u+i*(n+1)+n)/(*(u+i*(n+1)+i));

}

}

实习题4

1、

源代码：

#include

float Lagrange(float x[],float y[],float xx,int n) //n为（n+1）次插值；{

int i,j;

float *a,yy=0;

a=new float[n];

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

{

a[i]=y[i];

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

if(j!=i)a[i]*=(xx-x[j])/(x[i]-x[j]);

yy+=a[i];

}

delete a;

return yy;

}

void main()

{

float x[5]={-3.0,-1.0,1.0,2.0,3.0};

float y[5]={1.0,1.5,2.0,2.0,1.0};

float xx1=-2,xx2=0,xx3=2.75,yy1,yy2,yy3;

yy1=Lagrange(x,y,xx1,3);

yy2=Lagrange(x,y,xx2,3);

yy3=Lagrange(x,y,xx3,3);

printf("x1=%-20f,y1=%f\n",xx1,yy1);

printf("x2=%-20f,y2=%f\n",xx2,yy2);

printf("x3=%-20f,y3=%f\n",xx3,yy3);

}

结果：

2、

源代码：

#include

float Lagrange(float x[],float y[],float xx,int n) //n为（n+1）次插值；{

int i,j;

float *a,yy=0;

a=new float[n];

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

{

a[i]=y[i];

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

if(j!=i)a[i]*=(xx-x[j])/(x[i]-x[j]);

yy+=a[i];

}

delete a;

return yy;

}

void main()

{

float x[6]={0.30,0.42,0.50,0.58,0.66,0.72};

float y[6]={1.04403,1.08462,1.11803,1.15603,1.19817,1.23223};

float xx1=0.46,xx2=0.55,xx3=0.60,yy1,yy2,yy3;

yy1=Lagrange(x,y,xx1,6);

yy2=Lagrange(x,y,xx2,6);

yy3=Lagrange(x,y,xx3,6);

printf("x1=%-20f,y1=%f\n",xx1,yy1);

printf("x2=%-20f,y2=%f\n",xx2,yy2);

printf("x3=%-20f,y3=%f\n",xx3,yy3);

}

结果：

源代码：

#include

#define N 3

void Difference(float y[],float f[4][4],int n)

{

int k,i;

f[0][0]=y[0];f[1][0]=y[1];f[2][0]=y[2];f[3][0]=y[3];

for(k=1;k<=n;k++)

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

f[i][k]=f[i+1][k-1]-f[i][k-1];

return;

}

void main()

{

int i,k=1;

float a,b=1,m=21.4,t=1.4,f[4][4]={0};

float x[5]={20,21,22,23,24};

float y[5]={1.30103,1.32222,1.34242,1.36173,1.38021};

Difference(y,f,N);

a=f[0][0];

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

{

k=k*i;

b=b*(t-i+1);

a=a+b*f[0][i]/k;

}

printf("x(k)\n");

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

printf( "%-20f",x[i]);

printf("\ny(k)\n");

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

printf("%-20f",y[i]);

for(k=1;k<=3;k++)

{

printf("\nF(%d)\n ",k);

for(i=0;i<=(3-k);i++)

{

printf("%-20f",f[i][k]);

}

}

printf ("\n");

printf("f(%f)=%-20f",m,a);

printf ("\n");

结果：

实习题5

2、

源代码：

#include

#include

void main()

{

int i,n;

float a[2];

float x[15]={1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8},z[15];

float

y[15]={33.4,79.50,122.65,159.05,189.15,214.15,238.65,252.50,267.55,280.50,296.6 5,301.40,310.40,318.15,325.15};

for(n=0;n<=14;n++) //增加了数组z；

{

z[n]=log(y[n]/x[n]);

}

void Approx(float[],float[],int,int,float[]);

Approx(x,z,15,1,a); //变成一次拟合；

//for(i=0;i<=1;i++)

//printf("a[%d]=%f\n",i,a[i]);

printf("a=exp(a[0])=%f\n",exp(a[0]));

printf("b=-a[1]=%f\n",-a[1]); }

void Approx(float x[],float y[],int m,int n,float a[])

{

int i,j,t;

float *c=new float[(n+1)*(n+2)];

float power(int,float);

void ColPivot(float *,int,float[]);

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

{

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

{

*(c+i*(n+2)+j)=0;

for(t=0;t<=m-1;t++)

*(c+i*(n+2)+j)+=power(i+j,x[t]);

}

*(c+i*(n+2)+n+1)=0;

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

*(c+i*(n+2)+n+1)+=y[j]*power(i,x[j]);

}

ColPivot(c,n+1,a);

delete c;

}

void ColPivot(float *c,int n,float x[])

{

int i,j,t,k;

float p;

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

{

k=i;

for(j=i+1;j<=n-1;j++)

if(fabs(*(c+j*(n+1)+i))>(fabs(*(c+k*(n+1)+i))))

k=j;

if(k!=i)

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

{

p=*(c+i*(n+1)+j);

*(c+i*(n+1)+j)=*(c+k*(n+1)+j);

*(c+k*(n+1)+j)=p;

}

for(j=i+1;j<=n-1;j++)

{

p=(*(c+j*(n+1)+i))/(*(c+i*(n+1)+i));

for(t=i;t<=n;t++)*(c+j*(n+1)+t)-=p*(*(c+i*(n+1)+t));

}

}

for(i=n-1;i>=0;i--)

{

for(j=n-1;j>=i+1;j--)

(*(c+i*(n+1)+n))-=x[j]*(*(c+i*(n+1)+j));

x[i]=*(c+i*(n+1)+n)/(*(c+i*(n+1)+i));

}

float power(int i,float v)

{

float a=1;

while(i--)a*=v;

return a;

}

结果：

实习题6

1、

源代码：

（1）

#include

#include

float Cotes(float(*f)(float),float a,float b,int n)

{

int k;

float c,c1=0,c2,c3,c4;

float h=(b-a)/n;

c2=(*f)(a+h/4);

c3=(*f)(a+h/2);

c4=(*f)(a+3*h/4);

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

{

c1+=(*f)(a+k*h);

c2+=(*f)(a+k*h+h/4);

c3+=(*f)(a+k*h+h/2);

c4+=(*f)(a+k*h+3*h/4);

}

c=h/90*(7*((*f)(a)+(*f)(b))+14*c1+32*c2+12*c3+32*c4);

return c;

}

float f(float x)

{

return 1/sqrt(1+x*x*x);

}

void main()

int i,n=4;

float c;

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

{

c=Cotes(f,0,1,n);

printf("C(%d)=%f\n",n,c);

n*=2;

}

}

（2）

#include

#include

float Cotes(float(*f)(float),float a,float b,int n)

{

int k;

float c,c1=0,c2,c3,c4;

float h=(b-a)/n;

c2=(*f)(a+h/4);

c3=(*f)(a+h/2);

c4=(*f)(a+3*h/4);

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

{

c1+=(*f)(a+k*h);

c2+=(*f)(a+k*h+h/4);

c3+=(*f)(a+k*h+h/2);

c4+=(*f)(a+k*h+3*h/4);

}

c=h/90*(7*((*f)(a)+(*f)(b))+14*c1+32*c2+12*c3+32*c4);

return c;

}

float f(float x)

{ // return 1/sqrt(1+x*x*x);

if (x==0)

return 1;

else return sin(x)/x;

}

void main()

{

int i,n=4;

float c;

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

{

// c=Cotes(f,0,1,n);

c=Cotes(f,0,5,n);

printf("C(%d)=%f\n",n,c);

n*=2;

}

}

结果：

（1）

（2）

实习题7 一、改进欧拉法

1、

#include

#include

double Adams (

double (*f)(double x, double y),

double x0,

double y0,

double h,

int N

) {

for(int n=0; n

double x1=x0+h;

double yp=y0+h*f(x0,y0);

double y1=y0+h*f(x1,yp);

y1=(yp+y1)/2.0;

printf("ty=%f\n",y1);

x0=x1;

y0=y1;

}

}

int main(void)

{

double f(double x, double y); double x0=0,y0=0;

double x,y,step;

long i;

step=0.1;

Adams(f,x0,y0,0.1,10);

}

double f(double x, double y) {

double r;

r=x*x+y*y;

return(r);

}

2、

int main(void)

{

double f(double x, double y); double x0=0,y0=0;

double x,y,step;

long i;

step=0.1;

Adams(f,x0,y0,0.1,10);

}

double f(double x, double y) {

double r;

r=1/(1+y*y);

return(r);

}

3、

int main(void)

{

double f(double x, double y); double x0=0,y0=1;

double x,y,step;

long i;

step=0.1;

Adams(f,x0,y0,0.1,10);

}

double f(double x, double y) {

double r;

r=y-2*x/y;