排序大集锦（二）：N叉堆排序及优先级队列

首先需要明白何谓“堆”。逻辑上它是一种特殊的树形结构，而其物理实现仅仅是一维数组，如下所示：

上图中节点中的数表示元素的索引，将索引从1开始而不是0，可以使得整棵树从当前节点获得其孩子节点的操作对于根节点来说也是适用的，而额外“浪费”的这个空间可以用做哨兵使得比较操作的次数更少。可简单得出两点结论：
①除顶层与底层之外，每一层的节点数与上层节点数的比值固定。
②最底层的叶子节点尽可能向左分布。
注意在任意一种树形结构中并不一定保证同时满足上述两种性质，它仅仅适用于堆这种特殊的树，而往往也由于上述两种性质的成立导致在该数据结构中实施插入和删除操作时，在极为高效的前提下同样也能保证高效的空间利用。由结论①可知，在具有固定元素个数的集合中，如果我们改变每一层与上层节点数之间的比值，一般来说将存在如下关系：即当比值越大则堆的层数越小，因此从树的顶层到达底层将花费更少的时间。由结论②可知，因为叶子节点尽可能向左分布，因此在一维数组中内节点与叶子节点之间存在明确的边界。

根据上述分析，在给定一个节点的索引i之后，可以得到该节点的父、子节点。

/*make the number of node's sons generalization*/declare s as the number of node's sonsPARENT(i)  /*current node's parent*/return (i+s-2)/s//------------------------------------------------------FARLEFT(i) /*far left child*/ return ((i-1)*s+2)//------------------------------------------------------FARRIGHT(i)  /*far right child*/return (i*s+1)//------------------------------------------------------INTERNODE(length[Arr])return (length[Arr]+s-2)/s

证明过程省略，事实上观察到节点i的s倍为其所有孩子节点中的倒数第二个节点即可完成证明。而上述INTERNODE(length[Arr])子过程，则是根据当前数组的长度计算出内节点的最大索引项。
除此之外最重要的一点是，只有满足在任意节点中存储的关键字大于/小于其孩子节点中的关键字的数据结构方能称之为“堆”。

节点中存储的关键字大于其孩子节点中的关键字的堆称为大堆，反之则称为小堆。由于大堆及小堆的构造方式相同，所以我们以大堆为例剖析其构造过程。考虑一个简单集合A={2,1,3}，若需要将其变成大堆，只需将1和3简单交换即可，图示如下：

由上述操作可以发现如下事实：即父、子节点形成一个独立的组，该组中的最大值位于父节点处，而兄弟节点之间则不具备明显的关系。由此可以从堆的最后一个节点开始依次向前执行上述操作，最终在索引[1]处将得到整个集合的最大值。注意，当某对父子节点发生交换操作时，有可能使得交换后的子节点与其孩子节点之间的“堆性质”遭到破坏，因此需要使下层节点重新保持堆性质。

HEAPIFY(A,i,n)  /*n appointed to length[A] usually in the procedure of <BUILDHEAP>*/temp←A[i]maxnode←0while TRUEdo farleft←FARLEFT(i)   farright←FARRIGHT(i)   if farleft > n     then break   if farright > n  /*the number of the last internal node's children maybe less than n*/     then farright←n   t←farleft   for j←farleft+1 to farright  /*get the node has the max value*/     do t←max(t , j)   if A[t] > temp     then A[i]←A[t]          maxnode←i←t     else breakif maxnode≠0  then A[maxnode]←temp//---------------------------------------------------------------------------max(a,b)return ( A[a]>A[b] ? a : b )

假定单步操作的时间开销为O(1)，因为while循环的一次迭代过程所执行的比较次数与节点的孩子数相关，因此一次迭代操作所用时间为O(s)（s为节点的孩子数），并且该过程在最坏情况下将从根节点下降至堆的叶子节点处，此时执行㏒n次迭代，运行时间为O(s㏒n)。

从最后一个内节点到根节点，对每个节点使用以上子过程，使得堆中的任意一对父子节点均满足“堆性质”便完成堆的构造：

BUILDHEAP(A)start←INTERNODE(length[A])end←1for i←start to end  do HEAPIFY(A,i,length[A])

根据上述过程，简单分析其运行时间如下：

