Решить задачу модифицированным методом Ньютона

program MetodN;
uses crt;
 const
 Sigma0=5.67e-8;   {vvod konstant}
 Pi=3.14;
var T,R,P,a,en:real;       {Globalnie peremennie}
  function F(Ts:real):real;   {Uravnenie teplovogo balansa}
    begin
    F:=a*(T-Ts)+en*Sigma0*((T*T*T*T)-(Ts*Ts*Ts*Ts))-P/(4*Pi*R*R);
    end;
  function F1(Ts:real):real;     {Proizvodnaya po Ts}
    begin
    F1:=-a-4*en*Sigma0*Ts*Ts*Ts;
    end;
procedure Newton(x1,eps:real; var x2:real; var k:integer);   {Modif metod Newtona}
 var b:real;
  begin
  x2:=x1;
  eps:=0.1;
  k:=0;
     repeat
     b:=x2;
     k:=k+1;
     x2:=b-F(b)/F1(x1);
     until abs(x2-b)<=eps;
  end;
var x0,x,eps:real;
    k:integer;
begin
clrscr;
  T:=300;
  R:=0.05;
  P:=100;
  a:=5;
  en:=0.7;
    repeat
     write('Vvedite nachalnoe priblizenie  x0=');+
     readln(x0);
    until(x0>=-100000)and(x0<=100000);  {Interval poiska korney uravneniya}
    eps:=0.1;
Newton(x0,eps,x,k);
writeln('Ts=',x:0:2,'  iteracij=',k);
readln;
end.