Нахождение координат точек пересечения линий, ограничивающих закрашенные области методом проб (простого сканирования)

program CriticalPoints;
 
 
var
  r, a, b: real;
  step: real;
  analyticalX: array[1..4] of real;
  analyticalY: array[1..4] of real;
  i: integer;
  found: boolean;
  x_inter, y_inter, step_inter: real;
  is_in_circle1, is_in_ellipse1: boolean;
  is_in_circle2, is_in_ellipse2: boolean;
 
procedure CheckIntersection(x, y: real; var is_in_circle: boolean;
  var is_in_ellipse: boolean);
begin
  // Проверка попадания в круг
  is_in_circle := (x * x + y * y <= r * r);
  // Проверка попадания в эллипс
  is_in_ellipse := (x * x / (a * a) + y * y / (b * b) <= 1);
end;
 
function CalculateYFromEllipse(x: real): real;
begin
  CalculateYFromEllipse := sqrt(1225 * (1 - x * x / 441));
end;
 
begin
 
  r := 32;
  a := 21;
  b := 35;
  step := 0.001;
 
  writeln('Критические точки пересечения (круг и эллипс):');
  writeln('-------------------------------------------------');
  writeln('Имя       Aналитический метод     Вычислительный метод');
  writeln('точки      X          Y          X          Y');
  writeln('-------------------------------------------------------------------');
 
  // Аналитические точки (как было вычислено ранее)
 
   analyticalX[1] := 10.63;
   analyticalY[1] := 30.18;
   analyticalX[2] := 10.63;
   analyticalY[2] := -30.18;
   analyticalX[3] := -10.63;
   analyticalY[3] := 30.18;
   analyticalX[4] := -10.63;
   analyticalY[4] := -30.18;
 
 
  for i := 1 to 4 do
  begin
    if (abs(analyticalX[i]) > 0.00001) and (abs(analyticalY[i]) > 0.00001) then
    begin
      x_inter := 0;
      found := false;
      step_inter := step;
 
      while (not found) and (abs(x_inter) < abs(analyticalX[i] * 2)) do
      begin
        y_inter := CalculateYFromEllipse(x_inter);
        if analyticalY[i] < 0 then
          y_inter := -y_inter;
 
        CheckIntersection(x_inter, y_inter, is_in_circle1, is_in_ellipse1);
        CheckIntersection(x_inter + step_inter, y_inter, is_in_circle2, is_in_ellipse2);
 
        if (is_in_circle1 <> is_in_circle2) then
        begin
          found := true;
        end;
 
        if analyticalX[i] > 0 then
          x_inter := x_inter + step_inter
        else
          x_inter := x_inter - step_inter;
 
        if abs(x_inter - analyticalX[i]) < 0.001 then
          found := true;
      end;
 
      //Второй цикл для y
      found := false;
      step_inter := step;
 
      while (not found) and (abs(y_inter) < abs(analyticalY[i] * 2)) do
      begin
          CheckIntersection(x_inter, y_inter, is_in_circle1, is_in_ellipse1);
          CheckIntersection(x_inter, y_inter + step_inter, is_in_circle2, is_in_ellipse2);
 
 
        if (is_in_circle1 <> is_in_circle2) then
        begin
          found := true;
        end;
 
        if analyticalY[i] > 0 then
          y_inter := y_inter + step_inter
        else
          y_inter := y_inter - step_inter;
 
        if abs(y_inter - analyticalY[i]) < 0.001 then
          found := true;
      end;
      if i = 1 then
        write('A')
      else if i = 2 then
        write('B')
      else if i = 3 then
        write('C')
      else
        write('D');
      writeln('  ', analyticalX[i]: 10:4, '  ', analyticalY[i]: 10:4, '  ', x_inter: 10:4, '  ', y_inter: 10:4);
    end;
  end;
 
end.

program AreaCalculation;
 
 
var
r, a, b: real;
N: integer;
areaA_rect, areaA_trap, areaB_rect, areaB_trap: real;
xaA, xbA, xaB, xbB: real;
repeatCalculation: char;
 
 
function IsInCircle(x, y: real): boolean;
begin
IsInCircle := (x * x + y * y <= r * r);
end;
 
function IsInEllipse(x, y: real): boolean;
begin
IsInEllipse := (x * x / (a * a) + y * y / (b * b) <= 1);
end;
 
 
function CalculateYFromEllipse(x: real): real;
begin
CalculateYFromEllipse := sqrt(1024 * (1 - x * x / 441));
end;
 
function CalculateRectangleArea(xa, xb: real; N: integer): real;
var
dx, x, y, area: real;
i: integer;
begin
dx := abs(xb - xa) / N;
area := 0;
x := xa;
 
for i := 1 to N do
begin
y := CalculateYFromEllipse(x + dx/2);
area := area + y * dx;
x := x + dx;
end;
CalculateRectangleArea := area;
end;
 
function CalculateTrapezoidArea(xa, xb: real; N: integer): real;
var
dx, x, y1, y2, area: real;
i: integer;
begin
dx := abs(xb - xa) / N;
area := 0;
x := xa;
 
for i := 1 to N do
begin
y1 := CalculateYFromEllipse(x);
y2 := CalculateYFromEllipse(x+dx);
area := area + (y1 + y2) / 2 * dx;
x := x + dx;
end;
CalculateTrapezoidArea := area;
end;
 
 
begin
 
r := 32;
a := 21;
b := 35;
 
writeln('Вычислите площади закрашенных областей:');
writeln('---------------------------------------');
writeln;
 
 
xaA := 0;
xbA := a;
xaB := -a;
xbB := 0;
 
 
repeat
write('Введите число разбиений N: ');
readln(N);
 
areaA_rect := CalculateRectangleArea(xaA, xbA, N);
areaA_trap := CalculateTrapezoidArea(xaA, xbA, N);
 
areaB_rect := CalculateRectangleArea(xaB, xbB, N);
areaB_trap := CalculateTrapezoidArea(xaB, xbB, N);
 
 
writeln;
writeln('Число разбиений N Область A Область B');
writeln(' Прямоугольник Трапеция Прямоугольник Трапеция');
writeln('--------------------------------------------------------------------------');
writeln(' ', N, ' ', areaA_rect:10:4, ' ', areaA_trap:10:4, ' ', areaB_rect:10:4, ' ', areaB_trap:10:4);
writeln;
write('Хотите повторить вычисления с новым N? (y/n): ');
readln(repeatCalculation);
until (repeatCalculation = 'n') or (repeatCalculation = 'N');
 
end.