Chapter 16. Convex Hulls
Rovenski Vladimir, Haifa
The gift wrapping method.
> restart:
> convex:=proc(P) local i,Y,r,j,C,B,a,T; global L,m,n; n:=nops(P); L:=array(1..n+1); for i from 1 to n do L[i]:=P[i] od; Y:=op(2,L[1]); a:=1; for i from 2 to n do if op(2,L[i])<=Y then a:=i; Y:=op(2,L[i]) fi od; L[n+1]:=L[a]; B:=[1,0]; for j from 1 to n+1 do C:=-1; T:=L[a]; L[a]:=L[j]; L[j]:=T; if a=n+1 or j=n+1 then m:=j; j:=n+1 fi; for i from j+1 to n+1 do r:=L[i]-L[j]; if r[1]^2+r[2]^2>0 then Y:=evalf((r[1]*B[1]+r[2]*B[2])/(sqrt(r[1]^2+r[2]^2)*sqrt(B[1]^2+B[2]^2))); if Y>=C then C:=Y; a:=i fi fi od; B:=L[a]-L[j] od; plot([seq(L[i], i=1..m)], thickness=2) end:
> N:=37: P:=[seq([rand(1..40)(), rand(1..40)()], j=1..N)]:
> PP:=plot(P, style=point, symbol=circle, color=blue): plots[display]([PP, convex(P)]);
The set of layers of the finite planar set P (use the procedure convex(P) )
> G||0:=plot(P, style=point,symbol=circle,color=blue): G||1:=convex(P): for t from 2 while m<=n do G||t:=convex([seq(L[i], i=m+1..n+1)]); k:=t od: plots[display]([seq(G||i, i=0..k)]);