Obszar samo-przecinającego się wielokąta

32

Weźmy pod uwagę potencjalnie przecinający się wielokąt, zdefiniowany przez listę wierzchołków w przestrzeni 2D. Na przykład

{{0, 0}, {5, 0}, {5, 4}, {1, 4}, {1, 2}, {3, 2}, {3, 3}, {2, 3}, {2, 1}, {4, 1}, {4, 5}, {0, 5}}

Istnieje kilka sposobów definiowania obszaru takiego wielokąta, ale najciekawszym jest zasada parzysto-nieparzysta. Biorąc dowolny punkt na płaszczyźnie, narysuj linię od punktu do nieskończoności (w dowolnym kierunku). Jeśli linia ta przecina wielokąt nieparzystą liczbę razy, punkt jest częścią obszaru wielokąta, jeśli przecina wielokąt parzystą liczbę razy, punkt nie jest częścią wielokąta. W powyższym przykładzie wielokąta przedstawiono zarówno jego kontur, jak i jego parzysty nieparzysty obszar:

ZarysPowierzchnia

Wielobok na ogół nie będzie ortogonalny. Wybrałem tylko taki prosty przykład, aby ułatwić policzenie obszaru.

Obszar tego przykładu to 17(nie 24lub 33jak mogą dać inne definicje lub obszar).

Zauważ, że zgodnie z tą definicją obszar wielokąta jest niezależny od jego kolejności nawijania.

Wyzwanie

Biorąc pod uwagę listę wierzchołków o współrzędnych całkowitych definiujących wielokąt, określ jego obszar według reguły parzystej i nieparzystej.

Możesz napisać funkcję lub program, przyjmując dane wejściowe przez STDIN lub najbliższą alternatywę, argument wiersza poleceń lub argument funkcji i albo zwróć wynik, albo wydrukuj go do STDOUT lub najbliższej alternatywy.

Możesz pobierać dane wejściowe w dowolnym dogodnym formacie listy lub ciągu, o ile nie są one wstępnie przetwarzane.

Wynik powinien być albo liczbą zmiennoprzecinkową, z dokładnością do 6 cyfr znaczących (dziesiętnych), albo wynikiem racjonalnym, którego reprezentacja zmiennoprzecinkowa jest dokładna do 6 cyfr znaczących. (Jeśli uzyskasz racjonalne wyniki, prawdopodobnie będą one dokładne, ale nie mogę tego wymagać, ponieważ nie mam dokładnych wyników w celach informacyjnych).

Musisz być w stanie rozwiązać każdy z poniższych przypadków testowych w ciągu 10 sekund na rozsądnej maszynie stacjonarnej. (W tej regule jest pewna swoboda, więc dokonaj najlepszego osądu. Jeśli na moim laptopie zajmie to 20 sekund, dam ci wątpliwości, jeśli zajmie to chwilę, nie zrobię tego). Myślę, że ten limit powinno być bardzo hojne, ale powinno wykluczać podejścia, w których po prostu dyskrecjonujesz wielokąt na odpowiednio drobnej siatce i liczysz, lub używasz podejść probabilistycznych, takich jak Monte Carlo. Bądź dobrym sportowcem i nie próbuj optymalizować tych podejść, aby i tak dotrzymać terminu. ;)

Nie wolno używać żadnych istniejących funkcji związanych bezpośrednio z wielokątami.

To jest kod golfowy, więc wygrywa najkrótsze przesłanie (w bajtach).

Założenia

  • Wszystkie współrzędne są liczbami całkowitymi w zakresie 0 ≤ x ≤ 100, 0 ≤ y ≤ 100.
  • Będą przynajmniej 3i co najwyżej 50wierzchołki.
  • Nie będzie powtarzanych wierzchołków. Żaden wierzchołek nie będzie leżał na innej krawędzi. ( Jednak na liście mogą znajdować się punkty współliniowe.)

Przypadki testowe

{{0, 0}, {5, 0}, {5, 4}, {1, 4}, {1, 2}, {3, 2}, {3, 3}, {2, 3}, {2, 1}, {4, 1}, {4, 5}, {0, 5}}
17.0000

{{22, 87}, {6, 3}, {98, 77}, {20, 56}, {96, 52}, {79, 34}, {46, 78}, {52, 73}, {81, 85}, {90, 43}}
2788.39

{{90, 43}, {81, 85}, {52, 73}, {46, 78}, {79, 34}, {96, 52}, {20, 56}, {98, 77}, {6, 3}, {22, 87}}
2788.39

