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

西西

<>本人紧急求助，希望各位大虾们伸出你们友情之手，帮我一下。</P>
<>1。牛顿迭代法</P>
<>1。全主消元法</P>
<>3。改进平方根法</P>
<>4。牛顿向前，向后插值法</P>
+ M) i$ X* p: W+ H% m+ y0 {" B8 _6 t<>5。加权最小二乘法</P>
  w. V8 Y' Y& F<>6。龙贝格积分法</P>
5 G- g. N. n: S/ d: @% U<>拜托了。。。。</P>

loveyuan

<>1。牛顿迭代法</P><>
1 K3 I4 Y$ y1 K0 y  c! u0 u- u8 j  #include "stdio.h"
  #include "math.h"
  int dnewt(x,eps,js)
  int js;
. o: q4 J$ v$ M( b8 X  double *x,eps;
. `' ]! ~0 B9 V. `" `* V6 t& S- l  { 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;
( b- n" J+ i4 H9 Y! x5 P$ U    while ((d&gt;=eps)&amp;&amp;(l!=0))
/ T$ _* P( w% B$ W! ?      { if (fabs(y[1])+1.0==1.0)
5 O- J: F1 j7 @5 U          { printf("err\n"); return(-1);}
2 V  P. V( K- C. F* u        x1=x0-y[0]/y[1];
8 M' ^7 t& |: Y' k8 \        dnewtf(x1,y);
        d=fabs(x1-x0); p=fabs(y[0]);
' w' H2 f7 E+ }        if (p&gt;d) d=p;
4 d; ~$ B  y* Q& P# @+ u        x0=x1; l=l-1;
% j7 A4 v/ u9 q! K/ i( z" j      }
    *x=x1;
    k=js-l;
    return(k);
" t9 q: a* F' l  }</P><>全主消元法</P><>#include "stdlib.h"
  #include "stdio.h"
  int acgas(ar,ai,n,br,bi)
! L' e; i, V4 D' q  int n;
  double ar[],ai[],br[],bi[];
  { int *js,l,k,i,j,is,u,v;
1 X  K5 ^4 \! x! K  x9 j  g    double p,q,s,d;
( X; k& L2 H. @% c    js=malloc(n*sizeof(int));
    for (k=0;k&lt;=n-2;k++)
1 O: X/ F% u7 b; `2 Q      { d=0.0;
7 P  V9 Y' }4 q$ e9 m        for (i=k;i&lt;=n-1;i++)
        for (j=k;j&lt;=n-1;j++)
  a& ?/ O  f3 A6 X# e8 c" @0 r          { u=i*n+j;
            p=ar*ar+ai*ai;
            if (p&gt;d) {d=p;js[k]=j;is=i;}
          }
% v" I2 x: c- t) `7 Y$ T        if (d+1.0==1.0)
          { free(js); printf("err**fail\n");
1 e% l  H9 x) {' t7 {            return(0);
# }5 F8 k8 z- M+ N% ^% Y9 k          }
        if (is!=k)
          { for (j=k;j&lt;=n-1;j++)
7 |: W0 G7 ]8 b+ F  n& L0 m              { u=k*n+j; v=is*n+j;
4 g) Q" e. x  R8 J) F- a$ g                p=ar; ar=ar[v]; ar[v]=p;
- r- [5 B0 ]7 k: N/ C1 q                p=ai; ai=ai[v]; ai[v]=p;
0 k( y9 t' l, P4 L' t; t              }
* z3 m8 t0 E" t/ S! `) O" K            p=br[k]; br[k]=br[is]; br[is]=p;
            p=bi[k]; bi[k]=bi[is]; bi[is]=p;
          }
# R0 e7 m3 o% m- Y        if (js[k]!=k)
          for (i=0;i&lt;=n-1;i++)
            { u=i*n+k; v=i*n+js[k];
              p=ar; ar=ar[v]; ar[v]=p;
( H0 |/ i' X1 L4 C) ^5 Z' R              p=ai; ai=ai[v]; ai[v]=p;
            }
        v=k*n+k;
        for (j=k+1;j&lt;=n-1;j++)
. }7 [. I* c% ^          { u=k*n+j;
, o3 P, f3 Z" c) d) A  B" s            p=ar*ar[v]; q=-ai*ai[v];
            s=(ar[v]-ai[v])*(ar+ai);
% m) q5 W" n, u: Q4 }& o+ ~            ar=(p-q)/d; ai=(s-p-q)/d;
          }
- j; L# @. x& g2 L* x! |        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&lt;=n-1;i++)
          { u=i*n+k;
! p" l% ~3 z) P& s0 _            for (j=k+1;j&lt;=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]);
2 x1 ]( y2 j% ]. U                ar[l]=ar[l]-p+q;
                ai[l]=ai[l]-s+p+q;
6 t  W) R- O4 Q% z/ X  G" w- J              }
$ L: w" g; u4 I; k7 V            p=ar*br[k]; q=ai*bi[k];
. w5 R  l7 v0 i( J( h: n1 d: ~: p            s=(ar+ai)*(br[k]+bi[k]);
/ i( \5 }7 j& b, C6 S* D            br=br-p+q; bi=bi-s+p+q;
          }
( z/ t/ t- q8 k3 N" V( v      }
' M$ @, k% A! }* @    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);
  I  @/ T* h1 J4 W  f      }
    p=ar*br[n-1]; q=-ai*bi[n-1];
    s=(ar-ai)*(br[n-1]+bi[n-1]);
    br[n-1]=(p-q)/d; bi[n-1]=(s-p-q)/d;
2 \# P1 \  g0 f  Y5 H- ~    for (i=n-2;i&gt;=0;i--)
    for (j=i+1;j&lt;=n-1;j++)
' r) }. b0 e) r      { u=i*n+j;
  {+ ~$ T5 j" e* m9 y/ O9 I        p=ar*br[j]; q=ai*bi[j];
: s( l$ a4 n7 H+ r, b        s=(ar+ai)*(br[j]+bi[j]);
        br=br-p+q;
1 w/ ^0 }# e7 _! U. e8 W/ w        bi=bi-s+p+q;
      }
4 H$ g, ^4 [, `- h( J3 T    js[n-1]=n-1;
2 L6 O, H' e( H1 U    for (k=n-1;k&gt;=0;k--)
# E" E0 }+ o6 S      if (js[k]!=k)
        { p=br[k]; br[k]=br[js[k]]; br[js[k]]=p;
/ w  b( l4 y/ M1 V5 \3 k          p=bi[k]; bi[k]=bi[js[k]]; bi[js[k]]=p;
+ v1 d) @; B; y# `4 q        }
1 }6 b% J& F3 \& S1 U5 ^+ n    free(js);
5 T' V9 |/ s# U9 f0 b    return(1);
# @# k( a4 s$ g# @- J+ N/ e6 E# h  }</P><>平方根法</P><>#include "math.h"
4 f3 @# @$ C5 m& y( }" X  #include "stdio.h"
  int achol(a,n,m,d)
# \5 j0 F8 A, Q! P  int n,m;
5 d( X! r/ h0 j# z  double a[],d[];
  { int i,j,k,u,v;
    if ((a[0]+1.0==1.0)||(a[0]&lt;0.0))
      { printf("fail\n"); return(-2);}
    a[0]=sqrt(a[0]);
    for (j=1; j&lt;=n-1; j++) a[j]=a[j]/a[0];
4 ^; B2 {- K/ ]3 t# H9 k; [; n    for (i=1; i&lt;=n-1; i++)
      { u=i*n+i;
: I  r) n' b! ?2 h        for (j=1; j&lt;=i; j++)
( t4 u7 I3 K! ^$ ?          { v=(j-1)*n+i;
            a=a-a[v]*a[v];
" ^" d" G* |" S8 H' W! x          }
( Y: X* J, k% O# M, K9 _# b        if ((a+1.0==1.0)||(a&lt;0.0))
          { printf("fail\n"); return(-2);}
        a=sqrt(a);
        if (i!=(n-1))
! h5 F/ R- d2 t          { for (j=i+1; j&lt;=n-1; j++)
# N4 A+ P1 u2 e- g# t; n              { v=i*n+j;
                for (k=1; k&lt;=i; k++)
                  a[v]=a[v]-a[(k-1)*n+i]*a[(k-1)*n+j];
                a[v]=a[v]/a;
              }
, R0 y. f7 j2 y$ Z6 z          }
& d+ x) j# f, D, _& K+ @      }
    for (j=0; j&lt;=m-1; j++)
      { d[j]=d[j]/a[0];
        for (i=1; i&lt;=n-1; i++)
          { u=i*n+i; v=i*m+j;
            for (k=1; k&lt;=i; k++)
, _% X# i% N$ E0 d% b. J" I              d[v]=d[v]-a[(k-1)*n+i]*d[(k-1)*m+j];
" M3 |9 g; L9 `            d[v]=d[v]/a;
: H5 O/ x. z( Y& a! x4 w2 w          }
      }
2 M# Z) a: Z' M. l1 A    for (j=0; j&lt;=m-1; j++)
      { u=(n-1)*m+j;
9 e! K, F) f$ i! L* p: L6 p        d=d/a[n*n-1];
, N& K- }2 c' y8 \  u2 o        for (k=n-1; k&gt;=1; k--)
$ P* u) }1 U: ^7 s          { u=(k-1)*m+j;
            for (i=k; i&lt;=n-1; i++)
* P0 v6 s( B3 t$ ~              { v=(k-1)*n+i;
                d=d-a[v]*d[i*m+j];
              }
5 P3 x# G4 d! ~/ S- j            v=(k-1)*n+k-1;
            d=d/a[v];
          }
' Q5 Y$ b' ]) B  T3 O      }
    return(2);
: q9 D6 ~- b/ y  }</P><>牛顿向前，向后插值法</P><>double eelgr(x0,h,n,y,t)
; s+ `# j& p1 j& v  int n;
  double x0,h,t,y[];
  { int i,j,k,m;
$ [9 ~9 |. F6 s# y8 v: V    double z,s,xi,xj;
    float p,q;
' B' t* z2 J% Q    z=0.0;
    if (n&lt;1) return(z);
    if (n==1) { z=y[0]; return(z);}
% {2 G# q3 U! A; J9 ^5 B    if (n==2)
- b! {* B4 D% f; q' F3 d4 F0 c# J      { z=(y[1]*(t-x0)-y[0]*(t-x0-h))/h;
2 \' n) p8 ]- f9 K" G        return(z);
  L$ k2 [6 i5 ?( }" w$ g) u0 o5 X4 g      }
, I" C3 s: X7 \% z4 z5 g    if (t&gt;x0)
2 _. w/ W1 b8 E( A; X9 Z5 W      { p=(t-x0)/h; i=(int)p; q=(float)i;
        if (p&gt;q) i=i+1;
      }
    else i=0;
    k=i-4;
    if (k&lt;0) k=0;
& A5 ^2 J* m: X# Z# B    m=i+3;
    if (m&gt;n-1) m=n-1;
    for (i=k;i&lt;=m;i++)
9 Q; L" m& u# c( Q9 L1 [2 }' U      { s=1.0; xi=x0+i*h;
4 ?3 H3 n2 ~# k3 Z        for (j=k; j&lt;=m; j++)
          if (j!=i)
            { xj=x0+j*h;
- v! C# F$ {# ^% V2 P0 X: S              s=s*(t-xj)/(xi-xj);
            }
        z=z+s*y;
6 L1 T& i: b4 W/ }% u8 d      }
9 u& f2 _/ K% {    return(z);
  }
向前，向后是一样的思想！！</P><>加权最小二乘法</P><>#include "math.h"
  T8 k4 p; Y& I* {  void hpir1(x,y,n,a,m,dt)
0 m. y3 r/ S7 I  T6 W* s" X  int n,m;
# M- V& T9 [2 L- `6 F  double x[],y[],a[],dt[];
& N3 e& H: k3 P$ j% J9 A9 H  { int i,j,k;
    double z,p,c,g,q,d1,d2,s[20],t[20],b[20];
$ s3 r1 N5 b/ D    for (i=0; i&lt;=m-1; i++) a=0.0;
- b' [9 i# c" B& U% ~9 x    if (m&gt;n) m=n;
    if (m&gt;20) m=20;
# O0 U1 n; m; c7 u7 X7 o' y8 J8 C* X+ m* g    z=0.0;
    for (i=0; i&lt;=n-1; i++) z=z+x/(1.0*n);
    b[0]=1.0; d1=1.0*n; p=0.0; c=0.0;
8 g; V6 z- y  f6 w/ h+ E( z/ G    for (i=0; i&lt;=n-1; i++)
+ Q/ m6 }5 h/ }" u& B2 c7 g( p      { p=p+(x-z); c=c+y;}
    c=c/d1; p=p/d1;
  `7 U- q% q2 C    a[0]=c*b[0];
    if (m&gt;1)
      { t[1]=1.0; t[0]=-p;
        d2=0.0; c=0.0; g=0.0;
        for (i=0; i&lt;=n-1; i++)
0 s# c, C& C/ U! K" V; m) |1 G          { q=x-z-p; d2=d2+q*q;
            c=c+y*q;
, g+ L# L( V' H0 `* p+ f; N  ?3 U5 B            g=g+(x-z)*q*q;
          }
! C8 p* U" p  L- u- X% c+ j        c=c/d2; p=g/d2; q=d2/d1;
        d1=d2;
        a[1]=c*t[1]; a[0]=c*t[0]+a[0];
      }
    for (j=2; j&lt;=m-1; j++)
      { s[j]=t[j-1];
! p. C! R% Q' i        s[j-1]=-p*t[j-1]+t[j-2];
2 ]2 F! ?: V" I        if (j&gt;=3)
          for (k=j-2; k&gt;=1; k--)
            s[k]=-p*t[k]+t[k-1]-q*b[k];
) z* W9 v1 _6 q1 w7 i        s[0]=-p*t[0]-q*b[0];
        d2=0.0; c=0.0; g=0.0;
- P. ^0 }7 g/ B5 c+ s. y2 s6 c        for (i=0; i&lt;=n-1; i++)
9 ]6 p8 d! X* x/ @3 q! \0 ]          { q=s[j];
( C* r3 r$ d) [            for (k=j-1; k&gt;=0; k--)
. ^/ C+ P3 K4 E1 V1 ~              q=q*(x-z)+s[k];
            d2=d2+q*q; c=c+y*q;
            g=g+(x-z)*q*q;
          }
* e7 p: ^: S1 k- W7 K+ y        c=c/d2; p=g/d2; q=d2/d1;
* U' @; u) d4 r        d1=d2;
) q# N* P9 a" q1 X' m) {$ q, r& ?        a[j]=c*s[j]; t[j]=s[j];
        for (k=j-1; k&gt;=0; k--)
