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

西西

<>本人紧急求助，希望各位大虾们伸出你们友情之手，帮我一下。</P>
<>1。牛顿迭代法</P>
) Y2 \+ \; g% p! K* ?<>1。全主消元法</P>
<>3。改进平方根法</P>
<>4。牛顿向前，向后插值法</P>
<>5。加权最小二乘法</P>
<>6。龙贝格积分法</P>
. i: g- e9 W5 m% V; E& g9 |( I<>拜托了。。。。</P>

loveyuan

<>1。牛顿迭代法</P><>
  #include "stdio.h"
" k# m2 }+ O7 R" n9 T  #include "math.h"
/ j% \6 g9 x& w3 @8 b/ g$ p% O8 Z  int dnewt(x,eps,js)
0 _6 e  j! S. j6 i" m  N( A( a  int js;
  double *x,eps;
( {" M0 l1 g  z9 e  { extern void dnewtf();
& W9 r5 P, H  X# A    int k,l;
    double y[2],d,p,x0,x1;
+ P* c/ j: C1 }/ w, B' d1 A    l=js; k=1; x0=*x;
    dnewtf(x0,y);
    d=eps+1.0;
    while ((d&gt;=eps)&amp;&amp;(l!=0))
4 N, F6 }! `- y- S& R: P      { if (fabs(y[1])+1.0==1.0)
          { printf("err\n"); return(-1);}
        x1=x0-y[0]/y[1];
        dnewtf(x1,y);
# s! V9 L& X! K5 @6 H        d=fabs(x1-x0); p=fabs(y[0]);
        if (p&gt;d) d=p;
        x0=x1; l=l-1;
      }
$ T2 H# ~7 D8 S    *x=x1;
% p, Z6 ~0 j6 V    k=js-l;
    return(k);
1 w5 P+ O  C" F# X! P2 b2 r  }</P><>全主消元法</P><>#include "stdlib.h"
  #include "stdio.h"
  int acgas(ar,ai,n,br,bi)
. t; a% C7 c' P  int n;
  double ar[],ai[],br[],bi[];
' ^% {8 Z' j+ f* P2 g  { int *js,l,k,i,j,is,u,v;
    double p,q,s,d;
    js=malloc(n*sizeof(int));
2 O( K% C: r5 k1 ?1 s* O5 v; ~9 J    for (k=0;k&lt;=n-2;k++)
" M: f7 s$ y; D, W& c0 z      { d=0.0;
        for (i=k;i&lt;=n-1;i++)
        for (j=k;j&lt;=n-1;j++)
          { u=i*n+j;
+ n: R) q& h$ V. G, v9 O1 e            p=ar*ar+ai*ai;
+ E  ~. }+ E$ `" d0 |            if (p&gt;d) {d=p;js[k]=j;is=i;}
          }