{{70, 33}, {53, 89}, {76, 35}, {14, 56}, {14, 47}, {59, 49}, {12, 32}, {22, 66}, {85, 2}, {2, 81},
 {61, 39}, {1, 49}, {91, 62}, {67, 7}, {19, 55}, {47, 44}, {8, 24}, {46, 18}, {63, 64}, {23, 30}}
2037.98

{{42, 65}, {14, 59}, {97, 10}, {13, 1}, {2, 8}, {88, 80}, {24, 36}, {95, 94}, {18, 9}, {66, 64},
 {91, 5}, {99, 25}, {6, 66}, {48, 55}, {83, 54}, {15, 65}, {10, 60}, {35, 86}, {44, 19}, {48, 43},
 {47, 86}, {29, 5}, {15, 45}, {75, 41}, {9, 9}, {23, 100}, {22, 82}, {34, 21}, {7, 34}, {54, 83}}
3382.46

{{68, 35}, {43, 63}, {66, 98}, {60, 56}, {57, 44}, {90, 52}, {36, 26}, {23, 55}, {66, 1}, {25, 6},
 {84, 65}, {38, 16}, {47, 31}, {44, 90}, {2, 30}, {87, 40}, {19, 51}, {75, 5}, {31, 94}, {85, 56},
 {95, 81}, {79, 80}, {82, 45}, {95, 10}, {27, 15}, {18, 70}, {24, 6}, {12, 73}, {10, 31}, {4, 29},
 {79, 93}, {45, 85}, {12, 10}, {89, 70}, {46, 5}, {56, 67}, {58, 59}, {92, 19}, {83, 49}, {22,77}}
3337.62

{{15, 22}, {71, 65}, {12, 35}, {30, 92}, {12, 92}, {97, 31}, {4, 32}, {39, 43}, {11, 40}, 
 {20, 15}, {71, 100}, {84, 76}, {51, 98}, {35, 94}, {46, 54}, {89, 49}, {28, 35}, {65, 42}, 
 {31, 41}, {48, 34}, {57, 46}, {14, 20}, {45, 28}, {82, 65}, {88, 78}, {55, 30}, {30, 27}, 
 {26, 47}, {51, 93}, {9, 95}, {56, 82}, {86, 56}, {46, 28}, {62, 70}, {98, 10}, {3, 39}, 
 {11, 34}, {17, 64}, {36, 42}, {52, 100}, {38, 11}, {83, 14}, {5, 17}, {72, 70}, {3, 97}, 
 {8, 94}, {64, 60}, {47, 25}, {99, 26}, {99, 69}}
3514.46
Martin Ender
źródło
1
W szczególności chciałbym zastąpić ograniczniki w taki sposób, aby lista była prawidłową ścieżką użytkownika PostScript, dzięki czemu mogę przeanalizować całość za pomocą jednego upathoperatora. (W rzeczywistości jest to bardzo prosta konwersja 1: 1 między separatorami. Po }, {prostu staje się lineto, przecinek między xiy jest usuwany, a nawiasy otwierające i zamykające są zastępowane statycznym nagłówkiem i stopką ...)
AJMansfield
1
@AJMansfield Zwykle nie mam nic przeciwko korzystaniu z wygodnych, natywnych reprezentacji list, ale używanie upathi linetobrzmi, jakbyś faktycznie przetwarzał dane wejściowe. Tj. Nie bierzesz listy współrzędnych, ale rzeczywisty wielokąt.
Martin Ender,
1
@MattNoonan Oh, to dobra uwaga. Tak, możesz to założyć.
Martin Ender
2
@ Promień Kierunek może wpływać na liczbę przejazdów, ale zawsze zwiększa się lub zmniejsza o 2, zachowując parzystość. Spróbuję znaleźć referencję. Na początek SVG używa tej samej definicji.
Martin Ender
1
Mathematica 12.0 ma nową wbudowaną funkcję dla tego: CrossingPolygon.
alephalpha

Odpowiedzi:

14

Mathematica, 247 225 222