* b/ M) Z9 r) s% \2 ?- ~          { a[k]=c*s[k]+a[k];
5 i: W) s3 n/ _" K$ o+ m. j2 [2 [3 g            b[k]=t[k]; t[k]=s[k];
          }
. C+ h0 K0 \# L  h      }
    dt[0]=0.0; dt[1]=0.0; dt[2]=0.0;
: {" m7 P& H/ X5 A+ B    for (i=0; i&lt;=n-1; i++)
      { q=a[m-1];
        for (k=m-2; k&gt;=0; k--)
          q=a[k]+q*(x-z);
. S6 t3 Z$ I) i$ O( n$ L* u: v        p=q-y;
, F4 M' R/ k: @/ r) a        if (fabs(p)&gt;dt[2]) dt[2]=fabs(p);
        dt[0]=dt[0]+p*p;
        dt[1]=dt[1]+fabs(p);
      }
    return;
' O# Q/ H; N1 a, ~; z; P! m+ t$ H6 E! S  }</P><>龙贝格积分法</P><>#include "math.h"
  double fromb(a,b,eps)
8 O( ?/ p5 r* |& y9 K- h  double a,b,eps;
  { extern double frombf();
" V# g9 r! t& x; h& ~0 U    int m,n,i,k;
* k% p+ {0 i6 L$ |' m8 Y    double y[10],h,ep,p,x,s,q;
    h=b-a;
3 X5 L! E- ?# y+ [4 w    y[0]=h*(frombf(a)+frombf(b))/2.0;
0 m: [% y2 F& E$ \: J    m=1; n=1; ep=eps+1.0;
    while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
  a) f$ l; _% T1 y! ]) [/ Z( Z, `3 z      { p=0.0;
7 @2 y* m2 U# F( e: f        for (i=0;i&lt;=n-1;i++)
          { x=a+(i+0.5)*h;
            p=p+frombf(x);
, K+ ~7 N( F- ?3 A          }
        p=(y[0]+h*p)/2.0;
  c% N$ ~5 t; m3 J8 l        s=1.0;
        for (k=1;k&lt;=m;k++)
          { s=4.0*s;
1 i) r/ g+ q# n$ }2 x            q=(s*p-y[k-1])/(s-1.0);
# x) i' R+ ?* K4 r0 ]/ b4 d& h            y[k-1]=p; p=q;
  `! k9 Z, b+ p- i% w8 X# Y. G          }
        ep=fabs(q-y[m-1]);
        m=m+1; y[m-1]=q; n=n+n; h=h/2.0;
      }
3 w+ F! o' S( D4 t/ D& w    return(q);
  }</P><>呵呵 希望对你有用！！</P>

西西

<>谢谢你哦，呵呵</P>

陈了

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