堆排序
整个集合最大的元素位于整个堆的根节点处。由于按照升序排序，整个集合的最大元素应位于最后一个节点处，所以交换索引[1]与索引[n]处的元素之后，最大元素已经处于正确的位置，而发生交换操作后的根节点的“堆性质”遭到破坏，因此需要使其重新具备“堆性质”，这一过程是通过HEAPIFY完成。接着重复上述操作直至索引[2]处的节点处，便完成了整个数组的排序操作：

HEAPSORT(A)BUILDHEAP(A)for i←length[A] to 2  do swap(A[1] , A[i])     HEAPIFY(A,1,i-1)

简单分析其运行时间：由于BUILDHEAP(A)的时间代价为O(n)，同时调用n-1次HEAPIFY子过程，因此T(n)=O(n+(n-1)s㏒n)=O(n㏒n)，其中s为常量因子。

优先级队列
考虑在堆中执行插入操作。因为插入元素有可能大于其父节点，这导致父节点的堆性质遭到破坏，因此我们将插入元素与其父节点比较，若小于则说明元素的插入对整个堆不造成任何影响，反之则将该元素与其父节点交换。而交换之后的插入元素仍有可能大于其父节点，因此我们不断对该元素进行比较-交换操作直至其到达正确的位置，此时所有节点又全部具备了“堆性质”。

HEAPINSERT(t)length[A]←length[A]+1k←length[A]while A[PARENT(k)] < t   do A[k]←A[PARENT(k)]      k←PARENT(k)A[k]←t

接着考虑堆的删除操作，这可以通过将根节点与最后一个叶子节点进行交换，此后对根节点调用HEAPIFY过程使根节点重新具备“堆性质”。

EXTRACTMAX(A)t←A[1]swap(A[1], A[length[A]])length[A]←length[A]-1HEAPIFY(A,1,length[A])return t

最后是对某个节点的修改操作。当修改后的值小于原先的关键字时执行向下的筛选操作，反之则向上筛选。

CHANGEKEY(A,i,val)/*sift down*/if A[i] > val  /*decrease the node*/  then A[i]←val       HEAPIFY(A,i,length[A])       return/*sift up*/t←length[A]length[A]←i-1HEAPINSERT(val)length[A]←t

因为这两种算法的运行时间取决于BUILDHEAP和HEAPIFY子过程，而这两个过程的运行时间由节点的孩子数s而定，所以最终，问题被转换为确定s的值为多少可以使得算法的性能达到最优，建立如下测试：

