program Cube_Root_2; { This program is computing the cube root of 2.   (1-x)^(-m) = 1 + mx + m(m+1)/2! * x^2 + m(m+1)(m+2)/3! * x^3 +...   when m=1/3   1/(1-x)^(1/3) = 1 + 1/3 * x + 1*4/(3*6) * x^2 + 1*4*7/(3*6*9) * x^3 +...   Let x = 375/16000   At last, Answer = sum*25/20   Note: Last some digits may be Wrong. } label ext; const      dn=336; var     i,k:longint;     sum,a:array[1..dn] of integer;     ip:integer; procedure outp; var   i:integer; procedure writep(num:integer); begin     write(num div 1000);     write(num div 100 mod 10);     write(num div 10 mod 10);     write(num mod 10,' '); end; begin    writeln('Cube_root(2)=');    writeln(sum[1],'.');    for i:=2 to dn do         writep(sum[i]) ;    writeln;    writeln('Programmed by j.t.chang'); end; procedure m_div(k:longint); var    i:integer;    r1,c:longint; begin    c:=0;    for i:=ip to dn do       begin          r1:=c*10000+a[i];          a[i]:=r1 div k;          c:=r1 mod k;       end; end; procedure sm_div(k:longint); var    i:integer;    r1,c:longint; begin    c:=0;    for i:=1 to dn do       begin          r1:=c*10000+sum[i];          sum[i]:=r1 div k;          c:=r1 mod k;       end; end; procedure m_mul(k:longint); var    i:integer;    r1,c:longint; begin    c:=0;    for i:=dn downto ip do       begin          r1:=k*a[i]+c;          a[i]:=r1 mod 10000;          c:=r1 div 10000;       end;      if ip=1 then         begin            a[1]:=c;            exit;         end;      i:=ip-1;      while c>0 do       begin          if i=1 then             begin                a[1]:=c;                exit             end;          a[i]:=c mod 10000;          c:=c div 10000;          i:=i-1;       end; end; procedure sm_mul(k:longint); var    i:integer;    r1,c:longint; begin    c:=0;    for i:=dn downto 2 do       begin          r1:=k*sum[i]+c;          sum[i]:=r1 mod 10000;          c:=r1 div 10000;       end;       r1:=1;       r1:=r1*k*sum[1]+c;       sum[1]:=r1 ; end; procedure m_add; var    i:integer;    c:longint; begin    c:=0;    for i:=dn downto 1 do      begin         sum[i]:=sum[i]+a[i]+c;         if sum[i]>=10000 then           begin              c:=sum[i] div 10000;              sum[i]:=sum[i] mod 10000;           end         else c:=0;      end; end; begin     for i:=1 to dn do a[i]:=0;     a[1]:=1;     sum:=a;     k:=1;     ip:=1;     repeat        i:=ip;        while (a[i]=0) do          begin             i:=i+1;             if i>dn then goto ext;          end;        ip:=i;        m_div(3*k);        m_mul(3*k-2);        m_div(16000);        m_mul(375);        m_add;        k:=k+1;     until false; ext:     sm_div(20);     sm_mul(25);     outp; end.

评论