5 R* O' L. f" `% w        if (d+1.0==1.0)
          { free(js); printf("err**fail\n");
& M/ i7 F$ p2 d5 {- I# a  v/ f            return(0);
7 e; f3 N9 G5 f! G( N- `4 E          }
5 P' ?$ T( t) A7 E, d. ^        if (is!=k)
          { for (j=k;j&lt;=n-1;j++)
" A9 s' o+ j% n4 P8 N              { u=k*n+j; v=is*n+j;
                p=ar; ar=ar[v]; ar[v]=p;
2 l  `' W8 x: R                p=ai; ai=ai[v]; ai[v]=p;
" v& v- Q$ u  M) @! ~) Q# ?; N" I              }
            p=br[k]; br[k]=br[is]; br[is]=p;
1 [2 k4 t8 e6 B" x            p=bi[k]; bi[k]=bi[is]; bi[is]=p;
          }
        if (js[k]!=k)
          for (i=0;i&lt;=n-1;i++)
            { u=i*n+k; v=i*n+js[k];
- U5 \, _6 @4 i1 N1 r6 G5 E6 [" P              p=ar; ar=ar[v]; ar[v]=p;
: V6 C5 W6 z3 j; u' U1 x              p=ai; ai=ai[v]; ai[v]=p;
            }
        v=k*n+k;
1 Z! i/ I' t. a: x( _        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;
9 K4 M/ W! P2 I4 ^/ h. L. J          }
        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;
            for (j=k+1;j&lt;=n-1;j++)
              { v=k*n+j; l=i*n+j;
1 Z% ^3 {2 b2 a$ R0 t" c                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;
              }
% d+ `5 o, B0 T            p=ar*br[k]; q=ai*bi[k];
' d5 L( s: Q" T            s=(ar+ai)*(br[k]+bi[k]);
) n; \0 }- L# K9 J! R+ Y            br=br-p+q; bi=bi-s+p+q;
1 N- D% c, x7 I/ k  [- T& ^, v) W          }
      }
    u=(n-1)*n+n-1;
$ W$ L7 N2 z' M, m9 N    d=ar*ar+ai*ai;
    if (d+1.0==1.0)
$ x: R4 O8 Z: u/ n      { free(js); printf("err**fail\n");
) X& {) ~9 H6 s5 _& ^        return(0);
      }
0 \/ t% p2 U. U! G3 s5 A/ j    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;
    for (i=n-2;i&gt;=0;i--)
    for (j=i+1;j&lt;=n-1;j++)
, o' h8 V9 i, C7 c& a; m      { u=i*n+j;
        p=ar*br[j]; q=ai*bi[j];
  s0 h0 p7 R6 L. p        s=(ar+ai)*(br[j]+bi[j]);
        br=br-p+q;
0 T/ Y9 }2 J) q4 j9 O/ l        bi=bi-s+p+q;
      }
    js[n-1]=n-1;
9 Z- Y! C& Q! D4 f' E* l2 Y# `# _    for (k=n-1;k&gt;=0;k--)
      if (js[k]!=k)
0 O& Y' P  m/ ~. Y* z        { p=br[k]; br[k]=br[js[k]]; br[js[k]]=p;
! Q, y" }& ?0 j7 g* w- e          p=bi[k]; bi[k]=bi[js[k]]; bi[js[k]]=p;
        }
    free(js);
    return(1);
/ p4 A  a. ]& x6 i  }</P><>平方根法</P><>#include "math.h"
  #include "stdio.h"
  int achol(a,n,m,d)
! c) r7 E% F1 {( C! I7 u+ E  int n,m;
- Q% K- S. P$ Z8 v  double a[],d[];
: |/ [( r; A7 g! }5 M# D) l( T5 P  { int i,j,k,u,v;
% O) z! c& g2 h4 g    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];
    for (i=1; i&lt;=n-1; i++)
      { u=i*n+i;
( k7 \7 }- c' b; N9 G8 s6 ^        for (j=1; j&lt;=i; j++)
          { v=(j-1)*n+i;
1 W# z& [' a" s7 m            a=a-a[v]*a[v];
          }
/ |) }# t2 J6 ]# D- O) j        if ((a+1.0==1.0)||(a&lt;0.0))
          { printf("fail\n"); return(-2);}
        a=sqrt(a);
3 T3 B) p: ]8 x; f. k        if (i!=(n-1))
9 H3 s$ N% @& k' J          { for (j=i+1; j&lt;=n-1; j++)
  h; f" @1 @  q              { v=i*n+j;
                for (k=1; k&lt;=i; k++)