#include <iostream>#include <fstream>#include <string>#include <ctime>#include <exception>#include <Windows.h>using namespace std;#pragma once  /*xxoo~*/#define NumOfTestCount 20#define PreStoreMemPool 80000000/*can be changed, but fix it just to test the heap*/#define mintomillsec    (60*1000)#define sectomillsec    1000#define Start 1#define End 2#define HSStartHeapS 100000#define HSIntervalS 100000#define PQStartHeapS 10000000#define PQIntervalS 500000typedef enum {heapsort , pri_insert , pri_delete } FuncName;typedef size_t node;typedef size_t point;class Heap{private:    size_t NodeSon;        /*the number of sons A node has*/    node * temp;            /*temp space & just recover for other node*/    string str;    static size_t StartHeapSize;    static size_t IntervalSize;private:    Heap(const Heap&) { }    Heap& operator=(const Heap&) { }    node FarLeft(node i) const { return ((i-1)*NodeSon+2); }    node FarRight(node i) const { return (i*NodeSon+1); }    node InterNode() const { return (ArrLength+NodeSon-2)/NodeSon; }    void GeneRand(size_t n_Length);    node MaxInTwoNode(const node&, const node&) const;    size_t GetRunTime()    { Timer(Start); BuildHeap(); Timer(End);       return (e_millsec - s_millsec); }public:    virtual ~Heap();    explicit Heap(size_t);    Heap(node*, node*);    void BuildHeap();    void SetNodeSon(size_t ns) { NodeSon = ns; }    void SeqDisp() const;    void RecoverArr();    void ExtendHeapSize(size_t);    void GetStartAndIntervalSize();protected:    size_t ArrLength;        /*the size of Array*/    node * root;                /*a point to the array dynamic allocate*/    SYSTEMTIME tm;    DWORD s_millsec,e_millsec;    size_t RunTime[NumOfTestCount];  /*time buffer*/    ofstream fout;    streambuf * outbuf;    void Exit() { system("pause"); exit(0); }    node Parent(node i) { return (i+NodeSon-2)/NodeSon; }    void Heapify(node , size_t);    void swap(node& left, node& right)    { node temp = left; left = right; right = temp; }    void Timer(point);    void WriteHead();    void WriteTail();};class HeapSort : public Heap{private :     size_t GetRunTime()    { Timer(Start); Sort(); Timer(End);    return (e_millsec-s_millsec); }public :     explicit HeapSort(size_t n_Length) : Heap(n_Length) { }    void Sort();    void OpenFile()    { fout.open("f:\\HeapSort.m",ios::trunc | ios::out); }    void Go();    void WriteHSBdFile();    ~HeapSort() { }};class PriQueue : public Heap{private:    void TestPQ();public:    explicit PriQueue(size_t n_Length) : Heap(n_Length)     { root[0] = RAND_MAX; } /*set a guard*/    void PQInsert(node);    void OPNode()    { PQInsert( rand() ); ExtractMax(); }    node ExtractMax();    void ChangeKey(size_t, node);    void OpenFile()    { fout.open("f:\\PriQueue.m",ios::trunc | ios::out); }    void Go();    size_t GetOPTime()    { Timer(Start); OPNode(); Timer(End);     return (e_millsec - s_millsec); }    void WritePriBdFile();    ~PriQueue() { }};size_t Heap::StartHeapSize = 0;size_t Heap::IntervalSize = 0;int main(int argc, char *argv[]){    srand(static_cast<unsigned>( time(NULL) ));    HeapSort hs(rand()%PreStoreMemPool);    hs.Go();    PriQueue pq(rand()%PreStoreMemPool);    pq.Go();    system("pause");    return EXIT_SUCCESS;}void PriQueue::Go(){    /*eight times insertion*/    OpenFile();    cout<<"\nPriQueue already, go..."<<endl;    WriteHead();    WritePriBdFile();    WriteTail();}void HeapSort::Go(){    OpenFile();    cout<<"HeapSort already, go..."<<endl;    WriteHead();    WriteHSBdFile();    WriteTail();}void HeapSort::WriteHSBdFile(){    cout<<"clear;\nclc;\nord_x=["<<endl;    for(size_t index = 0; index != NumOfTestCount; ++index)    {        cout<<HSStartHeapS+index*HSIntervalS<<' ';        if((index+1)%10 == 0) cout<<"..."<<endl;    }    cout<<"];"<<endl;    for(size_t ns = 2; ns != 10; ++ns)    {        SetNodeSon(ns);  /*startup*/        cout.rdbuf(outbuf);        cout<<"NodeSon is "<<ns<<"\tStart Test:"<<endl;        ArrLength = HSStartHeapS;        for(size_t i = 0; i != NumOfTestCount; ++i)        {            RunTime[i] = GetRunTime();            cout<<"ArrLength = "<<ArrLength                <<"\t\tRun Time is "<<RunTime[i]<<endl;            ExtendHeapSize(HSStartHeapS+(i+1)*HSIntervalS);        }        outbuf = cout.rdbuf( fout.rdbuf() );        cout<<"NodeSon"<<ns<<"=["<<endl;        for(size_t j = 0; j != NumOfTestCount; ++j)        {            cout<<RunTime[j]<<' ';            if((j+1)%10 == 0) cout<<"..."<<endl;        }        cout<<"];\n"<<endl;    }}void PriQueue::WritePriBdFile(){    cout<<"clear;\nclc;\nord_x=["<<endl;    for(size_t index = 0; index != NumOfTestCount; ++index)    {        cout<<PQStartHeapS+index*PQIntervalS<<' ';        if((index+1)%10 == 0) cout<<"..."<<endl;    }    cout<<"];"<<endl;    for(size_t ns = 2; ns != 10; ++ns)    {        SetNodeSon(ns);  /*startup*/        cout.rdbuf(outbuf);        cout<<"NodeSon is "<<ns<<"\tStart Test:"<<endl;        ArrLength = PQStartHeapS;        for(size_t i = 0; i != NumOfTestCount; ++i)        {            RunTime[i] = GetOPTime();            cout<<"ArrLength = "<<ArrLength                <<"\t\tRun Time is "<<RunTime[i]<<endl;            ExtendHeapSize(PQStartHeapS+(i+1)*PQIntervalS);        }        outbuf = cout.rdbuf( fout.rdbuf() );        cout<<"NodeSon"<<ns<<"=["<<endl;        for(size_t j = 0; j != NumOfTestCount; ++j)        {            cout<<RunTime[j]<<' ';            if((j+1)%10 == 0) cout<<"..."<<endl;        }        cout<<"];\n"<<endl;    }}void Heap::WriteHead(){    outbuf = cout.rdbuf( fout.rdbuf() );    cout<<"%clear the var in current working space;"<<endl;}void Heap::WriteTail(){    cout<<"title('performance of each NodeSon');"<<endl;    cout<<"xlabel('x=[100000 : 2000000]');"<<endl;    cout<<"ylabel('the time of sort use(millisecond)');"<<endl;    cout<<"hold on;"<<endl;    for(size_t index = 2; index != 10; ++index)        cout<<"plot(ord_x,NodeSon"<<index<<",\'"<<str[index]<<"*-\');"<<endl;    cout<<"legend('NodeSon2'";    for(size_t index = 3; index != 10; ++index)        cout<<",\'NodeSon"<<index<<'\'';    cout<<");";    cout.rdbuf(outbuf); fout.close();}void Heap::Timer(point t_Time){    GetLocalTime(&tm);    point m_Temp;    m_Temp=tm.wMinute*mintomillsec+tm.wSecond*sectomillsec+\        tm.wMilliseconds;    if(t_Time == Start)        s_millsec = m_Temp;    else        e_millsec = m_Temp;}void PriQueue::ChangeKey(size_t index, node NewValue){    try {        if( index < 1 || index >ArrLength )            throw exception("Out of Range!");    } catch (exception &e) {        cout<<e.what()<<endl;        Exit();    }    /*sift down*/    if(root[index] > NewValue)  /*decrease the node*/    {        root[index] = NewValue;        Heapify(index,ArrLength);        return;    }    /*sift up*/    root[index] = NewValue;    size_t t_Length = ArrLength;    ArrLength = index-1;    PQInsert(NewValue);    ArrLength = t_Length;}node PriQueue::ExtractMax(){    try{        if(ArrLength  == 0)            throw exception("There has no more nodes!");    } catch(exception &e) {        cout<<e.what()<<endl;        Exit();    }    node m_MaxVal = root[1];    swap(root[1],root[ArrLength]);    --ArrLength;    Heapify(1,ArrLength);    return m_MaxVal;}void PriQueue::PQInsert(node Insertion){    ArrLength = ArrLength+1;    size_t m_Term = ArrLength;    while( root[Parent(m_Term)] < Insertion)    {        root[m_Term] = root[Parent(m_Term)];        m_Term = Parent(m_Term);    }    root[m_Term] = Insertion;}void HeapSort::Sort(){    BuildHeap();    for(size_t index = ArrLength; index != 1; --index)    {        swap(root[1], root[index]);        Heapify(1,index-1);    }}void Heap::GetStartAndIntervalSize(){    /*be sure the min heap size that can be measured &    a appropriate increment, here the define of measurement    is run_time greater than 50 millisec*/    size_t t_Shz, t_Is = t_Shz = 10000, t_RunTime;    bool ShzNeedToBeChanged ,         IntersNeedToBeChanged =ShzNeedToBeChanged = true;    ArrLength = t_Shz;    cout<<"Now Test StartHeapSize & IntervalSize:"<<endl;    cout<<"----------------------------------------------------------"<<endl;    for( ; IntersNeedToBeChanged || ShzNeedToBeChanged ; )    {        if( ShzNeedToBeChanged )        {            t_RunTime = GetRunTime();            cout<<"StartHeapSize = "<<t_Shz                <<"\t\tRunTime = "<<t_RunTime<<endl;            if(t_RunTime < 50)            {                t_Shz += 5000;                ExtendHeapSize(t_Shz);            }            else {                ShzNeedToBeChanged = false;                cout<<"StartHeapSize has been ensured!\n"                    " And then Test IntervalSize : "<<endl;            }        }        if(!ShzNeedToBeChanged && IntersNeedToBeChanged)        {            size_t m_Size = t_Shz+t_Is, m_RunTime, m_TimeDiff;            ExtendHeapSize(m_Size);            m_RunTime = GetRunTime();            m_TimeDiff = m_RunTime - t_RunTime;            cout<<"IntervalSize = "<<t_Is                <<"\tRunTime = "<<m_RunTime                <<"\tThe Difference = "<<m_TimeDiff<<endl;            if(m_RunTime < t_RunTime || m_TimeDiff < 80)                t_Is += 20000;            else                IntersNeedToBeChanged = false;        }    }    cout<<endl<<endl;    StartHeapSize = t_Shz;    IntervalSize = t_Is;}void Heap::BuildHeap(){    try {        if (NodeSon < 2)            throw exception("Please set the NodeSon!");    } catch(exception &e) {        cout<<e.what()<<endl;        Exit();    }    size_t m_Start = InterNode() , m_End = 0;    for(node index = m_Start; index != m_End; --index)        Heapify(index,ArrLength);}void Heap::Heapify(node i, size_t HeapSize){    node MaxNode = 0, t_ExChange = root[i];    for( ; ; )    {        node m_FarLeft = FarLeft(i);        node m_FarRight = FarRight(i);        if(m_FarLeft > HeapSize) break;        if(m_FarRight > HeapSize) m_FarRight = HeapSize;        node m_MaxNode = m_FarLeft;        for(node index = m_FarLeft+1; index != m_FarRight+1; ++index)            m_MaxNode = MaxInTwoNode(m_MaxNode , index);        if(root[m_MaxNode] > t_ExChange)        {            root[i] = root[m_MaxNode];            MaxNode = i = m_MaxNode;        }        else break;    }    if(MaxNode) root[MaxNode] = t_ExChange;}inline node Heap::MaxInTwoNode(const node& Left, const node& Right) const{    return root[Left] > root[Right] ? Left : Right;}void Heap::SeqDisp() const{    for(size_t i = 1; i != ArrLength+1; ++i)        cout<<root[i]<<' ';    cout<<endl;}inline void Heap::RecoverArr(){    node * t_Root = root+1, * t_Temp = temp+1;    memcpy(t_Root,t_Temp,ArrLength*sizeof(int));}void Heap::ExtendHeapSize(size_t n_Length){    RecoverArr();    ArrLength = n_Length;}void Heap::GeneRand(size_t n_Length){    //pay attention that first node start at index 1!    for(size_t i = 1; i != n_Length+1; ++i)        root[i] = static_cast<node>(rand());}Heap::Heap(size_t n_Length) : str("00ymcrgbwk"){    root = temp = NULL;    root = new node[PreStoreMemPool+1];    temp = new node[PreStoreMemPool+1];    root[0] = NodeSon = 0;    GeneRand(PreStoreMemPool);    memcpy(temp,root,(PreStoreMemPool+1)*sizeof(int));    ArrLength = n_Length;    memset(&tm,0,sizeof(SYSTEMTIME));    s_millsec = e_millsec = 0;}Heap::Heap(node * first, node * end) : str("00ymcrgbwk"){    root = temp = NULL;    root = new node[end - first+1];    temp = new node[end - first+1];    root[0] = NodeSon = 0;    for(size_t index = 1; index != end-first+1; ++index)        root[index] = *(first+index-1);    memcpy(temp,root,(end-first+1)*sizeof(int));    ArrLength = end-first;    memset(&tm,0,sizeof(SYSTEMTIME));    s_millsec = e_millsec = 0;}Heap::~Heap(){    if(root != NULL) delete [] root;    if(temp != NULL) delete [] temp;}

