Jump to content

Dijkstra's algorithm: Difference between revisions

Line 758:
[https://en.wikipedia.org/w/index.php?title=Dijkstra%27s_algorithm&oldid=1117533081#Using_a_priority_queue a Wikipedia article about Dijkstra's algorithm], although I use a different method to determine whether a queue entry is obsolete.
 
I prove that the algorithm terminates and that the priority queue has at least enough storage. The priority queue is a binary heap implemented as an array.
 
<syntaxhighlight lang="ats">
Line 783:
macdef infinity = $extval (flt, "INFINITY")
 
prfn
(*------------------------------------------------------------------*)
mul_compare_lte
(* A very inefficient implementation of a priority queue, for the sake
{i, j, n : nat | i <= j}
of demonstrating Dijkstra's algorithm. The queue is represented as
an association list. * ()
:<prf> [i * n <= j * n] void =
 
mul_gte_gte_gte {j - i, n} ()
typedef pqueue (priority_t : t@ype+,
value_t : t@ype+,
n : int) =
list (@(priority_t, value_t), n)
 
prfn
mul_compare_lt
lemma_pqueue_param
{ni, j, n : int | 0 <= i; i < j; 1 <= n}
{priority_t, value_t : t@ype}
(pq : pqueue (priority_t, value_t, n))
:<prf> [0 <= n]
void =
lemma_list_param pq
 
extern fn {priority_t : t@ype}
pqueue$cmp :
(priority_t, priority_t) -<> int
 
implement pqueue$cmp<double> (x, y) = compare (x, y)
 
fn
pqueue_make_empty
{priority_t : t@ype}
{value_t : t@ype}
()
:<!wrtprf> pqueue[i (priority_t,* value_t,n 0)< j * n] void =
mul_compare_lte {i + 1, j, n} ()
NIL
 
fn {}
pqueue_is_empty
{n : int}
{priority_t : t@ype}
{value_t : t@ype}
(pq : pqueue (priority_t, value_t, n))
:<> bool (n == 0) =
list_is_nil pq
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_insert
{n : int}
(pq : &pqueue (priority_t, value_t, n)
>> pqueue (priority_t, value_t, n + 1),
priority : priority_t,
value : value_t)
:<!wrt> void =
let
prval () = lemma_pqueue_param pq
in
pq := @(priority, value) :: pq
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_extract_min
{n : pos}
(pq : &pqueue (priority_t, value_t, n)
>> pqueue (priority_t, value_t, n - 1))
:<!wrt> @(priority_t, value_t) =
let
fun
loop {m1, m2 : nat | m1 + m2 + 1 == n}
.<m1>.
(current_min : @(priority_t, value_t),
to_be_scanned : list (@(priority_t, value_t), m1),
already_scanned : list (@(priority_t, value_t), m2))
:<> @(@(priority_t, value_t),
list (@(priority_t, value_t), n - 1)) =
case+ to_be_scanned of
| NIL => @(current_min, already_scanned)
| head :: tail =>
if pqueue$cmp<priority_t> (head.0, current_min.0) < 0 then
loop (head, tail, current_min :: already_scanned)
else
loop (current_min, tail, head :: already_scanned)
 
val+ head :: tail = pq
val @(retval, pq1) = loop (head, tail, NIL)
in
pq := pq1;
retval
end
 
(* A little unit testing of the priority queue. *)
var pq = pqueue_make_empty ()
val- true = pqueue_is_empty pq
val () = pqueue_insert (pq, 3.0, 3)
val () = pqueue_insert (pq, 5.0, 5)
val () = pqueue_insert (pq, 1.0, 1)
val () = pqueue_insert (pq, 2.0, 2)
val () = pqueue_insert (pq, 4.0, 4)
val- false = pqueue_is_empty pq
val- @(1.0, 1) = pqueue_extract_min<double> pq
val- @(2.0, 2) = pqueue_extract_min<double> pq
val- @(3.0, 3) = pqueue_extract_min<double> pq
val- @(4.0, 4) = pqueue_extract_min<double> pq
val- @(5.0, 5) = pqueue_extract_min<double> pq
val- true = pqueue_is_empty pq
 
(*------------------------------------------------------------------*)
Line 1,047 ⟶ 957:
end
end
 
(*------------------------------------------------------------------*)
(* A binary-heap priority queue, similar to the Pascal in Robert
Sedgewick, "Algorithms", 2nd ed. (reprinted with corrections),
1989. Note that Sedgewick does an extract-max, whereas we do an
extract-min.
Niklaus Wirth, within the heapsort implementation of "Algorithms +
Data Structures = Programs", has, I will note, some Pascal code
that is practically the same as Sedgewick's. Can we trace that code
back farther to Algol?
We do not have "goto" for Sedgewick's "downheap" (or Wirth's
"sift"), but do have mutual tail call as an obvious alternative to
the "goto". Nevertheless, because the code "jumped to" is small, I
simply use a macro to duplicate it. *)
 
typedef pqueue (priority_t : t@ype+,
value_t : t@ype+,
n : int,
n_max : int) =
[n <= n_max]
@{
arr = arrayref (@(priority_t, value_t), n_max + 1),
n = size_t n,
n_max = size_t n_max
}
 
prfn
lemma_pqueue_param
{n_max : int}
{n : int}
{priority_t, value_t : t@ype}
(pq : pqueue (priority_t, value_t, n, n_max))
:<prf> [0 <= n; n <= n_max] void =
lemma_g1uint_param (pq.n)
 
extern praxi
lemma_pqueue_size
{n_max : int}
{n : int}
{priority_t, value_t : t@ype}
(pq : pqueue (priority_t, value_t, n, n_max))
:<prf> [n1 : int | n1 == n] void
 
extern fn {priority_t : t@ype}
pqueue$cmp :
(priority_t, priority_t) -<> int
 
extern fn {priority_t : t@ype}
pqueue$priority_min :
() -<> priority_t
 
implement pqueue$cmp<double> (x, y) = compare (x, y)
implement pqueue$priority_min<double> () = neg infinity
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_make_empty
{n_max : int}
(n_max : size_t n_max,
arbitrary_entry : @(priority_t, value_t))
:<!wrt> pqueue (priority_t, value_t, 0, n_max) =
let
(* Currently an array is allocated whose size is the proven
bound. It might be better to use a smaller array and allow
reallocation up to this maximum size, or to break the array
into pieces. *)
 
prval () = lemma_g1uint_param n_max
val arr =
arrayref_make_elt<@(priority_t, value_t)>
(succ n_max, arbitrary_entry)
in
@{arr = arr,
n = i2sz 0,
n_max = n_max}
end
 
fn {}
pqueue_clear
{n_max : int}
{n : int}
{priority_t : t@ype}
{value_t : t@ype}
(pq : &pqueue (priority_t, value_t, n, n_max)
>> pqueue (priority_t, value_t, 0, n_max))
:<!wrt> void =
let
prval () = lemma_g1uint_param (pq.n_max)
in
pq := @{arr = pq.arr,
n = i2sz 0,
n_max = pq.n_max}
end
 
fn {}
pqueue_is_empty
{n_max : int}
{n : int}
{priority_t : t@ype}
{value_t : t@ype}
(pq : pqueue (priority_t, value_t, n, n_max))
:<> bool (n == 0) =
(pq.n) = i2sz 0
 
fn {}
pqueue_size
{n_max : int}
{n : int}
{priority_t : t@ype}
{value_t : t@ype}
(pq : pqueue (priority_t, value_t, n, n_max))
:<> size_t n =
pq.n
 
fn {priority_t : t@ype}
{value_t : t@ype}
_upheap {n_max : pos}
{n : int | n <= n_max}
{k0 : nat | k0 <= n}
(arr : arrayref (@(priority_t, value_t), n_max + 1),
k0 : size_t k0)
:<!refwrt> void =
let
macdef lt (x, y) = (pqueue$cmp<priority_t> (,(x), ,(y)) < 0)
macdef prio x = ,(x).0
 
val entry = arr[k0]
 
fun
loop {k : nat | k <= n}
.<k>.
(k : size_t k)
:<!refwrt> void =
if k = i2sz 0 then
arr[k] := entry
else
let
val kh = half k
in
if (prio entry) \lt (prio arr[kh]) then
begin
arr[k] := arr[kh];
loop kh
end
else
arr[k] := entry
end
in
arr[0] := @(pqueue$priority_min<priority_t> (), arr[0].1);
loop k0
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_insert
{n_max : int}
{n : int | n < n_max}
(pq : &pqueue (priority_t, value_t, n, n_max)
>> pqueue (priority_t, value_t, n + 1, n_max),
entry : @(priority_t, value_t))
:<!refwrt> void =
let
prval () = lemma_g1uint_param (pq.n)
val arr = pq.arr
and n1 = succ (pq.n)
in
arr[n1] := entry;
_upheap {n_max} {n + 1} (arr, n1);
pq := @{arr = arr,
n = n1,
n_max = pq.n_max}
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
_downheap {n_max : pos}
{n : pos | n <= n_max}
(arr : arrayref (@(priority_t, value_t), n_max + 1),
n : size_t n)
:<!refwrt> void =
let
macdef lt (x, y) = (pqueue$cmp<priority_t> (,(x), ,(y)) < 0)
macdef prio x = ,(x).0
 
val entry = arr[1]
and nh = half n
 
fun
loop {k : pos | k <= n}
.<n - k>.
(k : size_t k)
:<!refwrt> void =
let
macdef move_data i =
if (prio entry) \lt (prio arr[,(i)]) then
arr[k] := entry
else
begin
arr[k] := arr[,(i)];
loop ,(i)
end
in
if nh < k then
arr[k] := entry
else
let
stadef j = 2 * k
prval () = prop_verify {j <= n} ()
val j : size_t j = k + k
in
if j < n then
let
stadef j1 = j + 1
prval () = prop_verify {j1 <= n} ()
val j1 : size_t j1 = succ j
in
if ~((prio arr[j]) \lt (prio arr[j1])) then
move_data j1
else
move_data j
end
else
move_data j
end
end
in
loop (i2sz 1)
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_peek
{n_max : int}
{n : pos | n <= n_max}
(pq : pqueue (priority_t, value_t, n, n_max))
:<!ref> @(priority_t, value_t) =
let
val arr = pq.arr
in
arr[1]
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_delete
{n_max : int}
{n : pos | n <= n_max}
(pq : &pqueue (priority_t, value_t, n, n_max)
>> pqueue (priority_t, value_t, n - 1, n_max))
:<!refwrt> void =
let
val @{arr = arr,
n = n,
n_max = n_max} = pq
in
if i2sz 0 < pred n then
begin
arr[1] := arr[n];
_downheap {n_max} {n - 1} (arr, pred n)
end;
pq := @{arr = arr,
n = pred n,
n_max = n_max}
end
 
fn {priority_t : t@ype}
{value_t : t@ype}
pqueue_extract
{n_max : int}
{n : pos | n <= n_max}
(pq : &pqueue (priority_t, value_t, n, n_max)
>> pqueue (priority_t, value_t, n - 1, n_max))
:<!refwrt> @(priority_t, value_t) =
let
val retval = pqueue_peek<priority_t><value_t> {n_max} {n} pq
in
pqueue_delete<priority_t><value_t> {n_max} {n} pq;
retval
end
 
local (* A little unit testing of the priority queue
implementation. *)
#define NMAX 10
in
var pq = pqueue_make_empty<double><int> (i2sz NMAX, @(0.0, 0))
val- true = pqueue_is_empty pq
val- true = (pqueue_size pq = i2sz 0)
val () = pqueue_insert (pq, @(3.0, 3))
val () = pqueue_insert (pq, @(5.0, 5))
val () = pqueue_insert (pq, @(1.0, 1))
val () = pqueue_insert (pq, @(2.0, 2))
val () = pqueue_insert (pq, @(4.0, 4))
val- false = pqueue_is_empty pq
val- true = (pqueue_size pq = i2sz 5)
val- @(1.0, 1) = pqueue_extract<double> pq
val- @(2.0, 2) = pqueue_extract<double> pq
val- @(3.0, 3) = pqueue_extract<double> pq
val- @(4.0, 4) = pqueue_extract<double> pq
val- @(5.0, 5) = pqueue_extract<double> pq
val- true = pqueue_is_empty pq
val- true = (pqueue_size pq = i2sz 0)
end
 
(*------------------------------------------------------------------*)
Line 1,060 ⟶ 1,274:
(* Returns total-costs and previous-hops arrays. *)
:<!refwrt> @(arrayref (flt, n),
arrayref ([i : nat | i <= n] size_t i, n)) =
let
prval () = lemma_matrixref_param adj_matrix
Line 1,067 ⟶ 1,281:
typedef defined_index_t = [i : nat | i < n] size_t i
val index_t_undefined : size_t n = n
 
val arbitrary_pq_entry : @(flt, defined_index_t) =
@(0.0, i2sz 0)
 
val prev = arrayref_make_elt<index_t> (n, index_t_undefined)
and cost = arrayref_make_elt<flt> (n, infinity)
val () = cost[source] := zero
 
(* The priority queue never gets larger than m_max. There is code
below that proves this; thus there is no risk of overrunning
the queue's storage (unless the queue implementation itself is
made unsafe). FIXME: Is it possible to prove a tighter bound on
the size of the priority queue? *)
stadef m_max = (n * n) + n + n
prval () = mul_pos_pos_pos (mul_make {n, n} ())
prval () = prop_verify {n + n < m_max} ()
val m_max : size_t m_max = (n * n) + n + n
 
