查看: 6000|回复: 3

求助！！数值分析程序（c语言）

字体大小: 正常放大

西西

1 主题	0 听众	20 积分

升级 15.79%

该用户从未签到

电梯直达

1^#

发表于 2004-5-10 08:46 |只看该作者 |正序浏览

|招呼Ta 关注Ta

<

>本人紧急求助，希望各位大虾们伸出你们友情之手，帮我一下。
<

>1。牛顿迭代法
<

>1。全主消元法
<

>3。改进平方根法
<

>4。牛顿向前，向后插值法
<

>5。加权最小二乘法
<

>6。龙贝格积分法
<

>拜托了。。。。

zan

转播0 淘帖0 分享0 收藏0 支持0 反对0 微信

使用道具举报

陈了

15 主题	3 听众	553 积分

升级 84.33%

TA的每日心情

	开心 2014-11-16 19:37

签到天数: 16 天

[LV.4]偶尔看看III

自我介绍: 数学，没有不可能

4^#

发表于 2010-3-22 11:13 |只看该作者 |招呼Ta 关注Ta |邮箱已经成功绑定

哇哇哇哇哇哇。。。好复杂呀！看不懂

我也说一句

发表

使用道具举报

西西

1 主题	0 听众	20 积分

升级 15.79%

该用户从未签到

3^#

发表于 2004-5-14 10:54 |只看该作者 |招呼Ta 关注Ta

？

<

>谢谢你哦，呵呵

我也说一句

发表

使用道具举报

loveyuan

21 主题	2 听众	104 积分

该用户从未签到

2^#

发表于 2004-5-11 08:32 |只看该作者 |招呼Ta 关注Ta

<

>1。牛顿迭代法<

>
 #include "stdio.h"
 #include "math.h"
 int dnewt(x,eps,js)
 int js;
 double *x,eps;
 { extern void dnewtf();
int k,l;
double y[2],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])+1.0==1.0)
 { printf("err\n"); return(-1);}
 x1=x0-y[0]/y[1];
 dnewtf(x1,y);
 d=fabs(x1-x0); p=fabs(y[0]);
 if (p>d) d=p;
 x0=x1; l=l-1;
 }
*x=x1;
k=js-l;
return(k);
 }<

>全主消元法<

>#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[k]=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[v]; ar[v]=p;
 p=ai; ai=ai[v]; ai[v]=p;
 }
 p=br[k]; br[k]=br[is]; br[is]=p;
 p=bi[k]; bi[k]=bi[is]; bi[is]=p;
 }
 if (js[k]!=k)
 for (i=0;i<=n-1;i++)
 { u=i*n+k; v=i*n+js[k];
 p=ar; ar=ar[v]; ar[v]=p;
 p=ai; ai=ai[v]; ai[v]=p;
 }
 v=k*n+k;
 for (j=k+1;j<=n-1;j++)
 { u=k*n+j;
 p=ar*ar[v]; q=-ai*ai[v];
 s=(ar[v]-ai[v])*(ar+ai);
 ar=(p-q)/d; ai=(s-p-q)/d;
 }
 p=br[k]*ar[v]; q=-bi[k]*ai[v];
 s=(ar[v]-ai[v])*(br[k]+bi[k]);
 br[k]=(p-q)/d; bi[k]=(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[v]; q=ai*ai[v];
 s=(ar+ai)*(ar[v]+ai[v]);
 ar[l]=ar[l]-p+q;
 ai[l]=ai[l]-s+p+q;
 }
 p=ar*br[k]; q=ai*bi[k];
 s=(ar+ai)*(br[k]+bi[k]);
 br=br-p+q; bi=bi-s+p+q;' d* L6 ^ w1 Y1 D" e
 }' t7 L& ~7 {" l: v3 T
 }
