计算方法程序:
1、 秦九韶算法 2、二分法
3、拉格朗日插值 4、埃特金算法 5、复化梯形法
6、复化辛甫生算法 7、二阶龙格库塔方法
8、四阶龙格库塔方法 9、改进的欧拉方法 10、迭代法 11、埃特金加速方法:12、牛顿迭代法 13、追赶法 14、雅克比迭代 15、蛋白质设计:17高斯消去法:
18
格子法
实变函数程序
[10]可测函数
1.证明A=B;2.证明C=AB;3证明B为A的真子集;4 证明 C= A∩B,C= A∪B
5 证明A的余集为C;6 证明上下限集;7
证明B=S-A;8 证明D=A+B;9康托图
计算方法程序
1、秦九韶算法
P11.3利用秦九韶算法求多项式,在x=3时的值。
程序:
#include
#include
void main()
{float a[100],v,x;
int n,i,k;
scanf("%d%f",&n,&x);
for(i=0;i<=n;i++)
scanf("%f",&a[i]);
v=a[n];
k=1;
do
{v=x*v+a[n-k];
k=k+1;
}
while(k<=n);
printf("v=%f",v);}
运行结果:
2、二分法
P11.1用二分法求方程法x*x*x-x-1=0在[1,2]内的近似根,要求误差不超过
#include
#include
float fun(float);
void main()
{float a,b,c,x,y,y1;
scanf("%f%f%f",&a,&b,&c);
y1=fun(a);
do
{x=(a+b)/2;
y=fun(x);
{if(y*y1>0)
a=x;
else
b=x;}}
while((b-a)>=c);
printf("%f,%fn",x,y);
}
float fun(float m)
{float n;
n=m*m*m-m-1;
return(n);
}
运行结果:
3、拉格朗日插值
程序:
#include
main()
{float a,b,t,x[100],y[100];
int n,i,j,k;
scanf("%f%d",&a,&n);
for(i=0;i<=n;i++)
scanf("%f%f",&x[i],&y[i]);
k=0;
b=0;
for(k=0;k<=n;k++)
{t=1;
for(j=0;j<=n;j++)
{if(j!=k)
t=t*(a-x[j])/(x[k]-x[j]);}
b=b+t*y[k];
}
printf("%fn",b);
}
4、埃特金算法
程序:
#include
#include
main()
{float a,b,c,x[100],y[100];
int i,j,n,k;
scanf("%d%f",&n,&a);
for(i=0;i<=n;i++)
scanf("%f%f",&x[i],&y[i]);
for(k=1;k<=n;k++)
{for(i=k;i<=n;i++)
y[i]=y[k-1]+(y[i]-y[k-1])*(a-x[k-1])/(x[i]-x[k-1]);}
printf("%fn",y[n]);
}
5、复化梯形法
P95.9 设 ,用复化梯形法求积分 的近似值
程序:
#include
#include
double fun(double);
void main()
{
double a,b,h,s,x,y;
int n,k;
scanf("%lf%lf%d",&a,&b,&n);
h=(b-a)/n;
s=0;
x=a;
y=fun(x);
for(k=1;k<=n;k++)
{s=s+fun(x);
x=x+h;
s=s+fun(x);
}
s=(h/2)*s;
printf("s=%lfn",s);
}
double fun(double m)
{double n;
n=exp(-m)*sin(4*m)+1;
return(n);
}
运行结果:
6、复化辛甫生算法
P95.9 设 ,用复化辛甫生法求积分 的近似值
程序:
#include
#include
double fun(double);
void main()
{
double a,b,h,s,x,y;
int n,k;
scanf("%lf%lf%d",&a,&b,&n);
h=(b-a)/n;
s=0;
x=a;
y=fun(x);
for(k=1;k<=n;k++)
{s=s+fun(x);
x=x+h/2;
s=s+4*fun(x);
x=x+h/2;
s=s+fun(x);}
s=(h/6)*s;
printf("s=%lfn",s);
}
double fun(double m)
{double n;
n=exp(-m)*sin(4*m)+1;
return(n);
}
运行结果:
7、二阶龙格库塔方法
求解初值问题:
取h=0.2
程序:
#include
#include
float fun(float,float);
void main()
{float h,x0,y0,x1,y1,k1,k2;
int n,N;
scanf("%f%f%f%d",&x0,&y0,&h,&N);
n=1;
for(n=1;n<=N;n++)
{x1=x0+h;
k1=fun(x0,y0);
k2=fun(x0+h/2,y0+h/2*k1);
y1=y0+h*k2;
printf("%f,%fn",x1,y1);
x0=x1;y0=y1;}
}
float fun(float a,float b)
{float m;
m=b-2*a/b;
return(m);
}
运行结果:
8、四阶龙格库塔方法
求解初值问题:
取h=0.2
程序:
#include
#include
float fun(float,float);
void main()
{float h,x0,y0,x1,y1,k1,k2,k3,k4;
int n,N;
scanf("%f%f%f%d",&x0,&y0,&h,&N);
n=1;
for(n=1;n<=N;n++)
{x1=x0+h;
k1=fun(x0,y0);
k2=fun(x0+h/2,y0+h/2*k1);
k3=fun(x0+h/2,y0+h/2*k2);
k4=fun(x1,y0+h*k3);
y1=y0+h/6*(k1+2*k2+2*k3+k4);
printf("%f,%f",x1,y1);
x0=x1;y0=y1;}
}
{float m;
m=b-2*a/b;
return(m);
}
运行结果:
9、改进的欧拉方法
求解初值问题:
程序:
#include
#include
float fun(float,float);
void main()
{float x0,y0,h,x1,y1,yp,yc;
int n,N;
scanf("%f%f%f%d",&x0,&y0,&h,&N);
for(n=1;n<=N;n++)
{x1=x0+h;
yp=y0+h*fun(x0,y0);
yc=y0+h*fun(x1,yp);
y1=(yp+yc)/2;
printf("%f,%f",x1,y1);
x0=x1;
y0=y1;}
}
float fun(float a,float b)
{float m;
m=b-2*a/b;
return(m);
}
运行结果:
10、迭代法
P131 例2 用迭代法求方程 在 附近的一个根 ,要求精度为
程序:
#include
#include
float fun(float);
void main()
{float x0,x1,c;
int k,N;
scanf("%f%f%d",&x0,&c,&N);
for(k=1;k<=N;k++)
{x1=fun(x0);
printf("%fn",x1);
if(fabs(x1-x0)
break;
x0=x1;
}
if(k-1==N)
printf("Failure!n");
}
float fun(float m)
{float n;
n=exp(-m);
return(n);
}
运行结果:
11、埃特金加速方法:
程序:
#include
#include
float fun(float);
void main()
{float x0,x1,x2,c;
int k,N;
scanf("%f%f%d",&x0,&c,&N);
for(k=1;k<=N;k++)
{x1=fun(x0);
x2=fun(x1);
x2=x2-(x2-x1)*(x2-x1)/(x2-2*x1+x0);
if(fabs(x2-x0)
{printf("%fn",x2);
break;}
x0=x2;
}
if(k-1==N)
printf("Failure!n");
}
float fun(float m)
{float n;
n=exp(-m);
return(n);
}
运行结果:
12、牛顿迭代法:
例5、用牛顿法解方程
牛顿公式为: ,取x=0.5
程序:
#include
#include
float fun(float);
float ff(float);
void main()
{float x0,x1,c;
int k,N;
scanf("%f%f%d",&x0,&c,&N);
for(k=1;k<=N;k++)
{if(ff(x0)!=0)
{x1=x0-fun(x0)/ff(x0);
printf("%fn",x1);
if(fabs(x1-x0)
break;
x0=x1;
}
else
printf("****");
}
if(k==N)
printf("Failure!n");
}
float fun(float m)
{float n;
n=m-exp(-m);
return(n);
}
float ff(float m)
{float n;
n=1+m;
return(n);
}
运行结果:
13、追赶法:
P198.3.用追赶法求解下列方程组:
程序:
#include
#include
void main()
{float a[100],b[100],c[100],d[100],t;
int i,n;
scanf("%d",&n);
for(i=2;i<=n;i++)
scanf("%f",&a[i]);
for(i=1;i<=n;i++)
scanf("%f",&b[i]);
for(i=1;i<=n-1;i++)
scanf("%f",&c[i]);
for(i=1;i<=n;i++)
scanf("%f",&d[i]);
c[1]=c[1]/b[1];
d[1]=d[1]/b[1];
for(i=2;i
{t=b[i]-c[i-1]*a[i];
c[i]=c[i]/t;
d[i]=(d[i]-d[i-1]*a[i])/t;}
d[n]=(d[n]-d[n-1]*a[n])/(b[n]-c[n-1]*a[n]);
printf("%fn",d[n]);
for(i=n-1;i>=1;i--)
{d[i]=d[i]-c[i]*d[i+1];
printf("%fn",d[i]);}
}
运行结果:
14、雅克比迭代
程序:
#include
#include
#define N
50
#define M 4
void main()
{double x[M],y[M],a[M][M],b[M],d[M],c,t;
double ff(double [],int);
int k,i,j,n;
n=M-1;
scanf("%lf",&c);
for(i=0;i<=n;i++)
{scanf("%lf",&x[i]);
scanf("%lf",&b[i]);}
for(i=0;i<=n;i++)
for(j=0;j<=n;j++)
scanf("%lf",&a[0][i*M+j]);
for(k=1;k<=N;k++)
{for(i=1;i<=n;i++)
{for(j=1,t=0;j<=n;j++)
{if(j==i)
continue;
else
t=t+a[i][j]*x[j];}
y[i]=(b[i]-t)/a[i][i];
}
for(i=1;i<=n;i++)
d[i]=fabs(x[i]-y[i]);
t=ff(d,n+1);
if(t
else
for(i=1;i<=n;i++)
x[i]=y[i];
}
if(k==N)
printf("Failure!n");
if(k
{printf("k=%dn",k);
for(i=1;i<=n;i++)
printf("y[i]=%fn",y[i]);}
}
double ff(double a[],int n)
{double p;
int t;
p=a[1];
for(t=2;t
{if(p
p=a[t];}
return(p);
}
运行结果:
15、蛋白质设计:
程序:
#include
#include
#include
void main()
{printf("横坐标 纵坐标 竖坐标n");
FILE *fp;
char c[100],x[3000][7],y[3000][7],z[3000][7];
float a[3000],b[3000],d[3000],r[3000];
int i,k=0,j,n=0,m,p;
float X[3000],Y[3000],Z[3000],s1=0,s2=0,s3=0,av_x,av_y,av_z;
fp=fopen("1a1c.pdb","r");
for(i=0;i<=10000;i++)
{
fgets(c,81,fp);
if(c[0]=='A'&&c[1]=='T'&&c[2]=='O'&&c[3]=='M')
{
for(j=0;j<6;j++)
{x[k][j]=c[32+j];
printf("%c",x[k][j]);
}
printf(" ");
for(j=0;j<6;j++)
{y[k][j]=c[40+j];
printf("%c",y[k][j]);}
printf(" ");
for(j=0;j<6;j++)
{z[k][j]=c[48+j];
printf("%c",z[k][j]);
}
printf(" ");
X[k]=atof(x[k]);
s1+=X[k];
Y[k]=atof(y[k]);
s2+=Y[k];
Z[k]=atof(z[k]);
s3+=Z[k];
printf("n");
k++;
if(c[77]=='C') n++;
}
}
av_x=s1/k;
av_y=s2/k;
av_z=s3/k;
printf("共有:%dn",k+1);
printf("av_x=%f,av_y=%f,av_z=%fn",av_x,av_y,av_z);
printf("C原子个数为:%dn",n);
printf("平移原点后的坐标:n");
for(m=0;m
{a[m]=X[m]-av_x;printf("%f ",a[m]);
b[m]=Y[m]-av_y;printf("%f ",b[m]);
d[m]=Z[m]-av_z;printf("%f ",d[m]);
r[m]=sqrt(a[m]*a[m]+b[m]*b[m]+d[m]*d[m]);printf("%f ",r[m]);
printf("n");
}
}
17高斯消去法:
#include "stdio.h"
#include"math.h"
main()
{double a[3][3]={1,1,1,0,4,-1,2,-2,1},b[3]={6,5,1},x[10]={0};
int i,j,k,n=3;
for(k=0;k
{
for(i=k+1;i
{
for(j=k+1;j
{
a[i][j]=a[i][j]-a[k][j]*a[i][k]/a[k][k];
}
b[i]=b[i]-b[k]*a[i][k]/a[k][k];
}
}
x[n-1]=b[n-1]/a[n-1][n-1];
for(i=2;i<=n;i++)
{ k=n-i;
for(j=k+1;j
{
x[k]+=a[k][j]*x[j];
}
x[k]=(b[k]-x[k])/a[k][k];
}
for(k=0;k
printf("x[%d]=%f",k,x[k]);
}
18 格子法
求方程式组的极小解:x-y(x+y)-1=0; y(x+y)-y-2=0;
clear all;
x0=0;y0=0;h=0.8;ep=0.01;
f0=999999999;
for ii=1:999;f1=f0;
for
i=-1:1;for j=-1:1;
x=x0+i*h;y=y0+j*h;
f=(x-y*(x+y)-1)^2+(y*(x+y)-y-2)^2;
if(f
end;end;
h=h*2;
if(f1==f0);h=h/4;end;
if(f1
if(f0
end; % ii=1:
a(1)=x0;a(2)=y;a(3)=f0;a(4)=ii;a(5)=h
实变函数程序
[10]可测函数
s(y)=mE(f
h=0.02;
x=[-1:h:1]; y=[0:h:1];
for j=1:length(y);w=0;
for
i=1:length(x);f=x(i)*x(i);
if(f
end;
%i=1
s(j)=w;
end;%j=1
plot(y,s);
[1]//证明A=B
#include"stdio.h" main(){
double A[7]={4,5,6,7},B[5]={7,6,5,4};
int na,nb,i,j,k,nk,nc;
na=4;nb=4;
k=0;
for(i=0;i
for(j=0;j
if(A[j]==B[i]){k=k+1;}//计算A与B中相等元素个数
}
}
if(k==na)printf("A=B");
if(k !=na)printf("A is not Eque. B");
}
(2 )//证明C=AB
#include"stdio.h" main(){
double
A[7]={1,2,3,4,5,6,7},B[5]={5,6,7,8,9},C[3]={5,16,7},C1[5];
int na,nb,i,j,k,nk,nc;
//计算C1=AB;
na=7;nb=5;nc=3;
k=-1;
for(i=0;i
for(j=0;j
if(A[i]==B[j]){k=k+1;C1[k]=A[i];}//计算C1=AB
}
}
//判断C1=C
nk=k+1;
k=0;
for(i=0;i
for(j=0;j
if(C[j]==C1[i]){k=k+1;}//计算C1与C中相等元素个数
}
}
if(k==nc)printf("C=AB");
if(k !=nc)printf("C is not Eque. AB");
}
(3)证明B为A的真子集
#include
main()
{int A[5]={9,8,7,6,5},B[2]={6,7},i,j,t=0;
for(i=0;i<5;i++)
for(j=0;j<2;j++)
if(A[i]==B[j])
t++;
if(t==2)
printf("B为A的真子集n");
else
printf("B不是A的真子集n");
}
[4]
1、 已知A={4,2,3},B={5,3,6},证明:①A∩B={3};②A∪B={2,3,4,5,6}
#include
main()
{int A[3]={4,2,3},B[3]={5,3,6},S[5],i,j,k,t=0;
for(i=0;i<3;i++)
for(j=0;j<3;j++)
if(A[i]==B[j])
{t++;
k=A[i];
}
if(t==1)
{S[0]=k;
printf("A交B为%dn",S[0]);
j=1;
for(i=0;i<3;i++)
if(A[i]!=S[0])
{S[j]=A[i];
j++;
}
for(i=0;i<3;i++)
if(B[i]!=S[0])
{S[j]=B[i];
j++;
}
for(i=0;i<5;i++)
printf("%-3d",S[i]);
}
else
printf("无法证明题设n");
}
[5
] 若s={2,3,…,10},A={3,7,8,6},证明A的余集为{2,4,5,9,10}
#include
main()
{int s[9],A[4]={3,7,8,6},B[5]={2,4,5,9,10},i,j,k=0,t=0;
for(i=0;i<9;i++)
s[i]=i+2;
for(i=0;i<9;i++)
for(j=0;j<5;j++)
if(s[i]==B[j])
k++;
for(i=0;i<4;i++)
for(j=0;j<5;j++)
if(A[i]==B[j])
t++;
if(k==5&&t==0)
printf("A的余集为{2,4,5,9,10}n");
}
6 证明上下限集
#include
main()
{double x,an,w1,w2;
int i,j,n,m,N;
N=100;
n=N+1;
an=-1;
w1=an;
if((an-0)*(an-0)<0.01)
printf(“0属于An的上限集n”);
n=N+2;
an=1+1;
w2=an;
if((an-2)*(an-2)<0)
printf(“2属于An的上限集n”);
if((w1-w2)*(w1-w2)>0..1)
printf(“An的下限集为空集n”);
}
[7] 证明B=S-A
#include
#include
main()
{int S[6]={1,3,5,7,9,11},A[3]={1,3,6},B[100];
int i,j,k=0,m=0;
for(j=0;j<3;j++){
for(i=0;i<6;i++){
if(S[i] ==A[j]) {S[i]=-1;}
}
}
for(i=0;i<6;i++){if(S[i]>0){B[m]=S[i];m=m+1;
printf("%d %dn",m-1,B[m-1]);
}}
printf("S-A=B={%d,%d,%d,%d}n",B[0],B[1],B[2],B[3]);
}
[8] 证明D=A+B
#include
main()
{int
i,j,k=0,t=0,a[3]={1,2,3},b[4]={4,5,6,3},d[6]={1,2,3,4,5,6},h[10];
int g[9],p[9];
for(i=0;i<3;i++)h[i]=a[i];
for(i=0;i<4;i++)h[3+i]=b[i];
for(i=0;i<7;i++)g[i]=1;
for(j=0;j<6;j++)p[j]=1;
for(i=0;i<7;i++)
for(j=0;j<6;j++)
if(h[i]==d[j]){g[i]=0;p[j]=0;
}
k=0;
for(j=0;j<6;j++)k=k+p[j];
for(i=0;i<7;i++)k=k+g[i];
if(k==0)printf("C=D,%d",k);
if(k>0)printf("wrong,%d",k);
}
#include
main()
{int
i,j,k=0,t=0,a[3]={1,2,3},b[4]={4,5,6,3},d[6]={1,2,3,4,5,6},h[10];
for(i=0;i<3;i++)
h[i]=a[i];
for(i=0;i<4;i++)
h[3+i]=b[i];
for(i=0;i<7;i++)
for(j=0;j<6;j++)
if(h[i]==d[j])k++;
for(j=0;j<6;j++)
for(i=0;i<7;i++)
if(h[i]==d[j])
{t++;break;}
if(k==7&&t==6)
printf("C=D");
else printf("wrong");
}
9 康托图
clear all; hold off;
h(1)=1;n(1)=1;a(1,1)=0;
m=6; % the lay number
kc=2; %分叉数
r=0.2; %删除比例
for k=2:m;h(k)=h(k-1)/(kc+r*(kc-1));n(k)=n(k-1)*kc;end;
% Basic figure
z=[0.5 1
0.5 0
0 0
1 0]; plot(z(:,1),z(:,2)); axis([-0.1 1.1 0.7-m
1.1]);
% right end point a
for k=2:m;for i=1:n(k-1);for
j1=1:kc;j=kc-j1;a(k,kc*i-j)=a(k-1,i)+(1+r)*(j1-1)*h(k); q(:,1)=z(:,1)*h(k);q(:,2)=z(:,2);%伸缩变换
wx=a(k,kc*i-j)-q(3,1);wy=1-k-q(3,2);
vx=q(:,1)+wx;vy=q(:,2)+wy;%坐标平移
hold on; plot(vx,vy);
end;end;end; %k=2
二维康托图
clear all; hold off;
h(1)=1;n(1)=1; a(1,1)=0; h1(1)=1;
m=3; % 层数
kc=4; %分叉数
r=0.3; %删除比例
for k=2:m;h(k)=h(k-1)/(kc*(1+r)-r); n(k)=n(k-1)*kc;end;
% right end point a
for k=2:m; for i=1:n(k-1);for
j1=1:kc;j=kc-j1;a(k,kc*i-j)=a(k-1,i)+(j1-1)*h(k)*(1+r);end;end;end;
%k=2
% Basica Figure
z=[0 0
0 1
1 1
1 0
0 0
];
plot(z(:,1),z(:,2));axis([-0.1 1.1 -0.1
1.1]);
for k=2:m;
for i=1:n(k);for j=1:n(k);q=z*h(k);%伸缩变换
wx=a(k,i)-q(1,1);wy=a(k,j)-q(1,2);
vx=q(:,1)+wx;vy=q(:,2)+wy;%坐标平移
hold
on; plot(vx,vy);
end;end;
%i=,j=
end; %k=