- h) S. _/ x9 h# z2 [                  a[v]=a[v]-a[(k-1)*n+i]*a[(k-1)*n+j];
1 `3 s: Q4 q1 M! {) J" L. f7 a                a[v]=a[v]/a;
              }
          }
      }
    for (j=0; j&lt;=m-1; j++)
      { d[j]=d[j]/a[0];
1 W& L; n  [* Z# R/ V        for (i=1; i&lt;=n-1; i++)
0 z8 ]! b- U3 j# t7 a          { u=i*n+i; v=i*m+j;
            for (k=1; k&lt;=i; k++)
* F; o3 W" i. y& j4 `4 {/ [              d[v]=d[v]-a[(k-1)*n+i]*d[(k-1)*m+j];
% p, {3 ?( H* {3 w2 s* c' H# B6 S" [            d[v]=d[v]/a;
7 v* w7 P9 C+ n6 Z          }
      }
2 m- N- l2 Y2 X0 q& d; h5 X: V/ ^    for (j=0; j&lt;=m-1; j++)
      { u=(n-1)*m+j;
        d=d/a[n*n-1];
& ^: s) J" N* H) s$ S, y# G6 W1 V) u        for (k=n-1; k&gt;=1; k--)
          { u=(k-1)*m+j;
            for (i=k; i&lt;=n-1; i++)
              { v=(k-1)*n+i;
                d=d-a[v]*d[i*m+j];
              }
/ k3 P2 F9 ^# Y3 o, y6 {0 i- M            v=(k-1)*n+k-1;
            d=d/a[v];
9 s5 {* ~# Z9 q# p" a# K          }
+ h& F( G  D1 F+ g/ T      }
! O4 Y1 g4 @8 q# C& l: E, ]) G2 ]    return(2);
5 r$ N  R! W5 ^; i  }</P><>牛顿向前，向后插值法</P><>double eelgr(x0,h,n,y,t)
4 T# x, `& |" t7 q& D2 D  int n;
5 ?) f- @9 g0 d* H* F  double x0,h,t,y[];
. h* Y% z. X  t" J7 s) O7 ?9 M  { int i,j,k,m;
    double z,s,xi,xj;
% j4 a7 \" n1 f) ~2 z+ z    float p,q;
    z=0.0;
; Y- ]/ k" F/ _' x5 Q2 E    if (n&lt;1) return(z);
! q  [- I6 i; t4 S4 I, D0 a    if (n==1) { z=y[0]; return(z);}
- X: X  o+ w& e1 _! y    if (n==2)
      { z=(y[1]*(t-x0)-y[0]*(t-x0-h))/h;
, Z0 D* m) B; ~. I* P0 k        return(z);
6 b; p: H+ C7 N( n, V: v, ]      }
' U, L: U3 T' I4 @8 i+ L# U    if (t&gt;x0)
6 X3 r/ B$ s% U- F$ o6 A      { p=(t-x0)/h; i=(int)p; q=(float)i;
6 z. q, G! D  i; U% {. @# u, O4 O        if (p&gt;q) i=i+1;
7 {% u: j  K6 p7 E      }
    else i=0;
7 s; I4 T+ |* C: V. f: A    k=i-4;
    if (k&lt;0) k=0;
    m=i+3;
  g7 T+ v, E# i& X    if (m&gt;n-1) m=n-1;
    for (i=k;i&lt;=m;i++)
      { s=1.0; xi=x0+i*h;
/ U) g$ x1 Z) @4 z        for (j=k; j&lt;=m; j++)
+ U0 _3 O7 U5 }  T# [$ k9 a          if (j!=i)
! L( f+ K% J! `$ Q' G5 B3 A4 S  P            { xj=x0+j*h;
              s=s*(t-xj)/(xi-xj);
* n4 u3 u' f# f, C  q            }
        z=z+s*y;
3 @: ?# y7 s% q' [      }
2 e4 K3 m2 w) `$ r1 L- }; F    return(z);
  }
( a* U, b5 n1 }% L  y# X& s5 j! @2 W向前，向后是一样的思想！！</P><>加权最小二乘法</P><>#include "math.h"
  S% k/ E. j6 @# [, Q  void hpir1(x,y,n,a,m,dt)
2 W8 y  e1 K) A2 F6 r& Z$ H# u" u) e  int n,m;
# h( u% |4 M# Q- i' q  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;
3 N+ v2 P& e& _% K) c" u! w: m) x    if (m&gt;n) m=n;
    if (m&gt;20) m=20;
