Project Euler 252 - inamori’s diary

Q252.
与えられた疑似乱数で500点を発生させる。この500点のうち何点かで構成されて、かつ構成点以外を含まない凸多角形のうち最も大きい面積を求めよ。

最初、どうしていいか分からなかったが、地道にやってもなんとか解ける。
内部に他の点を含まない三角形を生成する。その三角形と三角形を併せて四角形を作る。四角形と三角形を併せて五角形を作る。とやっていけばよい。最初の三角形を生成するのがとても遅い。

やっと最後にたどり着いた。まだ20問も解けずに残っている問題があるが、ここにはもう新しい問題しか載せないことにする。



from itertools import imap
def gen_random():

    s = 290797

    while True:

        s = int(s * s % 50515093)

        yield s % 2000 - 1000
def make_points(n):

    a = [ ]

    g = gen_random()

    for k in range(n):

        a.append( (g.next(), g.next()) )

    return a
def has_points_in_triangle(i, j, k, a11, a12, a21, a22, points, sr):

    pt1 = points[i]

    det = a11 * a22 - a21 * a12

    if det == 0:

        return True

    

    for l in xrange(N):

        if l != i and l != j and l != k:

            pt = points[l]

            dx = pt[0] - pt1[0]

            dy = pt[1] - pt1[1]

            s =  a22 * dx - a12 * dy

            t = -a21 * dx + a11 * dy

            if det > 0:

                if s > 0 and t > 0 and s + t < det:

                    return True

                elif s < 0 and s + t > det:

                    sr.add(l)

            else:

                if s < 0 and t < 0 and s + t > det:

                    return True

                elif s > 0 and s + t < det:

                    sr.add(l)

    return False
def make_triangles(points):

    a = [ ]

    for i in xrange(N):

        pt1 = points[i]

        for j in xrange(i + 1, N):

            pt2 = points[j]

            a11 = pt2[0] - pt1[0]

            a21 = pt2[1] - pt1[1]

            s_remove = set()

            for k in xrange(j + 1, N):

                if k in s_remove:

                    continue

                pt3 = points[k]

                a12 = pt3[0] - pt1[0]

                a22 = pt3[1] - pt1[1]

                if not has_points_in_triangle(i, j, k,

                            a11, a12, a21, a22, points, s_remove):

                    area = a11 * a22 - a12 * a21

                    if area > 0:

                        t = (i, j, k, area)

                    else:

                        t = (i, k, j, -area)

                    a.append(t)

    return a
def gen_triangle_with_edge(v1, v2, triangles):

    edge = (v2, v1)

    if edge in dic_triangles:

        for t in dic_triangles[edge]:

            if t[0] == v2:

                yield (t[2], t[3])

            elif t[1] == v2:

                yield (t[0], t[3])

            else:

                yield (t[1], t[3])
def is_counterclockwise(pt0, pt1, pt2):

    ax = pt1[0] - pt0[0]

    ay = pt1[1] - pt0[1]

    bx = pt2[0] - pt1[0]

    by = pt2[1] - pt1[1]

    return ax * by - ay * bx > 0
def is_convex(p):

    n = len(p) - 1

    pt0 = points[p[-4] ]

    pt1 = points[p[-3] ]

    pt2 = points[p[-2] ]

    pt3 = points[p[0] ]

    pt4 = points[p[1] ]

    return (is_counterclockwise(pt0, pt1, pt2)

        and is_counterclockwise(pt2, pt3, pt4))
def make_polygons(polygons, triangles):

    p = [ ]

    for t in polygons:

        for v, area in gen_triangle_with_edge(t[-2], t[0], triangles):

            if v > t[0]:

                q = t[:-1] + (v, ) + (t[-1] + area, )

                if is_convex(q):

                    p.append(q)

    return p
def edge_to_triangles(triangles):

    d = { }

    def add(k, v):

        if k in d:

            d[k].append(v)

        else:

            d[k] = [ v ]

    

    for t in triangles:

        add( (t[0], t[1]), t )

        add( (t[1], t[2]), t )

        add( (t[2], t[0]), t )

    

    return d
N = 500

area = 0

points = make_points(N)

triangles = make_triangles(points)

dic_triangles = edge_to_triangles(triangles)

print len(triangles)

polygons = triangles

while len(polygons):

    area = reduce(max, imap(lambda e: e[-1], polygons), area)

    polygons = make_polygons(polygons, triangles)

    print len(polygons)

print area / 2.0