为什么在 rcpp 中实现的快速排序工作缓慢？

Question

我已经在 Rcpp 中实现了快速排序算法，但对于大型数组，它的运行速度明显慢于 sort(array, method="quick")。为什么？

这是我的 Rcpp 代码

// partition using hoare's scheme

#include <Rcpp.h>
using namespace Rcpp;

int partition(NumericVector a,int start,int end)
{
    double pivot = a[end];
    int i = start - 1;
    int j = end + 1;
    //Rcout << a <<"\n";
    while(1)
    {
        do {
            i++;
        } while (a[i] < pivot);

        do {
            j--;
        } while (pivot < a[j]);

        if(i > j)
            return j;

        //special.Swap(a, i, j);
        std::swap(a[i], a[j]);
    }
}

void qsort(NumericVector a,int start,int end)
{
  //Rcout << start <<"," << end <<"\n";
  if(start < end)
  {
    int P_index = partition(a, start, end);
    //Rcout << P_index << "\n";
    qsort(a, start, P_index);
    qsort(a, P_index + 1, end);
  }
}


// [[Rcpp::export]]
NumericVector QuickSortH_WC(NumericVector arr)
{
  int len = arr.size();
  qsort(arr, 0, len-1);
  //Rcout << arr <<"\n";
  return 1;
}

同样对于具有浮点值的数组，算法更差。我想和 hoare's 和 lomuto 分区方案做个比较，但我不知道这个实现是否有任何缺陷，因为哪个算法工作得更慢。

Answer 1

您的代码效率低下的主要原因似乎是混合了您要比较的两种分区方案。你声称使用Hoare分区方案，代码看起来很像，但是pivot是按照Lomuto分区方案计算的。此外，如果i >= j，则应返回j，而不是i > j。修复这两个问题并用稍快的++i 替换i++ 我得到：

// partition using hoare's scheme

#include <Rcpp.h>
using namespace Rcpp;

int partition(NumericVector a,int start,int end)
{
    double pivot = a[(start + end) / 2];
    int i = start - 1;
    int j = end + 1;
    //Rcout << a <<"\n";
    while(1)
    {
        do {
            ++i;
        } while (a[i] < pivot);

        do {
            --j;
        } while (pivot < a[j]);

        if(i >= j)
            return j;

        //special.Swap(a, i, j);
        std::swap(a[i], a[j]);
    }
}

void qsort(NumericVector a,int start,int end)
{
    //Rcout << start <<"," << end <<"\n";
    if(start < end)
    {
        int P_index = partition(a, start, end);
        //Rcout << P_index << "\n";
        qsort(a, start, P_index);
        qsort(a, P_index + 1, end);
    }
}


// [[Rcpp::export]]
NumericVector QuickSortH_WC(NumericVector arr)
{
    int len = arr.size();
    qsort(arr, 0, len-1);
    //Rcout << arr <<"\n";
    return arr;
}

/*** R
set.seed(42)
dat <- runif(1e6)
bench::mark(QuickSortH_WC(dat), sort(dat, method="quick"))
*/

输出

> bench::mark(QuickSortH_WC(dat), sort(dat, method="quick"))
# A tibble: 2 x 13
  expression                     min  median `itr/sec` mem_alloc `gc/sec` n_itr
  <bch:expr>                  <bch:> <bch:t>     <dbl> <bch:byt>    <dbl> <int>
1 QuickSortH_WC(dat)          95.7ms 100.5ms      8.63    2.49KB     43.2     5
2 sort(dat, method = "quick")   15ms  16.5ms     53.1    11.44MB     23.6    27
# … with 6 more variables: n_gc <dbl>, total_time <bch:tm>, result <list>,
#   memory <list>, time <list>, gc <list>
Warning message:
Some expressions had a GC in every iteration; so filtering is disabled. 

因此，虽然此方法比R's sort 慢约 7 倍，但它的运行时间至少具有可比的数量级。 （感谢@JosephWood 挖掘链接）。 Wikipedia 列出了对这两种模式的更多改进。

顺便说一句，我还更改了包装函数以返回更改后的数组。这允许我使用bench::mark 的默认行为，即比较返回的结果。我觉得这很有用...

Answer 2

Rcpp 不好地应用递归函数。 我建议迭代快速排序实现：

void _Quick_sorti( double _array[],int _l,int _h){
    int *_stack=new int [_h-_l+1]; double _tmp;int _i,_p,_top=-1;
    _stack[++_top]=_l;_stack[++_top]=_h;
    while(_top>=0){
      _h=_stack[_top--];_l=_stack[_top--];
      _tmp=_array[_h];
      _i=_l-1;
      for(int _j=_l;_j<=_h-1;_j++){
        if(_array[_j]<=_tmp){_i++;std::swap(_array[_i],_array[_j]);}
      }
      _p=_i+1;
      std::swap(_array[_p],_array[_h]);
      if(_p-1>_l){_stack[++_top]=_l;_stack[++_top]=_p-1;}
      if(_p+1<_h){_stack[++_top]=_p+1;_stack[++_top]=_h;}
    }
  delete _stack;
  }


     // [[Rcpp::export]]
SEXP Quick_sorti(SEXP &unsorted) { //run
  SEXP temp=clone(unsorted);// or Rf_duplicate
  double *z=REAL(temp);
  int N=LENGTH(temp)-1;
  int k=0;
  _Quick_sorti(z,k,N); // note that we have provide lvalue (if we put 0 it will not works int place of N)
  return temp;} 

代码改编自包含'_'前缀的宏，看起来很丑，而且它使用R 内部。添加堆栈意味着需要更多的内存。