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

西西

<>本人紧急求助，希望各位大虾们伸出你们友情之手，帮我一下。</P>
& d/ x" K. C" D/ d  C6 w<>1。牛顿迭代法</P>
<>1。全主消元法</P>
. C6 T2 [) K6 N5 k0 ~: [<>3。改进平方根法</P>
<>4。牛顿向前，向后插值法</P>
<>5。加权最小二乘法</P>
! d+ j- k2 U& b; e+ |1 o7 c! Y<>6。龙贝格积分法</P>
<>拜托了。。。。</P>

loveyuan

<>1。牛顿迭代法</P><>
0 C) n8 w3 w" `# c- z  #include "stdio.h"
  #include "math.h"
  int dnewt(x,eps,js)
  int js;
$ h/ w4 H2 ~: \8 c# U  double *x,eps;
  { extern void dnewtf();
: x" t3 M. z! ?5 X( [; Q* `3 n7 K    int k,l;
    double y[2],d,p,x0,x1;
# B" a3 \- e5 X. }5 g    l=js; k=1; x0=*x;
    dnewtf(x0,y);
    d=eps+1.0;
    while ((d&gt;=eps)&amp;&amp;(l!=0))
+ u$ c0 B. u5 r# T, G0 @: d# n      { if (fabs(y[1])+1.0==1.0)
7 N) A# h& x4 j& k' v! o          { printf("err\n"); return(-1);}
* L" h9 I7 Z4 O5 u" C& W        x1=x0-y[0]/y[1];
3 h* }+ L7 B) X, L! {! L! U0 L8 u        dnewtf(x1,y);
        d=fabs(x1-x0); p=fabs(y[0]);
        if (p&gt;d) d=p;
+ E- s* U0 X0 Z8 B/ k        x0=x1; l=l-1;
7 P9 s$ l% k, x/ v6 G/ A; a      }
) s! ^/ s  f# i2 D2 M  j9 [& @    *x=x1;
    k=js-l;
    return(k);
( i( u1 Z6 j) |, H$ s% |* C) c+ v  }</P><>全主消元法</P><>#include "stdlib.h"
  #include "stdio.h"
$ r# @" y5 T! R8 D' |$ F  int acgas(ar,ai,n,br,bi)
& v! y1 u5 E2 h# ], ~  int n;
; ?; d" n3 t: i" H) }- Q! d  double ar[],ai[],br[],bi[];
  { int *js,l,k,i,j,is,u,v;
% [$ w7 U" E! B3 p  d& ]2 _" K    double p,q,s,d;
    js=malloc(n*sizeof(int));
& r& t/ ^- e, f/ V9 X6 ]' {' j6 p    for (k=0;k&lt;=n-2;k++)
      { d=0.0;
        for (i=k;i&lt;=n-1;i++)
8 `0 ]$ g$ Y# \1 x" p        for (j=k;j&lt;=n-1;j++)
1 D! Q0 x& k1 u# }          { u=i*n+j;
* ?. g. R- y# v, U! t2 N            p=ar*ar+ai*ai;
            if (p&gt;d) {d=p;js[k]=j;is=i;}
. C2 N0 Z6 ^. h' K. @' ^% E          }
        if (d+1.0==1.0)
- [- {4 H* N( f! @+ ?( M  I5 T          { free(js); printf("err**fail\n");
            return(0);
9 u0 ~3 }& Z" d4 [; S          }
        if (is!=k)
          { for (j=k;j&lt;=n-1;j++)
              { u=k*n+j; v=is*n+j;
, l8 y: w; o/ z3 B) _) W                p=ar; ar=ar[v]; ar[v]=p;
                p=ai; ai=ai[v]; ai[v]=p;
8 [+ |+ n: o# P  E6 T              }
            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&lt;=n-1;i++)
4 d/ x: q( i$ Z            { u=i*n+k; v=i*n+js[k];
3 {. O8 X" {3 Q! i$ ?) T              p=ar; ar=ar[v]; ar[v]=p;
              p=ai; ai=ai[v]; ai[v]=p;
  x5 c2 {2 V# N& j1 O9 Y% p% w            }
        v=k*n+k;
- q& P* p& `. H5 ?- j1 \4 X        for (j=k+1;j&lt;=n-1;j++)
2 P# u/ I2 n0 S& h; ]4 c- M9 ~          { u=k*n+j;
            p=ar*ar[v]; q=-ai*ai[v];
            s=(ar[v]-ai[v])*(ar+ai);
( F0 E  H: d2 [% k6 _2 w* _% H/ U            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]);
! f' h: y* G6 \& e2 H$ s/ b7 F        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++)
& y# _. T2 H  M              { v=k*n+j; l=i*n+j;
" b, o' M( V" Y# s2 `                p=ar*ar[v]; q=ai*ai[v];
. i) o0 C8 W, ^. e5 c% o- S                s=(ar+ai)*(ar[v]+ai[v]);
0 {+ F# k, T- v4 f8 m) r7 v                ar[l]=ar[l]-p+q;
                ai[l]=ai[l]-s+p+q;
              }
8 K* l( ]/ _: a/ s; U- I            p=ar*br[k]; q=ai*bi[k];
            s=(ar+ai)*(br[k]+bi[k]);
            br=br-p+q; bi=bi-s+p+q;
: A! J3 z# b+ z  L% k  q          }
      }
    u=(n-1)*n+n-1;