typedef pqueue_t (m : int) =
[0 <= m]; pqueuem (flt,<= defined_index_t, m)m_max]
pqueue (flt, defined_index_t, m, m_max)
typedef pqueue_t =
[m : int] pqueue_t m
var pq : pqueue_t = pqueue_make_empty ()
 
fn
pq_make_empty ()
:<!wrt> pqueue_t 0 =
(* Create a priority queue, whose maximum size is our proven
upper bound on the queue size. *)
pqueue_make_empty<flt><defined_index_t>
(m_max, arbitrary_pq_entry)
 
var pq = pq_make_empty ()
val active = arrayref_make_elt<bool> (n, true)
var num_active : [i : nat | i <= n] size_t i = n
Line 1,084 ⟶ 1,320:
fill_pq {i : nat | i <= n}
.<n - i>.
(pq : &pqueue_t i >> _pqueue_t n,
i : size_t i)
:<!refwrt> void =
if i <> n then
begin
pqueue_insert {m_max} {i} (pq, @(cost[i], i));
fill_pq {i + 1} (pq, succ i)
end
 
fun
extract_min
{m0 : pos | m0 + n <= : natm_max}
.<m0>.
(pq : &pqueue_t m0 >> pqueue_t m1)
:<!refwrt> #[m1 : nat | m1 <= m0]
Option @(flt, defined_index_t) =
let
if pqueue_is_empty pq then
Noneval @(priority, vertex) =
pqueue_extract<flt><defined_index_t> {m_max} {m0} pq
else
letin
if val @(priority, active[vertex)] =then
@(priority, vertex)
pqueue_extract_min<flt><defined_index_t> pq
else if pqueue_is_empty {m_max} pq then
in
if ~active[vertex] thenarbitrary_pq_entry
extract_min pqelse
elseextract_min pq
end
Some @(priority, vertex)
end
 