$ o. U; h0 @& f0 M u=(n-1)*n+n-1;3 z# } i5 Q, z( f
d=ar*ar+ai*ai;' l5 }. P9 t8 G7 \' m0 d% J* q
if (d+1.0==1.0)
0 @( a: v8 E8 U; u1 Z { free(js); printf("err**fail\n");( |6 T& ]0 {( C7 H8 z* l2 y
 return(0);. u' P5 e, U" K1 s# {; Z! }
 }3 t5 n L. M! k9 m% R+ F& p
p=ar*br[n-1]; q=-ai*bi[n-1];
( [3 |% v- N0 I0 Q. b2 m s=(ar-ai)*(br[n-1]+bi[n-1]);
5 Z0 ]8 F: Y. B6 J$ X br[n-1]=(p-q)/d; bi[n-1]=(s-p-q)/d;* ]# U* a! z; ^
for (i=n-2;i>=0;i--)* F6 y) `) i" o) {+ d
for (j=i+1;j<=n-1;j++)
: _" Z- q1 W7 W1 t2 c { u=i*n+j;! X, v# ?! R" C. ~0 J
 p=ar*br[j]; q=ai*bi[j];$ ?5 s1 u- u, A# J- V. ~; |- x! m
 s=(ar+ai)*(br[j]+bi[j]);
9 E0 q q) d7 ~ br=br-p+q;# [2 j: m7 z" {* v! c6 q
 bi=bi-s+p+q;: t. x' W0 S4 R/ t. v0 ^) _
 }8 C. D+ s1 u4 q4 s6 K. @6 Y" \
js[n-1]=n-1;
+ ]' E' `+ m: z4 E5 ~, o for (k=n-1;k>=0;k--)
5 R, b* j0 R( B/ b: p if (js[k]!=k)+ l& B1 [' \0 H/ M, W$ n
 { p=br[k]; br[k]=br[js[k]]; br[js[k]]=p;; Y. K' t. L, S1 p
 p=bi[k]; bi[k]=bi[js[k]]; bi[js[k]]=p;+ Q" @- [% [$ W2 C
 }
4 J* `! O7 r$ ]8 l* n, }5 N free(js);$ j0 T! w& H; @, j4 m
return(1);; d, C' n5 U; A4 O
 }<>平方根法<>#include "math.h"+ _- |7 p5 i* `9 H0 `
 #include "stdio.h"
* E0 o5 ?2 w) }. C8 `' J! W int achol(a,n,m,d)
, c9 |7 R) m! ~1 R7 B2 ? int n,m;4 ?8 O5 x3 ~& s `3 d- n
 double a[],d[];, `5 k+ A" k8 ?2 ^) T
 { int i,j,k,u,v;
8 a2 h( q! ^# f- U2 L' `3 ^& b if ((a[0]+1.0==1.0)||(a[0]<0.0))0 g; {% m" }" u* z5 X- W' @ c
 { printf("fail\n"); return(-2);}
" s7 {) R: {- _6 f: F a[0]=sqrt(a[0]);
, H3 L7 `) o/ r% I for (j=1; j<=n-1; j++) a[j]=a[j]/a[0];4 N2 k4 [6 Y& @3 d
for (i=1; i<=n-1; i++)3 b1 d6 `& f4 t: K& H1 i& i& ^5 c
 { u=i*n+i;
- p- m4 T9 [6 N4 h1 \9 S; o for (j=1; j<=i; j++)8 g6 H, t5 e) [2 ]" a
 { v=(j-1)*n+i;
2 \/ M3 h6 z6 ?' P8 u3 P) Q a=a-a[v]*a[v];
) E4 g8 l& ]5 t# z. p- x# c9 b }
( Y- n9 H- Y" Y. A, B2 J) Z if ((a+1.0==1.0)||(a<0.0))
8 l: ~' ?+ p3 @- \- B { printf("fail\n"); return(-2);}4 @# n6 g( w1 A6 P2 B# q+ H4 B4 K
 a=sqrt(a);8 J1 q$ b8 j2 o1 f
 if (i!=(n-1))
& r( z( z; Z8 C, h' v+ ^ { for (j=i+1; j<=n-1; j++)
 [! H4 P7 Z! F! f+ U* f { v=i*n+j;( [9 P' M# H2 n6 O4 _0 d* E
 for (k=1; k<=i; k++)
! P: U; J i$ f& J, i/ Z$ a0 g: i a[v]=a[v]-a[(k-1)*n+i]*a[(k-1)*n+j];. h1 z3 s: @9 m5 e
 a[v]=a[v]/a;! i. ]0 U) Q; h- d) n. y/ @' G
 }" x! W6 M: w0 f) u
 }
' q. r/ n! {7 I3 T6 x x }
. q5 f1 u2 ^- s3 o for (j=0; j<=m-1; j++)& L8 f3 C. U; w: K" B) ^8 ~
 { d[j]=d[j]/a[0];
' k6 O6 \1 H7 z. y1 [, O( G for (i=1; i<=n-1; i++)
8 K$ d b" m! q/ z5 K; C { u=i*n+i; v=i*m+j;3 `- B' o+ P/ N \7 O
 for (k=1; k<=i; k++)" m R7 ~# y x0 ]$ P! r; Z2 G- V
 d[v]=d[v]-a[(k-1)*n+i]*d[(k-1)*m+j];% f) s4 r1 K5 Z7 H
 d[v]=d[v]/a;7 X0 {3 X* [! b, m
 }1 o: L- t/ a0 p) K" r4 w
 }
8 ^. O/ m: b% [& d/ A for (j=0; j<=m-1; j++)3 J \7 J% T' Y# N
 { u=(n-1)*m+j;0 L% a4 @* `5 v/ g" _5 b
 d=d/a[n*n-1];! g/ c. F: u; G/ W! F; x
 for (k=n-1; k>=1; k--)
/ y1 \" s2 r4 y8 B$ @, \, |9 C) O { u=(k-1)*m+j;3 s7 s! |% c; J
 for (i=k; i<=n-1; i++)1 z/ f7 S+ K4 t/ h" _
 { v=(k-1)*n+i;
- _* v6 F, C' b- Z' t6 L: m6 @ d=d-a[v]*d[i*m+j];
- D! S" P, Z2 `4 w( `1 z- |& G }
$ C- t# m/ [' Y7 T, C, K6 R( G v=(k-1)*n+k-1;
/ r! y0 R: S$ f d=d/a[v];
0 w; V9 Q( t- e& B. L }' t7 r6 s+ y8 G" C0 y9 t
 }5 v' C& j! a/ V" B
return(2);& p+ k' T: Q7 C7 D$ r
 }<>牛顿向前，向后插值法<>double eelgr(x0,h,n,y,t)
; w% ^& V" L7 n d- P int n;$ ]9 O0 k. D5 N: e. z
 double x0,h,t,y[];
6 r, C9 T9 F% _# a7 Y# P9 } { int i,j,k,m;
" Y# u* ^: j/ Y. E' F/ `6 b2 L double z,s,xi,xj;7 x! X8 \; X, Z( S2 ^
float p,q;* Y$ ^% _0 {! S U
z=0.0;. `) n; e( ~; C r5 x" T# N; t* C
if (n<1) return(z);$ n I7 F" P$ x3 M9 X, Z6 B2 s* A
if (n==1) { z=y[0]; return(z);}
" D9 }& W/ s5 q0 y) L1 x2 ? if (n==2)
( R. y E" d% B B+ }" M; B/ w9 O { z=(y[1]*(t-x0)-y[0]*(t-x0-h))/h;" v. X: P# R& ~6 d
 return(z);9 a0 ^& B7 o, J" R' |7 ^9 ~
 }
% }0 n- U: r' ^) N+ B if (t>x0)' E1 O2 f+ @; ^& n+ J. B% V' `0 n7 Y/ b
 { p=(t-x0)/h; i=(int)p; q=(float)i;
7 `0 ~; R2 h& C if (p>q) i=i+1;
0 f1 z' l8 ~( k y' O6 D1 l }1 |6 \( ^( V9 e0 f, o6 m
else i=0;
4 ^/ b8 E8 d* W k=i-4;
 S: |5 G( m2 Z0 x; Z% ? if (k<0) k=0;- r* K# R/ |3 ~- W! m
m=i+3;/ x: b/ p. a; @' _- q
if (m>n-1) m=n-1;
' L0 \$ u7 v0 Y5 g, \ for (i=k;i<=m;i++)
9 g6 m$ r% s+ r* a# S1 S { s=1.0; xi=x0+i*h; }' V' h" G) i8 ]! x. V
 for (j=k; j<=m; j++)7 u% B8 \' f- ~; ?. [, u
 if (j!=i)
3 E$ T( @/ W4 i4 u: U; c6 X7 n { xj=x0+j*h;/ o( H/ _& {2 E: a# o. Y
 s=s*(t-xj)/(xi-xj);! H9 K& s3 Q3 G/ v% ^
 }# o2 ]5 w' r: ], | A/ D7 E+ S
 z=z+s*y;
# R. H- S2 O- w/ {& K" E }
* Q( q( _8 B& S return(z);1 [; c3 \% P+ T2 M- \9 F, @
 }( ?6 E7 b6 ~( U+ g( A4 H
向前，向后是一样的思想！！<>加权最小二乘法<>#include "math.h"
6 t7 N2 ^1 I a* T5 ^ void hpir1(x,y,n,a,m,dt)
8 I6 D7 z$ N5 E: r. J% h5 d int n,m;6 {7 y/ ^9 Y- _' e
 double x[],y[],a[],dt[];
9 P: U5 y4 `& }9 w$ l" @' P { int i,j,k;( U/ R( e) u, o
double z,p,c,g,q,d1,d2,s[20],t[20],b[20];
/ {: T# d% i. @0 G for (i=0; i<=m-1; i++) a=0.0;
; A" Z3 k$ G5 p; t9 z6 @2 k0 W if (m>n) m=n;* A$ G. g( z+ c6 g) M
if (m>20) m=20;
, U& ?3 @1 Z$ d+ d% N z=0.0;9 ?/ E! @, d P
for (i=0; i<=n-1; i++) z=z+x/(1.0*n);1 ~8 s; e. C1 v% r# _. s
b[0]=1.0; d1=1.0*n; p=0.0; c=0.0;" K5 [/ a, ?! U( P6 {! O4 Q
for (i=0; i<=n-1; i++)/ F1 f: c* \0 e: p1 [* W
 { p=p+(x-z); c=c+y;}" d7 A" z) \: a4 U8 h
c=c/d1; p=p/d1;3 X5 o# N0 n. S6 p! @, i& ~! h
a[0]=c*b[0];1 s1 {% J3 p8 t; \% W
if (m>1)9 x; p/ |2 M" ^9 F3 ^8 h! H6 u
 { t[1]=1.0; t[0]=-p;
& v/ ]% W& U' X' U! q d2=0.0; c=0.0; g=0.0;
 x' r+ S2 ~$ {1 O3 s for (i=0; i<=n-1; i++), L e3 g$ L/ m( N; z& c/ K/ I
 { q=x-z-p; d2=d2+q*q;' N7 r3 C7 ^2 [! z9 [7 a
 c=c+y*q;5 P+ H7 W6 a8 D. G M+ Y
 g=g+(x-z)*q*q;
6 l2 m& U5 D: J0 k- S" ^7 s/ i }7 z8 H! Z7 `) r% p
 c=c/d2; p=g/d2; q=d2/d1;0 d; E( J& c0 l+ I3 a
 d1=d2;" G7 i4 s9 }0 g/ w8 ~% o# ~
 a[1]=c*t[1]; a[0]=c*t[0]+a[0];4 ?6 j2 a; [* W4 e5 |$ S
 } R2 }% p: U" A
for (j=2; j<=m-1; j++)9 X5 K, d! V5 W5 j4 |# m1 u
 { s[j]=t[j-1];! M1 ?6 K5 \; Z3 X4 H W
 s[j-1]=-p*t[j-1]+t[j-2];" f4 K" d1 B* @
 if (j>=3)
; `1 W& |/ O2 i5 h; p5 ]) u for (k=j-2; k>=1; k--)
: D/ o3 e9 Y( y. D7 L W. q s[k]=-p*t[k]+t[k-1]-q*b[k];
. r$ K% G9 V+ z1 R& m6 q9 f s[0]=-p*t[0]-q*b[0];) t# O1 n1 a3 u8 G
 d2=0.0; c=0.0; g=0.0;
1 Y: f; M+ Y m ^' Z for (i=0; i<=n-1; i++)
( P1 Q1 w( t7 l/ c$ O' {2 E" C { q=s[j];
% E# B$ J0 y4 S9 q6 ?* l for (k=j-1; k>=0; k--)9 i- n- r0 c. r- f0 _
 q=q*(x-z)+s[k];
 d2=d2+q*q; c=c+y*q; S% J* M( [7 g+ P& D4 p
 g=g+(x-z)*q*q;
+ o% A. b2 n( \2 H }8 ~( A8 P: l: ]# M! u$ ~, \; @
 c=c/d2; p=g/d2; q=d2/d1;
2 R0 j( |7 [5 ~! d3 z d1=d2;0 d+ U. H( `5 G, K3 @% }
 a[j]=c*s[j]; t[j]=s[j];
4 ~* N: d+ `. H* h for (k=j-1; k>=0; k--)' O# `$ \) O& \' B) }7 E
 { a[k]=c*s[k]+a[k];
: \) {; p6 H, e% Y1 _' @" L b[k]=t[k]; t[k]=s[k];: u8 ^ u; ?) E( g$ k9 G+ ^
 }" ~4 N" Z) f* U! k4 o- z
 }
: T! D+ d8 m E L+ G dt[0]=0.0; dt[1]=0.0; dt[2]=0.0;6 v2 R$ V2 c* U# ~
for (i=0; i<=n-1; i++), q5 }$ z6 L g% e& p; `( Q$ b* X
 { q=a[m-1];% T, ~) N+ u7 t
 for (k=m-2; k>=0; k--), e( P+ k& k3 b- P
 q=a[k]+q*(x-z);' \+ L8 w/ U+ A8 P* z1 k
 p=q-y;! s; @4 t9 y' s; A% d) u& ^. _
 if (fabs(p)>dt[2]) dt[2]=fabs(p);
- h" r, {1 y5 s) _ dt[0]=dt[0]+p*p;
- g# e3 Q6 |+ ? q% G dt[1]=dt[1]+fabs(p);
/ N* Q- b, |2 p" R" f, t }" `, P2 X9 h# N) x% _
return;
 E: T* V1 g, U% `" h }<>龙贝格积分法<>#include "math.h"' n- U/ p3 [3 S- x/ Q9 Q
 double fromb(a,b,eps)1 ]: U5 \+ Z' `3 t% J% ]- N7 k
 double a,b,eps;
. q, X; W8 e8 K* r { extern double frombf();
" p' M3 n- s) g$ i* L+ | int m,n,i,k;# E9 c+ {3 n& T
double y[10],h,ep,p,x,s,q;! R2 G( S9 _' ~* ?, m
h=b-a;4 T: f$ l) f% p$ N+ C
y[0]=h*(frombf(a)+frombf(b))/2.0;! K1 {2 F/ d' x1 b; y. ?
m=1; n=1; ep=eps+1.0;9 H+ D9 R2 w: Q9 h! W- z- C7 g
while ((ep>=eps)&&(m<=9))# l; Z+ C6 m/ e2 t
 { p=0.0;
' _/ }% @8 k7 j; c8 N for (i=0;i<=n-1;i++)
1 L1 a2 {5 f7 r4 d { x=a+(i+0.5)*h;5 ^- k Z& C/ |1 I @0 Q
 p=p+frombf(x);
6 H5 T, |* Y3 ^8 | }
4 g; c( F5 j1 E, ]' z4 e! ` p=(y[0]+h*p)/2.0;
( I& Y4 c% |, s1 C m s=1.0;2 h) U: J; ^3 y3 ]
 for (k=1;k<=m;k++)
/ N( j& T9 s a { s=4.0*s;5 V2 L3 N5 h) }( c, u7 g
 q=(s*p-y[k-1])/(s-1.0);
 y[k-1]=p; p=q;0 p8 p9 |2 V6 V$ B+ f* R( \
 }
. o# g+ y' ?5 i ep=fabs(q-y[m-1]);
6 N9 E8 ~8 v) B8 X% ]8 V( a m=m+1; y[m-1]=q; n=n+n; h=h/2.0;. @: w' Z4 [- F* ^
 }
) W* v N& ?" C# s( g4 _ return(q);! ]( v ]# |6 A9 b
 }<>呵呵希望对你有用！！

我也说一句

发表

May God grant me the wisdom to discover the right and strength to make it endure 因为陌生，所以勇敢，因为距离，所以美丽。

使用道具举报

返回列表

帐号		密码		只需要一步，快速开始		注册地址	找回密码

求助！！数值分析程序（c语言）

？

浏览过的版块

新人进步奖

元老勋章

QQ

电话咨询

关于我们| 联系我们| 诚征英才| 对外合作| 产品服务|

求助！！ 数值分析程序（c语言）

？

浏览过的版块

新人进步奖

元老勋章

QQ

电话咨询

关于我们| 联系我们| 诚征英才| 对外合作| 产品服务|

求助！！数值分析程序（c语言）