program assignment4_2;

uses math;

var I, alpha, converge, periodratio:real;
    j:integer;



function integrand(theta:real):real; far;
begin

  if (abs(theta-alpha)<1E-9) then integrand:=0
  else integrand:=1/sqrt(cos(theta)-cos(alpha))-1/sqrt((alpha-theta)*sin(alpha));

end;



type funcparam = function(x:real):real;

procedure integrate(func:funcparam; lowerlim, upperlim, converge:real; var I:real);
var j, n:integer;
    yodd, yeven, Iprev, Icurr:real;

begin

  n:=8;

  yodd:=0;
  yeven:=0;
  I:=func(lowerlim)+func(upperlim);

  for j:=1 to (n-1) do begin

    if odd(j) then	yodd := yodd+func(j*(upperlim-lowerlim)/n+lowerlim)
    else		yeven:=yeven+func(j*(upperlim-lowerlim)/n+lowerlim);

  end;

  Icurr:=(I+2*yeven+4*yodd)*(upperlim-lowerlim)/(3*n);

  repeat begin

    Iprev:=Icurr;

    yeven:=yeven+yodd;
    yodd:=0;

    for j:=1 to n do yodd := yodd+func((2*j-1)*(upperlim-lowerlim)/(2*n)+lowerlim);

    n:=n*2;

    Icurr:=(I+2*yeven+4*yodd)*(upperlim-lowerlim)/(3*n);

  end;
  until (abs(Icurr-Iprev)<abs(converge));

  I:=Icurr;

end;



begin

  converge:=1E-6;

  writeln('  alpha  |   T/T0   ');
  writeln('---------+----------');

  for j:=1 to 19 do begin

    alpha:=j*Pi()/20;

    integrate(@integrand, 0, alpha, converge, I);

    periodratio:=(sqrt(2)*I + 2*sqrt(2*alpha/sin(alpha)))/Pi();

    writeln(' ',j:2,'pi/20 | ',periodratio:8:5);

  end;

  writeln();
  write('Press any key to continue.');
  readln;

end.