fun
main_loop {num_active0 : nat | num_active0 <= n}
{qsize0 : nat}
{qlimit0 : int | 0 <= qlimit0;
qsize0 <= qlimit0 + n}
.<num_active0>.
(pq* The pf_qlimit0 prop helps us prove a bound :on &pqueue_t >> pqueue_t,the
size of the priority queue. We need it because the
proof capabilities built into ATS have very limited
ability to handle multiplication. *)
(pf_qlimit0 : MUL (n - num_active0, n, qlimit0) |
pq : &pqueue_t qsize0 >> pqueue_t 0,
num_active : &size_t num_active0 >> size_t num_active1)
:<!refwrt> #[num_active1 : nat | num_active1 <= num_active0]
void =
if num_active <>= i2sz 0 then
case+ extract_minpqueue_clear pq of
else if |pqueue_is_empty None{m_max} (){qsize0} =>pq then
let (* This branch should nevernot behappen. reached, but the*)
val- false = true
corresponding sanity check is done further below. Thus
in
the proof of termination of this loop does not depend on
an exception. *)end
()else
| Some @(priority, u) =>let
letprval () = mul_elim pf_qlimit0
valprval () = active[u] := false
andprop_verify {qsize0 <= ((n - num_active0) =* num_activen) :=+ predn} num_active()
prval () = mul_compare_lt {n - num_active0, n, n} ()
prval () = prop_verify {qsize0 < m_max} ()
 
