Last active
June 20, 2023 14:41
-
-
Save AsyaMatveeva/cd7833d41e4d71bdef3756e6c0c09646 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import math | |
import random | |
def cos(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round((pe1[0] * pe2[0] + pe1[1] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
def sin(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round(abs(pe1[1] * pe2[0] - pe1[0] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
def delaunay_correct(edge, p1, p2): | |
sin1 = sin(edge, p1) | |
cos1 = cos(edge, p1) | |
sin2 = sin(edge, p2) | |
cos2 = cos(edge, p2) | |
return sin1 * cos2 + sin2 * cos1 >= 0 | |
def line(p1, p2): | |
a = p1[1] - p2[1] | |
b = p2[0] - p1[0] | |
c = p1[0] * p2[1] - p2[0] * p1[1] | |
return [a, b, c] | |
def naive(points): | |
points = list(set(points)) | |
points.sort() | |
p0 = points[0] | |
edges = [] | |
triangulation = {} | |
max_abs_slope = [0] | |
for point in points[1:]: | |
if point[0] == p0[0]: | |
slope = float('inf') | |
else: | |
slope = (point[1] - p0[1])/(point[0]-p0[0]) | |
if abs(slope) > abs(max_abs_slope[0]): | |
max_abs_slope = [slope, point] | |
elif slope == max_abs_slope[0]: | |
max_abs_slope.append(point) | |
max_abs_slope[0] = p0 | |
for i in range(1, len(max_abs_slope)): | |
new_edge = tuple([max_abs_slope[i-1], max_abs_slope[i]]) | |
edges.append(new_edge) | |
triangulation[new_edge] = [] | |
while len(edges) > 0: | |
edge = edges.pop(0) | |
if len(triangulation[edge]) == 2: | |
continue | |
coeff = line(edge[0], edge[1]) | |
a, b, c = coeff[0], coeff[1], coeff[2] | |
if len(triangulation[edge]) == 0: | |
fixed_point = [edge[0][0] - 1, edge[0][1]] | |
else: | |
fixed_point = triangulation[edge][0] | |
fixed_dist = a * fixed_point[0] + b * fixed_point[1] + c | |
vertex = [2, []] | |
for point in points: | |
if (a * point[0] + b * point[1] + c) * fixed_dist < 0: | |
cosine = cos(edge, point) | |
angle = cos((edge[1], point), edge[0]) | |
if cosine < vertex[0]: | |
vertex[0] = cosine | |
vertex[1] = [[angle, point]] | |
elif cosine == vertex[0]: | |
vertex[1].append([angle, point]) | |
if vertex[0] != 2: | |
on_circle = vertex[1][::] | |
on_circle.sort() | |
on_circle.append([1, edge[1]]) | |
triangulation[edge].append(on_circle[0][1]) | |
new_edge = tuple(sorted([edge[0], on_circle[0][1]])) | |
if new_edge not in triangulation.keys(): | |
triangulation[new_edge] = [on_circle[1][1]] | |
edges.append(new_edge) | |
elif on_circle[1][1] not in triangulation[new_edge]: | |
triangulation[new_edge].append(on_circle[1][1]) | |
for i in range(len(on_circle)-1): | |
new_edge = tuple(sorted([on_circle[i][1], on_circle[i+1][1]])) | |
if new_edge not in triangulation.keys(): | |
triangulation[new_edge] = [edge[0]] | |
edges.append(new_edge) | |
elif edge[0] not in triangulation[new_edge]: | |
triangulation[new_edge].append(edge[0]) | |
for i in range(1, len(on_circle)-1): | |
new_edge = tuple(sorted([edge[0], on_circle[i][1]])) | |
triangulation[new_edge] = [on_circle[i-1][1], on_circle[i+1][1]] | |
return triangulation | |
flips = 0 | |
all_points = [] | |
flipped = [] | |
k = 5 | |
for x in range(k): | |
for y in range(k): | |
all_points.append((x, y)) | |
random.shuffle(all_points) | |
for i1 in range(k**2-5): | |
for i2 in range(i1+1, k**2-4): | |
for i3 in range(i2+1, k**2-3): | |
for i4 in range(i3+1, k**2-2): | |
for i5 in range(i4+1, k**2-1): | |
for i6 in range(i5+1, k**2): | |
points = [all_points[i1], all_points[i2], all_points[i3], all_points[i4], all_points[i5], all_points[i6]] | |
triangulation = naive(points) | |
keys = list(triangulation.keys()) | |
for edge in keys: | |
val = triangulation[edge] | |
if len(val) == 2 and not delaunay_correct(edge, val[0], val[1]): | |
flips += 1 | |
flipped.append(points) | |
print(flips) | |
for set_of_points in flipped: | |
print(set_of_points) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import math | |
import random | |
from turtle import * | |
#returns the cosine of an angle at which an edge is visible from a given point | |
def cos(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round((pe1[0] * pe2[0] + pe1[1] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
#returns the sine of an angle at which an edge is visible from a given point | |
def sin(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round(abs(pe1[1] * pe2[0] - pe1[0] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
#returns the coefficients of a line (a * x + b * y + c = 0) through two given points | |
def line(p1, p2): | |
a = p1[1] - p2[1] | |
b = p2[0] - p1[0] | |
c = p1[0] * p2[1] - p2[0] * p1[1] | |
return [a, b, c] | |
points = [] | |
for i in range(7): | |
x = random.randint(-5, 5) | |
y = random.randint(-5, 5) | |
points.append((x, y)) | |
points = list(set(points)) | |
points.sort() | |
p0 = points[0] | |
edges = [] #queue with edges | |
triangulation = {} #dictionary with edges as keys, and points connected to them as values | |
#First step is to find a line throw p0 and another point from points with the greatest absolute value of slope possible. Every two adjacent points from points[] on the line give a new edge in triangulation. For n points there will be (n-1) edges. | |
max_abs_slope = [0] | |
for point in points[1:]: | |
if point[0] == p0[0]: #Since the slope of a vertical line is undefined as the denominator is equal to zero, it is represented as an infinite value. | |
slope = float('inf') | |
else: | |
slope = (point[1] - p0[1])/(point[0]-p0[0]) | |
if abs(slope) > abs(max_abs_slope[0]): | |
max_abs_slope = [slope, point] | |
elif slope == max_abs_slope[0]: | |
max_abs_slope.append(point) | |
max_abs_slope[0] = p0 | |
for i in range(1, len(max_abs_slope)): | |
new_edge = tuple([max_abs_slope[i-1], max_abs_slope[i]]) | |
edges.append(new_edge) | |
triangulation[new_edge] = [] | |
#Second step is to find the point from which an edge is visible at the greatest angle for every edge in edges[]. Since every edge in edges[] is already connected to a point (except for the edges added in the beginning of the algorithm), the new point is searched in a different half-plane. | |
while len(edges) > 0: | |
edge = edges.pop(0) | |
if len(triangulation[edge]) == 2: | |
continue | |
coeff = line(edge[0], edge[1]) | |
a, b, c = coeff[0], coeff[1], coeff[2] | |
if len(triangulation[edge]) == 0: | |
fixed_point = [edge[0][0] - 1, edge[0][1]] | |
else: | |
fixed_point = triangulation[edge][0] | |
fixed_dist = a * fixed_point[0] + b * fixed_point[1] + c | |
vertex = [2, []] | |
for point in points: | |
if (a * point[0] + b * point[1] + c) * fixed_dist < 0: #fixed_point and new point are in different half-planes in relation to the edge | |
cosine = cos(edge, point) | |
angle = cos((edge[1], point), edge[0]) | |
if cosine < vertex[0]: | |
vertex[0] = cosine | |
vertex[1] = [[angle, point]] | |
elif cosine == vertex[0]: | |
vertex[1].append([angle, point]) | |
if vertex[0] != 2: #vertex[0] == 2 in case of the absence of points in the half-plane | |
on_circle = vertex[1][::] #this list contains the vertices of the inscribed polygon except edge[0] and edge[1]. | |
on_circle.sort() | |
on_circle.append([1, edge[1]]) | |
triangulation[edge].append(on_circle[0][1]) | |
new_edge = tuple(sorted([edge[0], on_circle[0][1]])) | |
if new_edge not in triangulation.keys(): #new_edge is not connected to any points yet. | |
triangulation[new_edge] = [on_circle[1][1]] | |
edges.append(new_edge) #new_edge is added to the queue because one more point could still be connected to this edge. | |
elif on_circle[1][1] not in triangulation[new_edge]: #connect new_edge to on_circle[1][1], if this has not been done already. | |
triangulation[new_edge].append(on_circle[1][1]) #new_edge is not added to the queue because it is already connected to two points. | |
for i in range(len(on_circle)-1): #building outer edges of the inscribed polygon. | |
new_edge = tuple(sorted([on_circle[i][1], on_circle[i+1][1]])) | |
if new_edge not in triangulation.keys(): #the same check as for the previous new_edge | |
triangulation[new_edge] = [edge[0]] | |
edges.append(new_edge) | |
elif edge[0] not in triangulation[new_edge]: | |
triangulation[new_edge].append(edge[0]) | |
for i in range(1, len(on_circle)-1): #building inner edges of the inscribed polygon. | |
new_edge = tuple(sorted([edge[0], on_circle[i][1]])) | |
triangulation[new_edge] = [on_circle[i-1][1], on_circle[i+1][1]] | |
for key, val in triangulation.items(): | |
print(key, ":", val) | |
def rendering(points, triangulation, k): | |
speed(1000) | |
up() | |
color("red") | |
for point in points: | |
goto(point[0] * k, point[1] * k) | |
dot(5) | |
up() | |
color("black") | |
for key in triangulation.keys(): | |
goto(key[0][0] * k, key[0][1] * k) | |
down() | |
goto(key[1][0] * k, key[1][1] * k) | |
up() | |
exitonclick() | |
rendering(points, triangulation, 30) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import math | |
import random | |
import time | |
import matplotlib.pyplot as plt | |
def cos(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round((pe1[0] * pe2[0] + pe1[1] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
def sin(edge, point): | |
pe1 = [edge[0][0] - point[0], edge[0][1] - point[1]] | |
pe2 = [edge[1][0] - point[0], edge[1][1] - point[1]] | |
return round(abs(pe1[1] * pe2[0] - pe1[0] * pe2[1])/math.dist(edge[0], point)/math.dist(edge[1], point), 14) | |
def delaunay_correct(edge, p1, p2): | |
sin1 = sin(edge, p1) | |
cos1 = cos(edge, p1) | |
sin2 = sin(edge, p2) | |
cos2 = cos(edge, p2) | |
return sin1 * cos2 + sin2 * cos1 >= 0 | |
def line(p1, p2): | |
a = p1[1] - p2[1] | |
b = p2[0] - p1[0] | |
c = p1[0] * p2[1] - p2[0] * p1[1] | |
return [a, b, c] | |
def naive(points): | |
points = list(set(points)) | |
points.sort() | |
p0 = points[0] | |
edges = [] | |
triangulation = {} | |
max_abs_slope = [0] | |
for point in points[1:]: | |
if point[0] == p0[0]: | |
slope = float('inf') | |
else: | |
slope = (point[1] - p0[1])/(point[0]-p0[0]) | |
if abs(slope) > abs(max_abs_slope[0]): | |
max_abs_slope = [slope, point] | |
elif slope == max_abs_slope[0]: | |
max_abs_slope.append(point) | |
max_abs_slope[0] = p0 | |
for i in range(1, len(max_abs_slope)): | |
new_edge = tuple([max_abs_slope[i-1], max_abs_slope[i]]) | |
edges.append(new_edge) | |
triangulation[new_edge] = [] | |
while len(edges) > 0: | |
edge = edges.pop(0) | |
if len(triangulation[edge]) == 2: | |
continue | |
coeff = line(edge[0], edge[1]) | |
a, b, c = coeff[0], coeff[1], coeff[2] | |
if len(triangulation[edge]) == 0: | |
fixed_point = [edge[0][0] - 1, edge[0][1]] | |
else: | |
fixed_point = triangulation[edge][0] | |
fixed_dist = a * fixed_point[0] + b * fixed_point[1] + c | |
vertex = [2, []] | |
for point in points: | |
if (a * point[0] + b * point[1] + c) * fixed_dist < 0: | |
cosine = cos(edge, point) | |
angle = cos((edge[1], point), edge[0]) | |
if cosine < vertex[0]: | |
vertex[0] = cosine | |
vertex[1] = [[angle, point]] | |
elif cosine == vertex[0]: | |
vertex[1].append([angle, point]) | |
if vertex[0] != 2: | |
on_circle = vertex[1][::] | |
on_circle.sort() | |
on_circle.append([1, edge[1]]) | |
triangulation[edge].append(on_circle[0][1]) | |
new_edge = tuple(sorted([edge[0], on_circle[0][1]])) | |
if new_edge not in triangulation.keys(): | |
triangulation[new_edge] = [on_circle[1][1]] | |
edges.append(new_edge) | |
elif on_circle[1][1] not in triangulation[new_edge]: | |
triangulation[new_edge].append(on_circle[1][1]) | |
for i in range(len(on_circle)-1): | |
new_edge = tuple(sorted([on_circle[i][1], on_circle[i+1][1]])) | |
if new_edge not in triangulation.keys(): | |
triangulation[new_edge] = [edge[0]] | |
edges.append(new_edge) | |
elif edge[0] not in triangulation[new_edge]: | |
triangulation[new_edge].append(edge[0]) | |
for i in range(1, len(on_circle)-1): | |
new_edge = tuple(sorted([edge[0], on_circle[i][1]])) | |
triangulation[new_edge] = [on_circle[i-1][1], on_circle[i+1][1]] | |
return triangulation | |
Time = [] | |
for n in range(10, 501, 10): | |
Time.append(0) | |
for k in range(10): | |
points = [(random.randint(-1000, 1000), random.randint(-1000, 1000)) for i in range(n)] | |
start = time.time() | |
triangulation = naive(points) | |
end = time.time() | |
Time[-1] += end - start | |
fig, ax = plt.subplots() | |
result = ax.plot([n for n in range(10, 501, 10)], Time, label='naive') | |
ax.legend() | |
plt.xlabel('the number of points') | |
plt.ylabel('time (sec)') | |
plt.show() |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment