求助!! 数值分析程序(c语言)
<P>本人紧急求助,希望各位大虾们伸出你们友情之手,帮我一下。</P><P>1。牛顿迭代法</P>
<P>1。全主消元法</P>
<P>3。改进平方根法</P>
<P>4。牛顿向前,向后插值法</P>
<P>5。加权最小二乘法</P>
<P>6。龙贝格积分法</P>
<P>拜托了。。。。</P> <P>1。牛顿迭代法</P><P>
#include "stdio.h"
#include "math.h"
int dnewt(x,eps,js)
int js;
double *x,eps;
{ extern void dnewtf();
int k,l;
double y,d,p,x0,x1;
l=js; k=1; x0=*x;
dnewtf(x0,y);
d=eps+1.0;
while ((d>=eps)&&(l!=0))
{ if (fabs(y)+1.0==1.0)
{ printf("err\n"); return(-1);}
x1=x0-y/y;
dnewtf(x1,y);
d=fabs(x1-x0); p=fabs(y);
if (p>d) d=p;
x0=x1; l=l-1;
}
*x=x1;
k=js-l;
return(k);
}</P><P>全主消元法</P><P>#include "stdlib.h"
#include "stdio.h"
int acgas(ar,ai,n,br,bi)
int n;
double ar[],ai[],br[],bi[];
{ int *js,l,k,i,j,is,u,v;
double p,q,s,d;
js=malloc(n*sizeof(int));
for (k=0;k<=n-2;k++)
{ d=0.0;
for (i=k;i<=n-1;i++)
for (j=k;j<=n-1;j++)
{ u=i*n+j;
p=ar*ar+ai*ai;
if (p>d) {d=p;js=j;is=i;}
}
if (d+1.0==1.0)
{ free(js); printf("err**fail\n");
return(0);
}
if (is!=k)
{ for (j=k;j<=n-1;j++)
{ u=k*n+j; v=is*n+j;
p=ar; ar=ar; ar=p;
p=ai; ai=ai; ai=p;
}
p=br; br=br; br=p;
p=bi; bi=bi; bi=p;
}
if (js!=k)
for (i=0;i<=n-1;i++)
{ u=i*n+k; v=i*n+js;
p=ar; ar=ar; ar=p;
p=ai; ai=ai; ai=p;
}
v=k*n+k;
for (j=k+1;j<=n-1;j++)
{ u=k*n+j;
p=ar*ar; q=-ai*ai;
s=(ar-ai)*(ar+ai);
ar=(p-q)/d; ai=(s-p-q)/d;
}
p=br*ar; q=-bi*ai;
s=(ar-ai)*(br+bi);
br=(p-q)/d; bi=(s-p-q)/d;
for (i=k+1;i<=n-1;i++)
{ u=i*n+k;
for (j=k+1;j<=n-1;j++)
{ v=k*n+j; l=i*n+j;
p=ar*ar; q=ai*ai;
s=(ar+ai)*(ar+ai);
ar=ar-p+q;
ai=ai-s+p+q;
}
p=ar*br; q=ai*bi;
s=(ar+ai)*(br+bi);
br=br-p+q; bi=bi-s+p+q;
}
}
u=(n-1)*n+n-1;
d=ar*ar+ai*ai;
if (d+1.0==1.0)
{ free(js); printf("err**fail\n");
return(0);
}
p=ar*br; q=-ai*bi;
s=(ar-ai)*(br+bi);
br=(p-q)/d; bi=(s-p-q)/d;
for (i=n-2;i>=0;i--)
for (j=i+1;j<=n-1;j++)
{ u=i*n+j;
p=ar*br; q=ai*bi;
s=(ar+ai)*(br+bi);
br=br-p+q;
bi=bi-s+p+q;
}
js=n-1;
for (k=n-1;k>=0;k--)
if (js!=k)
{ p=br; br=br]; br]=p;
p=bi; bi=bi]; bi]=p;
}
free(js);
return(1);
}</P><P>平方根法</P><P>#include "math.h"
#include "stdio.h"
int achol(a,n,m,d)
int n,m;
double a[],d[];
{ int i,j,k,u,v;
if ((a+1.0==1.0)||(a<0.0))
{ printf("fail\n"); return(-2);}
a=sqrt(a);
for (j=1; j<=n-1; j++) a=a/a;
for (i=1; i<=n-1; i++)
{ u=i*n+i;
for (j=1; j<=i; j++)
{ v=(j-1)*n+i;
a=a-a*a;
}
if ((a+1.0==1.0)||(a<0.0))
{ printf("fail\n"); return(-2);}
a=sqrt(a);
if (i!=(n-1))
{ for (j=i+1; j<=n-1; j++)
{ v=i*n+j;
for (k=1; k<=i; k++)
a=a-a[(k-1)*n+i]*a[(k-1)*n+j];
a=a/a;
}
}
}
for (j=0; j<=m-1; j++)
{ d=d/a;
for (i=1; i<=n-1; i++)
{ u=i*n+i; v=i*m+j;
for (k=1; k<=i; k++)
d=d-a[(k-1)*n+i]*d[(k-1)*m+j];
d=d/a;
}
}
for (j=0; j<=m-1; j++)
{ u=(n-1)*m+j;
d=d/a;
for (k=n-1; k>=1; k--)
{ u=(k-1)*m+j;
for (i=k; i<=n-1; i++)
{ v=(k-1)*n+i;
d=d-a*d;
}
v=(k-1)*n+k-1;
d=d/a;
}
}
return(2);
}</P><P>牛顿向前,向后插值法</P><P>double eelgr(x0,h,n,y,t)
int n;
double x0,h,t,y[];
{ int i,j,k,m;
double z,s,xi,xj;
float p,q;
z=0.0;
if (n<1) return(z);
if (n==1) { z=y; return(z);}
if (n==2)
{ z=(y*(t-x0)-y*(t-x0-h))/h;
return(z);
}
if (t>x0)
{ p=(t-x0)/h; i=(int)p; q=(float)i;
if (p>q) i=i+1;
}
else i=0;
k=i-4;
if (k<0) k=0;
m=i+3;
if (m>n-1) m=n-1;
for (i=k;i<=m;i++)
{ s=1.0; xi=x0+i*h;
for (j=k; j<=m; j++)
if (j!=i)
{ xj=x0+j*h;
s=s*(t-xj)/(xi-xj);
}
z=z+s*y;
}
return(z);
}
向前,向后是一样的思想!!</P><P>加权最小二乘法</P><P>#include "math.h"
void hpir1(x,y,n,a,m,dt)
int n,m;
double x[],y[],a[],dt[];
{ int i,j,k;
double z,p,c,g,q,d1,d2,s,t,b;
for (i=0; i<=m-1; i++) a=0.0;
if (m>n) m=n;
if (m>20) m=20;
z=0.0;
for (i=0; i<=n-1; i++) z=z+x/(1.0*n);
b=1.0; d1=1.0*n; p=0.0; c=0.0;
for (i=0; i<=n-1; i++)
{ p=p+(x-z); c=c+y;}
c=c/d1; p=p/d1;
a=c*b;
if (m>1)
{ t=1.0; t=-p;
d2=0.0; c=0.0; g=0.0;
for (i=0; i<=n-1; i++)
{ q=x-z-p; d2=d2+q*q;
c=c+y*q;
g=g+(x-z)*q*q;
}
c=c/d2; p=g/d2; q=d2/d1;
d1=d2;
a=c*t; a=c*t+a;
}
for (j=2; j<=m-1; j++)
{ s=t;
s=-p*t+t;
if (j>=3)
for (k=j-2; k>=1; k--)
s=-p*t+t-q*b;
s=-p*t-q*b;
d2=0.0; c=0.0; g=0.0;
for (i=0; i<=n-1; i++)
{ q=s;
for (k=j-1; k>=0; k--)
q=q*(x-z)+s;
d2=d2+q*q; c=c+y*q;
g=g+(x-z)*q*q;
}
c=c/d2; p=g/d2; q=d2/d1;
d1=d2;
a=c*s; t=s;
for (k=j-1; k>=0; k--)
{ a=c*s+a;
b=t; t=s;
}
}
dt=0.0; dt=0.0; dt=0.0;
for (i=0; i<=n-1; i++)
{ q=a;
for (k=m-2; k>=0; k--)
q=a+q*(x-z);
p=q-y;
if (fabs(p)>dt) dt=fabs(p);
dt=dt+p*p;
dt=dt+fabs(p);
}
return;
}</P><P>龙贝格积分法</P><P>#include "math.h"
double fromb(a,b,eps)
double a,b,eps;
{ extern double frombf();
int m,n,i,k;
double y,h,ep,p,x,s,q;
h=b-a;
y=h*(frombf(a)+frombf(b))/2.0;
m=1; n=1; ep=eps+1.0;
while ((ep>=eps)&&(m<=9))
{ p=0.0;
for (i=0;i<=n-1;i++)
{ x=a+(i+0.5)*h;
p=p+frombf(x);
}
p=(y+h*p)/2.0;
s=1.0;
for (k=1;k<=m;k++)
{ s=4.0*s;
q=(s*p-y)/(s-1.0);
y=p; p=q;
}
ep=fabs(q-y);
m=m+1; y=q; n=n+n; h=h/2.0;
}
return(q);
}</P><P>呵呵 希望对你有用!!</P>
?
<P>谢谢你哦,呵呵</P> 哇哇哇哇哇哇。。。好复杂呀!看不懂
页:
[1]