val @(priority, funu) = extract_min pq
prval [qsize : int] () = lemma_pqueue_size {m_max} pq
loop_over_neighbors
prval () = lemma_pqueue_param {m_max} {vqsize} : nat | v <= n}pq
prval () = prop_verify {qsize .<n -qsize0} v>.()
prval () = prop_verify {qsize < m_max} (pq : &pqueue_t >> pqueue_t,)
v : size_t v)
:<!refwrt> void =
if v = n then
()
else if ~active[v] then
loop_over_neighbors (pq, succ v)
else
let
val alternative = cost[u] + adj_matrix[u, n, v]
in
if alternative < cost[v] then
begin
cost[v] := alternative;
prev[v] := u;
 
val (* Rather than lower the) priority= ofactive[u] v,:= thisfalse
and () = num_active := pred num_active
implementation inserts a new entry for v and
ignores obsolete queue entries. Queue entries
are obsolete if the vertex's entry in the
"active" array is false. *)
 
fun
pqueue_insert<flt><defined_index_t>
(pq, alternative, v)loop_over_vertices
end;{v : nat | v <= n}
loop_over_neighbors (pq, succ{m0 : nat | qsize <= m0; m0 <= qsize + v)}
end .<n - v>.
(pq : &pqueue_t m0 >> pqueue_t m1,
in
loop_over_neighbors (pq, i2sz 0 v : size_t v);
:<!refwrt> #[m1 : int | qsize <= m1; m1 <= qsize + n]
main_loop {num_active0 - 1} (pq, num_active)
end void =
if v = n then
()
else if ~active[v] then
loop_over_vertices {v + 1} {m0} (pq, succ v)
else
let
val alternative = cost[u] + adj_matrix[u, n, v]
in
if alternative < cost[v] then
let
prval () = prop_verify {m0 < m_max} ()
in
cost[v] := alternative;
prev[v] := u;
 
(* Rather than lower the priority of v, this
val () = fill_pq (pq, i2sz 0)
implementation inserts a new entry for v and
val () = main_loop {n} (pq, num_active)
ignores obsolete queue entries. Queue entries
are obsolete if the vertex's entry in the
"active" array is false. *)
 
pqueue_insert<flt><defined_index_t>
(* Sanity check. *)
{m_max} {m0}
val- true = iseqz num_active
(pq, @(alternative, v));
 
loop_over_vertices {v + 1} {m0 + 1}
(pq, succ v)
end
else
loop_over_vertices {v + 1} {m0} (pq, succ v)
end
 
val () = loop_over_vertices {0} {qsize} (pq, i2sz 0)
in
main_loop {num_active0 - 1}
(MULind pf_qlimit0 | pq, num_active)
end
in
fill_pq {0} (pq, i2sz 0);
main_loop {n} {n} (MULbas () | pq, num_active);
@(cost, prev)
end
Line 1,276 ⟶ 1,537:
 
{{out}}
<pre>$ patscc -O3 -DATS_MEMALLOC_GCBDW dijkstra-algorithm.dats -lgc && ./a.out
The requested paths:
a -> c -> d -> e (cost = 26.000000)
1,448

edits

Cookies help us deliver our services. By using our services, you agree to our use of cookies.