( U* J9 @2 K/ \  ^    z=0.0;
    for (i=0; i&lt;=n-1; i++) z=z+x/(1.0*n);
9 i9 P+ U: {0 o% ^* b    b[0]=1.0; d1=1.0*n; p=0.0; c=0.0;
6 `- m" e% e; z# ^* `  ^    for (i=0; i&lt;=n-1; i++)
      { p=p+(x-z); c=c+y;}
! S; U) n1 \/ {/ u. d    c=c/d1; p=p/d1;
    a[0]=c*b[0];
5 A2 |' f" y0 B" m    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++)
          { q=x-z-p; d2=d2+q*q;
0 B  {) Q% x3 G0 t" P, ?            c=c+y*q;
            g=g+(x-z)*q*q;
          }
  {# f/ f4 g* L$ |        c=c/d2; p=g/d2; q=d2/d1;
        d1=d2;
  P6 _& y' {" p2 m        a[1]=c*t[1]; a[0]=c*t[0]+a[0];
4 E/ Y* t9 c0 m( Y; i, Q      }
    for (j=2; j&lt;=m-1; j++)
      { s[j]=t[j-1];
        s[j-1]=-p*t[j-1]+t[j-2];
0 U' A+ A$ }% Q2 z        if (j&gt;=3)
          for (k=j-2; k&gt;=1; k--)
( G0 M+ R+ g/ _5 K; `; A            s[k]=-p*t[k]+t[k-1]-q*b[k];
        s[0]=-p*t[0]-q*b[0];
- o: t9 H+ q: m: H' y9 u6 [9 u        d2=0.0; c=0.0; g=0.0;
        for (i=0; i&lt;=n-1; i++)
          { q=s[j];
* t% I* R+ f& N5 e            for (k=j-1; k&gt;=0; k--)
              q=q*(x-z)+s[k];
+ \2 r; s! F7 V1 l; W" ?& {5 A            d2=d2+q*q; c=c+y*q;
- I% S5 H+ }7 N3 ]6 Y* b( _            g=g+(x-z)*q*q;
          }
( S& y/ ^0 h) _8 `" P        c=c/d2; p=g/d2; q=d2/d1;
        d1=d2;
        a[j]=c*s[j]; t[j]=s[j];
        for (k=j-1; k&gt;=0; k--)
          { a[k]=c*s[k]+a[k];
" ~) O* K0 T$ m, J$ x            b[k]=t[k]; t[k]=s[k];
4 |( G- x& C: I" a  t          }
      }
    dt[0]=0.0; dt[1]=0.0; dt[2]=0.0;
2 w) M; H7 i9 u9 R2 L9 [2 q2 h3 ]& v    for (i=0; i&lt;=n-1; i++)
* |# `! A' z# _0 u, ], O& D) H      { q=a[m-1];
" I& B% j: e2 q+ H" a        for (k=m-2; k&gt;=0; k--)
. Q  c' ]( j7 \0 K3 T, T          q=a[k]+q*(x-z);
        p=q-y;
        if (fabs(p)&gt;dt[2]) dt[2]=fabs(p);
        dt[0]=dt[0]+p*p;
! b1 U8 ~( E% K+ e: z4 I& r        dt[1]=dt[1]+fabs(p);
      }
    return;
  }</P><>龙贝格积分法</P><>#include "math.h"
) V' d. D; ]- i( O1 L  double fromb(a,b,eps)
  double a,b,eps;
  { extern double frombf();
    int m,n,i,k;
0 I9 I  F" [7 k% R/ F% S# ^    double y[10],h,ep,p,x,s,q;
4 Y8 V' K$ J% \& I    h=b-a;
+ O- x# O8 D9 N! S9 M9 C    y[0]=h*(frombf(a)+frombf(b))/2.0;
& Z" r3 j+ J' }  r    m=1; n=1; ep=eps+1.0;
    while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
) W/ f8 R4 e% u' \      { p=0.0;
  m, ^- c$ P3 K/ I        for (i=0;i&lt;=n-1;i++)
. z4 {- _4 N$ d% H' g          { x=a+(i+0.5)*h;
            p=p+frombf(x);
5 V) ]0 \# }9 u. S# Q( F          }
        p=(y[0]+h*p)/2.0;
1 m0 P0 w3 o  I9 ^8 h& g        s=1.0;
        for (k=1;k&lt;=m;k++)
          { s=4.0*s;
            q=(s*p-y[k-1])/(s-1.0);
            y[k-1]=p; p=q;
% r$ y' ^% J, V2 S' p$ ?# G& F          }
        ep=fabs(q-y[m-1]);
        m=m+1; y[m-1]=q; n=n+n; h=h/2.0;
      }
    return(q);
  }</P><>呵呵 希望对你有用！！</P>

西西

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

陈了

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