Convex hull: Difference between revisions
Added Easylang
(Added Easylang) |
|||
| (34 intermediate revisions by 4 users not shown) | |||
Line 13:
{{trans|Nim}}
<
V val = (q.y - p.y) * (r.x - q.x) - (q.x - p.x) * (r.y - q.y)
I val == 0
Line 71:
V hull = calculateConvexHull(points)
L(i) hull
print(i)</
{{out}}
Line 85:
=={{header|Action!}}==
<
TYPE PointI=[INT x,y]
Line 186:
PrintE("Convex hull:")
PrintPoints(result,rLen)
RETURN</
{{out}}
[https://gitlab.com/amarok8bit/action-rosetta-code/-/raw/master/images/Convex_hull.png Screenshot from Atari 8-bit computer]
Line 199:
=={{header|Ada}}==
{{trans|D}}
<
with Ada.Containers.Vectors;
Line 286:
Show (Find_Convex_Hull (Vec));
Ada.Text_IO.New_Line;
end Convex_Hull;</
{{out}}
Line 296:
{{trans|Rust}}
<
orientation: function [P Q R][
Line 363:
hull: calculateConvexHull points
loop hull 'i ->
print i</
{{out}}
Line 374:
[x:5 y:19]
[x:-3 y:15]</pre>
=={{header|ATS}}==
{{trans|Standard ML}}
I purposely demonstrate various features, although I personally do not favor the syntax of the '''.x''' and '''.y''' overloads. (I have sometimes had this syntax get confused with record field syntax.)
<syntaxhighlight lang="ats">//
// Convex hulls by Andrew's monotone chain algorithm.
//
// For a description of the algorithm, see
// https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
//
//
// The implementation is designed for use with a garbage collector.
// In other words, some of the memory is allowed to leak.
//
// Sometimes, where it is slightly easier if one does so, we do try to
// free memory. Boehm GC sometimes cannot free allocated memory
// itself, so perhaps it is just as well if an ATS program explicitly
// frees some of its memory.
//
#include "share/atspre_staload.hats"
#define NIL list_nil ()
#define :: list_cons
(* Make the floating point type big, so use of a boxed representation
for points makes some sense. *)
typedef floatpt = ldouble
extern castfn int2floatpt : int -<> floatpt
overload i2fp with int2floatpt
(*------------------------------------------------------------------*)
//
// Enough plane geometry for our purpose.
//
(* Let us make the point type boxed. Here is one way to do that. The
name of the type will be "plane_point_t". The constructor will be
named "plane_point". *)
datatype plane_point_t =
| plane_point of (floatpt, floatpt)
fn {}
plane_point_x (p : plane_point_t) :<> floatpt =
let
val+ plane_point (x, _) = p
in
x
end
fn {}
plane_point_y (p : plane_point_t) :<> floatpt =
let
val+ plane_point (_, y) = p
in
y
end
fn {}
plane_point_eq (p : plane_point_t,
q : plane_point_t) :<> bool =
let
val+ plane_point (xp, yp) = p
val+ plane_point (xq, yq) = q
in
xp = xq && yp = yq
end
(* Impose a total order on points, making it one that will work for
Andrew's monotone chain algorithm. *)
fn {}
plane_point_order (p : plane_point_t,
q : plane_point_t) :<> bool =
let
val+ plane_point (xp, yp) = p
val+ plane_point (xq, yq) = q
in
xp < xq || (xp = xq && yp < yq)
end
(* Subtraction is really a vector or multivector operation. *)
fn {}
plane_point_sub (p : plane_point_t,
q : plane_point_t) :<> plane_point_t =
let
val+ plane_point (xp, yp) = p
val+ plane_point (xq, yq) = q
in
plane_point (xp - xq, yp - yq)
end
(* Cross product is really a multivector operation. *)
fn {}
plane_point_cross (p : plane_point_t,
q : plane_point_t) :<> floatpt =
let
val+ plane_point (xp, yp) = p
val+ plane_point (xq, yq) = q
in
(xp * yq) - (yp * xq)
end
overload .x with plane_point_x
overload .y with plane_point_y
overload = with plane_point_eq
overload order with plane_point_order
overload < with plane_point_order
overload - with plane_point_sub
(* Make "cross" a left-associative infix operator with the same
precedence as "*". *)
infixl ( * ) cross
overload cross with plane_point_cross
(*------------------------------------------------------------------*)
//
// Sorting an array of points.
//
fn
plane_point_array_sort
{n : int}
(arr : &(@[plane_point_t][n]) >> _,
n : size_t n) :<!wrt> void =
let
(* The comparison will be inlined by template expansion. *)
implement
array_quicksort$cmp<plane_point_t> (p, q) =
if order (p, q) then (* An overload for plane_point_order. *)
~1
else if q < p then (* Another overload for plane_point_order. *)
1
else
0
in
(* The following sort routine is in the ATS2 prelude. *)
array_quicksort<plane_point_t> (arr, n)
end
(*------------------------------------------------------------------*)
//
// Removing duplicates from a sorted array. Returns a new array.
//
extern fun {a : t@ype}
array_delete_neighbor_dups$eq (p : a, q : a) :<> bool
fn {a : t@ype}
array_delete_neighbor_dups
{n : int}
(arr : &(@[a][n]),
n : size_t n)
:<!wrt> [m : nat | m <= n]
[parr1 : addr]
@(@[a][m] @ parr1,
mfree_gc_v parr1 |
ptr parr1,
size_t m) =
let
macdef eq = array_delete_neighbor_dups$eq<a>
fn
nondups_list {n : int}
(arr : &(@[a][n]),
n : size_t n)
:<> [m : nat | m <= n] list (a, m) =
let
fun
loop {i : nat | i < n}
{k : nat | k < n - i}
.<i>. (* <-- proof of termination. *)
(arr : &(@[a][n]),
i : size_t i,
lst : list (a, k))
:<> [m : nat | m <= n] list (a, m) =
(* Cons a list of non-duplicates, going backwards through
the array so the list will be in forwards order. (The
order does not really matter in ATS, though, because in
the prelude there are both array_initize_list *and*
array_initize_rlist. *)
if i = i2sz 0 then
arr[i] :: lst
(* The "\" in the following line makes eq temporarily
infix. *)
else if arr[i - 1] \eq arr[i] then
loop (arr, pred i, lst)
else
loop (arr, pred i, arr[i] :: lst)
prval () = lemma_array_param arr
in
if n = i2sz 0 then
NIL
else
loop (arr, pred n, NIL)
end
val lst = nondups_list (arr, n)
prval () = lemma_list_param lst
val m = i2sz (length lst)
val @(pfarr1, pfgc1 | parr1) = array_ptr_alloc<a> (m)
val () = array_initize_list<a> (!parr1, sz2i m, lst)
in
@(pfarr1, pfgc1 | parr1, m)
end
(*------------------------------------------------------------------*)
//
// Removing duplicates from a sorted plane_point_t array. Returns a
// new array.
//
fn
plane_point_array_delete_neighbor_dups
{n : int}
(arr : &(@[plane_point_t][n]),
n : size_t n)
:<!wrt> [m : nat | m <= n]
[parr1 : addr]
@(@[plane_point_t][m] @ parr1,
mfree_gc_v parr1 |
ptr parr1,
size_t m) =
let
implement
array_delete_neighbor_dups$eq<plane_point_t> (p, q) =
p = q (* Here = is an overload for plane_point_eq. *)
in
array_delete_neighbor_dups<plane_point_t> (arr, n)
end
(*------------------------------------------------------------------*)
//
// The convex hull algorithm.
//
fn
cross_test {m, k : int | 1 <= k; k < m}
(pt_i : plane_point_t,
hull : &(@[plane_point_t][m]),
k : size_t k) :<> bool =
let
val hull_k = hull[k]
and hull_k1 = hull[k - 1]
in
i2fp 0 < (hull_k - hull_k1) cross (pt_i - hull_k1)
end
fn
construct_lower_hull {n : int | 2 <= n}
(pt : &(@[plane_point_t][n]),
n : size_t n)
:<!wrt> [m : int | 2 <= m; m <= n]
[phull : addr]
@(@[plane_point_t][n] @ phull,
mfree_gc_v phull |
ptr phull,
size_t m) =
let
val @(pfhull, pfgc | phull) = array_ptr_alloc<plane_point_t> n
(* It is easier to work with an array if it is fully
initialized. (Yes, there are also ways to cheat and so merely
pretend the array has been initialized.) *)
val arbitrary_point = plane_point (i2fp 0, i2fp 0)
val () = array_initize_elt<plane_point_t> (!phull, n,
arbitrary_point)
(* The macro "hull" can be used with index notation "[]". *)
macdef hull = !phull
val () = hull[0] := pt[0]
val () = hull[1] := pt[1]
fun
outer_loop {i : int | 0 <= i; i <= n}
{j : int | 1 <= j; j < n}
.<n - i>. (* <-- proof of termination. *)
(pt : &(@[plane_point_t][n]),
hull : &(@[plane_point_t][n]),
i : size_t i,
j : size_t j)
:<!wrt> [m : int | 2 <= m; m <= n] size_t m =
if i = n then
succ j
else
let
val pt_i = pt[i]
fun
inner_loop {k : int | 0 <= k; k < n - 1}
.<k>. (* <-- proof of termination. *)
(hull : &(@[plane_point_t][n]),
k : size_t k)
:<!wrt> [j : int | 1 <= j; j < n] size_t j =
if k = i2sz 0 then
begin
hull[succ k] := pt_i;
succ k
end
else if cross_test (pt_i, hull, k) then
begin
hull[succ k] := pt_i;
succ k
end
else
inner_loop (hull, pred k)
(* I do not know how to write a proof of the following, so
let us test it at runtime. *)
val () = $effmask_exn assertloc (j < n - 1)
in
outer_loop (pt, hull, succ i, inner_loop (hull, j))
end
val hull_size = outer_loop (pt, hull, i2sz 2, i2sz 1)
in
@(pfhull, pfgc | phull, hull_size)
end
fn
construct_upper_hull {n : int | 2 <= n}
(pt : &(@[plane_point_t][n]),
n : size_t n)
:<!wrt> [m : int | 2 <= m; m <= n]
[phull : addr]
@(@[plane_point_t][n] @ phull,
mfree_gc_v phull |
ptr phull,
size_t m) =
let
val @(pfhull, pfgc | phull) = array_ptr_alloc<plane_point_t> n
(* It is easier to work with an array if it is fully
initialized. (Yes, there are also ways to cheat and so merely
pretend the array has been initialized.) *)
val arbitrary_point = plane_point (i2fp 0, i2fp 0)
val () = array_initize_elt<plane_point_t> (!phull, n,
arbitrary_point)
(* The macro "hull" can be used with index notation "[]". *)
macdef hull = !phull
val () = hull[0] := pt[n - 1]
val () = hull[1] := pt[n - 2]
fun
outer_loop {i1 : int | 0 <= i1; i1 <= n - 2}
{j : int | 1 <= j; j < n}
.<i1>. (* <-- proof of termination. *)
(pt : &(@[plane_point_t][n]),
hull : &(@[plane_point_t][n]),
i1 : size_t i1,
j : size_t j)
:<!wrt> [m : int | 2 <= m; m <= n] size_t m =
if i1 = i2sz 0 then
succ j
else
let
val i = pred i1
val pt_i = pt[i]
fun
inner_loop {k : int | 0 <= k; k < n - 1}
.<k>. (* <-- proof of termination. *)
(hull : &(@[plane_point_t][n]),
k : size_t k)
:<!wrt> [j : int | 1 <= j; j < n] size_t j =
if k = i2sz 0 then
begin
hull[succ k] := pt_i;
succ k
end
else if cross_test (pt_i, hull, k) then
begin
hull[succ k] := pt_i;
succ k
end
else
inner_loop (hull, pred k)
(* I do not know how to write a proof of the following, so
let us test it at runtime. *)
val () = $effmask_exn assertloc (j < n - 1)
in
outer_loop (pt, hull, pred i1, inner_loop (hull, j))
end
val hull_size = outer_loop (pt, hull, n - i2sz 2, i2sz 1)
in
@(pfhull, pfgc | phull, hull_size)
end
fn
construct_hull {n : int | 2 <= n}
(pt : &(@[plane_point_t][n]),
n : size_t n)
:<!wrt> [hull_size : int | 2 <= hull_size]
[phull : addr]
@(@[plane_point_t][hull_size] @ phull,
mfree_gc_v phull |
ptr phull,
size_t hull_size) =
(* The following implementation demonstrates slightly complicated
"safe" initialization. By "safe" I mean so one never *reads* from
an uninitialized entry. Elsewhere in the program I did this more
simply, with prelude routines.
I also demonstrate freeing of a linear object. If you remove the
calls to array_ptr_free, you cannot compile the program. *)
let
(* Side note: Construction of the lower and upper hulls can be
done in parallel. *)
val [lower_hull_size : int]
[p_lower_hull : addr]
@(pf_lower_hull, pfgc_lower | p_lower_hull, lower_hull_size) =
construct_lower_hull (pt, n)
and [upper_hull_size : int]
[p_upper_hull : addr]
@(pf_upper_hull, pfgc_upper | p_upper_hull, upper_hull_size) =
construct_upper_hull (pt, n)
stadef hull_size = lower_hull_size + upper_hull_size - 2
val hull_size : size_t hull_size =
lower_hull_size + upper_hull_size - i2sz 2
val [phull : addr] @(pfhull, pfgc_hull | phull) =
array_ptr_alloc<plane_point_t> hull_size
(* Split off the part of each partial hull's view that we actually
will use. *)
prval @(pf_lower, pf_lower_rest) =
array_v_split {plane_point_t} {p_lower_hull} {n}
{lower_hull_size - 1}
pf_lower_hull
prval @(pf_upper, pf_upper_rest) =
array_v_split {plane_point_t} {p_upper_hull} {n}
{upper_hull_size - 1}
pf_upper_hull
(* Split the new array's view in two. The question mark means
that the array elements are uninitialized. *)
prval @(pfleft, pfright) =
array_v_split {plane_point_t?} {phull} {hull_size}
{lower_hull_size - 1}
pfhull
(* Copy the lower hull, thus initializing the left side of the new
array. *)
val () = array_copy<plane_point_t> (!phull, !p_lower_hull,
pred lower_hull_size)
(* Copy the upper hull, thus initializing the left side of the new
array. *)
val phull_right =
ptr_add<plane_point_t> (phull, pred lower_hull_size)
val () = array_copy<plane_point_t> (!phull_right, !p_upper_hull,
pred upper_hull_size)
(* Join the views of the initialized halves. *)
prval pfhull = array_v_unsplit (pfleft, pfright)
(* Restore the views of the partial hulls. *)
prval () = pf_lower_hull :=
array_v_unsplit (pf_lower, pf_lower_rest)
prval () = pf_upper_hull :=
array_v_unsplit (pf_upper, pf_upper_rest)
(* We do not need the lower and upper hulls anymore, and so can
free them. (Of course, if there is a garbage collector, you
could make freeing be a no-operation.) *)
val () = array_ptr_free (pf_lower_hull, pfgc_lower | p_lower_hull)
val () = array_ptr_free (pf_upper_hull, pfgc_upper | p_upper_hull)
in
@(pfhull, pfgc_hull | phull, hull_size)
end
fn
plane_convex_hull {num_points : int}
(points_lst : list (plane_point_t, num_points))
:<!wrt> [hull_size : int | 0 <= hull_size]
[phull : addr]
@(@[plane_point_t][hull_size] @ phull,
mfree_gc_v phull |
ptr phull,
size_t hull_size) =
(* Takes an arbitrary list of points, which may be in any order and
may contain duplicates. Returns an ordered array of points that
make up the convex hull. If the initial list of points is empty,
the returned array is empty. *)
let
prval () = lemma_list_param points_lst
val num_points = i2sz (length points_lst)
(* Copy the list to an array. *)
val @(pf_points, pfgc_points | p_points) =
array_ptr_alloc<plane_point_t> num_points
val () =
array_initize_list<plane_point_t> (!p_points, sz2i num_points,
points_lst)
(* Sort the array. *)
val () = plane_point_array_sort (!p_points, num_points)
(* Create a new sorted array that has the duplicates removed. *)
val @(pf_pt, pfgc_pt | p_pt, n) =
plane_point_array_delete_neighbor_dups (!p_points, num_points)
(* The original array no longer is needed. *)
val () = array_ptr_free (pf_points, pfgc_points | p_points)
in
if n <= 2 then
@(pf_pt, pfgc_pt | p_pt, n)
else
let
val @(pfhull, pfgc_hull | phull, hull_size) =
construct_hull (!p_pt, n)
val () = array_ptr_free (pf_pt, pfgc_pt | p_pt)
in
@(pfhull, pfgc_hull | phull, hull_size)
end
end
(*------------------------------------------------------------------*)
implement
main0 () =
{
val example_points =
$list (plane_point (i2fp 16, i2fp 3),
plane_point (i2fp 12, i2fp 17),
plane_point (i2fp 0, i2fp 6),
plane_point (i2fp ~4, i2fp ~6),
plane_point (i2fp 16, i2fp 6),
plane_point (i2fp 16, i2fp ~7),
plane_point (i2fp 16, i2fp ~3),
plane_point (i2fp 17, i2fp ~4),
plane_point (i2fp 5, i2fp 19),
plane_point (i2fp 19, i2fp ~8),
plane_point (i2fp 3, i2fp 16),
plane_point (i2fp 12, i2fp 13),
plane_point (i2fp 3, i2fp ~4),
plane_point (i2fp 17, i2fp 5),
plane_point (i2fp ~3, i2fp 15),
plane_point (i2fp ~3, i2fp ~9),
plane_point (i2fp 0, i2fp 11),
plane_point (i2fp ~9, i2fp ~3),
plane_point (i2fp ~4, i2fp ~2),
plane_point (i2fp 12, i2fp 10))
val (pf_hull, pfgc_hull | p_hull, hull_size) =
plane_convex_hull example_points
macdef hull = !p_hull
val () =
let
var i : [i : nat] size_t i
in
for (i := i2sz 0; i < hull_size; i := succ i)
println! ("(", hull[i].x(), " ", hull[i].y(), ")")
end
val () = array_ptr_free (pf_hull, pfgc_hull | p_hull)
}
(*------------------------------------------------------------------*)</syntaxhighlight>
{{out}}
<pre>$ patscc -O3 -DATS_MEMALLOC_GCBDW convex_hull_task.dats -lgc && ./a.out
(-9.000000 -3.000000)
(-3.000000 -9.000000)
(19.000000 -8.000000)
(17.000000 5.000000)
(12.000000 17.000000)
(5.000000 19.000000)
(-3.000000 15.000000)</pre>
=={{header|C}}==
{{trans|C++}}
<
#include <stdbool.h>
#include <stdio.h>
Line 494 ⟶ 1,079:
return 0;
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
=={{header|C sharp|C#}}==
<
using System.Collections.Generic;
Line 594 ⟶ 1,179:
}
}
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
Line 600 ⟶ 1,185:
=={{header|C++}}==
{{trans|D}}
<
#include <iostream>
#include <ostream>
Line 687 ⟶ 1,272:
return 0;
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
=={{header|Common Lisp}}==
{{trans|Scheme}}
The program is wrapped in a [https://roswell.github.io/ Roswell] script, which is one of the popular ways to get around differences between Common Lisp implementations.
(Side note: It occurs to me that the "complete with tail recursions" comment in the code is not meant to include the Scheme '''do''' constructs, which are secretly really tail recursions. How Common Lisps implement the '''loop''' macro I do not know; it might be tail recursion, or it might not. Note also how the Scheme uses "named '''let'''" syntactic sugar, where the Common Lisp has '''defun''' fully written out.)
<syntaxhighlight lang="lisp">#!/bin/sh
#|-*- mode:lisp -*-|#
#|
exec ros -Q -- $0 "$@"
|#
(progn ;;init forms
(ros:ensure-asdf)
#+quicklisp(ql:quickload '() :silent t)
)
(defpackage :ros.script.convex-hull-task.3861520611
(:use :cl))
(in-package :ros.script.convex-hull-task.3861520611)
;;;
;;; Convex hulls by Andrew's monotone chain algorithm.
;;;
;;; For a description of the algorithm, see
;;; https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
;;;
;;; This program is translated rather faithfully from the Scheme,
;;; complete with tail recursions.
;;;
;; x and y coordinates of a "point". A "point" is represented by a
;; list of length 2.
(defun x@ (u) (car u))
(defun y@ (u) (cadr u))
(defun cross (u v)
;; Cross product (as a signed scalar).
(- (* (x@ u) (y@ v)) (* (y@ u) (x@ v))))
(defun point- (u v)
(list (- (x@ u) (x@ v)) (- (y@ u) (y@ v))))
(defun sort-points-vector (points-vector)
;; Ascending sort on x-coordinates, followed by ascending sort
;; on y-coordinates.
(sort points-vector #'(lambda (u v)
(or (< (x@ u) (x@ v))
(and (= (x@ u) (x@ v))
(< (y@ u) (y@ v)))))))
(defun construct-lower-hull (sorted-points-vector)
(let* ((pt sorted-points-vector)
(n (length pt))
(hull (make-array n))
(j 1))
(setf (aref hull 0) (aref pt 0))
(setf (aref hull 1) (aref pt 1))
(loop for i from 2 to (1- n)
do (progn
(defun inner-loop ()
(if (or (zerop j)
(plusp
(cross (point- (aref hull j)
(aref hull (1- j)))
(point- (aref pt i)
(aref hull (1- j))))))
(progn
(setf j (1+ j))
(setf (aref hull j) (aref pt i)))
(progn
(setf j (1- j))
(inner-loop))))
(inner-loop)))
(values (+ j 1) hull))) ; Hull size, hull points.
(defun construct-upper-hull (sorted-points-vector)
(let* ((pt sorted-points-vector)
(n (length pt))
(hull (make-array n))
(j 1))
(setf (aref hull 0) (aref pt (- n 1)))
(setf (aref hull 1) (aref pt (- n 2)))
(loop for i from (- n 3) downto 0
do (progn
(defun inner-loop ()
(if (or (zerop j)
(plusp
(cross (point- (aref hull j)
(aref hull (1- j)))
(point- (aref pt i)
(aref hull (1- j))))))
(progn
(setf j (1+ j))
(setf (aref hull j) (aref pt i)))
(progn
(setf j (1- j))
(inner-loop))))
(inner-loop)))
(values (+ j 1) hull))) ; Hull size, hull points.
(defun construct-hull (sorted-points-vector)
;; Notice that the lower and upper hulls could be constructed in
;; parallel. (The Scheme "let-values" macro made this apparent,
;; despite not actually doing the computation in parallel. The
;; coding here makes it less obvious.)
(multiple-value-bind (lower-hull-size lower-hull)
(construct-lower-hull sorted-points-vector)
(multiple-value-bind (upper-hull-size upper-hull)
(construct-upper-hull sorted-points-vector)
(let* ((hull-size (+ lower-hull-size upper-hull-size -2))
(hull (make-array hull-size)))
(loop for i from 0 to (- lower-hull-size 2)
do (setf (aref hull i) (aref lower-hull i)))
(loop for i from 0 to (- upper-hull-size 2)
do (setf (aref hull (+ i (1- lower-hull-size)))
(aref upper-hull i)))
hull))))
(defun vector-delete-neighbor-dups (elt= v)
;; A partial clone of the SRFI-132 procedure of the same name. This
;; implementation is similar to the reference implementation for
;; SRFI-132, and may use a bunch of stack space. That reference
;; implementation is by Olin Shivers and rests here:
;; https://github.com/scheme-requests-for-implementation/srfi-132/blob/master/sorting/delndups.scm
;; The license is:
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;; This code is
;;; Copyright (c) 1998 by Olin Shivers.
;;; The terms are: You may do as you please with this code, as long as
;;; you do not delete this notice or hold me responsible for any outcome
;;; related to its use.
;;;
;;; Blah blah blah. Don't you think source files should contain more lines
;;; of code than copyright notice?
;;;
(let ((start 0)
(end (length v)))
(let ((x (aref v start)))
(defun recur (x i j)
(if (< i end)
(let ((y (aref v i))
(nexti (1+ i)))
(if (funcall elt= x y)
(recur x nexti j)
(let ((ansvec (recur y nexti (1+ j))))
(setf (aref ansvec j) y)
ansvec)))
(make-array j)))
(let ((ans (recur x start 1)))
(setf (aref ans 0) x)
ans))))
(defun vector-convex-hull (points)
(let* ((points-vector (coerce points 'vector))
(sorted-points-vector
(vector-delete-neighbor-dups
#'equalp
(sort-points-vector points-vector))))
(if (<= (length sorted-points-vector) 2)
sorted-points-vector
(construct-hull sorted-points-vector))))
(defun list-convex-hull (points)
(coerce (vector-convex-hull points) 'list))
(defconstant example-points
'((16 3) (12 17) (0 6) (-4 -6) (16 6)
(16 -7) (16 -3) (17 -4) (5 19) (19 -8)
(3 16) (12 13) (3 -4) (17 5) (-3 15)
(-3 -9) (0 11) (-9 -3) (-4 -2) (12 10)))
(defun main (&rest argv)
(declare (ignorable argv))
(write (list-convex-hull example-points))
(terpri))
;;; vim: set ft=lisp lisp:</syntaxhighlight>
{{out}}
<pre>$ ./convex-hull-task.ros
((-9 -3) (-3 -9) (19 -8) (17 5) (12 17) (5 19) (-3 15))</pre>
=={{header|D}}==
{{trans|Kotlin}}
<
import std.stdio;
Line 757 ⟶ 1,527:
auto hull = convexHull(points);
writeln("Convex Hull: ", hull);
}</
{{out}}
<pre>Convex Hull: [(-9,-3), (-3,-9), (19,-8), (17,5), (12,17), (5,19), (-3,15)]</pre>
Line 768 ⟶ 1,538:
{{libheader| System.Generics.Defaults}}
{{libheader| System.Generics.Collections}}
<syntaxhighlight lang="delphi">
program ConvexHulls;
Line 878 ⟶ 1,648:
end.
</syntaxhighlight>
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
=={{header|EasyLang}}==
{{trans|Nim}}
<syntaxhighlight>
func orientation p[] q[] r[] .
return (q[2] - p[2]) * (r[1] - q[1]) - (q[1] - p[1]) * (r[2] - q[2])
.
proc calcConvexHull . pts[][] res[][] .
res[][] = [ ]
indMinX = 1
for i to len pts[][]
if pts[i][1] < pts[indMinX][1]
indMinX = i
.
.
p = indMinX
repeat
res[][] &= pts[p][]
q = (p + 1) mod1 len pts[][]
for i to len pts[][]
if orientation pts[p][] pts[i][] pts[q][] < 0
q = i
.
.
p = q
until p = indMinX
.
.
#
pts[][] = [ [ 16 3 ] [ 12 17 ] [ 0 6 ] [ -4 -6 ] [ 16 6 ] [ 16 -7 ] [ 16 -3 ] [ 17 -4 ] [ 5 19 ] [ 19 -8 ] [ 3 16 ] [ 12 13 ] [ 3 -4 ] [ 17 5 ] [ -3 15 ] [ -3 -9 ] [ 0 11 ] [ -9 -3 ] [ -4 -2 ] [ 12 10 ] ]
calcConvexHull pts[][] res[][]
print res[][]
</syntaxhighlight>
{{out}}
<pre>
[
[ -9 -3 ]
[ -3 -9 ]
[ 19 -8 ]
[ 17 5 ]
[ 12 17 ]
[ 5 19 ]
[ -3 15 ]
]
</pre>
=={{header|F_Sharp|F#}}==
{{trans|C}}
<
type Point =
Line 931 ⟶ 1,746:
affiche (convexHull (List.sortBy (fun (x : Point) -> x.X, x.Y) poly))
</syntaxhighlight>
{{out}}
<pre>Convex Hull: [(-9,-3), (-3,-9), (19,-8), (17,5), (12,17), (5,19), (-3,15)]</pre>
Line 940 ⟶ 1,755:
<
!
! Convex hulls by Andrew's monotone chain algorithm.
Line 1,021 ⟶ 1,836:
end do
end subroutine sort_points
subroutine delete_neighbor_duplicates (n, pt)
integer, intent(inout) :: n
complex, intent(inout) :: pt(0:*)
call delete_trailing_duplicates
call delete_nontrailing_duplicates
contains
subroutine delete_trailing_duplicates
integer :: i
logical :: done
i = n - 1
done = .false.
do while (.not. done)
if (i == 0) then
n = 1
done = .true.
else if (pt(i - 1) /= pt(i)) then
n = i + 1
done = .true.
else
i = i - 1
end if
end do
end subroutine delete_trailing_duplicates
subroutine delete_nontrailing_duplicates
integer :: i, j, num_deleted
logical :: done
i = 0
do while (i < n - 1)
j = i + 1
done = .false.
do while (.not. done)
if (j == n) then
done = .true.
else if (pt(j) /= pt(i)) then
done = .true.
else
j = j + 1
end if
end do
if (j /= i + 1) then
num_deleted = j - i - 1
do while (j /= n)
pt(j - num_deleted) = pt(j)
j = j + 1
end do
n = n - num_deleted
end if
i = i + 1
end do
end subroutine delete_nontrailing_duplicates
end subroutine delete_neighbor_duplicates
subroutine construct_lower_hull (n, pt, hull_size, hull)
Line 1,124 ⟶ 1,998:
complex :: pt(0 : n - 1)
integer :: ipoints0, ihull0, numpt
ipoints0 = lbound (points, 1)
ihull0 = lbound (hull, 1)
numpt = n
call sort_points (numpt, pt)
call delete_neighbor_duplicates (numpt, pt)
if (numpt == 0) then
hull_size = 0
else if (
hull_size =
hull(:ihull0 + numpt - 1) =
else
end if
end subroutine find_convex_hull
Line 1,167 ⟶ 2,042:
write (*, fmt) (points(i), i = 1, n)
end program convex_hull_task</
{{out}}
Line 1,175 ⟶ 2,050:
=={{header|FreeBASIC}}==
<
#include "crt.bi"
Screen 20
Line 1,267 ⟶ 2,142:
Next
Sleep
</syntaxhighlight>
{{out}}
<pre>
Line 1,281 ⟶ 2,156:
=={{header|Go}}==
<
import (
Line 1,343 ⟶ 2,218:
hull := pts.ConvexHull()
fmt.Println("Convex Hull:", hull)
}</
{{out}}
<pre>
Line 1,351 ⟶ 2,226:
=={{header|Groovy}}==
{{trans|Java}}
<
private static class Point implements Comparable<Point> {
private int x, y
Line 1,436 ⟶ 2,311:
println("Convex Hull: $hull")
}
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
=={{header|Haskell}}==
<
import Data.Ord (comparing)
Line 1,501 ⟶ 2,376:
, [-4, -2]
, [12, 10]
]</
{{Out}}
<pre>[-3.0,-9.0]
Line 1,513 ⟶ 2,388:
=={{header|Haxe}}==
{{trans|Nim}}
<
class Main {
Line 1,605 ⟶ 2,480:
Sys.println(']');
}
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, 15), (5, 19), (12, 17), (17, 5), (19, -8), (-3, -9)]</pre>
=={{header|Icon}}==
{{trans|ObjectIcon}}
<syntaxhighlight lang="icon">#
# Convex hulls by Andrew's monotone chain algorithm.
#
# For a description of the algorithm, see
# https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
#
record PlanePoint (x, y)
######################################################################
#
# Merge sort adapted from the Object Icon IPL (public domain code).
#
# A merge sort implementation. This returns a sorted copy, leaving the
# original unchanged.
#
# :Parameters :
# : `l` - the list to sort
# : `cmp` - a comparator function
#
procedure mergesort (l, cmp)
return mergesort1 (l, cmp, 1, *l)
end
procedure mergesort1 (l, cmp, first, last)
local l1, l2, l3, m, v1
if last <= first then
return l[first:last + 1]
m := (first + last) / 2
l1 := mergesort1 (l, cmp, first, m)
l2 := mergesort1 (l, cmp, m + 1, last)
l3 := []
every v1 := !l1 do {
while cmp (v1, l2[1]) > 0 do
put (l3, get(l2))
put (l3, v1)
}
every put(l3, !l2)
return l3
end
######################################################################
procedure point_equals (p, q)
if p.x = q.x & p.y = q.y then return else fail
end
# Impose a total order on points, making it one that will work for
# Andrew's monotone chain algorithm. *)
procedure point_comes_before (p, q)
if (p.x < q.x) | (p.x = q.x & p.y < q.y) then return else fail
end
# Subtraction is really a vector or multivector operation.
procedure point_subtract (p, q)
return PlanePoint (p.x - q.x, p.y - q.y)
end
# Cross product is really a multivector operation.
procedure point_cross (p, q)
return (p.x * q.y) - (p.y * q.x)
end
procedure point_to_string (p)
return "(" || string (p.x) || " " || string (p.y) || ")"
end
######################################################################
# Comparison like C's strcmp(3).
procedure compare_points (p, q)
local cmp
if point_comes_before (p, q) then
cmp := -1
else if point_comes_before (q, p) then
cmp := 1
else
cmp := 0
return cmp
end
procedure sort_points (points)
# Non-destructive sort.
return mergesort (points, compare_points)
end
procedure delete_neighbor_dups (arr, equals)
local arr1, i
if *arr = 0 then {
arr1 := []
} else {
arr1 := [arr[1]]
i := 2
while i <= *arr do {
if not (equals (arr[i], arr1[-1])) then
put (arr1, arr[i])
i +:= 1
}
}
return arr1
end
procedure construct_lower_hull (pt)
local hull, i, j
hull := list (*pt)
hull[1] := pt[1]
hull[2] := pt[2]
j := 2
every i := 3 to *pt do {
while (j ~= 1 &
point_cross (point_subtract (hull[j], hull[j - 1]),
point_subtract (pt[i], hull[j - 1])) <= 0) do j -:= 1
j +:= 1
hull[j] := pt[i]
}
return hull[1 : j + 1]
end
procedure construct_upper_hull (pt)
local hull, i, j
hull := list (*pt)
hull[1] := pt[-1]
hull[2] := pt[-2]
j := 2
every i := 3 to *pt do {
while (j ~= 1 &
point_cross (point_subtract (hull[j], hull[j - 1]),
point_subtract (pt[-i], hull[j - 1])) <= 0) do j -:= 1
j +:= 1
hull[j] := pt[-i]
}
return hull[1 : j + 1]
end
procedure construct_hull (pt)
local lower_hull, upper_hull
lower_hull := construct_lower_hull (pt)
upper_hull := construct_upper_hull (pt)
return lower_hull[1 : -1] ||| upper_hull [1 : -1]
end
procedure find_convex_hull (points)
local pt, hull
if *points = 0 then {
hull := []
} else {
pt := delete_neighbor_dups (sort_points (points), point_equals)
if *pt <= 2 then {
hull := pt
} else {
hull := construct_hull (pt)
}
}
return hull
end
procedure main ()
local example_points, hull
example_points :=
[PlanePoint (16.0, 3.0),
PlanePoint (12.0, 17.0),
PlanePoint (0.0, 6.0),
PlanePoint (-4.0, -6.0),
PlanePoint (16.0, 6.0),
PlanePoint (16.0, -7.0),
PlanePoint (16.0, -3.0),
PlanePoint (17.0, -4.0),
PlanePoint (5.0, 19.0),
PlanePoint (19.0, -8.0),
PlanePoint (3.0, 16.0),
PlanePoint (12.0, 13.0),
PlanePoint (3.0, -4.0),
PlanePoint (17.0, 5.0),
PlanePoint (-3.0, 15.0),
PlanePoint (-3.0, -9.0),
PlanePoint (0.0, 11.0),
PlanePoint (-9.0, -3.0),
PlanePoint (-4.0, -2.0),
PlanePoint (12.0, 10.0)]
hull := find_convex_hull (example_points)
every write (point_to_string (!hull))
end
######################################################################</syntaxhighlight>
{{out}}
<pre>$ icon convex_hull_task-Icon.icn
(-9.0 -3.0)
(-3.0 -9.0)
(19.0 -8.0)
(17.0 5.0)
(12.0 17.0)
(5.0 19.0)
(-3.0 15.0)</pre>
=={{header|J}}==
Line 1,614 ⟶ 2,698:
Restated from the implementation at [https://web.archive.org/web/20150919194242/http://kukuruku.co/hub/funcprog/introduction-to-j-programming-language-2004] which in turn is Kukuruku Hub's [https://web.archive.org/web/20150908060956/http://kukuruku.co/] translation of Dr. K. L. Metlov's original Russian article [http://dr-klm.livejournal.com/42312.html].
<
crossproduct =: 11 o. [: (* +)/ }. - {.
removeinner =: #~ (1 , (0 > 3 crossproduct\ ]) , 1:)
hull =: [: removeinner^:_ counterclockwise</
Example use:
<
_9j_3 _3j_9 19j_8 17j5 12j17 5j19 _3j15
Line 1,654 ⟶ 2,738:
1 7
0 4
</syntaxhighlight>
=={{header|Java}}==
{{trans|Kotlin}}
<
import java.util.Arrays;
import java.util.List;
Line 1,745 ⟶ 2,829:
System.out.printf("Convex Hull: %s\n", hull);
}
}</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
=={{header|Javascript}}==
<syntaxhighlight lang="javascript">
function convexHull(points) {
points.sort(comparison);
Line 1,779 ⟶ 2,863:
return (a.x - o.x) * (b.y - o.y) - (a.y - o.y) * (b.x - o.x);
}
</syntaxhighlight>
'''For usage''':
<nowiki>convexhull.js</nowiki>
<syntaxhighlight lang="javascript">
var points = [];
var hull = [];
Line 1,853 ⟶ 2,937:
return (a.x - o.x) * (b.y - o.y) - (a.y - o.y) * (b.x - o.x);
}
</syntaxhighlight>
<nowiki>convexhull.html</nowiki>
<
<html>
Line 1,875 ⟶ 2,959:
</html>
</syntaxhighlight>
=={{header|jq}}==
Line 1,881 ⟶ 2,965:
{{works with|jq}}
'''Works with gojq, the Go implementation of jq'''
<
def ccw(a; b; c):
a as [$ax, $ay]
Line 1,902 ⟶ 2,986:
| . + [$pt])
| .[:-1]
end ;</
'''The task'''
<
[16, 3], [12, 17], [ 0, 6], [-4, -6], [16, 6],
[16, -7], [16, -3], [17, -4], [ 5, 19], [19, -8],
Line 1,911 ⟶ 2,995:
];
"Convex Hull: \(pts|convexHull)"</
{{out}}
<pre>
Line 1,919 ⟶ 3,003:
=={{header|Julia}}==
<
# https://github.com/JuliaPolyhedra/Polyhedra.jl/blob/master/examples/operations.ipynb
using Polyhedra, CDDLib
Line 1,927 ⟶ 3,011:
Pch = convexhull(P, P)
removevredundancy!(Pch)
println("$Pch")</
{{out}}
<pre>convexhull([5.0, 19.0], [19.0, -8.0], [17.0, 5.0], [-3.0, 15.0], [-9.0, -3.0], [12.0, 17.0], [-3.0, -9.0])</pre>
Line 1,933 ⟶ 3,017:
=={{header|Kotlin}}==
{{trans|Go}}
<
class Point(val x: Int, val y: Int) : Comparable<Point> {
Line 1,982 ⟶ 3,066:
val hull = convexHull(points)
println("Convex Hull: $hull")
}</
{{out}}
Line 2,017 ⟶ 3,101:
=={{header|Lua}}==
{{trans|C++}}
<
io.write("("..p.x..", "..p.y..")")
return nil
Line 2,086 ⟶ 3,170:
io.write("Convex Hull: ")
print_points(hull)
print()</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
Line 2,094 ⟶ 3,178:
Function are available in the geometry, ComputationalGeometry and simplex packages.
<
[3,16],[12,13],[3,-4],[17,5],[-3,15],[-3,-9],[0,11],[-9,-3],[-4,-2],[12,10]]:
Line 2,106 ⟶ 3,190:
simplex:-convexhull(pts);
# [[-9, -3], [-3, -9], [19, -8], [17, 5], [12, 17], [5, 19], [-3, 15]]</
=={{header|Mathematica}}/{{header|Wolfram Language}}==
<
hullPoints[{{16, 3}, {12, 17}, {0, 6}, {-4, -6}, {16, 6}, {16, -7}, {16, -3}, {17, -4}, {5, 19}, {19, -8}, {3, 16}, {12, 13}, {3, -4}, {17, 5}, {-3, 15}, {-3, -9}, {0, 11}, {-9, -3}, {-4, -2}, {12, 10}}]</
{{out}}{{12., 17.}, {5., 19.}, {19., -8.}, {17., 5.}, {-3., 15.}, {-3., -9.}, {-9., -3.}}
=={{header|Mercury}}==
{{trans|Standard ML}}
{{works with|Mercury|20.06.1}}
<syntaxhighlight lang="mercury">:- module convex_hull_task.
:- interface.
:- import_module io.
:- pred main(io, io).
:- mode main(di, uo) is det.
:- implementation.
:- import_module exception.
:- import_module float.
:- import_module int.
:- import_module list.
:- import_module pair.
:- import_module string.
:- import_module version_array.
%%--------------------------------------------------------------------
%% fetch_items/3 for version_array, similar to the library function
%% for regular array.
:- func fetch_items(version_array(T), int, int) = list(T).
fetch_items(Arr, I, J) = fetch_items_(Arr, I, J, []).
:- func fetch_items_(version_array(T), int, int, list(T)) = list(T).
fetch_items_(Arr, I, J, Lst0) = Lst :-
if (J < I) then (Lst = Lst0)
else (J1 = J - 1,
Lst = fetch_items_(Arr, I, J1, [Arr^elem(J) | Lst0])).
%%--------------------------------------------------------------------
:- type point == pair(float).
:- type point_list == list(point).
:- type point_array == version_array(point).
:- pred point_comes_before(point, point).
:- mode point_comes_before(in, in) is semidet.
point_comes_before(P, Q) :-
(fst(P) < fst(Q); fst(P) - fst(Q) = (0.0),
snd(P) < snd(Q)).
:- pred point_compare(point, point, comparison_result).
:- mode point_compare(in, in, out) is det.
point_compare(P, Q, Cmp) :-
if (point_comes_before(P, Q)) then (Cmp = (<))
else if (point_comes_before(Q, P)) then (Cmp = (>))
else (Cmp = (=)).
:- func point_subtract(point, point) = point.
point_subtract(P, Q) = pair(fst(P) - fst(Q),
snd(P) - snd(Q)).
:- func point_cross_product(point, point) = float.
point_cross_product(P, Q) = (fst(P) * snd(Q)) - (snd(P) * fst(Q)).
:- func point_to_string(point) = string.
point_to_string(P) = ("(" ++ from_float(fst(P)) ++
" " ++ from_float(snd(P)) ++ ")").
%%--------------------------------------------------------------------
:- func convex_hull(point_list) = point_list.
convex_hull(Pt) = Hull :-
Pt1 = unique_points_sorted(Pt),
(if (Pt1 = []; Pt1 = [_]; Pt1 = [_, _]) then (Hull = Pt1)
else (Hull = construct_hull(Pt1))).
:- func unique_points_sorted(point_list) = point_list.
unique_points_sorted(Pt0) = Pt :-
sort_and_remove_dups(point_compare, Pt0, Pt).
:- func construct_hull(point_list) = point_list.
construct_hull(Pt) = (construct_lower_hull(Pt) ++
construct_upper_hull(Pt)).
:- func construct_lower_hull(point_list) = point_list.
construct_lower_hull(Pt) = Hull :-
if (Pt = [P0, P1 | Rest])
then (N = length(Pt),
Arr0 = (version_array.init(N, P0)),
Arr1 = (Arr0^elem(1) := P1),
hull_construction(Rest, Arr1, Arr2, 1, N_Hull),
%% In the fetch_items/3 call, we leave out the last item. It
%% is redundant with the first item of the upper hull.
N_Hull2 = N_Hull - 2,
Hull = fetch_items(Arr2, 0, N_Hull2))
else throw("construct_lower_hull expects list of length >= 3").
:- func construct_upper_hull(point_list) = point_list.
%% An upper hull is merely a lower hull for points going the other
%% way.
construct_upper_hull(Pt) = construct_lower_hull(reverse(Pt)).
:- pred hull_construction(point_list, point_array, point_array,
int, int).
:- mode hull_construction(in, in, out, in, out) is det.
hull_construction([], Arr0, Arr, J, N_Hull) :-
Arr = Arr0,
N_Hull = J + 1.
hull_construction([P | Rest], Arr0, Arr, J, N_Hull) :-
if cross_test(P, Arr0, J)
then (J1 = J + 1,
Arr1 = (Arr0^elem(J1) := P),
hull_construction(Rest, Arr1, Arr, J1, N_Hull))
else (J1 = J - 1,
hull_construction([P | Rest], Arr0, Arr, J1, N_Hull)).
:- pred cross_test(point, point_array, int).
:- mode cross_test(in, in, in) is semidet.
cross_test(P, Arr, J) :-
if (J = 0) then true
else (Elem_J = Arr^elem(J),
J1 = J - 1,
Elem_J1 = Arr^elem(J1),
0.0 < point_cross_product(point_subtract(Elem_J, Elem_J1),
point_subtract(P, Elem_J1))).
%%--------------------------------------------------------------------
main(!IO) :-
Example_points = [pair(16.0, 3.0),
pair(12.0, 17.0),
pair(0.0, 6.0),
pair(-4.0, -6.0),
pair(16.0, 6.0),
pair(16.0, -7.0),
pair(16.0, -3.0),
pair(17.0, -4.0),
pair(5.0, 19.0),
pair(19.0, -8.0),
pair(3.0, 16.0),
pair(12.0, 13.0),
pair(3.0, -4.0),
pair(17.0, 5.0),
pair(-3.0, 15.0),
pair(-3.0, -9.0),
pair(0.0, 11.0),
pair(-9.0, -3.0),
pair(-4.0, -2.0),
pair(12.0, 10.0)],
Hull = convex_hull(Example_points),
HullStr = join_list(" ", map(point_to_string, Hull)),
write_string(HullStr, !IO),
nl(!IO).
%%%-------------------------------------------------------------------
%%% local variables:
%%% mode: mercury
%%% prolog-indent-width: 2
%%% end:</syntaxhighlight>
{{out}}
<pre>$ mmc convex_hull_task.m && ./convex_hull_task
(-9.0 -3.0) (-3.0 -9.0) (19.0 -8.0) (17.0 5.0) (12.0 17.0) (5.0 19.0) (-3.0 15.0)</pre>
=={{header|Modula-2}}==
<
FROM FormatString IMPORT FormatString;
FROM Storage IMPORT ALLOCATE, DEALLOCATE;
Line 2,362 ⟶ 3,606:
DeleteNode(nodes);
ReadChar
END ConvexHull.</
{{out}}
<pre>[(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
Line 2,368 ⟶ 3,612:
=={{header|Nim}}==
{{trans|Rust}}
<
Point = object
x: float
Line 2,439 ⟶ 3,683:
let hull = calculateConvexHull(points)
for i in hull:
echo i</
{{out}}
<pre>
Line 2,450 ⟶ 3,694:
(x: -3.0, y: 15.0)
</pre>
=={{header|ObjectIcon}}==
{{trans|Fortran}}
{{trans|Standard ML}}
<syntaxhighlight lang="objecticon"># -*- ObjectIcon -*-
#
# Convex hulls by Andrew's monotone chain algorithm.
#
# For a description of the algorithm, see
# https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
#
import io
import ipl.sort
class PlanePoint () # Enough plane geometry for our purpose.
private readable x, y
public new (x, y)
self.x := x
self.y := y
return
end
public equals (other)
if self.x = other.x & self.y = other.y then
return
else
fail
end
# Impose a total order on points, making it one that will work for
# Andrew's monotone chain algorithm. *)
public comes_before (other)
if (self.x < other.x) | (self.x = other.x & self.y < other.y) then
return
else
fail
end
# Subtraction is really a vector or multivector operation.
public minus (other)
return PlanePoint (self.x - other.x, self.y - other.y)
end
# Cross product is really a multivector operation.
public cross (other)
return (self.x * other.y) - (self.y * other.x)
end
public to_string ()
return "(" || string (self.x) || " " || string (self.y) || ")"
end
end
# Comparison like C's strcmp(3).
procedure compare_points (p, q)
local cmp
if p.comes_before (q) then
cmp := -1
else if q.comes_before (p) then
cmp := 1
else
cmp := 0
return cmp
end
procedure sort_points (points)
# Non-destructive sort.
return mergesort (points, compare_points)
end
procedure delete_neighbor_dups (arr, equals)
local arr1, i
if *arr = 0 then {
arr1 := []
} else {
arr1 := [arr[1]]
i := 2
while i <= *arr do {
unless equals (arr[i], arr1[-1]) then
put (arr1, arr[i])
i +:= 1
}
}
return arr1
end
procedure construct_lower_hull (pt)
local hull, i, j
hull := list (*pt)
hull[1] := pt[1]
hull[2] := pt[2]
j := 2
every i := 3 to *pt do {
while (j ~= 1 &
(hull[j].minus (hull[j - 1])).cross (pt[i].minus (hull[j - 1])) <= 0) do
j -:= 1
j +:= 1
hull[j] := pt[i]
}
return hull[1 : j + 1]
end
procedure construct_upper_hull (pt)
local hull, i, j
hull := list (*pt)
hull[1] := pt[-1]
hull[2] := pt[-2]
j := 2
every i := 3 to *pt do {
while (j ~= 1 &
(hull[j].minus (hull[j - 1])).cross (pt[-i].minus (hull[j - 1])) <= 0) do
j -:= 1
j +:= 1
hull[j] := pt[-i]
}
return hull[1 : j + 1]
end
procedure construct_hull (pt)
local lower_hull, upper_hull
lower_hull := construct_lower_hull (pt)
upper_hull := construct_upper_hull (pt)
return lower_hull[1 : -1] ||| upper_hull [1 : -1]
end
procedure points_equal (p, q)
if p.equals (q) then
return
else
fail
end
procedure find_convex_hull (points)
local pt, hull
if *points = 0 then {
hull := []
} else {
pt := delete_neighbor_dups (sort_points (points), points_equal)
if *pt <= 2 then
hull := pt
else
hull := construct_hull (pt)
}
return hull
end
procedure main ()
local example_points, hull
example_points :=
[PlanePoint (16.0, 3.0),
PlanePoint (12.0, 17.0),
PlanePoint (0.0, 6.0),
PlanePoint (-4.0, -6.0),
PlanePoint (16.0, 6.0),
PlanePoint (16.0, -7.0),
PlanePoint (16.0, -3.0),
PlanePoint (17.0, -4.0),
PlanePoint (5.0, 19.0),
PlanePoint (19.0, -8.0),
PlanePoint (3.0, 16.0),
PlanePoint (12.0, 13.0),
PlanePoint (3.0, -4.0),
PlanePoint (17.0, 5.0),
PlanePoint (-3.0, 15.0),
PlanePoint (-3.0, -9.0),
PlanePoint (0.0, 11.0),
PlanePoint (-9.0, -3.0),
PlanePoint (-4.0, -2.0),
PlanePoint (12.0, 10.0)]
hull := find_convex_hull (example_points)
every write ((!hull).to_string ())
end</syntaxhighlight>
{{out}}
<pre>$ oiscript convex_hull_task-OI.icn
(-9.0 -3.0)
(-3.0 -9.0)
(19.0 -8.0)
(17.0 5.0)
(12.0 17.0)
(5.0 19.0)
(-3.0 15.0)</pre>
=={{header|OCaml}}==
{{trans|Standard ML}}
{{works with|OCaml|4.14.0}}
<syntaxhighlight lang="ocaml">(*
* Convex hulls by Andrew's monotone chain algorithm.
*
* For a description of the algorithm, see
* https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
*)
(*------------------------------------------------------------------*)
(* Just enough plane geometry for our purpose. *)
module Plane_point =
struct
type t = float * float
let make (xy : t) = xy
let to_tuple ((x, y) : t) = (x, y)
let x ((x, _) : t) = x
let y ((_, y) : t) = y
let equal (u : t) (v : t) = (x u = x v && y u = y v)
(* Impose a total order on points, making it one that will work
for Andrew's monotone chain algorithm. *)
let order (u : t) (v : t) =
(x u < x v) || (x u = x v && y u < y v)
(* The Array module's sort routines expect a "cmp" function. *)
let cmp u v =
if order u v then
(-1)
else if order v u then
(1)
else
(0)
(* Subtraction is really a vector or multivector operation. *)
let sub (u : t) (v : t) = make (x u -. x v, y u -. y v)
(* Cross product is really a multivector operation. *)
let cross (u : t) (v : t) = (x u *. y v) -. (y u *. x v)
let to_string ((x, y) : t) =
"(" ^ string_of_float x ^ " " ^ string_of_float y ^ ")"
end
;;
(*------------------------------------------------------------------*)
(* We want something akin to array_delete_neighbor_dups of Scheme's
SRFI-132. *)
let array_delete_neighbor_dups equal arr =
(* Returns a Seq.t rather than an array. *)
let rec loop i lst =
(* Cons a list of non-duplicates, going backwards through the
array so the list will be in forwards order. *)
if i = 0 then
arr.(i) :: lst
else if equal arr.(i - 1) arr.(i) then
loop (i - 1) lst
else
loop (i - 1) (arr.(i) :: lst)
in
let n = Array.length arr in
List.to_seq (if n = 0 then [] else loop (n - 1) [])
;;
(*------------------------------------------------------------------*)
(* The convex hull algorithm. *)
let cross_test pt_i hull j =
let hull_j = hull.(j)
and hull_j1 = hull.(j - 1) in
0.0 < Plane_point.(cross (sub hull_j hull_j1)
(sub pt_i hull_j1))
let construct_lower_hull n pt =
let hull = Array.make n (Plane_point.make (0.0, 0.0)) in
let () = hull.(0) <- pt.(0)
and () = hull.(1) <- pt.(1) in
let rec outer_loop i j =
if i = n then
j + 1
else
let pt_i = pt.(i) in
let rec inner_loop j =
if j = 0 || cross_test pt_i hull j then
begin
hull.(j + 1) <- pt_i;
j + 1
end
else
inner_loop (j - 1)
in
outer_loop (i + 1) (inner_loop j)
in
let hull_size = outer_loop 2 1 in
(hull_size, hull)
let construct_upper_hull n pt =
let hull = Array.make n (Plane_point.make (0.0, 0.0)) in
let () = hull.(0) <- pt.(n - 1)
and () = hull.(1) <- pt.(n - 2) in
let rec outer_loop i j =
if i = (-1) then
j + 1
else
let pt_i = pt.(i) in
let rec inner_loop j =
if j = 0 || cross_test pt_i hull j then
begin
hull.(j + 1) <- pt_i;
j + 1
end
else
inner_loop (j - 1)
in
outer_loop (i - 1) (inner_loop j)
in
let hull_size = outer_loop (n - 3) 1 in
(hull_size, hull)
let construct_hull n pt =
(* Side note: Construction of the lower and upper hulls can be done
in parallel. *)
let (lower_hull_size, lower_hull) = construct_lower_hull n pt
and (upper_hull_size, upper_hull) = construct_upper_hull n pt in
let hull_size = lower_hull_size + upper_hull_size - 2 in
let hull = Array.make hull_size (Plane_point.make (0.0, 0.0)) in
begin
Array.blit lower_hull 0 hull 0 (lower_hull_size - 1);
Array.blit upper_hull 0 hull (lower_hull_size - 1)
(upper_hull_size - 1);
hull
end
let plane_convex_hull points =
(* Takes an arbitrary sequence of points, which may be in any order
and may contain duplicates. Returns an ordered array of points
that make up the convex hull. If the sequence of points is empty,
the returned array is empty. *)
let pt = Array.of_seq points in
let () = Array.fast_sort Plane_point.cmp pt in
let pt = Array.of_seq
(array_delete_neighbor_dups Plane_point.equal pt) in
let n = Array.length pt in
if n <= 2 then
pt
else
construct_hull n pt
;;
(*------------------------------------------------------------------*)
let example_points =
[Plane_point.make (16.0, 3.0);
Plane_point.make (12.0, 17.0);
Plane_point.make (0.0, 6.0);
Plane_point.make (-4.0, -6.0);
Plane_point.make (16.0, 6.0);
Plane_point.make (16.0, -7.0);
Plane_point.make (16.0, -3.0);
Plane_point.make (17.0, -4.0);
Plane_point.make (5.0, 19.0);
Plane_point.make (19.0, -8.0);
Plane_point.make (3.0, 16.0);
Plane_point.make (12.0, 13.0);
Plane_point.make (3.0, -4.0);
Plane_point.make (17.0, 5.0);
Plane_point.make (-3.0, 15.0);
Plane_point.make (-3.0, -9.0);
Plane_point.make (0.0, 11.0);
Plane_point.make (-9.0, -3.0);
Plane_point.make (-4.0, -2.0);
Plane_point.make (12.0, 10.0)]
;;
let hull = plane_convex_hull (List.to_seq example_points)
;;
Array.iter
(fun p -> (print_string (Plane_point.to_string p);
print_string " "))
hull;
print_newline ()
;;
(*------------------------------------------------------------------*)</syntaxhighlight>
{{out}}
<pre>$ ocaml convex_hull_task.ml
(-9. -3.) (-3. -9.) (19. -8.) (17. 5.) (12. 17.) (5. 19.) (-3. 15.)</pre>
=={{header|Pascal}}==
{{trans|Fortran}}
<syntaxhighlight lang="pascal">{$mode ISO}
program convex_hull_task (output);
{ Convex hulls, by Andrew's monotone chain algorithm.
For a description of the algorithm, see
https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169 }
const max_points = 1000;
type points_range = 0 .. max_points - 1;
type point =
record
x, y : real
end;
type point_array = array [points_range] of point;
var ciura_gaps : array [1 .. 8] of integer;
var example_points : point_array;
var hull : point_array;
var hull_size : integer;
var index : integer;
function make_point (x, y : real) : point;
begin
make_point.x := x;
make_point.y := y;
end;
{ The cross product as a signed scalar. }
function cross (u, v : point) : real;
begin
cross := (u.x * v.y) - (u.y * v.x)
end;
function point_subtract (u, v : point) : point;
begin
point_subtract := make_point (u.x - v.x, u.y - v.y)
end;
function point_equal (u, v : point) : boolean;
begin
point_equal := (u.x = v.x) and (u.y = v.y)
end;
procedure sort_points (num_points : integer;
var points : point_array);
{ Sort first in ascending order by x-coordinates, then in
ascending order by y-coordinates. Any decent sort algorithm will
suffice; for the sake of interest, here is the Shell sort of
https://en.wikipedia.org/w/index.php?title=Shellsort&oldid=1084744510 }
var
i, j, k, gap, offset : integer;
temp : point;
done : boolean;
begin
for k := 1 to 8 do
begin
gap := ciura_gaps[k];
for offset := 0 to gap - 1 do
begin
i := offset;
while i <= num_points - 1 do
begin
temp := points[i];
j := i;
done := false;
while not done do
begin
if j < gap then
done := true
else if points[j - gap].x < temp.x then
done := true
else if ((points[j - gap].x = temp.x)
and (points[j - gap].y < temp.y)) then
done := true
else
begin
points[j] := points[j - gap];
j := j - gap
end
end;
points[j] := temp;
i := i + gap
end
end
end
end; { sort_points }
procedure delete_neighbor_duplicates (var n : integer;
var pt : point_array);
procedure delete_trailing_duplicates;
var
i : integer;
done : boolean;
begin
i := n - 1;
done := false;
while not done do
begin
if i = 0 then
begin
n := 1;
done := true
end
else if not point_equal (pt[i - 1], pt[i]) then
begin
n := i + 1;
done := true
end
else
i := i + 1
end
end;
procedure delete_nontrailing_duplicates;
var
i, j, num_deleted : integer;
done : boolean;
begin
i := 0;
while i < n - 1 do
begin
j := i + 1;
done := false;
while not done do
begin
if j = n then
done := true
else if not point_equal (pt[j], pt[i]) then
done := true
else
j := j + 1
end;
if j <> i + 1 then
begin
num_deleted := j - i - 1;
while j <> n do
begin
pt[j - num_deleted] := pt[j];
j := j + 1
end;
n := n - num_deleted
end;
i := i + 1
end
end;
begin
delete_trailing_duplicates;
delete_nontrailing_duplicates
end; { delete_neighbor_duplicates }
procedure construct_lower_hull (n : integer;
pt : point_array;
var hull_size : integer;
var hull : point_array);
var
i, j : integer;
done : boolean;
begin
j := 1;
hull[0] := pt[0];
hull[1] := pt[1];
for i := 2 to n - 1 do
begin
done := false;
while not done do
begin
if j = 0 then
begin
j := j + 1;
hull[j] := pt[i];
done := true
end
else if 0.0 < cross (point_subtract (hull[j],
hull[j - 1]),
point_subtract (pt[i],
hull[j - 1])) then
begin
j := j + 1;
hull[j] := pt[i];
done := true
end
else
j := j - 1
end
end;
hull_size := j + 1
end; { construct_lower_hull }
procedure construct_upper_hull (n : integer;
pt : point_array;
var hull_size : integer;
var hull : point_array);
var
i, j : integer;
done : boolean;
begin
j := 1;
hull[0] := pt[n - 1];
hull[1] := pt[n - 2];
for i := n - 3 downto 0 do
begin
done := false;
while not done do
begin
if j = 0 then
begin
j := j + 1;
hull[j] := pt[i];
done := true
end
else if 0.0 < cross (point_subtract (hull[j],
hull[j - 1]),
point_subtract (pt[i],
hull[j - 1])) then
begin
j := j + 1;
hull[j] := pt[i];
done := true
end
else
j := j - 1
end
end;
hull_size := j + 1
end; { construct_upper_hull }
procedure contruct_hull (n : integer;
pt : point_array;
var hull_size : integer;
var hull : point_array);
var
i : integer;
lower_hull_size, upper_hull_size : integer;
lower_hull, upper_hull : point_array;
begin
{ A side note: the calls to construct_lower_hull and
construct_upper_hull could be done in parallel. }
construct_lower_hull (n, pt, lower_hull_size, lower_hull);
construct_upper_hull (n, pt, upper_hull_size, upper_hull);
hull_size := lower_hull_size + upper_hull_size - 2;
for i := 0 to lower_hull_size - 2 do
hull[i] := lower_hull[i];
for i := 0 to upper_hull_size - 2 do
hull[lower_hull_size - 1 + i] := upper_hull[i]
end; { contruct_hull }
procedure find_convex_hull (n : integer;
points : point_array;
var hull_size : integer;
var hull : point_array);
var
pt : point_array;
numpt : integer;
i : integer;
begin
for i := 0 to n - 1 do
pt[i] := points[i];
numpt := n;
sort_points (numpt, pt);
delete_neighbor_duplicates (numpt, pt);
if numpt = 0 then
hull_size := 0
else if numpt <= 2 then
begin
hull_size := numpt;
for i := 0 to numpt - 1 do
hull[i] := pt[i]
end
else
contruct_hull (numpt, pt, hull_size, hull)
end; { find_convex_hull }
begin
ciura_gaps[1] := 701;
ciura_gaps[2] := 301;
ciura_gaps[3] := 132;
ciura_gaps[4] := 57;
ciura_gaps[5] := 23;
ciura_gaps[6] := 10;
ciura_gaps[7] := 4;
ciura_gaps[8] := 1;
example_points[0] := make_point (16, 3);
example_points[1] := make_point (12, 17);
example_points[2] := make_point (0, 6);
example_points[3] := make_point (-4, -6);
example_points[4] := make_point (16, 6);
example_points[5] := make_point (16, -7);
example_points[6] := make_point (16, -3);
example_points[7] := make_point (17, -4);
example_points[8] := make_point (5, 19);
example_points[9] := make_point (19, -8);
example_points[10] := make_point (3, 16);
example_points[11] := make_point (12, 13);
example_points[12] := make_point (3, -4);
example_points[13] := make_point (17, 5);
example_points[14] := make_point (-3, 15);
example_points[15] := make_point (-3, -9);
example_points[16] := make_point (0, 11);
example_points[17] := make_point (-9, -3);
example_points[18] := make_point (-4, -2);
example_points[19] := make_point (12, 10);
find_convex_hull (19, example_points, hull_size, hull);
for index := 0 to hull_size - 1 do
writeln (hull[index].x, ' ', hull[index].y)
end.
{--------------------------------------------------------------------}
{ The Emacs Pascal mode is intolerable.
Until I can find a substitute: }
{ local variables: }
{ mode: fundamental }
{ end: }</syntaxhighlight>
{{out}}
<pre>$ fpc convex_hull_task.pas && ./convex_hull_task
Free Pascal Compiler version 3.2.2 [2021/06/27] for x86_64
Copyright (c) 1993-2021 by Florian Klaempfl and others
Target OS: Linux for x86-64
Compiling convex_hull_task.pas
Linking convex_hull_task
324 lines compiled, 0.1 sec
-9.0000000000000000e+000 -3.0000000000000000e+000
-3.0000000000000000e+000 -9.0000000000000000e+000
1.9000000000000000e+001 -8.0000000000000000e+000
1.7000000000000000e+001 5.0000000000000000e+000
1.2000000000000000e+001 1.7000000000000000e+001
5.0000000000000000e+000 1.9000000000000000e+001
-3.0000000000000000e+000 1.5000000000000000e+001</pre>
=={{header|Perl}}==
{{trans|Raku}}
<
use warnings;
use feature 'say';
Line 2,524 ⟶ 4,511:
$list = join ' ', map { Point::print($_) } @hull_2;
say "Convex Hull (@{[scalar @hull_2]} points): [$list]";
</syntaxhighlight>
{{out}}
<pre>Convex Hull (7 points): [(-3, -9) (19, -8) (17, 5) (12, 17) (5, 19) (-3, 15) (-9, -3)]
Line 2,533 ⟶ 4,520:
{{libheader|Phix/online}}
You can run this online [http://phix.x10.mx/p2js/convexhull.htm here].
<!--<
<span style="color: #000080;font-style:italic;">--
-- demo/rosetta/Convex_hull.exw
Line 2,635 ⟶ 4,622:
<span style="color: #7060A8;">IupClose</span><span style="color: #0000FF;">()</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<!--</
{{out}}
<pre>
Line 2,656 ⟶ 4,643:
{{libheader|Shapely}}
An approach that uses the shapely library:
<
from shapely.geometry import MultiPoint
if __name__=="__main__":
pts = MultiPoint([(16,3), (12,17), (0,6), (-4,-6), (16,6), (16,-7), (16,-3), (17,-4), (5,19), (19,-8), (3,16), (12,13), (3,-4), (17,5), (-3,15), (-3,-9), (0,11), (-9,-3), (-4,-2), (12,10)])
print (pts.convex_hull)</
{{out}}
Line 2,669 ⟶ 4,656:
Also an implementation of https://en.wikibooks.org/wiki/Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain (therefore kinda {{trans|Go}}
<
(define-type Point (Pair Real Real))
(define-type Points (Listof Point))
Line 2,716 ⟶ 4,703:
'(5 . 19) '(19 . -8) '(3 . 16) '(12 . 13) '(3 . -4) '(17 . 5) '(-3 . 15) '(-3 . -9)
'(0 . 11) '(-9 . -3) '(-4 . -2) '(12 . 10)))
(list '(-9 . -3) '(-3 . -9) '(19 . -8) '(17 . 5) '(12 . 17) '(5 . 19) '(-3 . 15))))</
{{out}}
Line 2,728 ⟶ 4,715:
Modified the angle sort method as the original could fail if there were multiple points on the same y coordinate as the starting point. Sorts on tangent rather than triangle area. Inexpensive since it still doesn't do any trigonometric math, just calculates the ratio of opposite over adjacent.
<syntaxhighlight lang="raku"
has Real $.x is rw;
has Real $.y is rw;
Line 2,791 ⟶ 4,778:
);
say "Convex Hull ({+@hull} points): ", @hull;</
{{out}}
<pre>Convex Hull (7 points): [(-3, -9) (19, -8) (17, 5) (12, 17) (5, 19) (-3, 15) (-9, -3)]
Convex Hull (7 points): [(-3, -9) (20, -9) (17, 5) (12, 17) (5, 19) (-3, 15) (-9, -3)]</pre>
=={{header|RATFOR}}==
{{trans|Fortran}}
{{works with|ratfor77|[https://sourceforge.net/p/chemoelectric/ratfor77/ public domain 1.0]}}
{{works with|gfortran|11.3.0}}
{{works with|f2c|20100827}}
<syntaxhighlight lang="ratfor">#
# Convex hulls by Andrew's monotone chain algorithm.
#
# For a description of the algorithm, see
# https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
#
# As in the Fortran 2018 version upon which this code is based, I
# shall use the built-in "complex" type to represent "points" in the
# plane. This is merely for convenience, rather than to express a
# mathematical equivalence.
#
define(point, complex)
function x (u)
# Return the x-coordinate of "point" u.
implicit none
point u
real x
x = real (u)
end
function y (u)
# Return the y-coordinate of "point" u.
implicit none
point u
real y
y = aimag (u)
end
function cross (u, v)
# Return, as a signed scalar, the cross product of "points" u and v
# (regarded as "vectors" or multivectors).
implicit none
point u, v
real cross, x, y
cross = (x (u) * y (v)) - (y (u) * x (v))
end
subroutine sortpt (numpt, pt)
# Sort "points" in ascending order, first by the x-coordinates and
# then by the y-coordinates. Any decent sort algorithm will suffice;
# for the sake of interest, here is the Shell sort of
# https://en.wikipedia.org/w/index.php?title=Shellsort&oldid=1084744510
implicit none
integer numpt
point pt(0:*)
real x, y
integer i, j, k, gap, offset
point temp
logical done
integer gaps(1:8)
data gaps / 701, 301, 132, 57, 23, 10, 4, 1 /
for (k = 1; k <= 8; k = k + 1)
{
gap = gaps(k)
for (offset = 0; offset <= gap - 1; offset = offset + 1)
for (i = offset; i <= numpt - 1; i = i + gap)
{
temp = pt(i)
j = i
done = .false.
while (!done)
{
if (j < gap)
done = .true.
else if (x (pt(j - gap)) < x (temp))
done = .true.
else if (x (pt(j - gap)) == x (temp) _
&& (y (pt(j - gap)) <= y (temp)))
done = .true.
else
{
pt(j) = pt(j - gap)
j = j - gap
}
}
pt(j) = temp
}
}
end
subroutine deltrd (n, pt)
# Delete trailing neighbor duplicates.
implicit none
integer n
point pt(0:*)
integer i
logical done
i = n - 1
done = .false.
while (!done)
{
if (i == 0)
{
n = 1
done = .true.
}
else if (pt(i - 1) != pt(i))
{
n = i + 1
done = .true.
}
else
i = i - 1
}
end
subroutine delntd (n, pt)
# Delete non-trailing neighbor duplicates.
implicit none
integer n
point pt(0:*)
integer i, j, numdel
logical done
i = 0
while (i < n - 1)
{
j = i + 1
done = .false.
while (!done)
{
if (j == n)
done = .true.
else if (pt(j) != pt(i))
done = .true.
else
j = j + 1
}
if (j != i + 1)
{
numdel = j - i - 1
while (j != n)
{
pt(j - numdel) = pt(j)
j = j + 1
}
n = n - numdel
}
i = i + 1
}
end
subroutine deldup (n, pt)
# Delete neighbor duplicates.
implicit none
integer n
point pt(0:*)
call deltrd (n, pt)
call delntd (n, pt)
end
subroutine cxlhul (n, pt, hullsz, hull)
# Construct the lower hull.
implicit none
integer n # Number of points.
point pt(0:*)
integer hullsz # Output.
point hull(0:*) # Output.
real cross
integer i, j
logical done
j = 1
hull(0) = pt(0)
hull(1) = pt(1)
for (i = 2; i <= n - 1; i = i + 1)
{
done = .false.
while (!done)
{
if (j == 0)
{
j = j + 1
hull(j) = pt(i)
done = .true.
}
else if (0.0 < cross (hull(j) - hull(j - 1), _
pt(i) - hull(j - 1)))
{
j = j + 1
hull(j) = pt(i)
done = .true.
}
else
j = j - 1
}
}
hullsz = j + 1
end
subroutine cxuhul (n, pt, hullsz, hull)
# Construct the upper hull.
implicit none
integer n # Number of points.
point pt(0:*)
integer hullsz # Output.
point hull(0:*) # Output.
real cross
integer i, j
logical done
j = 1
hull(0) = pt(n - 1)
hull(1) = pt(n - 2)
for (i = n - 3; 0 <= i; i = i - 1)
{
done = .false.
while (!done)
{
if (j == 0)
{
j = j + 1
hull(j) = pt(i)
done = .true.
}
else if (0.0 < cross (hull(j) - hull(j - 1), _
pt(i) - hull(j - 1)))
{
j = j + 1
hull(j) = pt(i)
done = .true.
}
else
j = j - 1
}
}
hullsz = j + 1
end
subroutine cxhull (n, pt, hullsz, lhull, uhull)
# Construct the hull.
implicit none
integer n # Number of points.
point pt(0:*) # Overwritten with hull.
integer hullsz # Output.
point lhull(0:*) # Workspace.
point uhull(0:*) # Workspace
integer lhulsz, uhulsz, i
# A side note: the calls to construct_lower_hull and
# construct_upper_hull could be done in parallel.
call cxlhul (n, pt, lhulsz, lhull)
call cxuhul (n, pt, uhulsz, uhull)
hullsz = lhulsz + uhulsz - 2
for (i = 0; i <= lhulsz - 2; i = i + 1)
pt(i) = lhull(i)
for (i = 0; i <= uhulsz - 2; i = i + 1)
pt(lhulsz - 1 + i) = uhull(i)
end
subroutine cvxhul (n, pt, hullsz, lhull, uhull)
# Find a convex hull.
implicit none
integer n # Number of points.
point pt(0:*) # The contents of pt is replaced with the hull.
integer hullsz # Output.
point lhull(0:*) # Workspace.
point uhull(0:*) # Workspace
integer numpt
numpt = n
call sortpt (numpt, pt)
call deldup (numpt, pt)
if (numpt == 0)
hullsz = 0
else if (numpt <= 2)
hullsz = numpt
else
call cxhull (numpt, pt, hullsz, lhull, uhull)
end
program cvxtsk
# The Rosetta Code convex hull task.
implicit none
integer n, i
point points(0:100)
point lhull(0:100)
point uhull(0:100)
character*100 fmt
point exampl(0:19)
data exampl / (16, 3), (12, 17), (0, 6), (-4, -6), (16, 6), _
(16, -7), (16, -3), (17, -4), (5, 19), (19, -8), _
(3, 16), (12, 13), (3, -4), (17, 5), (-3, 15), _
(-3, -9), (0, 11), (-9, -3), (-4, -2), (12, 10) /
n = 20
for (i = 0; i <= n - 1; i = i + 1)
points(i) = exampl(i)
call cvxhul (n, points, n, lhull, uhull)
write (fmt, '("(", I20, ''("(", F3.0, 1X, F3.0, ") ")'', ")")') n
write (*, fmt) (points(i), i = 0, n - 1)
end</syntaxhighlight>
{{out}}
<pre>$ ratfor77 convex_hull_task.r > cvxtsk.f && f2c -C cvxtsk.f && cc cvxtsk.c -lf2c && ./a.out
(-9. -3.) (-3. -9.) (19. -8.) (17. 5.) (12. 17.) ( 5. 19.) (-3. 15.)</pre>
=={{header|REXX}}==
===version 1===
<
* Compute the Convex Hull for a set of points
* Format of the input file:
Line 3,262 ⟶ 5,611:
Trace ?R
Nop
Exit 12</
{{out}}
<pre> 1 x= -9 yl=-3
Line 3,280 ⟶ 5,629:
===version 2===
After learning from https://www.youtube.com/watch?v=wRTGDig3jx8
<
* Compute the Convex Hull for a set of points
* Format of the input file:
Line 3,554 ⟶ 5,903:
*--------------------------------------------------------------------*/
If arg(2)=1 Then say arg(1)
Return lineout(g.0oid,arg(1))</
{{out}}
<pre> 1 x= -9 yl=-3
Line 3,572 ⟶ 5,921:
=={{header|Ruby}}==
{{trans|D}}
<
include Comparable
attr :x, :y
Line 3,638 ⟶ 5,987:
end
main()</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
Line 3,644 ⟶ 5,993:
=={{header|Rust}}==
Calculates convex hull from list of points (f32, f32). This can be executed entirely in the Rust Playground.
<syntaxhighlight lang="rust">
#[derive(Debug, Clone)]
struct Point {
Line 3,711 ⟶ 6,060:
return Point {x:x as f32, y:y as f32};
}
</syntaxhighlight>
{{out}}
<pre>
Line 3,725 ⟶ 6,074:
=={{header|Scala}}==
Scala Implementation to find Convex hull of given points collection. Functional Paradigm followed
<syntaxhighlight lang="scala">
object convex_hull{
def get_hull(points:List[(Double,Double)], hull:List[(Double,Double)]):List[(Double,Double)] = points match{
Line 3,760 ⟶ 6,109:
}
}
</syntaxhighlight>
{{out}}
<pre>
Line 3,777 ⟶ 6,126:
The program is in R7RS Scheme. For CHICKEN, you need to install the '''r7rs''' and '''srfi-132''' eggs.
See also [[#Common_Lisp|Common Lisp]], [[#Standard_ML|Standard ML]], [[#OCaml|OCaml]], and [[#ATS|ATS]]. These implementations were based closely on the Scheme. The last includes proofs of various constraints and of termination of the recursions.
Interestingly, I also derived a [[#Fortran|Fortran]] implementation from this Scheme, which was itself based on Fortran code I wrote over a decade beforehand.
<syntaxhighlight lang="scheme">(define-library (convex-hulls)
(export vector-convex-hull)
Line 3,910 ⟶ 6,263:
(write (list-convex-hull example-points))
(newline)</
{{out}}
Line 3,919 ⟶ 6,272:
<pre>$ csc -O3 -R r7rs convex-hull-task.scm && ./convex-hull-task
((-9 -3) (-3 -9) (19 -8) (17 5) (12 17) (5 19) (-3 15))</pre>
=== A second implementation ===
{{trans|Mercury}}
From the first Scheme implementation, above, I derived an implementation in [[#Standard_ML|Standard ML]]. From that I derived one in [[#Mercury|Mercury]]. Now, starting from that [[#Mercury|Mercury implementation]], I can come back and write a Scheme program more concise than the original.
One could pass these changes along to other languages as well.
If you are using CHICKEN Scheme you will need the '''srfi-1''' egg, in addition to those you needed for the first Scheme implementation.
<syntaxhighlight lang="scheme">(define-library (convex-hulls)
(export convex-hull)
(import (scheme base))
(import (only (srfi 1) fold))
(import (only (srfi 1) append!))
(import (only (srfi 132) list-sort))
(import (only (srfi 132) list-delete-neighbor-dups))
(begin
;;
;; Andrew's monotone chain algorithm for the convex hull in a
;; plane.
;;
;; For a description of the algorithm, see
;; https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=4016938
;;
(define x@ car)
(define y@ cadr)
(define (cross u v)
;; Cross product (as a signed scalar).
(- (* (x@ u) (y@ v)) (* (y@ u) (x@ v))))
(define (point- u v)
(list (- (x@ u) (x@ v)) (- (y@ u) (y@ v))))
(define (point=? u v)
(and (= (x@ u) (x@ v)) (= (y@ u) (y@ v))))
(define (point<? u v)
(let ((xu (x@ u))
(xv (x@ v)))
(or (< xu xv) (and (= xu xv) (< (y@ u) (y@ v))))))
(define (convex-hull points-list)
(let* ((points (list-delete-neighbor-dups
point=? (list-sort point<? points-list)))
(n (length points)))
(cond
((<= n 2) points)
(else
(let ((half-hull (make-vector n)))
(define (cross-test pt j)
(or (zero? j)
(let ((elem-j (vector-ref half-hull j))
(elem-j1 (vector-ref half-hull (- j 1))))
(positive? (cross (point- elem-j elem-j1)
(point- pt elem-j1))))))
(define (construct-half-hull points)
(vector-set! half-hull 0 (car points))
(vector-set! half-hull 1 (cadr points))
(fold (lambda (pt j)
(let loop ((j j))
(if (cross-test pt j)
(let ((j1 (+ j 1)))
(vector-set! half-hull j1 pt)
j1)
(loop (- j 1)))))
1 (cddr points)))
(let* ((lower-hull
;; Leave out the last point, which is the same
;; as the first point of the upper hull.
(let ((j (construct-half-hull points)))
(vector->list half-hull 0 j)))
(upper-hull
;; Leave out the last point, which is the same
;; as the first point of the lower hull.
(let ((j (construct-half-hull (reverse points))))
(vector->list half-hull 0 j))))
(append! lower-hull upper-hull)))))))
)) ;; end of library convex-hulls.
;;
;; A demonstration.
;;
(import (scheme base))
(import (scheme write))
(import (convex-hulls))
(define example-points
'((16 3) (12 17) (0 6) (-4 -6) (16 6)
(16 -7) (16 -3) (17 -4) (5 19) (19 -8)
(3 16) (12 13) (3 -4) (17 5) (-3 15)
(-3 -9) (0 11) (-9 -3) (-4 -2) (12 10)))
(write (convex-hull example-points))
(newline)</syntaxhighlight>
{{out}}
<pre>$ gosh convex-hull-task-2.scm
((-9 -3) (-3 -9) (19 -8) (17 5) (12 17) (5 19) (-3 15))</pre>
Footnote: My Mercury implementation originally had a "fold" in it, the way this Scheme program does, but that got lost during debugging. (The bug turned out to be elsewhere.) The fold and inner loop became one loop that I think is easier to understand. However, here I wanted to show it could be done with a fold.
=={{header|Sidef}}==
{{trans|Raku}}
<
method to_s {
"(#{x}, #{y})"
Line 3,995 ⟶ 6,461:
[-3,15], [-3,-9], [ 0,11], [-9,-3], [-4,-2], [12,10], [14,-9], [1,-9])
say("Convex Hull (#{hull.len} points): ", hull.join(" "))</
{{out}}
<pre>
Line 4,001 ⟶ 6,467:
Convex Hull (9 points): (-3, -9) (1, -9) (14, -9) (19, -8) (17, 5) (12, 17) (5, 19) (-3, 15) (-9, -3)
</pre>
=={{header|Standard ML}}==
{{trans|Scheme}}
{{trans|Fortran}}
{{works with|MLton|20210117}}
{{works with|Poly/ML|5.9}}
<syntaxhighlight lang="sml">(*
* Convex hulls by Andrew's monotone chain algorithm.
*
* For a description of the algorithm, see
* https://en.wikibooks.org/w/index.php?title=Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain&stableid=40169
*)
(*------------------------------------------------------------------*)
(*
* Just enough plane geometry for our purpose.
*)
signature PLANE_POINT =
sig
type planePoint
val planePoint : real * real -> planePoint
val toTuple : planePoint -> real * real
val x : planePoint -> real
val y : planePoint -> real
val == : planePoint * planePoint -> bool
(* Impose a total order on points, making it one that will work for
Andrew's monotone chain algorithm. *)
val order : planePoint * planePoint -> bool
(* Subtraction is really a vector or multivector operation. *)
val subtract : planePoint * planePoint -> planePoint
(* Cross product is really a multivector operation. *)
val cross : planePoint * planePoint -> real
val toString : planePoint -> string
end
structure PlanePoint : PLANE_POINT =
struct
type planePoint = real * real
fun planePoint xy = xy
fun toTuple (x, y) = (x, y)
fun x (x, _) = x
fun y (_, y) = y
fun == ((a : real, b : real), (c : real, d : real)) =
Real.== (a, c) andalso Real.== (b, d)
fun order ((a : real, b : real), (c : real, d : real)) =
Real.< (a, c) orelse (Real.== (a, c) andalso Real.< (b, d))
fun subtract ((a : real, b : real), (c : real, d : real)) =
(a - c, b - d)
fun cross ((a : real, b : real), (c : real, d : real)) =
(a * d) - (b * c)
fun toString (x, y) =
"(" ^ Real.toString x ^ " " ^ Real.toString y ^ ")"
end
(*------------------------------------------------------------------*)
(*
* Rather than rely on compiler extensions for sorting, let us write
* our own.
*
* For no special reason, let us use the Shell sort of
* https://en.wikipedia.org/w/index.php?title=Shellsort&oldid=1084744510
*)
val ciura_gaps = Array.fromList [701, 301, 132, 57, 23, 10, 4, 1]
fun sort_in_place less arr =
let
open Array
fun span_gap gap =
let
fun iloop i =
if length arr <= i then
()
else
let
val temp = sub (arr, i)
fun jloop j =
if j < gap orelse
less (sub (arr, j - gap), temp) then
update (arr, j, temp)
else
(update (arr, j, sub (arr, j - gap));
jloop (j - gap))
in
jloop i;
iloop (i + gap)
end
in
iloop 0
end
in
app span_gap ciura_gaps
end
(*------------------------------------------------------------------*)
(*
* To go with our sort routine, we want something akin to
* array_delete_neighbor_dups of Scheme's SRFI-132.
*)
fun array_delete_neighbor_dups equal arr =
let
open Array
fun loop i lst =
(* Cons a list of non-duplicates, going backwards through
the array so the list will be in forwards order. *)
if i = 0 then
sub (arr, i) :: lst
else if equal (sub (arr, i - 1), sub (arr, i)) then
loop (i - 1) lst
else
loop (i - 1) (sub (arr, i) :: lst)
val n = length arr
in
fromList (if n = 0 then [] else loop (n - 1) [])
end
(*------------------------------------------------------------------*)
(*
* The convex hull algorithm.
*)
fun cross_test (pt_i, hull, j) =
let
open PlanePoint
val hull_j = Array.sub (hull, j)
val hull_j1 = Array.sub (hull, j - 1)
in
0.0 < cross (subtract (hull_j, hull_j1),
subtract (pt_i, hull_j1))
end
fun construct_lower_hull (n, pt) =
let
open PlanePoint
val hull = Array.array (n, planePoint (0.0, 0.0))
val () = Array.update (hull, 0, Array.sub (pt, 0))
val () = Array.update (hull, 1, Array.sub (pt, 1))
fun outer_loop i j =
if i = n then
j + 1
else
let
val pt_i = Array.sub (pt, i)
fun inner_loop j =
if j = 0 orelse cross_test (pt_i, hull, j) then
(Array.update (hull, j + 1, pt_i);
j + 1)
else
inner_loop (j - 1)
in
outer_loop (i + 1) (inner_loop j)
end
val hull_size = outer_loop 2 1
in
(hull_size, hull)
end
fun construct_upper_hull (n, pt) =
let
open PlanePoint
val hull = Array.array (n, planePoint (0.0, 0.0))
val () = Array.update (hull, 0, Array.sub (pt, n - 1))
val () = Array.update (hull, 1, Array.sub (pt, n - 2))
fun outer_loop i j =
if i = ~1 then
j + 1
else
let
val pt_i = Array.sub (pt, i)
fun inner_loop j =
if j = 0 orelse cross_test (pt_i, hull, j) then
(Array.update (hull, j + 1, pt_i);
j + 1)
else
inner_loop (j - 1)
in
outer_loop (i - 1) (inner_loop j)
end
val hull_size = outer_loop (n - 3) 1
in
(hull_size, hull)
end
fun construct_hull (n, pt) =
let
(* Side note: Construction of the lower and upper hulls can be
done in parallel. *)
val (lower_hull_size, lower_hull) = construct_lower_hull (n, pt)
and (upper_hull_size, upper_hull) = construct_upper_hull (n, pt)
val hull_size = lower_hull_size + upper_hull_size - 2
val hull =
Array.array (hull_size, PlanePoint.planePoint (0.0, 0.0))
fun copy_lower i =
if i = lower_hull_size - 1 then
()
else
(Array.update (hull, i, Array.sub (lower_hull, i));
copy_lower (i + 1))
fun copy_upper i =
if i = upper_hull_size - 1 then
()
else
(Array.update (hull, i + lower_hull_size - 1,
Array.sub (upper_hull, i));
copy_upper (i + 1))
in
copy_lower 0;
copy_upper 0;
hull
end
fun plane_convex_hull points_lst =
(* Takes an arbitrary list of points, which may be in any order
and may contain duplicates. Returns an ordered array of points
that make up the convex hull. If the list of points is empty,
the returned array is empty. *)
let
val pt = Array.fromList points_lst
val () = sort_in_place PlanePoint.order pt
val pt = array_delete_neighbor_dups (PlanePoint.==) pt
val n = Array.length pt
in
if n <= 2 then
pt
else
construct_hull (n, pt)
end
(*------------------------------------------------------------------*)
fun main () =
let
open PlanePoint
val example_points =
[planePoint (16.0, 3.0),
planePoint (12.0, 17.0),
planePoint (0.0, 6.0),
planePoint (~4.0, ~6.0),
planePoint (16.0, 6.0),
planePoint (16.0, ~7.0),
planePoint (16.0, ~3.0),
planePoint (17.0, ~4.0),
planePoint (5.0, 19.0),
planePoint (19.0, ~8.0),
planePoint (3.0, 16.0),
planePoint (12.0, 13.0),
planePoint (3.0, ~4.0),
planePoint (17.0, 5.0),
planePoint (~3.0, 15.0),
planePoint (~3.0, ~9.0),
planePoint (0.0, 11.0),
planePoint (~9.0, ~3.0),
planePoint (~4.0, ~2.0),
planePoint (12.0, 10.0)]
val hull = plane_convex_hull example_points
in
Array.app (fn p => (print (toString p);
print " "))
hull;
print "\n"
end;
main ();
(*------------------------------------------------------------------*)
(* local variables: *)
(* mode: sml *)
(* sml-indent-level: 2 *)
(* sml-indent-args: 2 *)
(* end: *)</syntaxhighlight>
{{out}}
Output with MLton:
<pre>$ mlton convex_hull_task.sml && ./convex_hull_task
(~9 ~3) (~3 ~9) (19 ~8) (17 5) (12 17) (5 19) (~3 15)</pre>
Output with Poly/ML (yes, the result is printed twice):
<pre>$ polyc convex_hull_task.sml && ./a.out
(~9.0 ~3.0) (~3.0 ~9.0) (19.0 ~8.0) (17.0 5.0) (12.0 17.0) (5.0 19.0) (~3.0 15.0)
(~9.0 ~3.0) (~3.0 ~9.0) (19.0 ~8.0) (17.0 5.0) (12.0 17.0) (5.0 19.0) (~3.0 15.0)</pre>
When you use Poly/ML, if you want the result printed only once, comment out the call to '''main()''' near the bottom of the source file. The call is redundant, because a program compiled with Poly/ML automatically calls '''main()'''.
=={{header|Swift}}==
Line 4,006 ⟶ 6,782:
{{trans|Rust}}
<
public var x: Double
public var y: Double
Line 4,082 ⟶ 6,858:
print("Input: \(points)")
print("Output: \(calculateConvexHull(fromPoints: points))")</
{{out}}
Line 4,093 ⟶ 6,869:
Translation of: https://en.wikibooks.org/wiki/Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain#Python
<
namespace eval test_convex_hull {
Line 4,145 ⟶ 6,921:
puts [convex_hull $tpoints]
}
</syntaxhighlight>
{{out}}
<pre>Test points:
Line 4,154 ⟶ 6,930:
=={{header|Visual Basic .NET}}==
{{trans|C#}}
<
Module Module1
Line 4,245 ⟶ 7,021:
End Sub
End Module</
{{out}}
<pre>Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]</pre>
Line 4,253 ⟶ 7,029:
{{libheader|Wren-sort}}
{{libheader|Wren-trait}}
{{libheader|Wren-iterate}}
<syntaxhighlight lang="wren">import "./
import "./trait" for Comparable
import "./iterate" for Stepped
class Point is Comparable {
Line 4,302 ⟶ 7,080:
Point.new(-3, -9), Point.new( 0, 11), Point.new(-9, -3), Point.new(-4, -2), Point.new(12, 10)
]
System.print("Convex Hull: %(convexHull.call(pts))")</
{{out}}
Line 4,308 ⟶ 7,086:
Convex Hull: [(-9, -3), (-3, -9), (19, -8), (17, 5), (12, 17), (5, 19), (-3, 15)]
</pre>
=={{header|XPL0}}==
<syntaxhighlight lang "XPL0">func real CosAng(A, B); \Return cosine of angle between vectors A and B
int A, B; \Cos(Ang) = DotProd / (|A|*|B|)
real DotProd, Magnitude;
[DotProd:= float( A(0)*B(0) + A(1)*B(1));
Magnitude:= sqrt(float( sq(B(0)-A(0)) + sq(B(1)-A(1)) ));
return DotProd / Magnitude;
];
proc ConvexHull(N, P);
int N, P;
int Min, I, HullI, EndI, A(2), B(2), J, SJ;
real Ang, MinAng;
[Min:= -1>>1; \find index of point with smallest X coordinate
for I:= 0 to N-1 do
if P(I,0) < Min then
[Min:= P(I,0); HullI:= I];
EndI:= HullI;
A(0):= 0; A(1):= -1;
repeat ChOut(0, ^();
IntOut(0, P(HullI,0)); ChOut(0, ^,); IntOut(0, P(HullI,1));
ChOut(0, ^));
MinAng:= -1e12; \find index of point with smallest diverting angle
for J:= 0 to N-1 do
[B(0):= P(J,0) - P(HullI,0); B(1):= P(J,1) - P(HullI,1);
Ang:= CosAng(A, B);
if Ang > MinAng and J # HullI then
[MinAng:= Ang; SJ:= J];
];
A(0):= P(SJ,0) - P(HullI,0); A(1):= P(SJ,1) - P(HullI,1);
HullI:= SJ;
if HullI # EndI then Text(0, ", ");
until HullI = EndI;
];
ConvexHull(20, [[16,3], [12,17], [0,6], [-4,-6], [16,6], [16,-7], [16,-3],
[17,-4], [5,19], [19,-8], [3,16], [12,13], [3,-4], [17,5],
[-3,15], [-3,-9], [0,11], [-9,-3], [-4,-2], [12,10]])
</syntaxhighlight>
{{out}}
<pre>
(-9,-3), (-3,-9), (19,-8), (17,5), (12,17), (5,19), (-3,15)</pre>
=={{header|zkl}}==
<
// https://en.wikipedia.org/wiki/Graham_scan, O(n log n)
// http://www.geeksforgeeks.org/convex-hull-set-2-graham-scan/
Line 4,341 ⟶ 7,163:
fcn ccw(a,b,c){ // a,b,c are points: (x,y)
((b[0] - a[0])*(c[1] - a[1])) - ((b[1] - a[1])*(c[0] - a[0]))
}</
<
T(16, -7), T(16,-3),T(17,-4), T(5,19), T(19,-8),
T(3,16), T(12,13), T(3,-4), T(17,5), T(-3,15),
Line 4,348 ⟶ 7,170:
.apply(fcn(xy){ xy.apply("toFloat") }).copy();
hull:=grahamScan(pts);
println("Convex Hull (%d points): %s".fmt(hull.len(),hull.toString(*)));</
{{out}}
<pre>
| |||