浙大-计算方法作业

来自14级机械专业的《计算方法》课程课本为黄皮小书只有当年作业布置的部分题目2.4.7.8.10.13.17

第2题

二分法

#include

#include

double calculate_y(double n)

{

using namespace std;

return exp(n)+10*n-2;

}

int main()

{

using namespace std;

int k;

double a,b,c1;

double y1,y2,x,y;

cout<<"请依次输入a（左端点），b（右端点），c（控制精度），中间用enter隔开"<

cin>>a>>b>>c1;

y1 = calculate_y(a);

y2 = calculate_y(b);

if(y1*y2<=0);

else{

do{

cout<<"["<

cout<<"实根不再此范围内，请重新输入！" <

cin>>a>>b>>c1;

y1 = calculate_y(a);

y2 = calculate_y(b);

}while(y1*y2>0);

}

k=0;

k=k+1;

x=(a+b)/2;

y=calculate_y(x);

if(y1*y<0) b=x;

else {a=x;y1=y;}

}while(b-a>=c1);

cout<<"k="<

cin.get();

cin.get();

}

迭代法

#include

#include

int main()

{

using namespace std;

int k=1,N;

double x0,c,x1;

cout<<"请依次输入N，x0，c"<

cin>>N>>x0>>c;

do{

k=k+1;

x1=(2-exp(x0))/10;

if(fabs(x1-x0)

{cout<<"k="<

cin.get();

cin.get();

return 0;

}

牛顿法

#include

#include

double fx(double n1)

{

using namespace std;

return exp(n1)+10*n1-2;

}

double f1x(double n2)

{

using namespace std;

return exp(n2)+10;

}

int main()

{

using namespace std;

int N,k=0,I;

double x0,x1,c,f_x0,f1_x0;

cout<<"请依次输入x0（迭代初值），c（精度），N（最大允许迭代次数）"<

cin>>x0>>c>>N;

do{

k++;

f_x0=fx(x0);

f1_x0=f1x(x0);

if(f1_x0==0){I=-1;break;}

else

x1=x0-f_x0/f1_x0;

}

if(fabs(x1-x0)

}

cout<

}

for(i=0;i

cout<

double t=0.0;

x[N-1]=b[N-1]/a[N-1][N-1];

for(i=N-2;i>=0;i--)

{

x[i]=b[i];

for(j=i+1;j

{

x[i]-=a[i][j]*x[j];

}

x[i]/=a[i][i];

}

cout<<"方程的解依次为"<

for(i=0;i

{

cout<<"x["<

}

}

列主元消元法

#include

#include

using namespace std;

#define N 4

void display(double (*a)[N],double *b)

{

for(int i=0;i

{

for(int j=0;j

{

cout<

}

cout<

}

}

int main()

{

double a[N][N]={{1.1348,3.8326,1.1651,3.4017}, {0.5301,1.7875,2.5330,1.5435},

{3.4129,4.9317,8.7643,1.3142},

{1.2371,4.9998,10.6721,0.0147}};

double b[N]={9.5342,6.3941,18.4231,16.9237};

double m[N][N];

double x[N];

int k,i,j,l;

double max;

double p;

for(k=0;k

{

max=a[k][k];

l=k;

for(i=k+1;i

{

if(fabs(a[i][k])>fabs(max))

{

max=a[i][k];

l=i;

}

}

if(fabs(max)*100000<1)

cout<<"det A=0"<

else{

if(l!=k)

{

p=b[l];b[l]=b[k];b[k]=p;

for(j=0;j

{

p=a[l][j];

a[l][j]=a[k][j];

a[k][j]=p;

}

display(a,b);

cout<

}

for(i=k+1;i

{

m[i][k]=a[i][k]/a[k][k];

for(j=k+1;j

{

a[i][j]-=m[i][k]*a[k][j];

}

b[i]-=m[i][k]*b[k];

}

display(a,b);

cout<

}

}

if(fabs(a[N-1][N-1])*100000<1)

cout<<"det A=0"<

else{

x[N-1]=b[N-1]/a[N-1][N-1];

for(i=N-2;i>=0;i--)

{

x[i]=b[i];

for(j=i+1;j

{

x[i]-=a[i][j]*x[j];

}

x[i]/=a[i][i];

}

}

cout<<"方程的解依次为"<

for(i=0;i

{

cout<<"x["<

}

第7题

雅可比迭代法

#include

#include

#include

using namespace std;

int main()

{

double s=0;

double max(double array[6]);

double a[6][6]={

{4,-1,0,-1,0,0},{-1,4,-1,0,-1,0},

{0,-1,4,-1,0,-1},{-1,0,-1,4,-1,0},

{0,-1,0,-1,4,-1},{0,0,-1,0,-1,4,},

};

double b[6]={0,5,-2,5,-2,6};

double c[6]={0,0,0,0,0,0};

double x[6]={0,0,0,0,0,0};

double x0[6]={0,0,0,0,0,0};

int i,k,j;

double r,sum=0;

cout<<"输入精度："<

cin>>s;

for(k=1;;k++)

{

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

{

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

{

sum=a[i][j]*x0[j]+sum;

}

x[i]=x0[i]+((b[i]-sum)/a[i][i]);

c[i]=fabs(x[i]-x0[i]);

sum=0;

}

r=max(c);

if(r

{for(i=0;i<6;i++)

cout<<"x"<

return 0;

}

else

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

x0[i]=x[i];

}

}

double max(double array[6])

{

double a=array[0];

int i;

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

{

if(a

a=array[i];

}

return a;

}

高斯-赛德尔迭代法

#include

#include

#include

using namespace std;

int main()

{

double s=0;

double max(double array[3]);

double a[6][6]={

{4,-1,0,-1,0,0},{-1,4,-1,0,-1,0},

{0,-1,4,-1,0,-1},{-1,0,-1,4,-1,0},

{0,-1,0,-1,4,-1},{0,0,-1,0,-1,4},

};

double b[6]={0,5,-2,5,-2,6};

double c[6]={0,0,0,0,0,0};

double x[6]={0,0,0,0,0,0};

double x0[10]={0,0,0,0,0,0};

int i,k,j;

double r,sum=0;

cout<<"输入精度："<

cin>>s;

for(k=1;;k++)

{

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

{

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

{

sum=a[i][j]*x0[j]+sum;

}

x[i]=x0[i]+((b[i]-sum)/a[i][i]);

c[i]=x[i]-x0[i];

if(c[i]<0)

c[i]=-c[i];

x0[i]=x[i];

sum=0;

}

r=max(c);

if(r

{for(i=0;i<6;i++)

cout<<"x"<

return 0;

}

}

}

double max(double array[6])

{

double a=array[0];

int i;

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

{

if(a

a=array[i];

}

return a;

}

SOR迭代法

#include

#include

#include

using namespace std;

int main()

{

double s=0,w=0;

double max(double array[3]);

double a[6][6]={

{4,-1,0,-1,0,0},{-1,4,-1,0,-1,0},

{0,-1,4,-1,0,-1},{-1,0,-1,4,-1,0},

{0,-1,0,-1,4,-1},{0,0,-1,0,-1,4}

};

double b[6]={0,5,-2,5,-2,6};

double c[6]={0,0,0};

double x[6]={0,0,0,0,0,0};

double x0[6]={0,0,0,0,0,0};

int i,k,j;

double r,sum=0;

cout<<"输入精度："<

cin>>s;

cout<<"松弛因子："<

cin>>w;

for(k=1;;k++)

{

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

{

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

{

sum=a[i][j]*x0[j]+sum;

}

x[i]=x0[i]+(w*(b[i]-sum)/a[i][i]);

c[i]=x[i]-x0[i];

if(c[i]<0)

c[i]=-c[i];

x0[i]=x[i];

sum=0;

}

r=max(c);

if(r

{for(i=0;i<6;i++)

cout<<"x"<

cout<<"迭代次数："<

return 0;

}

}

}

double max(double array[6])

{

double a=array[0];

int i;

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

{

if(a

a=array[i];

}

return a;

}

第8题

分段线性插值

#include

int main(void)

{

int i,j;

double a[6],b[6],x,y;

printf("Enter x:");

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

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

printf("Enter f(x):");

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

scanf("%lf",&b[i]);

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

printf("Enter x:");

scanf("%lf",&x);

for(i=0;x>a[i];i++);

y=b[i-1]*(x-a[i])/(a[i-1]-a[i])+b[i]*(x-a[i-1])/(a[i]-a[i-1]); printf("x=%lf , ",x);

printf("y=%lf\n",y);

}

return 0;

}

分段二次插值

#include

int main(void)

{

int i,j;

double a[6],b[6],x,y;

printf("Enter x:");

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

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

printf("Enter f(x):");

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

scanf("%lf",&b[i]);

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

printf("Enter x:");

scanf("%lf",&x);

for(i=0;x>a[i];i++);

if((a[i+1]-x)<(x-a[i-2]))

{

y=b[i-1]*(x-a[i])/(a[i-1]-a[i])*(x-a[i+1])/(a[i-1]-a[i+1]);

y=y+b[i]*(x-a[i-1])/(a[i]-a[i-1])*(x-a[i+1])/(a[i]-a[i+1]);

y=y+b[i+1]*(x-a[i-1])/(a[i+1]-a[i-1])*(x-a[i])/(a[i+1]-a[i]);

printf("x=%lf , ",x); printf("y=%lf\n",y);

}

else{

y=b[i-2]*(x-a[i-1])/(a[i-2]-a[i-1])*(x-a[i])/(a[i-2]-a[i]);

y=y+b[i-1]*(x-a[i-2])/(a[i-1]-a[i-2])*(x-a[i])/(a[i-1]-a[i]);

y=y+b[i]*(x-a[i-2])/(a[i]-a[i-2])*(x-a[i-1])/(a[i]-a[i-1]);

printf("x=%lf , ",x); printf("y=%lf\n",y);

}

}

return 0;

}

全区间上拉格朗日插值

#include

int main(void)

{

int i,j,k;

double a[6],b[6],x,y,p;

printf("Enter x:");

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

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

printf("Enter f(x):");

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

scanf("%lf",&b[i]);

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

{

y=0;

printf("Enter x:");

scanf("%lf",&x);

for(i=0;i<6;i++){

p=1;

for(k=0;k<6;k++){

if(k!=i)

p=p*(x-a[k])/(a[i]-a[k]);}

y=y+p*b[i];

}

printf("x=%lf , ",x); printf("y=%lf\n",y);

}

return 0;

}

第10题

最小二乘法二次回归

#include

#include

#include

float cf(int i,float f)

{ float g=1.0;

while(i--)

g*=f;

return g;

}

void ColPivot(float a[99][99], int n, float b[99][1],float x[99])

{

float max,m,temp;

int k,q=0,i,j,I;

for(k=0;k

{ max=a[k][k];

for(i=k+1;i

if(fabs(max)

{

max=a[i][k];

I=i;

}

if(max==0)

q=1;

else

{ if(I!=k)

{

temp=b[I][0]; b[I][0]=b[k][0]; b[k][0]=temp;

for(j=k;j

{

temp=a[I][j]; a[I][j]=a[k][j]; a[k][j]=temp;

}

}

for(i=k+1;i

{

m=a[i][k]/a[k][k];

b[i][0]=b[i][0]-b[k][0]*m;

for(j=0;j

a[i][j]=a[i][j]-a[k][j]*m;

}

}

}

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

{

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

b[i][0]=b[i][0]-a[i][j]*x[j];

x[i]=b[i][0]/a[i][i];

}

float x0,y0;

x0=-x[1]/(2*x[2]);

y0=x[0]-x[1]*x[1]/(4*x[2]);

printf("\nx0=%f",x0);

printf("\ny0=%f",y0);

}

int main()

{

float x[99],y[99],z[99],s[99][99],t[99][1];

int n,m,i,j,k;

printf("n(变量数)=");

scanf("%d",&n);

printf("m（拟合阶数）=");

scanf("%d",&m);

printf("x[%]=",n);

for(i=0;i

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

printf("y[%d]=",n);

for(i=0;i

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

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

{

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

{s[i][j]=0.0;

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

s[i][j]+=cf(i+j,y[k]);

}

t[i][0]=0.0;

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

t[i][0]+=z[j]*cf(i,y[j]);

}

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

s[i][m+1]=t[i][0];

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

{

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

printf("%f\t",s[i][j]);

printf("\n");

}

ColPivot(s,m+1,t,x);

printf("\n\n");

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

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

printf("\n\n");

}

第13题

复合梯形递推算法

#include

#include

/***************计算函数值**************/

double f(double x)//

{

using namespace std;

double s,e=exp(1);

s=sqrt((1-pow(e,-x)))/x;//预先输入的待积分函数

return(s);

}

/***************************************/

/***************计算2^k值***************/

double power(int K)//

{

using namespace std;

double P;

P=pow(2.0,double(K));

return P;

}

/***************************************/

/**************复化梯形公式的递推过程***************/

double DTGCH(double a,double b,double precision)//

{

using namespace std;

int i,k;

double fa,fb,t1,p,S,X,t;

fa=1;

fb=f(b);

k=1;

//h=b-a;

t1=(b-a)/2*(fa+fb);//公式T1=(b-a)/2[f(a)+f(b)]

p=precision+1;//误差值初始化

while(p>=precision)

{

double n;

S=0;

for(i=1;i<=power(k-1);i++)//计算公式中:∑（2^(k-1),i=1）f[a+(2i-1)(b-a)/2^k]

{

X=a+(2*i-1)*((b-a)/power(k));

S=S+f(X);

}

t=(t1)/2+(b-a)/power(k)*S;

p=fabs(t1-t);//误差值

n=power(k-1);

cout<<"步长2^(k-1)为"<

t1=t;

k++;

}

return(t);

}

/**************************************************/

int main()

{

using namespace std;

double result,a,b,precision;//result为计算出的积分值、a为积分下限、b为积分上限、precision为要求的精确度

cout<<"需要求解的积分式为f(x)=(1-e(-x))(1/2)/x"<

cout<<"请输入积分下线:a="<<'\t';

cin>>a;

cout<<"请输入积分上限:b="<<'\t';

cin>>b;

cout<<"请输入精确度"<<'\t';

cin>>precision;

result=DTGCH(a,b,precision);

cout<<"经过复化梯形公式的递推算法得方程结果为:"<

}

第17题

欧拉方法

#include

#include

double f(double x,double y)

{

double s;

s=2.0/3*(x*pow(y,-2));//修改这里改变常微分方程

return s;

}

main()

{

double x=0,y=1,xi=0,yi=1,h=0.1,e;//改变x，y，xi，yi，h来控制

int i,j;

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

{

yi=yi+h*f(x,yi);

xi=x;

x=x+h;

y=pow((1+pow(x,2)),1.0/3);//修改这里改变准确解

printf("xi=%lf,yi=%lf,y=%lf\n",xi,yi,y);

}

}

改进欧拉方法

#include

#include

double f(double x,double y)

{

double s;

s=2.0/3*(x*pow(y,-2));

return s;

}

main()

{

double x=0,y=1,xi=0,yi=1,h=0.1,T1,T2;

int i,j;

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

{

T1=yi+h*f(x,yi);

xi=x;

x=x+h;

T2=yi+h*f(x,T1);

yi=(T1+T2)/2;

y=pow((1+pow(x,2)),1.0/3);

printf("xi=%lf,yi=%lf,y=%lf\n",xi,yi,y);

}

}

经典R-K方法

#include

#include

double f(double x,double y)

{

double s;

s=2.0/3*(x*pow(y,-2));

return s;