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

西西

<>本人紧急求助，希望各位大虾们伸出你们友情之手，帮我一下。</P>
  z$ y  w3 O+ |. q<>1。牛顿迭代法</P>
<>1。全主消元法</P>
<>3。改进平方根法</P>
<>4。牛顿向前，向后插值法</P>
<>5。加权最小二乘法</P>
( _$ z4 \% y1 Y: }& J( g<>6。龙贝格积分法</P>
<>拜托了。。。。</P>

loveyuan

<>1。牛顿迭代法</P><>
  #include "stdio.h"
; _: m9 ?3 R0 R7 q% N9 i  #include "math.h"
  int dnewt(x,eps,js)
0 V& O3 y$ \+ i2 G+ ^# c  int js;
  double *x,eps;
. e% `" G3 T+ J# r' Q5 b# S  { extern void dnewtf();
- h3 }. ]. N% w/ [( O    int k,l;
+ o9 M- q$ A0 {1 I    double y[2],d,p,x0,x1;
    l=js; k=1; x0=*x;
    dnewtf(x0,y);
# n. e2 X. x8 H6 M, d2 u    d=eps+1.0;
3 |5 w4 O' @2 N1 l# [% g6 r, b0 B    while ((d&gt;=eps)&amp;&amp;(l!=0))
      { if (fabs(y[1])+1.0==1.0)
          { printf("err\n"); return(-1);}
        x1=x0-y[0]/y[1];
* F* B' {7 u# I        dnewtf(x1,y);
        d=fabs(x1-x0); p=fabs(y[0]);
        if (p&gt;d) d=p;
        x0=x1; l=l-1;
* n! X, ?: a' z. }) y      }
    *x=x1;
    k=js-l;
3 P' I* A& `6 |9 o: K0 }. X' [' X9 V    return(k);
8 [2 g* ?2 L" ?8 a  }</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;
5 c* Z, P, v: \) @    double p,q,s,d;
    js=malloc(n*sizeof(int));
    for (k=0;k&lt;=n-2;k++)
      { d=0.0;
7 G# W& I- f8 _; {; q* t        for (i=k;i&lt;=n-1;i++)
/ D% s3 O7 j% f) z        for (j=k;j&lt;=n-1;j++)
          { u=i*n+j;
* i: Q) I8 o9 e2 l            p=ar*ar+ai*ai;
            if (p&gt;d) {d=p;js[k]=j;is=i;}