5 p; v+ v( W; J; z* `' S8 L! ^  l, S    d=ar*ar+ai*ai;
    if (d+1.0==1.0)
      { free(js); printf("err**fail\n");
  A; s3 d- k% P& {        return(0);
6 m' @# e( F* K# y6 K. ?8 _      }
/ \% m( s' S6 u: f9 w  K    p=ar*br[n-1]; q=-ai*bi[n-1];
    s=(ar-ai)*(br[n-1]+bi[n-1]);
4 E0 |, B" d0 \) g7 L    br[n-1]=(p-q)/d; bi[n-1]=(s-p-q)/d;
, T) F7 x' ]) a! f1 B    for (i=n-2;i&gt;=0;i--)
7 t( A6 _8 p2 d% K' `- I2 r4 M0 M4 B    for (j=i+1;j&lt;=n-1;j++)
  n6 ~% |6 ~! ~" A+ \. E. V( o' a      { u=i*n+j;
        p=ar*br[j]; q=ai*bi[j];
* c: @" L. z9 S" Y+ f        s=(ar+ai)*(br[j]+bi[j]);
        br=br-p+q;
; p' h9 ]) L$ h3 P6 o# N7 r        bi=bi-s+p+q;
      }
    js[n-1]=n-1;
0 k) V, g' l& U# O  j    for (k=n-1;k&gt;=0;k--)
; p. o  u' L+ l8 g      if (js[k]!=k)
        { p=br[k]; br[k]=br[js[k]]; br[js[k]]=p;
          p=bi[k]; bi[k]=bi[js[k]]; bi[js[k]]=p;
$ e0 I4 G+ Q4 Q2 _+ s. g4 _. Q        }
    free(js);
+ K/ W! ]0 u; o' Y5 @. }* L5 ^    return(1);
8 `/ \$ U1 t/ }8 L  r7 C4 O  }</P><>平方根法</P><>#include "math.h"
9 F; b) z; ~0 X; F+ V, m  #include "stdio.h"
) n. o  {! d/ k- X  }5 s+ x# `3 m6 `  int achol(a,n,m,d)
/ n3 n% o& B. m/ e3 U) w: ?  int n,m;
  double a[],d[];
7 R+ W9 Q3 L# y7 B; H( o8 u7 L# Y  { 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]);
& r! D/ h- [' f" N' x. @- J6 M6 l    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;
        for (j=1; j&lt;=i; j++)
% T3 s4 {, X/ m/ u1 E0 S          { v=(j-1)*n+i;
  \3 ~8 Z$ J; E& ]/ N. `" J            a=a-a[v]*a[v];
3 X; N& [. v4 s4 D( H  }          }
8 K/ J7 v% h5 p3 l) ~$ U+ ?        if ((a+1.0==1.0)||(a&lt;0.0))
          { printf("fail\n"); return(-2);}
: C0 p  g; O# `- ^! x+ r* D& f        a=sqrt(a);
        if (i!=(n-1))
& d' R) Q2 V% R: q2 B          { for (j=i+1; j&lt;=n-1; j++)
              { v=i*n+j;
7 `  U3 w9 A1 n- f$ R                for (k=1; k&lt;=i; k++)
3 D2 Q1 {5 s- i( m; L2 Q                  a[v]=a[v]-a[(k-1)*n+i]*a[(k-1)*n+j];
                a[v]=a[v]/a;
              }
          }
      }
- v" }' T3 F  @    for (j=0; j&lt;=m-1; j++)
0 Q% @) i2 }  i. E  `- ~      { 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++)
              d[v]=d[v]-a[(k-1)*n+i]*d[(k-1)*m+j];
            d[v]=d[v]/a;
          }
      }