这里需要指出的一点是，虽然在上述的测试中设置了不同的节点来测量优先级队列插入和删除操作的运行时间，但几乎用光了全部的内存，各个节点下的插入和删除开销还是无法计算出来（windows API的测量精度只能达到毫秒级，如果想要更精确的测度需要用到汇编指令），因此这里只列出了每种节点下堆排序的运行时间：

HeapSort运行时间

可以发现各种节点下的运行时间及其增长幅度差别不大。不过这里还有一个有用的优化技巧，由上图可以发现当孩子节点数达到4时，运行时间系数降到最低，由4开始系数逐渐增加，这正是我们想要的，因为对2的n次方的随机访问下标的计算可以通过移位操作实现。另外，优先级队列中的插入和删除操作的开销在大多数情况下都花费在了HEAPIFY子过程上，因此如果增加孩子节点的个数，那么显然这两类操作的运行时间将变小，，这是因为同等规模下节点孩子数越多则整个堆的深度越小。综上，我们可以将堆的孩子节点的个数设置为4，从而使整个结构的性能达到最优：

/*make the number of node's sons generalization*//*the number of node's sons equ 4*/PARENT(i)  /*current node's parent*/return (i+2)>>2//------------------------------------------------------FARLEFT(i) /*far left child*/ return ((i-1)<<2+2)//------------------------------------------------------FARRIGHT(i)  /*far right child*/return (i<<2+1)//------------------------------------------------------INTERNODE(length[Arr])return (length[Arr]+2)>>2