Вычисления площади произвольного многоугольника
Работа из раздела: «
Программирование и комп-ры»
АННОТАЦИЯ
В курсовом проекте решается задача вычисления площади произвольного
многоугольника итерационным алгоритмом.
ЗАДАНИЕ.
Многоугольник (не обязательно выпуклый) задан на плоскости
пересечением координат вершин в порядке обхода его границ. Определить
площадь многоугольника.
СОДЕРЖАНИЕ
Аннотация
Задание на выполнение курсового проекта
Содержание
Введение
1 Разработка программной реализации
2 Проверка на контрольных примерах
3 Заключение
Приложение 1. Блок-схема.
Приложение 2. Программа.
ВВЕДЕНИЕ
Системы, подобные представленной, часто можно встретить в повседневной
жизни.
Данная задача не имеет аналитического решения. В геометрии существуют
формулы, позволяющие вычислять площади правильных многоугольников, но для
произвольных многоугольников таких формул нет. Решение задачи можно
получить численными методами. Рассмотрим два из них.
1. Площадь произвольной фигуры можно вычислить методом Монте-Карло.
Фигура вписывается в другую фигуру с известной площадью. Случайным
образом на последнюю ставятся произвольное количество точек. Площадь
определяется по формуле [pic], где Nф – количество точек попавших в
заданную фигуру, N – общее количество точек. Достоинство данного
метода заключается в простоте реализации, сложность состоит только в
определении попадания точки внутрь заданной фигуры. Очевидно, что
точность вычисленной площади зависит от количества точек. Приемлемая
точность может быть достигнута только при большом их количестве. В
этом заключается один из недостатков метода. Точность также сильно
зависит от качества генератора случайных чисел.
2. Из курса геометрии известно, что любой многоугольник можно разбить
на несколько треугольников, соединяя отрезками несмежные вершины.
Площадь многоугольника при этом будет равна сумме площадей
полученных треугольников. В этом заключается второй метод
определения площади. Площадь треугольника по заданным вершинам легко
определяется по аналитическим формулам, поэтому этот метод позволяет
получить большую точность при меньших затратах вычислительных
ресурсов.
РАЗРАБОТКА ПРОГРАММНОЙ МОДЕЛИ.
Решение задачи будем производить, разбивая одну большую и трудную
задачу на несколько небольших и несложных.
В черновом виде данный алгоритм можно представить в следующем виде:
1. Ввод вершин
2. Предварительная обработка
3. Пока количество вершин больше трех повторяем:
. Найти выпуклую вершин, т.е. вершину, внутренний угол которой меньше
1800. Например на рисунке вершины 1,3,4,5 являются выпуклыми.
. Отрезаем треугольник образованный этой вершиной и двумя смежными.
4. Площадь многоугольника будет равна сумме площадей отрезанных
треугольников и площади оставшегося (при выходе из цикла) треугольника.
Рассмотрим все пункты алгоритма.
1) Ввод данных. Данные будем хранить в текстовом файле ,каждая первая
строка которого содержит количество вершин, а последующие – пары
координат (X,Y), разделенных пробелом. Координаты вершин и внутренние
углы будем хранить в структуре типа:
sd: array[1..100] of
record
x,y: real;
angle: real;
end;
А количество вершин в глобальной переменной n.
Следующая процедура осуществляет ввод данных:
procedure input;
var f: text;
i: integer;
begin
Assign(f,'points.dat');
reset(f);
readln(f, n);
for i:=1 to n do readln(f, sd[i].x, sd[i].y);
end;
2) Предварительная обработка.
В данном пункте алгоритма осуществляется вычисление внутренних углов
многоугольника.
Рассмотрим часть произвольного многоугольника:
Пусть вектор A образует с ось OX угол (1, а вектор B – угол (2. Тогда
угол между ними (внутренний угол многоугольника) будет равен
180–(1–(2. Здесь нельзя использовать формулу угла между векторами
через скалярное произведение, т.к таким образом вычисляется
минимальный угол. Но при этом возможен такой случай:
Угол будет внешним.
Так вычислим либо все внутренние, либо все внешние углы
многоугольника. Чтобы выяснить какие углы мы нашли, рассмотрим
следующую теорему:
Сумма внешних углов произвольного многоугольника больше суммы
внутренних.
Доказательство проведем по индукции:
1) Очевидно, что теорема справедлива для треугольника
2) Предположим, что теорема справедлива для k-угольника
3) Докажем теперь, что теорема справедлива для (k+1)-угольника.
Пусть сумма внутренних углов k-угольника равна (1, а внешних (2.
Из п.2 следует, что (1<(2. k–угольник можно сделать (k+1)-
угольником 'нарастив' его на один треугольник:
Тогда сумма внутренних узлов (k+1)-угольника – (1+(+(+(, а
внешних (2+(360-()-(-(. Из геометрии известно, что сумма углов
треугольника равна 1800. Тогда:
сумма внутренних углов: (1+180
сумма внешних углов: (2+360-((+(+()=(2+180.
Но из п.2 следует, что (1<(2. Следовательно, сумма внешних углов
(k+1)-угольника больше суммы внутренних.(k+1)-угольника.
Теорема доказана.
Внутренние углы многоугольника будем вычислять следующим образом:
. для i-той вершины, имеющей координаты (Xi, Yi) найдем координаты
входящих и выходящих векторов:
A{Xi–Xi–1, Yi–Yi–1} – входящий вектор
B{Xi+1–Xi, Yi+1–Yi} – выходящий вектор
. Вычисляем углы, образованные этими векторами с осями координат
[pic]
[pic]
. Вычисляем угол i-той вершины (i=180–(1–(2.
. Находим сумму [pic]
. Находим сумму [pic]
. Если S1=0 then res:=pi+res
else res:=-pi-res
else
if dy<0 then res:=-res;
ArcCos:=res
end;
begin
dxp:=sd[1].x-sd[n].x;
dyp:=sd[1].y-sd[n].y;
a:=sqrt(dxp*dxp+dyp*dyp);
dxp:=dxp/a;
dyp:=dyp/a;
i:=1;
while i<=(n-1) do
begin
dx:=sd[i+1].x-sd[i].x;
dy:=sd[i+1].y-sd[i].y;
a:=sqrt(dx*dx+dy*dy);
dx:=dx/a;
dy:=dy/a;
if (dx=dxp) and (dy=dyp) then
begin
dec(n);
for j:=i to n do sd[j]:=sd[j+1];
end;
dxp:=dx; dyp:=dy;
inc(i)
end;
dx:=sd[1].x-sd[n].x;
dy:=sd[1].y-sd[n].y;
al1:=ArcCos(dx/sqrt(dx*dx+dy*dy));
for i:=1 to n-1 do
begin
dx:=sd[i+1].x-sd[i].x;
dy:=sd[i+1].y-sd[i].y;
al2:=ArcCos(dx/sqrt(dx*dx+dy*dy));
sd[i].angle:=pi-al1+al2;
if sd[i].angle>2*pi then sd[i].angle:=sd[i].angle-2*pi
else
if sd[i].angle<0 then sd[i].angle:=2*pi+sd[i].angle;
al1:=al2
end;
dx:=sd[1].x-sd[n].x;
dy:=sd[1].y-sd[n].y;
al2:=ArcCos(dx/sqrt(dx*dx+dy*dy));
sd[n].angle:=pi-al1+al2;
if sd[n].angle>2*pi then sd[n].angle:=sd[n].angle-2*pi
else
if sd[n].angle<0 then sd[n].angle:=2*pi+sd[n].angle;
s_in:=0;
s_out:=0;
for i:=1 to n do
begin
if sd[i].angle<0 then sd[i].angle:=2*pi+sd[i].angle;
S_in:=S_in+sd[i].angle;
S_out:=S_out+(2*pi-sd[i].angle);
end;
if S_out0) and (i<>c)
then exit;
cross:=false;
exit
end;
for i:=1 to n-1 do
begin
AA1:=sd[i+1].y-sd[i].y;
BB1:=-(sd[i+1].x-sd[i].x);
CC1:=sd[i].y*(sd[i+1].x-sd[i].x)-sd[i].x*(sd[i+1].y-sd[i].y);
if Mline(sd[i].x, sd[i].y)*Mline(sd[i+1].x,sd[i+1].y)<0 then
if Sline(sd[a].x, sd[a].y)*Sline(sd[b].x,sd[b].y)<0 then exit;
end;
AA1:=sd[1].y-sd[n].y;
BB1:=-(sd[1].x-sd[n].x);
CC1:=sd[n].y*(sd[1].x-sd[n].x)-sd[n].x*(sd[1].y-sd[n].y);
if Mline(sd[n].x, sd[n].y)*Mline(sd[1].x,sd[1].y)<0 then
if Sline(sd[a].x, sd[a].y)*Sline(sd[b].x,sd[b].y)<0 then exit;
cross:=false;
end;
4) Вычисление площадей отсеченных треугольников будем по формуле Герона
[pic]
где [pic].
function St(x1,y1, x2,y2, x3,y3: real): real;
var a, b, c, p: real;
begin
a:=sqrt(sqr(x1-x2)+sqr(y1-y2));
b:=sqrt(sqr(x1-x3)+sqr(y1-y3));
c:=sqrt(sqr(x3-x2)+sqr(y3-y2));
p:=(a+b+c)/2;
St:=sqrt(p*(p-a)*(p-b)*(p-c));
end;
5) Отсечение i-той вершины
dec(n);
for j:=i to n do sd[j]:=sd[j+1];
После отсечения какой-либо вершины необходимо заново рассчитать внутренние
углы многоугольника, т.е. вызвать процедуру Angles.
После вычисления площади, выводим ее на экран и ожидаем нажатия любой
клавиши.
Writeln('Площадь фигуры: ', S:3:3);
Readkey
Полный текст программы приведен в приложении 2.
Проверка на контрольных примерах.
Проверим работу программы на фигурах, площади которых можно вычислить по
формулам.
1) Треугольник
Содержимое файла points.dat
3
0 0
5 0
2 3
Площадь треугольника по формуле: [pic]
Результат работы программы:
Площадь фигуры: 7.500
2) Прямоугольник
Содержимое файла points.dat
4
0 0
5 0
5 3
0 3
Площадь прямоугольника по формуле [pic]
Результат работы программы
Площадь фигуры: 15.000
3) Невыпуклая фигура
Содержимое файла points.dat
4
0 0
3 2
2 0
3 -2
Эта фигура симметрична относительно оси OX. Ее площадь будет равна
[pic]
Результат работы программы:
Площадь фигуры: 4.000
Заключение.
Рассмотренный алгоритм является комбинацией аналитического и численного
методов. Поэтому он включает в себя достоинства обоих. Т.к. основной
операцией вычисления площади многоугольника является вычисление площади
элементарного треугольника, то на результат вычисления не будут влиять
методические погрешности, т.е. погрешности вносимые самим алгоритмом. Этим
приведенный алгоритм отличается от метода Монте-Карло, точность которого
зависит от количества точек. Погрешность будет вноситься лишь на этапе
вычислений и будет определяться конкретной ЭВМ, на которой ведется расчет.
Точность зависит от вещественного типа Real, в котором представляются
основные переменные. Данный тип представим в компьютере лишь с определенной
точностью, зависящей от внутреннего формата числа. Для персонального
компьютера типа IBM PC/AT тип Real имеет следующие параметры:
Длинна, байт 6
Количество значащих цифр 11…12
Диапазон десятичного порядка -39…+38
Эта точность вполне удовлетворительна для нашей задачи.
Приложение 1.
Блок-схема
Приложение 2. Текст программы.
Uses Crt;
const max=100;
var
n, i, j: integer;
sd: array[1..100] of
record
x,y: real;
angle: real;
end;
S: real;
procedure Angles;
var
al1,al2,
dx, dy, dxp, dyp,
s_in, s_out, a: real;
i,j: integer;
function ArcCos(a: real): real;
var res: real;
begin
if abs(a)<1.0E-30 then res:=pi/2
else res:=ArcTan(sqrt(1-a*a)/a);
if dx<0 then
if dy>=0 then res:=pi+res
else res:=-pi-res
else
if dy<0 then res:=-res;
ArcCos:=res
end;
begin
dxp:=sd[1].x-sd[n].x;
dyp:=sd[1].y-sd[n].y;
a:=sqrt(dxp*dxp+dyp*dyp);
dxp:=dxp/a;
dyp:=dyp/a;
i:=1;
while i<=(n-1) do
begin
dx:=sd[i+1].x-sd[i].x;
dy:=sd[i+1].y-sd[i].y;
a:=sqrt(dx*dx+dy*dy);
dx:=dx/a;
dy:=dy/a;
if (dx=dxp) and (dy=dyp) then
begin
dec(n);
for j:=i to n do sd[j]:=sd[j+1];
end;
dxp:=dx; dyp:=dy;
inc(i)
end;
dx:=sd[1].x-sd[n].x;
dy:=sd[1].y-sd[n].y;
al1:=ArcCos(dx/sqrt(dx*dx+dy*dy));
for i:=1 to n-1 do
begin
dx:=sd[i+1].x-sd[i].x;
dy:=sd[i+1].y-sd[i].y;
al2:=ArcCos(dx/sqrt(dx*dx+dy*dy));
sd[i].angle:=pi-al1+al2;
if sd[i].angle>2*pi then sd[i].angle:=sd[i].angle-2*pi
else
if sd[i].angle<0 then sd[i].angle:=2*pi+sd[i].angle;
al1:=al2
end;
dx:=sd[1].x-sd[n].x;
dy:=sd[1].y-sd[n].y;
al2:=ArcCos(dx/sqrt(dx*dx+dy*dy));
sd[n].angle:=pi-al1+al2;
if sd[n].angle>2*pi then sd[n].angle:=sd[n].angle-2*pi
else
if sd[n].angle<0 then sd[n].angle:=2*pi+sd[n].angle;
s_in:=0;
s_out:=0;
for i:=1 to n do
begin
if sd[i].angle<0 then sd[i].angle:=2*pi+sd[i].angle;
S_in:=S_in+sd[i].angle;
S_out:=S_out+(2*pi-sd[i].angle);
end;
if S_out0) and (i<>c)
then exit;
cross:=false;
exit
end;
for i:=1 to n-1 do
begin
AA1:=sd[i+1].y-sd[i].y;
BB1:=-(sd[i+1].x-sd[i].x);
CC1:=sd[i].y*(sd[i+1].x-sd[i].x)-sd[i].x*(sd[i+1].y-sd[i].y);
if Mline(sd[i].x, sd[i].y)*Mline(sd[i+1].x,sd[i+1].y)<0 then
if Sline(sd[a].x, sd[a].y)*Sline(sd[b].x,sd[b].y)<0 then exit;
end;
AA1:=sd[1].y-sd[n].y;
BB1:=-(sd[1].x-sd[n].x);
CC1:=sd[n].y*(sd[1].x-sd[n].x)-sd[n].x*(sd[1].y-sd[n].y);
if Mline(sd[n].x, sd[n].y)*Mline(sd[1].x,sd[1].y)<0 then
if Sline(sd[a].x, sd[a].y)*Sline(sd[b].x,sd[b].y)<0 then exit;
cross:=false;
end;
begin
ClrScr;
input;
Angles;
S:=0;
while n>3 do
begin
i:=1;
while (sd[i].angle>pi) or (cross(i)) do
inc(i);
if i=1 then
S:=S+St(sd[1].x,sd[1].y, sd[2].x,sd[2].y, sd[n].x,sd[n].y)
else
if i=n then
S:=S+St(sd[n].x,sd[n].y, sd[1].x,sd[1].y, sd[n-1].x,sd[n-
1].y)
else S:=S+St(sd[i].x,sd[i].y, sd[i-1].x,sd[i-1].y,
sd[i+1].x,sd[i+1].y);
dec(n);
for j:=i to n do sd[j]:=sd[j+1];
Angles
end;
S:=S+St(sd[1].x,sd[1].y, sd[2].x,sd[2].y, sd[3].x,sd[3].y);
Writeln('Площадь фигуры: ', S:3:3);
Readkey
end.
-----------------------
2
1
4
5
3
A
B
(
(
(
А
A1
A2
B1
B2
А