p=Partition[#,2,1,1]&;{a_,b_}~r~{c_,d_}=Det/@{{a-c,c-d},{a,c}-b}/Det@{a-b,c-d};f=Abs@Tr@MapIndexed[Det@#(-1)^Tr@#2&,p[Join@@MapThread[{1-#,#}&/@#.#2&,{Sort/@Cases[{s_,t_}/;0<=s<=1&&0<=t<=1:>s]/@Outer[r,#,#,1],#}]&@p@#]]/2&

Najpierw dodaj punkty przecięcia do wielokąta, następnie odwróć niektóre krawędzie, a następnie jego pole można obliczyć jak zwykły wielokąt.

wprowadź opis zdjęcia tutaj

Przykład:

In[2]:= f[{{15, 22}, {71, 65}, {12, 35}, {30, 92}, {12, 92}, {97, 31}, {4, 32}, {39, 43}, {11, 40}, 
 {20, 15}, {71, 100}, {84, 76}, {51, 98}, {35, 94}, {46, 54}, {89, 49}, {28, 35}, {65, 42}, 
 {31, 41}, {48, 34}, {57, 46}, {14, 20}, {45, 28}, {82, 65}, {88, 78}, {55, 30}, {30, 27}, 
 {26, 47}, {51, 93}, {9, 95}, {56, 82}, {86, 56}, {46, 28}, {62, 70}, {98, 10}, {3, 39}, 
 {11, 34}, {17, 64}, {36, 42}, {52, 100}, {38, 11}, {83, 14}, {5, 17}, {72, 70}, {3, 97}, 
 {8, 94}, {64, 60}, {47, 25}, {99, 26}, {99, 69}}]

Out[2]= 3387239559852305316061173112486233884246606945138074528363622677708164\
 6419838924305735780894917246019722157041758816629529815853144003636562\
 9161985438389053702901286180223793349646170997160308182712593965484705\
 3835036745220226127640955614326918918917441670126958689133216326862597\
 0109115619/\
 9638019709367685232385259132839493819254557312303005906194701440047547\
 1858644412915045826470099500628074171987058850811809594585138874868123\
 9385516082170539979030155851141050766098510400285425157652696115518756\
 3100504682294718279622934291498595327654955812053471272558217892957057\
 556160

In[3]:= N[%] (*The numerical value of the last output*)

Out[3]= 3514.46
alephalpha
źródło
Niestety nie jestem pewien, czy ta logika zadziała we wszystkich sytuacjach. Można spróbować {1,2},{4,4},{4,2},{2,4},{2,1},{5,3}? Powinieneś wyjść z 3.433333333333309. Patrzyłem na użycie podobnej logiki.
MickyT,
@MickyT Tak, działa. Zwrócił 103/30, a wartość liczbowa to 3.43333.
alephalpha
Przepraszam za to. Dobre rozwiązanie
MickyT,
44

Python 2, 323 319 bajtów

exec u"def I(s,a,b=1j):c,d=s;d-=c;c-=a;e=(d*bX;return e*(0<=(b*cX*e<=e*e)and[a+(d*cX*b/e]or[]\nE=lambda p:zip(p,p[1:]+p);S=sorted;P=E(input());print sum((t-b)*(r-l)/2Fl,r@E(S(i.realFa,b@PFe@PFi@I(e,a,b-a)))[:-1]Fb,t@E(S(((i+j)XFe@PFi@I(e,l)Fj@I(e,r)))[::2])".translate({70:u" for ",64:u" in ",88:u".conjugate()).imag"})

Pobiera listę wierzchołków przez STDIN jako liczby zespolone, w następującej formie

[  X + Yj,  X + Yj,  ...  ]

i zapisuje wynik w STDOUT.

Ten sam kod po zamianie łańcucha i niektórych odstępach:

def I(s, a, b = 1j):
    c, d = s; d -= c; c -= a;
    e = (d*b.conjugate()).imag;
    return e * (0 <= (b*c.conjugate()).imag * e <= e*e) and \
           [a + (d*c.conjugate()).imag * b/e] or []

E = lambda p: zip(p, p[1:] + p);
S = sorted;

P = E(input());

print sum(
    (t - b) * (r - l) / 2

    for l, r in E(S(
        i.real for a, b in P for e in P for i in I(e, a, b - a)
    ))[:-1]

    for b, t in E(S(
        ((i + j).conjugate()).imag for e in P for i in I(e, l) for j in I(e, r)
    ))[::2]
)

Wyjaśnienie

Dla każdego punktu przecięcia dwóch boków wielokąta wejściowego (w tym wierzchołków), przeprowadź linię pionową przez ten punkt.

Rycina 1

(W rzeczywistości, z powodu gry w golfa, program przechodzi kilka kolejnych linii; to nie ma znaczenia, o ile mijamy przynajmniej te linie.) Ciało wielokąta między dowolnymi dwiema kolejnymi liniami składa się z pionowych trapezoidów ( oraz trójkąty i segmenty linii, jako ich szczególne przypadki). Tak musi być, ponieważ jeśli którykolwiek z tych kształtów miałby dodatkowy wierzchołek między dwiema podstawami, między tymi dwoma liniami byłaby kolejna pionowa linia. Suma powierzchni wszystkich takich trapezoidów to powierzchnia wielokąta.

Oto jak odnajdujemy te trapezoidy: Dla każdej pary kolejnych pionowych linii znajdujemy segmenty każdej strony wielokąta, które (prawidłowo) leżą między tymi dwiema liniami (które mogą nie istnieć dla niektórych boków). Na powyższej ilustracji są to sześć czerwonych segmentów, biorąc pod uwagę dwie czerwone pionowe linie. Zauważ, że te segmenty nie przecinają się prawidłowo (tzn. Mogą się spotkać tylko w punktach końcowych, całkowicie się pokrywają lub wcale nie przecinają się, ponieważ, ponownie, jeśli prawidłowo się przecinają, pomiędzy nimi będzie inna linia pionowa;) dlatego warto rozmawiać o zamawianiu ich od góry do dołu, co robimy. Zgodnie z zasadą parzysto-nieparzystą, gdy przekroczymy pierwszy segment, jesteśmy wewnątrz wielokąta; kiedy przekroczymy drugi, jesteśmy poza; trzeci znowu; czwarte wyjście; i tak dalej...

Ogólnie jest to algorytm O ( n 3 log n ).

Łokieć
źródło
4
To jest genialne! Wiedziałem, że mogę na ciebie liczyć. ;) (Być może zechcesz odpowiedzieć na to pytanie przy przepełnieniu stosu.)
Martin Ender
@ MartinBüttner Keep 'em come :)
Ell
7
Świetna robota i świetne wytłumaczenie
MickyT,
1
To imponująca odpowiedź. Czy sam opracowałeś algorytm, czy jest już praca nad tym problemem? Jeśli istnieje praca, byłbym wdzięczny za wskazówkę, gdzie mogę ją znaleźć. Nie miałem pojęcia, jak sobie z tym poradzić.
Logic Knight
5
@CarpetPython Opracowałem go sam, ale byłbym bardzo zaskoczony, gdyby nie zostało to zrobione wcześniej.
Ell,
9

