ATS и smoothsort: продолжение
Feb. 2nd, 2012 02:50 pm![[personal profile]](https://www.dreamwidth.org/img/silk/identity/user.png){kind=small}
thedeemonАкт 5. Мемоизация доказательств
Для обеспечения коррекности структур и убеждения компилятора в том, что мы не делаем ничего предосудительного, у нас используются доказательства утверждений LEON(k,sz) о том, что sz - действительно k-е число Леонардо, L[k]=sz. В одном из прошлых постов этой серии мы видели, как построить функцию, вычисляющую нужное число Леонардо и поставляющую вместе с ним доказательство того, что это правильное число. Там было несколько вариантов, и лучший из них имел сложность O(k). Если такую функцию задействовать в smoothsort'e, его сложность резко возрастет и из O(n * log n) превратится во что-то вроде O(n * log n * log n). Нужна мемоизация, и оригинальный алгоритм подразумевает использование массива с уже вычисленными числами Леонардо, так получение нужного числа отнимает O(1) операций. Но такой массив имеет размер опять же O(log n), поэтому рассказы про потребление памяти О(1) - враки. Попытка добавить мемоизацию в функцию вычисления чисел Леонардо сталкивается с проблемой. Функция определена как
fun leon {n:nat} (n : int n) : [x:nat] (LEON(n,x) | int x) тип ее результата зависит от значения аргумента. Зависимые типы, понимаш. Если мы захотим сохранить результаты в массиве, то каждый его элемент будет иметь свой тип! А так нельзя - в массивах все элементы одного типа. Не знаю, как эту проблему решают во взрослых языках с зависимыми типами, а в ATS можно воспользоваться хаком. Сперва дадим типу результата функции leon свое имя:
typedef L(n:int) = [x:nat | x > 0] (LEON(n,x) | int x)
Определим интерфейс для функций сохранения вычисленных значений в массиве-кэше и извлечения оттуда.
extern fun leon_cache_get {n:nat} (n: int n): Option(L(n)) = "leon_cache_get"
extern fun leon_cache_set {n:nat; x:int} (n: int n, v: L(n)): void = "leon_cache_set"
А теперь опишем похожий тип, но без параметра. L0 - "существует натуральное число n, для которого L(n)". И опишем аналогичные функции с этим типом.
typedef L0 = [n:nat] L(n)
extern fun leon_cache_get0 {n:nat} (n: int n): Option(L0) = "leon_cache_get"
extern fun leon_cache_set0 {n,x:nat} (n: int n, v: L0): void = "leon_cache_set"
Теперь тип результата не зависит от значения, и мы можем сохранить результат в массиве:
#define LSIZE 43 val Lcache : array(Option(L0), LSIZE) = array_make_elt(LSIZE, None) implement leon_cache_get0(n) = if n < LSIZE then Lcache[n] else None implement leon_cache_set0(n, v) = if n < LSIZE then Lcache[n] := Some v else ()
И тут срабатывает хак: когда мы после декларации функции написали = "leon_cache_get", мы задали имя функции на Си, в которую та компилируется. Задав обоим вариантам функций работы с кэшем одинаковые Си-имена, мы сделали так, что реализация leon_cache_get0 является одновременно и реализацией leon_cache_get. Компилятор это устраивает, т.к. при трансляции в Си вся шелуха статических условий и деталей исчезает, остаются только реальные данные, а они в обоих вариантах совпадают - на входе инт и на выходе инт. В результате мы получили функции для сохранения и извлечения вычисленных значений зависимых типов за О(1), что и хотелось. С ними функция вычисления чисел Леонардо с мемоизацией выглядит так:
fun leon {n:nat} (n: int n): L(n) =
if n < 2 then (LEONbase | 1)
else
case+ leon_cache_get n of
| Some x => x
| None () => let
val (pa | a) = leon(n-2)
val (pb | b) = leon(n-1)
val res = (LEONind(pa, pb) | a + b + 1)
in
leon_cache_set(n, res);
res
end
Насколько я знаю, ни в одном мануале такой трюк не описан, но он был использован автором языка в одном из примеров на сайте.
Акт 6. Метаморфозы
Выше мы определили три типа структур, описывающих кучу Леонардо, отличающиеся тем, какие свойства для кучи выполнены (и доказаны). Если куча просто имеет правильную форму и положение, то это heap(r,k,sz), где r - позиция корня, k - индекс (номер числа Леонардо) и sz - размер. Если для нее доказано утверждение MAXRIGHT(r-sz+1, r), т.е. корневой элемент является максимальным во всей куче (выполнено свойство кучи), то это heap1. Когда же эта куча - звено правильного ряда куч, для которого выполнено свойство ряда, выраженное как MAXRIGHT(0,r), то это heap2. Ряд куч описан типом heaps(m), где m - суммарное число элементов в нем, и является списком из heap2. В первой части алгоритма мы в цикле поэлементно растим кучу от пустого списка до heaps(n), где n - размер массива, путем либо добавления новой кучи из одного элемента, либо объединения двух последних куч в одну с добавляемым элементом в качестве корня.
Маленькую кучу из одного элемента будет создавать функция small_heap:
fn small_heap {r:nat} (r:int r, prev_k : int):<> [k:two] heap(r,k,1) =
if prev_k = 1 then #[.. | @{root=r, k=0, sz=1, pf_sz=LEONbase, pf_r = LHpf, pf_gc = GCsmall}]
else #[.. | @{root=r, k=1, sz=1, pf_sz=LEONbase, pf_r = LHpf, pf_gc = GCsmall}]
Ее индекс - 0 или 1 в зависимости от индекса предыдущей кучи. Поскольку размер ее равен 1, доказательства утверждений про размер тривиальны (используются базы индукции). Синтаксическая конструкция #[.. | some data ] здесь превращает значение типа heap(r,k,1) в значение типа [k:two] heap(r,k,1), т.е. добавляет статической информации, в данном случае о том, что k:two, т.е. k:nat и k < 2. Эту конструкцию вы не найдете в документации на язык, это пример тайного знания, передаваемого лично от мастера к ученику.
Для кучи из одного элемента свойство кучи выполнено автоматически. Теперь, чтобы ее можно было добавить в список-ряд, нужно превратить ее в heap2, т.е. обеспечить выполненность свойства ряда, для чего может потребоваться переставить местами какие-то элементы. Свойства кучи и ряда формулируются через отношения между корнями куч, поэтому для выполнения этих свойств нужно обеспечить правильный порядок среди корней, т.е. меняться местами будут только корни каких-то куч. Если у большой кучи, для которой было выполнено свойство кучи, меняется значение в ее корне, то она теряет свое свойство, но две входящих в нее кучи поменьше не затронуты, они свои свойства сохраняют. Чтобы восстановить свойство кучи, нужно сравнить новое значение корня с корнями подкуч и при необходимости произвести перестановку, для этого нужно эти подкучи получить в виде нужных структур. Наши структуры для кучи не содержат подкучи в явном виде, иначе объем требуемой памяти под них был бы заметно больше, чем сам сортируемый массив. Поэтому нам нужна функция извлечения из большой кучи двух ее подкуч. Поскольку большие кучи мы всегда строим из правильных куч с выполненными свойствами, то при извлечении подкуч они должны быть с выполненными свойствами, т.е. они должны иметь тип heap1, содержащий доказательство MAXRIGHT(r-sz+1, r). Но откуда взять это доказательство? Из аналогичного MAXRIGHT для большой кучи оно не следует. Тут мы немного срежем угол и зададим нужное свойство аксиомой. При сборке большой кучи из маленьких мы потребуем выполненности свойств кучи для обеих подкуч. Для этого в определение кучи добавлено утверждение GOODCHILDREN. Для маленькой кучи оно выполнено бесплатно, это говорит конструктор GCsmall, определенный для k=0 и k=1 без пререквизитов.
// a heap is either a single element or contains proper sub-heaps
dataprop GOODCHILDREN(int) =
| {k:two} GCsmall(k)
| {k:nat | k > 1; a,b,c:nat} GCbig(k) of (MAXRIGHT(a,b), MAXRIGHT(b+1,c))
extern praxi good_children_mr {a,b,c,d:nat} {k:nat | k > 1}
(pf : GOODCHILDREN(k)): (MAXRIGHT(a,b), MAXRIGHT(c,d))
extern praxi good_child {d:nat | d < 3} {k:nat | k > 1} (pf : GOODCHILDREN(k)): GOODCHILDREN(k-d)
А при разделении кучи и извлечении из нее подкуч выполненность свойств кучи для них будет следовать из доказательства GOODCHILDREN для большой кучи аксиоматично. Т.е. все честно, и при разделении большой кучи доказанное свойство кучи для подкучи может быть получено только если оно было выполнено при создании этой большой кучи. Тогда функция разделения большой кучи и извлечения подкуч выглядит так:
fn split_heap {r,k,sz:nat | k > 1} (h: heap(r,k,sz)):
[cr1,csz1, cr2,csz2:nat | cr1 - csz1 == r - sz; cr2 == r-1; cr2 == cr1 + csz2; csz1 > 0; csz2 > 0 ]
(heap1(cr1, k-1, csz1), heap1(cr2, k-2, csz2)) = let
val (pfsz1 | child_size1) = leon(h.k - 1)
prval () = lh_use h.pf_r
prval () = leon_mono(pfsz1, h.pf_sz)
val child_size2 = h.sz - 1 - child_size1
prval pfsz2 = leon_size_dif(h.pf_sz, pfsz1)
prval pfgc1 = good_child {1} h.pf_gc
prval pfgc2 = good_child {2} h.pf_gc
val left_h = @{root = h.root - h.sz + child_size1, k = h.k - 1, sz = child_size1,
pf_sz = pfsz1, pf_r = LHpf, pf_gc = pfgc1}
val right_h = @{root = h.root - 1, k = h.k - 2, sz = child_size2,
pf_sz = pfsz2, pf_r = LHpf, pf_gc = pfgc2}
prval (pfmr1, pfmr2) = good_children_mr h.pf_gc
in #[.. | (@{hp=left_h, pf_mr=pfmr1}, @{hp=right_h, pf_mr=pfmr2})]
end
Она из большой кучи типа heap получает две ее подкучи типа heap1. Размер первой подкучи вычисляется через числа Леонардо, размер второй - вычитанием размера первого из размера большой кучи минус 1. Доказательство того, что это нужное число Леонардо нам дает доказанная выше теорема leon_size_dif. Доказательство того, что новый корень не окажется за пределами массива, дают теоремы lh_use (про связь корня и размера большой кучи) и leon_mono (про то, что размер подкучи меньше размера большой кучи, т.к. это предыдущее число Леонардо). Доказательства свойств кучи для подкуч дает аксиома good_children_mr. Помимо самих структур для подкуч функция поставляет статическую информацию о том, как связаны в них значения, в виде пяти равенств и неравенств.
Теперь, если наш ряд состоит из двух куч с соседними убывающими индексами, и мы добавляем новый элемент, то из этих двух куч и элемента мы легко можем сделать одну большую кучу типа heap, т.е. имеющую правильную форму, имеющую правильных детей (GOODCHILDREN), но пока еще без выполненности свойства кучи. Аналогично, если у нас была хорошая куча со всеми свойствами, но у нее поменялось значение в корне, то опять она стала простой heap, и нужно восстановить свойство кучи. Для восстановления свойства кучи и получения структуры типа heap1 нам потребуются такие функции:
datatype compare_res(int,int) =
| {i,j:nat} LeftGr(i,j) of (LTE(j,i) | int)
| {i,j:nat} RightGr(i,j) of (LTE(i,j) | int)
fn {a:t@ype} compare_elements {n,i,j:nat | j < n; i <= j}
(A: array(a,n), i : int i, j : int j, gt : gt(a)) : compare_res(i,j) =
if A[i] \gt A[j] then LeftGr (LTEcompared | 0)
else RightGr (LTEcompared | 1)
fun restore_heap_prop {r,k,sz:nat | r < n} (h: heap(r,k,sz)):<cloref1> heap1(r,k,sz) =
if h.k < 2 then
let prval () = leon_base_is1 h.pf_sz in @{hp = h, pf_mr = MRsingle} end
else let
val (left, right) = split_heap h in
restore_big_heap_prop(left, right, h)
end
and restore_big_heap_prop {r,k,sz, cr1,csz1, cr2,csz2:nat | r < n; k > 1;
cr1 - csz1 == r - sz; cr2 == r-1; cr2 == cr1 + csz2}
(left: heap1(cr1, k-1, csz1), right: heap1(cr2, k-2, csz2), h: heap(r,k,sz)):<cloref1> heap1(r,k,sz) = let
val start = left.hp.root - left.hp.sz + 1
prval () = lh_use left.hp.pf_r
prval () = leon_positive right.hp.pf_sz
in
case+ compare_elements(A, left.hp.root, right.hp.root, gt) of
| LeftGr(pf_lte_rl | _) => let
val (pf_lte_lh | swapped) = order_elements(A, left.hp.root, h.root, gt)
val rstart = right.hp.root - right.hp.sz + 1
prval pf_lte_rh = lte_trans(pf_lte_rl, pf_lte_lh)
prval pf_mrr = mr_grow_r(right.pf_mr, pf_lte_rh, rstart, right.hp.root)
in
if swapped then let
val left1 = restore_heap_prop left.hp
prval pf_mrheap = mr_join(left1.pf_mr, pf_mrr, pf_lte_lh, start, left1.hp.root)
in @{hp = h, pf_mr = pf_mrheap}
end else let
prval pf_mrheap = mr_join(left.pf_mr, pf_mrr, pf_lte_lh, start, left.hp.root)
in @{hp = h, pf_mr = pf_mrheap}
end
end
| RightGr(pf_lte_lr | _) => let
val (pf_lte_rh | swapped) = order_elements(A, right.hp.root, h.root, gt)
in
if swapped then let
val right1 = restore_heap_prop right.hp
prval pf_mr_leftright = mr_join(left.pf_mr, right1.pf_mr, pf_lte_lr, start, left.hp.root)
prval pf_mrheap = mr_grow_r( pf_mr_leftright, pf_lte_rh, start, right1.hp.root)
in @{hp = h, pf_mr = pf_mrheap}
end else let
prval pf_mr_leftright = mr_join(left.pf_mr, right.pf_mr, pf_lte_lr, start, left.hp.root)
prval pf_mrheap = mr_grow_r( pf_mr_leftright, pf_lte_rh, start, right.hp.root)
in @{hp = h, pf_mr = pf_mrheap}
end
end
end // end of restore_big_heap_prop
Функция compare_elements сравнивает два элемента и вместе с результатом поставляет факт-доказательство того, что один элемент не меньше другого. Маленькая куча (с индексом 0 или 1) обладает свойством кучи автоматически, т.к. имеет размер 1, что доказывает теорема leon_base_is1. В большой же куче нужно сравнить между собой корни ее подкуч, и тот, что больше, сравнить с корнем большой кучи. Если корень большой кучи больше максимального из корней подкуч, то он больше их обоих, и из этого факта (получаемого через аксиому о транзитивности lte_trans) и фактов-доказательств свойств кучи подкуч (MAXRIGHT) можно соорудить доказательство свойства кучи для всей большой кучи (используя теоремы mr_grow_r и mr_join). Если же корень большой кучи меньше, то он меняется местами с максимальным из корней подкуч, после чего та подкуча теряет свое свойство и его нужно восстановить, применив рекурсивно ту же самую функцию восстановления свойства. С т.з. доказательств случаи, когда корень левой подкучи больше правого, и когда правый больше левого, не одинаковы, там доказательство свойства кучи получаются немного по-разному, поэтому кода получилось так много. При реализации на языках без доказательств эти случаи неразличимы, и код заметно проще.
Умея восстанавливать свойство кучи, т.е. получать heap1 из heap, можно сделать восстановление свойства ряда. Когда мы к ряду хотим присоединить очередную кучу, у нас есть текущий ряд, для которого свойство ряда
уже выполнено (он состоит из куч типа heap2) и новая куча типа heap, из которой мы хотим получить heap2, которую можно добавить в наш ряд-список типа heaps. При этом могут потребоваться перестановки корней куч, последняя куча ряда может потерять свои свойства и превратиться из heap2 опять в heap, тогда мы получим аналогичную ситуацию уже с меньшим рядом. Восстановлением свойства ряда будут заниматься функции с такими сигнатурами:
fun restore_heapstring_prop {m,r,k,sz:nat | m + sz - 1 == r; r < n}
(hs : heaps(m), h:heap(r,k,sz)):<cloref1> heaps(r+1)
and restore_heapstring_prop_small {m,r,k,sz,r1,k1,sz1:nat | m + sz - 1 == r; r1 < n; r1 == r+1; sz1 == 1}
(hs: heaps(m), h: heap2(r,k,sz), small: heap(r1,k1,sz1)):<cloref1> heaps(r+2)
and restore_heapstring_prop_big {m,phr,phk,phsz,lr,lsz,rr,rsz,hr,hk,hsz:nat |
m + phsz - 1 == phr; hr == phr + hsz; hsz == lsz + rsz + 1; rr + 1 == hr; rr == lr + rsz;
lsz > 0; rsz > 0; hr < n; hk > 1}
(hs: heaps(m), prev: heap2(phr,phk,phsz),
left: heap1(lr,hk-1,lsz), right: heap1(rr,hk-2,rsz), h: heap(hr,hk,hsz)):<cloref1> heaps(hr+1)
Да, пять строк только на описание типа одной функции. Выкинь одно из статических условий, и код не скомпилируется. Тело функций не привожу, там много громоздкого кода. Смысл его в том, что ежели ряд hs пустой, то достаточно восстановить свойство кучи для h, превратив ее из heap в heap1, и свойство кучи для нее одновременно будет и свойством ряда, т.к. она в нем будет одна. Если же в текущем ряду есть хоть одна куча, то в зависимости от размера добавляемой кучи h вызывается либо restore_heapstring_prop_small, либо restore_heapstring_prop_big. Для маленькой кучи нужно ее корень сравнить с корнем предыдущей, при необходимости сделав перестановку и восстановив свойство кучи для предыдущей, в результате можно получить доказательство, что новый корень не меньше всех предыдущих элементов и свойство ряда будет доказано. Для большой же кучи корень последней кучи текущего ряда нужно сравнивать как с корнем добавляемой, так и с корнями ее подкуч, и делать перестановку только если он больше их всех. Иначе достаточно восстановить свойство кучи в добавляемой, тогда ее корень станет не меньше корня последней кучи текущего ряда, и свойство ряда для новой кучи будет выполнено. Полный код можно посмотреть здесь, там все эти случаи аккуратно рассматриваются и нужные доказательства кропотливо конструируются.
Акт 7. Сортировка
Основная функция сортировки будет состоять из "выращивания" ряда от heaps(0) (т.е. nil) до heaps(n), где n - размер массива, и его последующего "сжимания" обратно, попутно формируя доказательство отсортированности, используя доказательства свойств куч ряда. В процессе, глядишь, и массив отсортируется. Выращиваем мы ряд путем вызова grow_loop(0, nil), где
fun grow_loop {i:nat | i <= n} (i : int i, hs : heaps(i)) :<cloref1> heaps(n) =
if i = n then hs else grow_loop(i+1, grow hs)
Где функция grow растит ряд на один элемент:
fun grow {m:nat; m < n} (hs : heaps(m)):<cloref1> heaps(m+1) =
case+ hs of
| nil () => @{hp = small_heap(0, ~1), pf_mr = MRsingle, pf_totalmr = MRsingle} :: nil
| h :: nil () => let
val small = small_heap(h.hp.root + 1, h.hp.k)
in restore_heapstring_prop_small(nil, h, small)
end
| h0 :: (h0rest as h1 :: rest) =>
if h0.hp.k + 1 = h1.hp.k then let //join two top heaps into a bigger one
prval () = lh_use h1.hp.pf_r
prval () = leon_positive h0.hp.pf_sz //prove that child sizes > 0 to use LEONind
prval pfgc = GCbig(h1.pf_mr, h0.pf_mr)
val bighp = @{root = h0.hp.root + 1, k = h1.hp.k + 1, sz = h1.hp.sz + h0.hp.sz + 1,
pf_sz = LEONind(h0.hp.pf_sz, h1.hp.pf_sz), pf_r = LHpf, pf_gc = pfgc }
in
case+ rest of
| prev :: hps => restore_heapstring_prop_big(hps, prev, heap2to1 h1, heap2to1 h0, bighp)
| nil () => let
val bighp1 = restore_big_heap_prop (heap2to1 h1, heap2to1 h0, bighp)
in @{hp = bighp, pf_mr = bighp1.pf_mr, pf_totalmr = bighp1.pf_mr} :: nil
end
end else let //add a small heap
val small = small_heap(h0.hp.root + 1, h0.hp.k)
in restore_heapstring_prop_small(h0rest, h0, small)
end
Здесь паттерн-матчингом по списку куч мы смотрим на две последних (если они есть), и если их индексы - подряд идущие убывающие числа (h0.hp.k + 1 = h1.hp.k), то объединяем их с новым элементом в одну кучу побольше. Иначе добавляем в ряд маленькую кучу из одного элемента.
Кстати, тип функции grow_loop как раз является теоремой о том, что любое натуральное число можно представить суммой чисел Леонардо, а код ее служит доказательством.
Сжиманием непустого ряда на один элемент занимается функция shrink:
fun shrink {m:nat | m > 0; m <= n} (hs : heaps(m)) :<cloref1> heaps(m-1) =
case+ hs of h :: rest =>
if h.hp.k < 2 then let prval () = leon_base_is1 h.hp.pf_sz in rest end
else let
val (left,right) = split_heap h.hp
val hs1 = restore_heapstring_prop(rest, left.hp)
in restore_heapstring_prop(hs1, right.hp)
end
Если последняя куча имеет индекс меньше 2, то ее размер равен 1 (что доказывает теорема leon_base_is1) и ее можно просто выкинуть, вернув остаток списка. Если же это большая куча, то ее можно разбить на две ее подкучи и корень, две этих кучи добавить в список вместо исходной, а корень выбросить. Добавлением куч занимается restore_heapstring_prop, которая заодно сделает необходимые перестановки и обеспечит выполнение всех требуемых свойств.
Дальше самое главное - цикл скукоживания ряда, в котором мы формируем доказательство отсортированности. Имея непустой ряд из i элементов (hs : heaps(i)), можно взять его последнюю кучу htop типа heap2, она содержит доказательство свойства ряда "типа" MAXRIGHT(0,i-1). Из нее теоремой lte_from_maxright можно получить факт-доказательство LTE(i-2,i-1) того, что A[i-2] <= A[i-1] (когда i больше 1). Если у нас к этому моменту будет доказательство отсортированности элементов A[i-1 .. n-1], то из них мы можем построить доказательство отсортированности A[i-2 .. n-1], которое передать на следующую итерацию цикла. Когда же i = 1, переданное доказательство SORTED(i-1,n-1) и есть SORTED(0,n-1) - искомое доказательство отсортированности всего массива!
fn top {m:int | m > 0}(hs : heaps(m)):<> [r,k,sz:nat | r + 1 == m] heap2(r,k,sz) =
case+ hs of h :: rest => h
fun shrink_loop{i:nat | i < n; i > 0}
(pf_sorted : SORTED(i-1,n-1) | i : int i, hs : heaps(i)) :<cloref1> (SORTED(0,n-1) | void) =
if i = 1 then (pf_sorted | ())
else let
val htop = top hs
prval pf_lte = lte_from_maxright(htop.pf_totalmr, i-2, 0, i-1)
in
shrink_loop(SORTEDjoin(pf_lte, pf_sorted) | i-1, shrink hs)
end
Тогда в теле основной функции smoothsort нужно лишь вызвать grow_loop для выращивания ряда куч, потом сократить его на один элемент, получив необходимые параметры для цикла скукоживания shrink_loop, и вызвать его, получив доказательство отсортированности.
if n > 1 then let
val hps = grow_loop(0, nil) //hps : heaps(n)
val htop = top hps
prval pf_lte21 = lte_from_maxright(htop.pf_totalmr, n-2, 0, n-1) // LTE(n-2, n-1)
prval pf_sorted0 : SORTED(n-1, n-1) = SORTEDsingle
prval pf_sorted1 : SORTED(n-2, n-1) = SORTEDjoin(pf_lte21, pf_sorted0)
val hps1 = shrink hps
in shrink_loop(pf_sorted1 | n-1, hps1)
end else (SORTEDsingle | ())
Все это можно делать, когда в массиве есть хотя бы 2 элемента. А массив из одного элемента и так отсортирован, доказательство этого тривиально - база индукции в SORTED.
Вот так легко и просто можно написать на ATS сортировку smoothsort с доказательством ровно в 400 строк. А теперь обещанный психологический тест. Если вы осилили этот пост и поняли прочитанное, а не просто пролистали, то ваш IQ не меньше 111, и у вас слишком много свободного времени. :)
