辛甫生算法c语言,[转载]实变函数,计算方法程序,康托图

计算方法程序：

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=