Legend:
Page
Library
Module
Module type
Parameter
Class
Class type
Source
Source file array.ml
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869open!ImportincludeArray0type'at='aarray[@@deriving_inlinecompare,sexp,sexp_grammar]letcompare:'a.('a->'a->int)->'at->'at->int=compare_arraylett_of_sexp:'a.(Sexplib0.Sexp.t->'a)->Sexplib0.Sexp.t->'at=array_of_sexpletsexp_of_t:'a.('a->Sexplib0.Sexp.t)->'at->Sexplib0.Sexp.t=sexp_of_arraylet(t_sexp_grammar:'aSexplib0.Sexp_grammar.t->'atSexplib0.Sexp_grammar.t)=fun_'a_sexp_grammar->array_sexp_grammar_'a_sexp_grammar;;[@@@end](* This module implements a new in-place, constant heap sorting algorithm to replace the
one used by the standard libraries. Its only purpose is to be faster (hopefully
strictly faster) than the base sort and stable_sort.
At a high level the algorithm is:
- pick two pivot points by:
- pick 5 arbitrary elements from the array
- sort them within the array
- take the elements on either side of the middle element of the sort as the pivots
- sort the array with:
- all elements less than pivot1 to the left (range 1)
- all elements >= pivot1 and <= pivot2 in the middle (range 2)
- all elements > pivot2 to the right (range 3)
- if pivot1 and pivot2 are equal, then the middle range is sorted, so ignore it
- recurse into range 1, 2 (if pivot1 and pivot2 are unequal), and 3
- during recursion there are two inflection points:
- if the size of the current range is small, use insertion sort to sort it
- if the stack depth is large, sort the range with heap-sort to avoid n^2 worst-case
behavior
See the following for more information:
- "Dual-Pivot Quicksort" by Vladimir Yaroslavskiy.
Available at
http://www.kriche.com.ar/root/programming/spaceTimeComplexity/DualPivotQuicksort.pdf
- "Quicksort is Optimal" by Sedgewick and Bentley.
Slides at http://www.cs.princeton.edu/~rs/talks/QuicksortIsOptimal.pdf
- http://www.sorting-algorithms.com/quick-sort-3-way *)moduleSort=struct(* For the sake of speed we could use unsafe get/set throughout, but speed tests don't
show a significant improvement. *)letget=getletset=setletswaparrij=lettmp=getarriinsetarri(getarrj);setarrjtmp;;moduletypeSort=sigvalsort:'at->compare:('a->'a->int)->left:int(* leftmost index of sub-array to sort *)->right:int(* rightmost index of sub-array to sort *)->unitend(* http://en.wikipedia.org/wiki/Insertion_sort *)moduleInsertion_sort:Sort=structletsortarr~compare~left~right=(* loop invariant:
[arr] is sorted from [left] to [pos - 1], inclusive *)forpos=left+1torightdo(* loop invariants:
1. the subarray arr[left .. i-1] is sorted
2. the subarray arr[i+1 .. pos] is sorted and contains only elements > v
3. arr[i] may be thought of as containing v
Note that this does not allocate a closure, but is left in the for
loop for the readability of the documentation. *)letreclooparr~left~compareiv=leti_next=i-1inifi_next>=left&&compare(getarri_next)v>0then(setarri(getarri_next);looparr~left~comparei_nextv)elseiinletv=getarrposinletfinal_pos=looparr~left~compareposvinsetarrfinal_posvdone;;end(* http://en.wikipedia.org/wiki/Heapsort *)moduleHeap_sort:Sort=struct(* loop invariant:
root's children are both either roots of max-heaps or > right *)letrecheapifyarr~compareroot~left~right=letrelative_root=root-leftinletleft_child=(2*relative_root)+left+1inletright_child=(2*relative_root)+left+2inletlargest=ifleft_child<=right&&compare(getarrleft_child)(getarrroot)>0thenleft_childelserootinletlargest=ifright_child<=right&&compare(getarrright_child)(getarrlargest)>0thenright_childelselargestiniflargest<>rootthen(swaparrrootlargest;heapifyarr~comparelargest~left~right);;letbuild_heaparr~compare~left~right=(* Elements in the second half of the array are already heaps of size 1. We move
through the first half of the array from back to front examining the element at
hand, and the left and right children, fixing the heap property as we go. *)fori=(left+right)/2downtoleftdoheapifyarr~comparei~left~rightdone;;letsortarr~compare~left~right=build_heaparr~compare~left~right;(* loop invariants:
1. the subarray arr[left ... i] is a max-heap H
2. the subarray arr[i+1 ... right] is sorted (call it S)
3. every element of H is less than every element of S *)fori=rightdowntoleft+1doswaparrlefti;heapifyarr~compareleft~left~right:(i-1)done;;end(* http://en.wikipedia.org/wiki/Introsort *)moduleIntro_sort:sigincludeSortvalfive_element_sort:'at->compare:('a->'a->int)->int->int->int->int->int->unitend=structletfive_element_sortarr~comparem1m2m3m4m5=letcompare_and_swapij=ifcompare(getarri)(getarrj)>0thenswaparrijin(* Optimal 5-element sorting network:
{v
1--o-----o-----o--------------1
| | |
2--o-----|--o--|-----o--o-----2
| | | | |
3--------o--o--|--o--|--o-----3
| | |
4-----o--------o--o--|-----o--4
| | |
5-----o--------------o-----o--5
v} *)compare_and_swapm1m2;compare_and_swapm4m5;compare_and_swapm1m3;compare_and_swapm2m3;compare_and_swapm1m4;compare_and_swapm3m4;compare_and_swapm2m5;compare_and_swapm2m3;compare_and_swapm4m5;;(* choose pivots for the array by sorting 5 elements and examining the center three
elements. The goal is to choose two pivots that will either:
- break the range up into 3 even partitions
or
- eliminate a commonly appearing element by sorting it into the center partition
by itself
To this end we look at the center 3 elements of the 5 and return pairs of equal
elements or the widest range *)letchoose_pivotsarr~compare~left~right=letsixth=(right-left)/6inletm1=left+sixthinletm2=m1+sixthinletm3=m2+sixthinletm4=m3+sixthinletm5=m4+sixthinfive_element_sortarr~comparem1m2m3m4m5;letm2_val=getarrm2inletm3_val=getarrm3inletm4_val=getarrm4inifcomparem2_valm3_val=0thenm2_val,m3_val,trueelseifcomparem3_valm4_val=0thenm3_val,m4_val,trueelsem2_val,m4_val,false;;letdual_pivot_partitionarr~compare~left~right=letpivot1,pivot2,pivots_equal=choose_pivotsarr~compare~left~rightin(* loop invariants:
1. left <= l < r <= right
2. l <= p <= r
3. l <= x < p implies arr[x] >= pivot1
and arr[x] <= pivot2
4. left <= x < l implies arr[x] < pivot1
5. r < x <= right implies arr[x] > pivot2 *)letreclooplpr=letpv=getarrpinifcomparepvpivot1<0then(swaparrpl;cont(l+1)(p+1)r)elseifcomparepvpivot2>0then((* loop invariants: same as those of the outer loop *)letrecscan_backwardsr=ifr>p&&compare(getarrr)pivot2>0thenscan_backwards(r-1)elserinletr=scan_backwardsrinswaparrrp;contlp(r-1))elsecontl(p+1)randcontlpr=ifp>rthenl,relselooplprinletl,r=contleftleftrightinl,r,pivots_equal;;letrecintro_sortarr~max_depth~compare~left~right=letlen=right-left+1in(* This takes care of some edge cases, such as left > right or very short arrays,
since Insertion_sort.sort handles these cases properly. Thus we don't need to
make sure that left and right are valid in recursive calls. *)iflen<=32thenInsertion_sort.sortarr~compare~left~rightelseifmax_depth<0thenHeap_sort.sortarr~compare~left~rightelse(letmax_depth=max_depth-1inletl,r,middle_sorted=dual_pivot_partitionarr~compare~left~rightinintro_sortarr~max_depth~compare~left~right:(l-1);ifnotmiddle_sortedthenintro_sortarr~max_depth~compare~left:l~right:r;intro_sortarr~max_depth~compare~left:(r+1)~right);;letsortarr~compare~left~right=letheap_sort_switch_depth=(* We bail out to heap sort at a recursion depth of 32. GNU introsort uses 2lg(n).
The expected recursion depth for perfect 3-way splits is log_3(n).
Using 32 means a balanced 3-way split would work up to 3^32 elements (roughly
2^50 or 10^15). GNU reaches a depth of 32 at 65536 elements.
For small arrays, this makes us less likely to bail out to heap sort, but the
32*N cost before we do is not that much.
For large arrays, this means we are more likely to bail out to heap sort at
some point if we get some bad splits or if the array is huge. But that's only a
constant factor cost in the final stages of recursion.
All in all, this seems to be a small tradeoff and avoids paying a cost to
compute a logarithm at the start. *)32inintro_sortarr~max_depth:heap_sort_switch_depth~compare~left~right;;endendletsort?pos?lenarr~compare=letpos,len=Ordered_collection_common.get_pos_len_exn()?pos?len~total_length:(lengtharr)inSort.Intro_sort.sortarr~compare~left:pos~right:(pos+len-1);;letto_arrayt=tletis_emptyt=lengtht=0letis_sortedt~compare=leti=ref(lengtht-1)inletresult=reftrueinwhile!i>0&&!resultdoletelt_i=unsafe_gett!iinletelt_i_minus_1=unsafe_gett(!i-1)inifcompareelt_i_minus_1elt_i>0thenresult:=false;decridone;!result;;letis_sorted_strictlyt~compare=leti=ref(lengtht-1)inletresult=reftrueinwhile!i>0&&!resultdoletelt_i=unsafe_gett!iinletelt_i_minus_1=unsafe_gett(!i-1)inifcompareelt_i_minus_1elt_i>=0thenresult:=false;decridone;!result;;letmergea1a2~compare=letl1=Array.lengtha1inletl2=Array.lengtha2inifl1=0thencopya2elseifl2=0thencopya1elseifcompare(unsafe_geta20)(unsafe_geta1(l1-1))>=0thenappenda1a2elseifcompare(unsafe_geta10)(unsafe_geta2(l2-1))>0thenappenda2a1else(letlen=l1+l2inletmerged=create~len(unsafe_geta10)inleta1_index=ref0inleta2_index=ref0infori=0tolen-1doletuse_a1=ifl1=!a1_indexthenfalseelseifl2=!a2_indexthentrueelsecompare(unsafe_geta1!a1_index)(unsafe_geta2!a2_index)<=0inifuse_a1then(unsafe_setmergedi(unsafe_geta1!a1_index);a1_index:=!a1_index+1)else(unsafe_setmergedi(unsafe_geta2!a2_index);a2_index:=!a2_index+1)done;merged);;letcopy_matrix=map~f:copyletfolding_mapt~init~f=letacc=refinitinmapt~f:(funx->letnew_acc,y=f!accxinacc:=new_acc;y);;letfold_mapt~init~f=letacc=refinitinletresult=mapt~f:(funx->letnew_acc,y=f!accxinacc:=new_acc;y)in!acc,result;;letfold_resultt~init~f=Container.fold_result~fold~init~ftletfold_untilt~init~f=Container.fold_until~fold~init~ftletcountt~f=Container.count~foldt~fletsummt~f=Container.sum~foldmt~fletmin_eltt~compare=Container.min_elt~foldt~compareletmax_eltt~compare=Container.max_elt~foldt~compareletfoldit~init~f=letacc=refinitinfori=0tolengtht-1doacc:=fi!acc(unsafe_getti)done;!acc;;letfolding_mapit~init~f=letacc=refinitinmapit~f:(funix->letnew_acc,y=fi!accxinacc:=new_acc;y);;letfold_mapit~init~f=letacc=refinitinletresult=mapit~f:(funix->letnew_acc,y=fi!accxinacc:=new_acc;y)in!acc,result;;letcountit~f=foldit~init:0~f:(funidxcounta->iffidxathencount+1elsecount);;letconcat_mapt~f=concat(to_list(map~ft))letconcat_mapit~f=concat(to_list(mapi~ft))letrev_inplacet=leti=ref0inletj=ref(lengtht-1)inwhile!i<!jdoswapt!i!j;incri;decrjdone;;letrevt=lett=copytinrev_inplacet;t;;letof_list_revl=matchlwith|[]->[||]|a::l->letlen=1+List.lengthlinlett=create~lenainletr=reflin(* We start at [len - 2] because we already put [a] at [t.(len - 1)]. *)fori=len-2downto0domatch!rwith|[]->assertfalse|a::l->t.(i)<-a;r:=ldone;t;;(* [of_list_map] and [of_list_rev_map] are based on functions from the OCaml
distribution. *)letof_list_mapxs~f=matchxswith|[]->[||]|hd::tl->leta=create~len:(1+List.lengthtl)(fhd)inletrecfilli=function|[]->a|hd::tl->unsafe_setai(fhd);fill(i+1)tlinfill1tl;;letof_list_mapixs~f=matchxswith|[]->[||]|hd::tl->leta=create~len:(1+List.lengthtl)(f0hd)inletrecfillai=function|[]->a|hd::tl->unsafe_setai(fihd);filla(i+1)tlinfilla1tl;;letof_list_rev_mapxs~f=lett=of_list_mapxs~finrev_inplacet;t;;letof_list_rev_mapixs~f=lett=of_list_mapixs~finrev_inplacet;t;;letfilter_mapit~f=letr=ref[||]inletk=ref0infori=0tolengtht-1domatchfi(unsafe_getti)with|None->()|Somea->if!k=0thenr:=create~len:(lengtht)a;unsafe_set!r!ka;incrkdone;if!k=lengthtthen!relseif!k>0thensub~pos:0~len:!k!relse[||];;letfilter_mapt~f=filter_mapit~f:(fun_ia->fa)letfilter_optt=filter_mapt~f:Fn.idletraise_length_mismatchnamen1n2=invalid_argf"length mismatch in %s: %d <> %d"namen1n2()[@@cold][@@inlinenever][@@localnever][@@specialisenever];;letcheck_length2_exnnamet1t2=letn1=lengtht1inletn2=lengtht2inifn1<>n2thenraise_length_mismatchnamen1n2;;letiter2_exnt1t2~f=check_length2_exn"Array.iter2_exn"t1t2;iterit1~f:(funix1->fx1(unsafe_gett2i));;letmap2_exnt1t2~f=check_length2_exn"Array.map2_exn"t1t2;init(lengtht1)~f:(funi->f(unsafe_gett1i)(unsafe_gett2i));;letfold2_exnt1t2~init~f=check_length2_exn"Array.fold2_exn"t1t2;foldit1~init~f:(funiacx->facx(unsafe_gett2i));;letfiltert~f=filter_mapt~f:(funx->iffxthenSomexelseNone)letfilterit~f=filter_mapit~f:(funix->iffixthenSomexelseNone)letexistst~f=leti=ref(lengtht-1)inletresult=reffalseinwhile!i>=0&¬!resultdoiff(unsafe_gett!i)thenresult:=trueelsedecridone;!result;;letexistsit~f=leti=ref(lengtht-1)inletresult=reffalseinwhile!i>=0&¬!resultdoiff!i(unsafe_gett!i)thenresult:=trueelsedecridone;!result;;letmemta~equal=existst~f:(equala)letfor_allt~f=leti=ref(lengtht-1)inletresult=reftrueinwhile!i>=0&&!resultdoifnot(f(unsafe_gett!i))thenresult:=falseelsedecridone;!result;;letfor_allit~f=letlength=lengthtinleti=ref(length-1)inletresult=reftrueinwhile!i>=0&&!resultdoifnot(f!i(unsafe_gett!i))thenresult:=falseelsedecridone;!result;;letexists2_exnt1t2~f=check_length2_exn"Array.exists2_exn"t1t2;leti=ref(lengtht1-1)inletresult=reffalseinwhile!i>=0&¬!resultdoiff(unsafe_gett1!i)(unsafe_gett2!i)thenresult:=trueelsedecridone;!result;;letfor_all2_exnt1t2~f=check_length2_exn"Array.for_all2_exn"t1t2;leti=ref(lengtht1-1)inletresult=reftrueinwhile!i>=0&&!resultdoifnot(f(unsafe_gett1!i)(unsafe_gett2!i))thenresult:=falseelsedecridone;!result;;letequalequalt1t2=lengtht1=lengtht2&&for_all2_exnt1t2~f:equalletmap_inplacet~f=fori=0tolengtht-1dounsafe_setti(f(unsafe_getti))done;;let[@inlinealways]findi_internalt~f~if_found~if_not_found=letlength=lengthtiniflength=0thenif_not_found()else(leti=ref0inletfound=reffalseinletvalue_found=ref(unsafe_gett0)inwhile(not!found)&&!i<lengthdoletvalue=unsafe_gett!iiniff!ivaluethen(value_found:=value;found:=true)elseincridone;if!foundthenif_found~i:!i~value:!value_foundelseif_not_found());;letfindit~f=findi_internalt~f~if_found:(fun~i~value->Some(i,value))~if_not_found:(fun()->None);;letfindi_exnt~f=findi_internalt~f~if_found:(fun~i~value->i,value)~if_not_found:(fun()->raise(Not_found_s(Atom"Array.findi_exn: not found")));;letfind_exnt~f=findi_internalt~f:(fun_ix->fx)~if_found:(fun~i:_~value->value)~if_not_found:(fun()->raise(Not_found_s(Atom"Array.find_exn: not found")));;letfindt~f=Option.map(findit~f:(fun_ix->fx))~f:(fun(_i,x)->x)letfind_mapt~f=letlength=lengthtiniflength=0thenNoneelse(leti=ref0inletvalue_found=refNoneinwhileOption.is_none!value_found&&!i<lengthdoletvalue=unsafe_gett!iinvalue_found:=fvalue;incridone;!value_found);;letfind_map_exn=letnot_found=Not_found_s(Atom"Array.find_map_exn: not found")inletfind_map_exnt~f=matchfind_mapt~fwith|None->raisenot_found|Somex->xin(* named to preserve symbol in compiled binary *)find_map_exn;;letfind_mapit~f=letlength=lengthtiniflength=0thenNoneelse(leti=ref0inletvalue_found=refNoneinwhileOption.is_none!value_found&&!i<lengthdoletvalue=unsafe_gett!iinvalue_found:=f!ivalue;incridone;!value_found);;letfind_mapi_exn=letnot_found=Not_found_s(Atom"Array.find_mapi_exn: not found")inletfind_mapi_exnt~f=matchfind_mapit~fwith|None->raisenot_found|Somex->xin(* named to preserve symbol in compiled binary *)find_mapi_exn;;letfind_consecutive_duplicatet~equal=letn=lengthtinifn<=1thenNoneelse(letresult=refNoneinleti=ref1inletprev=ref(unsafe_gett0)inwhile!i<ndoletcur=unsafe_gett!iinifequalcur!prevthen(result:=Some(!prev,cur);i:=n)else(prev:=cur;incri)done;!result);;letreducet~f=iflengtht=0thenNoneelse(letr=ref(unsafe_gett0)infori=1tolengtht-1dor:=f!r(unsafe_getti)done;Some!r);;letreduce_exnt~f=matchreducet~fwith|None->invalid_arg"Array.reduce_exn"|Somev->v;;letpermute=Array_permute.permuteletrandom_element_exn?(random_state=Random.State.default)t=ifis_emptytthenfailwith"Array.random_element_exn: empty array"elset.(Random.State.intrandom_state(lengtht));;letrandom_element?(random_state=Random.State.default)t=trySome(random_element_exn~random_statet)with|_->None;;letzipt1t2=iflengtht1<>lengtht2thenNoneelseSome(map2_exnt1t2~f:(funx1x2->x1,x2));;letzip_exnt1t2=iflengtht1<>lengtht2thenfailwith"Array.zip_exn"elsemap2_exnt1t2~f:(funx1x2->x1,x2);;letunzipt=letn=lengthtinifn=0then[||],[||]else(letx,y=t.(0)inletres1=create~len:nxinletres2=create~len:nyinfori=1ton-1doletx,y=t.(i)inres1.(i)<-x;res2.(i)<-ydone;res1,res2);;letsorted_copyt~compare=lett1=copytinsortt1~compare;t1;;letpartitioni_tft~f=letboth=mapit~f:(funix->iffixthenEither.FirstxelseEither.Secondx)inlettrues=filter_mapboth~f:(function|Firstx->Somex|Second_->None)inletfalses=filter_mapboth~f:(function|First_->None|Secondx->Somex)intrues,falses;;letpartition_tft~f=partitioni_tft~f:(fun_ix->fx)letlastt=t.(lengtht-1)(* Convert to a sequence but does not attempt to protect against modification
in the array. *)letto_sequence_mutablet=Sequence.unfold_step~init:0~f:(funi->ifi>=lengthtthenSequence.Step.DoneelseSequence.Step.Yield(t.(i),i+1));;letto_sequencet=to_sequence_mutable(copyt)letcartesian_productt1t2=ifis_emptyt1||is_emptyt2then[||]else(letn1=lengtht1inletn2=lengtht2inlett=create~len:(n1*n2)(t1.(0),t2.(0))inletr=ref0infori1=0ton1-1dofori2=0ton2-1dot.(!r)<-t1.(i1),t2.(i2);incrrdonedone;t);;lettransposett=iflengthtt=0thenSome[||]else(letwidth=lengthttinletdepth=lengthtt.(0)inifexiststt~f:(funt->lengtht<>depth)thenNoneelseSome(initdepth~f:(fund->initwidth~f:(funw->tt.(w).(d)))));;lettranspose_exntt=matchtransposettwith|None->invalid_arg"Array.transpose_exn"|Somett'->tt';;includeBinary_searchable.Make1(structtypenonrec'at='atletget=getletlength=lengthend)includeBlit.Make1(structtypenonrec'at='atletlength=lengthletcreate_like~lent=iflen=0then[||]else(assert(lengtht>0);create~lent.(0));;letunsafe_blit=unsafe_blitend)letinvariantinvariant_at=itert~f:invariant_amodulePrivate=structmoduleSort=Sortend