Haskell, 549

Nie wygląda na to, bym mógł zagrać w tę grę wystarczająco daleko, ale koncepcja była inna niż dwie pozostałe odpowiedzi, więc pomyślałem, że i tak podzielę się tym. Wykonuje operacje wymierne O (N ^ 2) w celu obliczenia obszaru.

import Data.List
_%0=2;x%y=x/y
h=sort
z f w@(x:y)=zipWith f(y++[x])w
a=(%2).sum.z(#);(a,b)#(c,d)=b*c-a*d
(r,p)?(s,q)=[(0,p)|p==q]++[(t,v t p r)|u t,u$f r]where f x=(d q p#x)%(r#s);t=f s;u x=x^2<x
v t(x,y)(a,b)=(x+t*a,y+t*b);d=v(-1)
s x=zip(z d x)x
i y=h.(=<<y).(?)=<<y
[]!x=[x];x!_=x
e n(a@(x,p):y)|x>0=(n!y,a):(e(n!y)$tail$dropWhile((/=p).snd)y)|0<1=(n,a):e n y
c[p]k=w x[]where((_,q):x)=e[]p;w((n,y):z)b|q==y=(k,map snd(q:b)):c n(-k)|0<1=w z(y:b);c[]_=[]
b(s,p)=s*a p
u(_,x)(_,y)=h x==h y
f p=abs$sum$map b$nubBy u$take(length p^2)$c[cycle$i$s p]1

Przykład:

λ> f test''
33872395598523053160611731124862338842466069451380745283636226777081646419838924305735780894917246019722157041758816629529815853144003636562916198543838905370290128618022379334964617099716030818271259396548470538350367452202261276409556143269189189174416701269586891332163268625970109115619 % 9638019709367685232385259132839493819254557312303005906194701440047547185864441291504582647009950062807417198705885081180959458513887486812393855160821705399790301558511410507660985104002854251576526961155187563100504682294718279622934291498595327654955812053471272558217892957057556160
λ> fromRational (f test'')
3514.4559380388832

Pomysł polega na ponownym okablowaniu wielokąta przy każdym skrzyżowaniu, w wyniku czego powstanie połączenie wielokątów bez przecinających się krawędzi. Następnie możemy obliczyć (podpisany) obszar każdego z wielokątów za pomocą wzoru sznurowadła Gaussa ( http://en.wikipedia.org/wiki/Shoelace_formula ). Zasada parzystej nieparzystości wymaga, aby po przekształceniu skrzyżowania pole nowego wielokąta było liczone ujemnie w stosunku do starego wielokąta.

Rozważmy na przykład wielokąt w pierwotnym pytaniu. Skrzyżowanie w lewym górnym rogu jest konwertowane na dwie ścieżki, które spotykają się tylko w jednym punkcie; obie ścieżki są zorientowane zgodnie z ruchem wskazówek zegara, więc obszary każdej z nich byłyby dodatnie, chyba że zadeklarujemy, że ścieżka wewnętrzna jest ważona o -1 w stosunku do ścieżki zewnętrznej. Jest to równoważne odwróceniu ścieżki alphaalpha.

Wieloboki pochodzące z oryginalnego przykładu

Jako kolejny przykład rozważ wielokąt z komentarza MickyT:

Wieloboki pochodzą z komentarza MickyT

Tutaj niektóre wielokąty są zorientowane zgodnie z ruchem wskazówek zegara, a niektóre przeciwnie do ruchu wskazówek zegara. Reguła „przerzucanie na krzyż” gwarantuje, że regiony zorientowane zgodnie z ruchem wskazówek zegara podniosą dodatkowy współczynnik -1, powodując, że wniosą one pozytywną wartość do tego obszaru.

Oto jak działa program:

import Data.List  -- for sort and nubBy

-- Rational division, with the unusual convention that x/0 = 2
_%0=2;x%y=x/y

-- Golf
h=sort

-- Define a "cyclic zipWith" operation. Given a list [a,b,c,...z] and a binary
-- operation (@), z (@) [a,b,c,..z] computes the list [b@a, c@b, ..., z@y, a@z]
z f w@(x:y)=zipWith f(y++[x])w

-- The shoelace formula for the signed area of a polygon
a=(%2).sum.z(#)

-- The "cross-product" of two 2d vectors, resulting in a scalar.
(a,b)#(c,d)=b*c-a*d

-- Determine if the line segment from p to p+r intersects the segment from
-- q to q+s.  Evaluates to the singleton list [(t,x)] where p + tr = x is the
-- point of intersection, or the empty list if there is no intersection. 
(r,p)?(s,q)=[(0,p)|p==q]++[(t,v t p r)|u t,u$f r]where f x=(d q p#x)%(r#s);t=f s;u x=x^2<x

-- v computes an affine combination of two vectors; d computes the difference
-- of two vectors.
v t(x,y)(a,b)=(x+t*a,y+t*b);d=v(-1)

-- If x is a list of points describing a polygon, s x will be the list of
-- (displacement, point) pairs describing the edges.
s x=zip(z d x)x

-- Given a list of (displacement, point) pairs describing a polygon's edges,
-- create a new polygon which also has a vertex at every point of intersection.
-- Mercilessly golfed.
i y=h.(=<<y).(?)=<<y


-- Extract a simple polygon; when an intersection point is reached, fast-forward
-- through the polygon until we return to the same point, then continue.  This
-- implements the edge rewiring operation. Also keep track of the first
-- intersection point we saw, so that we can process that polygon next and with
-- opposite sign.
[]!x=[x];x!_=x
e n(a@(x,p):y)|x>0=(n!y,a):(e(n!y)$tail$dropWhile((/=p).snd)y)|0<1=(n,a):e n y

-- Traverse the polygon from some arbitrary starting point, using e to extract
-- simple polygons marked with +/-1 weights.
c[p]k=w x[]where((_,q):x)=e[]p;w((n,y):z)b|q==y=(k,map snd(q:b)):c n(-k)|0<1=w z(y:b);c[]_=[]

-- If the original polygon had N vertices, there could (very conservatively)
-- be up to N^2 points of intersection.  So extract N^2 polygons using c,
-- throwing away duplicates, and add up the weighted areas of each polygon.
b(s,p)=s*a p
u(_,x)(_,y)=h x==h y
f p=abs$sum$map b$nubBy u$take(length p^2)$c[cycle$i$s p]1
Matt Noonan
źródło