home *** CD-ROM | disk | FTP | other *** search
/ Programmer 7500 / MAX_PROGRAMMERS.iso / INFO / TURBOPAS / PAS_ENG.ZIP / ERF4.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1985-07-18  |  1.5 KB  |  80 lines
  1. program erfd4;  { -> 334 }
  2.  
  3. { evaluation of the gaussian error function }
  4.  
  5. var x,er,ec  : real;
  6.  done  : boolean;
  7.  
  8.  
  9. function erf(x: real): real;
  10. { infinite series expansion of the Gaussian error function }
  11.  
  12. const sqrtpi  = 1.7724538;
  13.  t2  = 0.66666667;
  14.  t3  = 0.66666667;
  15.  t4  = 0.07619048;
  16.  t5  = 0.01693122;
  17.  t6  = 3.078403E-3;
  18.  t7  = 4.736005E-4;
  19.  t8  = 6.314673E-5;
  20.  t9  = 7.429027E-6;
  21.  t10  = 7.820028E-7;
  22.  t11  = 7.447646E-8;
  23.  t12  = 6.476214E-9;
  24.  
  25. var x2,sum  : real;
  26.  i  : integer;
  27.  
  28. begin
  29.   x2:=x*x;
  30.   sum:=t5+x2*(t6+x2*(t7+x2*(t8+x2*(t9+x2*(t10+x2*(t11+x2*t12))))));
  31.   erf:=2.0*exp(-x2)/sqrtpi*(x*(1+x2*(t2+x2*(t3+x2*(t4+x2*sum)))))
  32. end; { function erf }
  33.  
  34. function erfc(x: real): real;
  35. { complement of error function }
  36. const sqrtpi  = 1.7724538;
  37.  
  38. var x2,v,sum : real;
  39.  
  40. begin
  41.   x2:=x*x;
  42.   v:=1.0/(2.0*x2);
  43.   sum:=v/(1+8*v/(1+9*v/(1+10*v/(1+11*v/(1+12*v)))));
  44.   sum:=v/(1+3*v/(1+4*v/(1+5*v/(1+6*v/(1+7*sum)))));
  45.   erfc:=1.0/(exp(x2)*x*sqrtpi*(1+v/(1+2*sum)))
  46. end;  { function ercf }
  47.  
  48. begin  { main }
  49.   clrscr;
  50.   done:=false;
  51.   writeln;
  52.   repeat
  53.     write('Arg? ');
  54.     readln(x);
  55.     if x<0.0 then done:=true
  56.     else
  57.       begin
  58.  if x=0.0 then
  59.    begin
  60.      er:=0.0;
  61.      ec:=1.0
  62.    end
  63.  else
  64.    begin
  65.      if x<1.5 then
  66.        begin
  67.   er:=erf(x);
  68.   ec:=1.0-er
  69.        end
  70.      else
  71.        begin
  72.   ec:=erfc(x);
  73.   er:=1.0-ec
  74.      end { if }
  75.  end;
  76.  writeln('X= ',x:6:2,', Erf= ',er:7:4,', Erfc= ',ec)
  77.       end { if }
  78.     until done
  79. end.
  80.