1 S+ L% p/ Z  g+ c( `6 Z3 }    for (j=0; j&lt;=m-1; j++)
      { u=(n-1)*m+j;
' U* n' {7 K0 u# S9 x4 N$ S        d=d/a[n*n-1];
" L7 Y( `, {: {1 R3 X1 g        for (k=n-1; k&gt;=1; k--)
          { u=(k-1)*m+j;
, g( J- ~, X, V) b            for (i=k; i&lt;=n-1; i++)
0 P6 |2 Z" U1 f9 k" A! y              { v=(k-1)*n+i;
; J4 U4 l) U* d+ [& B/ D                d=d-a[v]*d[i*m+j];
7 {5 U  q( Y' l8 N/ ?' b# n              }
            v=(k-1)*n+k-1;
            d=d/a[v];
5 S' m& p6 t1 Y; L) l1 Z8 H          }
/ j$ e+ d( v& I/ I& W8 W      }
    return(2);
1 m: y4 z3 A2 I  }</P><>牛顿向前，向后插值法</P><>double eelgr(x0,h,n,y,t)
& F- f( y/ ]; f5 Q) B! s  int n;
6 P" R* x( ^) w3 A8 @" s$ o, w  double x0,h,t,y[];
  { int i,j,k,m;
, |- O5 i7 l7 j2 B! L$ K! ]    double z,s,xi,xj;
6 ^& o7 p! F; @0 V; x. h    float p,q;
    z=0.0;
) ^5 s; F& J3 i' B. i    if (n&lt;1) return(z);
    if (n==1) { z=y[0]; return(z);}
    if (n==2)
( |' P( L/ p( _) m( U      { z=(y[1]*(t-x0)-y[0]*(t-x0-h))/h;
        return(z);
      }
( ?0 [- E$ m3 i/ A% k0 e$ A/ K    if (t&gt;x0)
      { p=(t-x0)/h; i=(int)p; q=(float)i;
        if (p&gt;q) i=i+1;
" h$ q$ ]. i" o+ `$ E9 i. u      }
    else i=0;
; y; p+ q6 B5 s# n1 w9 f    k=i-4;
    if (k&lt;0) k=0;
    m=i+3;
* t/ x5 G- Z; G# t; ]    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++)
          if (j!=i)
. [# N: S4 W8 q- H( c, a            { xj=x0+j*h;
              s=s*(t-xj)/(xi-xj);
& H8 J" A* `8 S' o% X            }
4 I! {# [* P& @5 x$ s9 \9 m        z=z+s*y;
      }
    return(z);
  }
2 I0 q+ d3 c. H/ C2 {; ]3 Z2 L向前，向后是一样的思想！！</P><>加权最小二乘法</P><>#include "math.h"
7 a( r% z8 N/ F( L6 R: i$ w  void hpir1(x,y,n,a,m,dt)
  int n,m;
# x$ m/ V2 J$ T. j) l# k/ E  double x[],y[],a[],dt[];
  { int i,j,k;
    double z,p,c,g,q,d1,d2,s[20],t[20],b[20];
( j; j4 a' J+ {    for (i=0; i&lt;=m-1; i++) a=0.0;
# B! ?8 B7 v6 J5 _    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);
1 n1 |9 a, g7 t( z" |    b[0]=1.0; d1=1.0*n; p=0.0; c=0.0;
    for (i=0; i&lt;=n-1; i++)
+ ]2 k& ~* e" D6 Z6 R+ k. B6 ]      { p=p+(x-z); c=c+y;}
    c=c/d1; p=p/d1;
    a[0]=c*b[0];
) B# a, ~+ b$ d& d    if (m&gt;1)
7 r3 ]+ D, e+ g5 h1 V      { t[1]=1.0; t[0]=-p;
! d/ I- R: c. N. l6 L        d2=0.0; c=0.0; g=0.0;
6 ]# r: L0 q' g4 c+ e        for (i=0; i&lt;=n-1; i++)
          { q=x-z-p; d2=d2+q*q;
9 K9 B+ u& N7 s! A6 `            c=c+y*q;
* ^4 ~' F' i# M3 U, k, r3 m+ ]9 J            g=g+(x-z)*q*q;
& q9 g  e- c( ]- U2 l          }
        c=c/d2; p=g/d2; q=d2/d1;
" Y6 H% B. q* R; }- j        d1=d2;
        a[1]=c*t[1]; a[0]=c*t[0]+a[0];
      }
1 C& U' N/ b; s    for (j=2; j&lt;=m-1; j++)
* {2 L9 I. m$ q- K3 w( `$ H      { s[j]=t[j-1];
! D3 Y$ r5 S) C) ~2 T; M        s[j-1]=-p*t[j-1]+t[j-2];
        if (j&gt;=3)
5 _7 p. c# r  c$ [3 F! r          for (k=j-2; k&gt;=1; k--)
            s[k]=-p*t[k]+t[k-1]-q*b[k];
' x, Y# P- A9 I! i: e        s[0]=-p*t[0]-q*b[0];
        d2=0.0; c=0.0; g=0.0;
        for (i=0; i&lt;=n-1; i++)
          { q=s[j];
            for (k=j-1; k&gt;=0; k--)
$ Q; ]- v) W$ L, f. \0 e$ W              q=q*(x-z)+s[k];
$ ~' c* i$ C/ R2 T1 a            d2=d2+q*q; c=c+y*q;
8 q! Q3 U! R% a6 z3 G            g=g+(x-z)*q*q;
          }
7 Y* n8 A# a6 Y" M, i1 t        c=c/d2; p=g/d2; q=d2/d1;
  o. B- z. _4 q9 h, ~/ m/ R        d1=d2;
        a[j]=c*s[j]; t[j]=s[j];
+ q# p; [& ]! G2 ]: ?        for (k=j-1; k&gt;=0; k--)
          { a[k]=c*s[k]+a[k];
            b[k]=t[k]; t[k]=s[k];
8 J2 l+ R- e4 T" E. V* L# B$ k0 V          }
      }
8 J$ @! v* m  u# p4 C    dt[0]=0.0; dt[1]=0.0; dt[2]=0.0;
0 M4 ]5 n$ X5 _6 r# c. H0 k2 r    for (i=0; i&lt;=n-1; i++)
7 ?( n( {- K* ^1 N4 f      { q=a[m-1];
+ i, e- @4 L9 K3 W( n% f        for (k=m-2; k&gt;=0; k--)
          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;
- X) D  H1 y0 |5 X# e. H( I        dt[1]=dt[1]+fabs(p);
  e5 ?6 Z- S" h& I      }
    return;
, c5 y) |- d; U7 p( @  }</P><>龙贝格积分法</P><>#include "math.h"
  double fromb(a,b,eps)
  double a,b,eps;
4 s) k' _2 l# O( R  { extern double frombf();
$ Q* T# R* [2 w* a# N) z( \    int m,n,i,k;
0 k1 F3 k/ Y6 c& I    double y[10],h,ep,p,x,s,q;
# A. a2 w/ a" {/ d6 ~/ a+ `+ W+ [: G    h=b-a;
    y[0]=h*(frombf(a)+frombf(b))/2.0;
7 t+ o' D) R& l) Y    m=1; n=1; ep=eps+1.0;
    while ((ep&gt;=eps)&amp;&amp;(m&lt;=9))
      { p=0.0;
        for (i=0;i&lt;=n-1;i++)
1 R4 q$ X$ _" L$ `          { x=a+(i+0.5)*h;
            p=p+frombf(x);
          }
        p=(y[0]+h*p)/2.0;
        s=1.0;
, s9 r% p8 C+ \4 H6 g$ g        for (k=1;k&lt;=m;k++)
% W6 X* H, C; s, {          { s=4.0*s;
            q=(s*p-y[k-1])/(s-1.0);
: M9 N7 L& D: n- O7 z, ?            y[k-1]=p; p=q;
          }
        ep=fabs(q-y[m-1]);
( x& ?7 w/ T6 f3 g! ]$ F        m=m+1; y[m-1]=q; n=n+n; h=h/2.0;
      }
2 M4 G, f1 Y+ G% N5 D    return(q);
, L2 s- U  Q( f4 G& `9 z2 z5 k/ X. k  }</P><>呵呵 希望对你有用！！</P>

西西

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

陈了

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