( ~% b2 f6 ~4 v- I" Z. B7 R          }
0 [6 |! i, N: L        if (d+1.0==1.0)
. f/ e* p) y* g) P( a4 Q/ a$ |          { free(js); printf("err**fail\n");
            return(0);
( z# A! B+ x+ k) l" u) |          }
        if (is!=k)
          { for (j=k;j&lt;=n-1;j++)
8 B5 {+ O4 p$ |  s/ a              { u=k*n+j; v=is*n+j;
4 d8 L4 b4 @9 b- a* q                p=ar; ar=ar[v]; ar[v]=p;
                p=ai; ai=ai[v]; ai[v]=p;
' w" _! R* k) X) @              }
            p=br[k]; br[k]=br[is]; br[is]=p;
            p=bi[k]; bi[k]=bi[is]; bi[is]=p;
          }
9 c9 E3 |0 n- J; a0 B! k9 E        if (js[k]!=k)
          for (i=0;i&lt;=n-1;i++)
; [* x) I' c" z' e            { u=i*n+k; v=i*n+js[k];
9 _% N& {4 D; \, U              p=ar; ar=ar[v]; ar[v]=p;
              p=ai; ai=ai[v]; ai[v]=p;
/ a  n, U- S+ L( o: V: R! f. }            }
        v=k*n+k;
0 X; _5 j; m7 R$ k; j        for (j=k+1;j&lt;=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;
          }
% t5 t8 H) B2 r5 S4 e7 ^        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;
- G1 G, @$ L# }( d7 d' F# z        for (i=k+1;i&lt;=n-1;i++)
# o; B) `& s4 r          { u=i*n+k;
            for (j=k+1;j&lt;=n-1;j++)
5 K& P6 J1 B% d. t/ O              { v=k*n+j; l=i*n+j;
6 g6 `: S: B. m+ k! _                p=ar*ar[v]; q=ai*ai[v];
                s=(ar+ai)*(ar[v]+ai[v]);
8 e1 K! A% j- U4 \9 v6 v                ar[l]=ar[l]-p+q;
: |" @) K8 x$ t3 |0 L+ }                ai[l]=ai[l]-s+p+q;
2 l5 S. o* X0 t1 O: Z! X              }
            p=ar*br[k]; q=ai*bi[k];
: k8 `% k' h3 Y6 U2 c: _/ U4 g            s=(ar+ai)*(br[k]+bi[k]);
            br=br-p+q; bi=bi-s+p+q;
          }
1 o6 O6 x6 ~- i      }
    u=(n-1)*n+n-1;
    d=ar*ar+ai*ai;
    if (d+1.0==1.0)
9 Y+ `# H. s  i5 `& O      { free(js); printf("err**fail\n");
$ W' m& j6 L- e1 h        return(0);
7 L" c' t$ I2 q9 h4 {1 ^      }
* k, e  m! Z% {    p=ar*br[n-1]; q=-ai*bi[n-1];
) G( W! y1 ?3 s8 S5 i    s=(ar-ai)*(br[n-1]+bi[n-1]);
" a, _. K. E6 D3 R    br[n-1]=(p-q)/d; bi[n-1]=(s-p-q)/d;
7 I9 C$ y* T" c  Z- E; }    for (i=n-2;i&gt;=0;i--)
    for (j=i+1;j&lt;=n-1;j++)
      { u=i*n+j;
        p=ar*br[j]; q=ai*bi[j];
/ S: S. {4 _, @; z        s=(ar+ai)*(br[j]+bi[j]);
# L& i$ l# A8 M  C1 g1 C. g8 t( d        br=br-p+q;
3 x5 l1 a( l. z/ }2 M( L# ^        bi=bi-s+p+q;
; t2 L% m3 o8 S! C% ]: T      }
    js[n-1]=n-1;
    for (k=n-1;k&gt;=0;k--)
* a* p& @1 e7 N6 r7 \      if (js[k]!=k)
        { p=br[k]; br[k]=br[js[k]]; br[js[k]]=p;
/ m& G$ L- v' Y% j) R- K2 A! H          p=bi[k]; bi[k]=bi[js[k]]; bi[js[k]]=p;
        }
    free(js);
, K7 Q% z0 i3 x0 O    return(1);
  }</P><>平方根法</P><>#include "math.h"
" T3 U- q& N8 z5 e# [  #include "stdio.h"
2 s& ?6 |0 |1 x6 n# C# ^. @  int achol(a,n,m,d)
/ h' L9 A7 F8 L# Z2 R  int n,m;
& h; ^3 P+ d$ U  double a[],d[];
  { int i,j,k,u,v;
    if ((a[0]+1.0==1.0)||(a[0]&lt;0.0))
, [3 c( o# V, @9 b/ t9 k, G      { printf("fail\n"); return(-2);}
    a[0]=sqrt(a[0]);
  X5 d" [3 |8 ^( X) A8 u) n' ^    for (j=1; j&lt;=n-1; j++) a[j]=a[j]/a[0];
) C. h% S) E' J: M! R2 t# h    for (i=1; i&lt;=n-1; i++)
7 H, y. S/ ^7 G5 y9 @      { u=i*n+i;
  {4 ]6 J4 Y9 D: O/ K4 F, d4 l        for (j=1; j&lt;=i; j++)
          { v=(j-1)*n+i;
, l% k0 V5 j% m2 h8 d$ Z0 {8 n  O            a=a-a[v]*a[v];
+ F2 B+ w- k8 |( r$ W/ V          }
& H6 `/ F6 i; J        if ((a+1.0==1.0)||(a&lt;0.0))
! t' Z3 k5 [; x" E* Q: q/ f          { printf("fail\n"); return(-2);}
        a=sqrt(a);
, n5 y6 K0 H8 M/ k        if (i!=(n-1))
          { for (j=i+1; j&lt;=n-1; j++)
              { 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;
# y; I; g) W! c  g5 h              }
. ]6 }% I1 L: t2 s8 e          }
      }
. i; i( S; ]. D, c- D- X8 i    for (j=0; j&lt;=m-1; j++)
* h$ \% t; t" M7 j$ N9 k9 v) [      { d[j]=d[j]/a[0];
        for (i=1; i&lt;=n-1; i++)
          { u=i*n+i; v=i*m+j;
- O' p' n. X9 a5 s9 ~            for (k=1; k&lt;=i; k++)
- B9 m3 ^6 S( X  t              d[v]=d[v]-a[(k-1)*n+i]*d[(k-1)*m+j];
( S, g& }- W4 E' ]' a2 n7 k0 b            d[v]=d[v]/a;
. ^; P+ ]' `6 l- S: t% h0 {          }
      }
1 z% t! e  T# G9 T* t  |" j    for (j=0; j&lt;=m-1; j++)
      { u=(n-1)*m+j;
5 e7 G* E" ?& Z! m* S        d=d/a[n*n-1];
$ ~* r0 y3 K) j5 `* c% \        for (k=n-1; k&gt;=1; k--)
9 B1 w  Y+ Q' q- K          { u=(k-1)*m+j;
8 g4 K6 T+ F( f) Q4 t# S3 `  r( a            for (i=k; i&lt;=n-1; i++)
              { v=(k-1)*n+i;
. |+ K' Y2 u: `5 d# D                d=d-a[v]*d[i*m+j];
              }
% Q5 G% M- o& I$ c" r( s* x            v=(k-1)*n+k-1;
            d=d/a[v];
          }
" w- X& u* O0 J3 d: U' G! E      }
    return(2);
1 N5 |9 D8 E# |# @" h, Q0 I  }</P><>牛顿向前，向后插值法</P><>double eelgr(x0,h,n,y,t)
  int n;
  double x0,h,t,y[];
  { int i,j,k,m;
' S; g, ]6 e- m3 q: a6 A    double z,s,xi,xj;
1 j2 G# ]2 K# N  }7 P0 b- Q    float p,q;
    z=0.0;
    if (n&lt;1) return(z);
    if (n==1) { z=y[0]; return(z);}
% W% c, e& j! [( R  ]& A    if (n==2)
% D$ A6 u9 l& g+ D' E! h) z# S      { z=(y[1]*(t-x0)-y[0]*(t-x0-h))/h;
        return(z);
  N( ]! _- Y" x$ y! _2 }      }
$ X# @! W7 o* K% z( w4 n    if (t&gt;x0)
      { p=(t-x0)/h; i=(int)p; q=(float)i;
% x# r7 o* V9 v! s; T        if (p&gt;q) i=i+1;
      }
    else i=0;
1 V8 S0 ^+ F& `/ o    k=i-4;
    if (k&lt;0) k=0;
    m=i+3;
  T/ v* j; \! V% S  `5 V; i) H" z    if (m&gt;n-1) m=n-1;
    for (i=k;i&lt;=m;i++)
      { s=1.0; xi=x0+i*h;
        for (j=k; j&lt;=m; j++)
2 F% N/ x. r- ^# Y. [/ v          if (j!=i)
            { xj=x0+j*h;
              s=s*(t-xj)/(xi-xj);
% a! w; `6 T7 [            }
        z=z+s*y;
$ P4 h; f8 ^# ~  F      }
! w; I' ?: l( |3 t. A! O    return(z);
# s$ M6 e! g1 d, t! h7 m  }
向前，向后是一样的思想！！</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[20],t[20],b[20];
    for (i=0; i&lt;=m-1; i++) a=0.0;
    if (m&gt;n) m=n;
    if (m&gt;20) m=20;
    z=0.0;
    for (i=0; i&lt;=n-1; i++) z=z+x/(1.0*n);
( X1 z' j( I( p$ ~! e    b[0]=1.0; d1=1.0*n; p=0.0; c=0.0;
, T, t8 I- B! n7 f    for (i=0; i&lt;=n-1; i++)
: X( s5 o5 S; G1 W6 X- {      { p=p+(x-z); c=c+y;}
2 j9 o2 U2 J5 y, g0 M    c=c/d1; p=p/d1;
    a[0]=c*b[0];
    if (m&gt;1)
8 y8 ?* y. n  J7 M& U      { t[1]=1.0; t[0]=-p;
9 ]9 V$ D$ {0 w3 H        d2=0.0; c=0.0; g=0.0;
        for (i=0; i&lt;=n-1; i++)
9 i) K* u  m4 H& l5 {# m/ [          { q=x-z-p; d2=d2+q*q;
            c=c+y*q;
7 [4 H* @$ `" \' h/ a            g=g+(x-z)*q*q;
          }
! G& N3 f) {0 k- u8 z7 w        c=c/d2; p=g/d2; q=d2/d1;
" x' o  ?$ K& n% u9 C1 G        d1=d2;
# _6 k/ m; V, K6 t/ I4 b        a[1]=c*t[1]; a[0]=c*t[0]+a[0];
6 |( ]( N9 X+ X' L, c3 ~  d3 ^      }
2 {( J) R. m3 U2 b. \% u+ r, ~4 N: {    for (j=2; j&lt;=m-1; j++)
      { s[j]=t[j-1];
/ l) A4 X' a% Y8 r# `" T- y        s[j-1]=-p*t[j-1]+t[j-2];
. O( A8 s  [6 l7 `+ l! x1 i        if (j&gt;=3)
          for (k=j-2; k&gt;=1; k--)
1 x- o3 g- j& [0 M, [& m6 p# n& I            s[k]=-p*t[k]+t[k-1]-q*b[k];
4 {) q7 k2 V, j- }        s[0]=-p*t[0]-q*b[0];
/ n& {% Q' E) [- Q3 D        d2=0.0; c=0.0; g=0.0;
        for (i=0; i&lt;=n-1; i++)
          { q=s[j];
4 m7 W% d( g2 N( A            for (k=j-1; k&gt;=0; k--)
4 D' R2 K. g" b; U* H; P) c" H              q=q*(x-z)+s[k];
  @% [) B  s, e8 z            d2=d2+q*q; c=c+y*q;
            g=g+(x-z)*q*q;
          }
        c=c/d2; p=g/d2; q=d2/d1;
3 {" p" D& g5 A  G, e        d1=d2;
5 l6 L% ]$ X: a" s2 \        a[j]=c*s[j]; t[j]=s[j];
* A6 i; T2 R. _7 h8 H        for (k=j-1; k&gt;=0; k--)
. N1 Z: C# x: ^          { a[k]=c*s[k]+a[k];
! i+ a' C. {2 J! j            b[k]=t[k]; t[k]=s[k];
          }
      }
- l) b0 X" F9 u" O    dt[0]=0.0; dt[1]=0.0; dt[2]=0.0;
    for (i=0; i&lt;=n-1; i++)
      { q=a[m-1];
        for (k=m-2; k&gt;=0; k--)
8 z; o- J: H, U" T" [! ~& w          q=a[k]+q*(x-z);
        p=q-y;
3 G- l& P9 m0 a; N2 [        if (fabs(p)&gt;dt[2]) dt[2]=fabs(p);
& b% D+ b. F9 z$ \! s' @& L/ ]8 p        dt[0]=dt[0]+p*p;
        dt[1]=dt[1]+fabs(p);
      }
    return;
  }</P><>龙贝格积分法</P><>#include "math.h"
  double fromb(a,b,eps)
  double a,b,eps;
  { extern double frombf();
2 x# _8 b# m2 Z  w    int m,n,i,k;
    double y[10],h,ep,p,x,s,q;
    h=b-a;
$ o" S0 n9 |; ~$ Y$ w% B) X    y[0]=h*(frombf(a)+frombf(b))/2.0;
    m=1; n=1; ep=eps+1.0;
- |! s' ?8 b: T' f( d    while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
      { p=0.0;
        for (i=0;i&lt;=n-1;i++)
          { x=a+(i+0.5)*h;
- v: S& _0 R, k5 z/ s1 C            p=p+frombf(x);
$ |0 w+ D' o5 j" t4 B$ c          }
0 ]& }8 J& K+ c7 U: Q4 `        p=(y[0]+h*p)/2.0;
        s=1.0;
& J/ j/ Y7 K% h( N        for (k=1;k&lt;=m;k++)
0 q+ V6 i+ R$ s' s+ e9 h! G          { s=4.0*s;
# a# u3 c4 `2 C8 M! z. e+ ~" t, I            q=(s*p-y[k-1])/(s-1.0);
( ~3 r. J( x$ r5 e% @1 Q            y[k-1]=p; p=q;
          }
/ }' _" p6 G6 {$ V3 m; K9 W        ep=fabs(q-y[m-1]);
        m=m+1; y[m-1]=q; n=n+n; h=h/2.0;
      }
    return(q);
* M' Y" f. R! w4 l1 c  }</P><>呵呵 希望对你有用！！</P>

西西

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

陈了

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