Быстрый выбор C Алгоритм быстрее, чем C Qsort - PullRequest
0 голосов
/ 25 августа 2018

Я попытался реализовать алгоритм C QuickSelect, как описано в этом посте ( 3-кратная быстрая сортировка (реализация C) ).Тем не менее, все, что я получаю, - это производительность в 5-10 раз меньше, чем стандартная qsort (даже с начальной перетасовкой).Я попытался покопаться в исходном исходном коде qsort, указанном здесь (https://github.com/lattera/glibc/blob/master/stdlib/qsort.c),, но он слишком сложный. У кого-нибудь есть более простой и лучший алгоритм? Любая идея приветствуется. Спасибо, примечание: моя первоначальная проблема состоит в том, чтобы попробоватьчтобы получить наименьшие K-ые значения массива для первых K-ых индексов. Поэтому я планировал вызвать quickselect K times EDIT 1: Вот код Cython, скопированный и адаптированный по ссылке выше

cdef void qswap(void* a, void* b, const size_t size) nogil:
    cdef char temp[size]# C99, use malloc otherwise
    #char serves as the type for "generic" byte arrays

    memcpy(temp, b,    size)
    memcpy(b,    a,    size)
    memcpy(a,    temp, size)

cdef void qshuffle(void* base, size_t num, size_t size) nogil: #implementation of Fisher
    cdef int i, j, tmp# create local variables to hold values for shuffle

    for i in range(num - 1, 0, -1): # for loop to shuffle
        j = c_rand() % (i + 1)#randomise j for shuffle with Fisher Yates
        qswap(base + i*size, base + j*size, size)

cdef void partition3(void* base,
                      size_t *low, size_t *high, size_t size,
                      QComparator compar) nogil:       
    # Modified median-of-three and pivot selection.                      
    cdef void *ptr = base
    cdef size_t lt = low[0]
    cdef size_t gt = high[0] # lt is the pivot
    cdef size_t i = lt + 1# (+1 !) we don't compare pivot with itself
    cdef int c = 0

    while (i <= gt):
        c = compar(ptr + i * size, ptr + lt * size)
        if (c < 0):# base[i] < base[lt] => swap(i++,lt++)
            qswap(ptr + lt * size, ptr + i * size, size)
            i += 1
            lt += 1
        elif (c > 0):#base[i] > base[gt] => swap(i, gt--)
            qswap(ptr + i * size, ptr + gt* size, size)
            gt -= 1
        else:#base[i] == base[gt]
            i += 1
    #base := [<<<<<lt=====gt>>>>>>]
    low[0] = lt                                          
    high[0] = gt 


cdef void qselectk3(void* base, size_t lo, size_t hi, 
   size_t size, size_t k, 
   QComparator compar) nogil:                                             
    cdef size_t low = lo                                          
    cdef size_t high = hi                                                      

    partition3(base, &low, &high,  size, compar)    

    if ((k - 1) < low): #k lies in the less-than-pivot partition           
        high = low - 1
        low = lo                      
    elif ((k - 1) >= low and  (k - 1) <= high): #k lies in the equals-to-pivot partition
        qswap(base, base + size*low, size)
        return                              
    else: # k > high => k lies in the greater-than-pivot partition                    
        low = high + 1
        high = hi 
    qselectk3(base, low, high, size, k, compar)

"""
 A selection algorithm to find the nth smallest elements in an unordered list. 
 these elements ARE placed at the nth positions of the input array                                                                         
"""
cdef void qselect(void* base, size_t num, size_t size,
                              size_t n,
                              QComparator compar) nogil:
    cdef int k
    qshuffle(base, num, size)
    for k in range(n):
        qselectk3(base + size*k, 0, num - k - 1, size, 1, compar)

, которую я используювремя Python для получения производительности как метода pyselect (с N = 50), так и pysort. Вот так

def testPySelect():
    A = np.random.randint(16, size=(10000), dtype=np.int32)
    pyselect(A, 50)
timeit.timeit(testPySelect, number=1)

def testPySort():
    A = np.random.randint(16, size=(10000), dtype=np.int32)
    pysort(A)
timeit.timeit(testPySort, number=1)

Ответы [ 2 ]

0 голосов
/ 25 августа 2018

Ответ @chqrlie - хороший и окончательный ответ, но чтобы завершить публикацию, я публикую версию Cython вместе с результатами бенчмаркинга.Короче говоря, предлагаемое решение в 2 раза быстрее, чем qsort на длинных векторах!


    cdef void qswap2(void *aptr, void *bptr, size_t size) nogil:
        cdef uint8_t* ac = <uint8_t*>aptr
        cdef uint8_t* bc = <uint8_t*>bptr
        cdef uint8_t t
        while (size > 0): t = ac[0]; ac[0] = bc[0]; bc[0] = t; ac += 1; bc += 1; size -= 1

    cdef struct qselect2_stack:
        uint8_t *base
        uint8_t *last

    cdef void qselect2(void *base, size_t nmemb, size_t size,
                      size_t k, QComparator compar) nogil:
        cdef qselect2_stack stack[64]
        cdef qselect2_stack *sp = &stack[0]

        cdef uint8_t *lb
        cdef uint8_t*ub
        cdef uint8_t *p
        cdef uint8_t *i
        cdef uint8_t *j
        cdef uint8_t *top

        if (nmemb < 2 or size <= 0):
            return

        top = <uint8_t *>base
        if(k < nmemb): 
            top += k*size 
        else: 
            top += nmemb*size

        sp.base = <uint8_t *>base
        sp.last = <uint8_t *>base + (nmemb - 1) * size
        sp += 1

        cdef size_t offset

        while (sp > stack):
            sp -= 1
            lb = sp.base
            ub = sp.last

            while (lb < ub and lb < top):
                #select middle element as pivot and exchange with 1st element
                offset = (ub - lb) >> 1
                p = lb + offset - offset % size
                qswap2(lb, p, size)

                #partition into two segments
                i = lb + size
                j = ub
                while 1:
                    while (i < j and compar(lb, i) > 0):
                        i += size
                    while (j >= i and compar(j, lb) > 0):
                        j -= size
                    if (i >= j):
                        break
                    qswap2(i, j, size)
                    i += size
                    j -= size

                # move pivot where it belongs
                qswap2(lb, j, size)

                # keep processing smallest segment, and stack largest
                if (j - lb <= ub - j):
                    sp.base = j + size
                    sp.last = ub
                    sp += 1
                    ub = j - size
                else:
                    sp.base = lb
                    sp.last = j - size
                    sp += 1
                    lb = j + size

    cdef int int_comp(void* a, void* b) nogil:
        cdef int ai = (<int*>a)[0] 
        cdef int bi = (<int*>b)[0]
        return (ai > bi ) - (ai < bi)

    def pyselect2(numpy.ndarray[int, ndim=1, mode="c"] na, int n):
        cdef int* a = <int*>&na[0]
        qselect2(a, len(na), sizeof(int), n, int_comp)

Вот результаты тестов (1000 тестов):

#of elements   K      #qsort (s)                     #qselect2 (s)
1,000          50     0.1261                         0.0895
1,000          100    0.1261                         0.0910

10,000         50     0.8113                         0.4157
10,000         100    0.8113                         0.4367
10,000         1,000  0.8113                         0.4746

100,000        100    7.5428                         3.8259
100,000        1,000  7,5428                         3.8325
100,000        10,000 7,5428                         4.5727

Для тех, ктоЛюбопытно, этот кусок кода является жемчужиной в области реконструкции поверхности с использованием нейронных сетей.Еще раз спасибо @chqrlie, ваш код уникален в Интернете.

0 голосов
/ 25 августа 2018

Вот краткая реализация для вашей цели: qsort_select - это простая реализация qsort с автоматическим сокращением ненужных диапазонов.

Без && lb < top он ведет себя как обычный qsort за исключением патологических случаев, когда более продвинутые версии имеют лучшую эвристику. Этот дополнительный тест предотвращает полную сортировку диапазонов, которые находятся за пределами цели 0 .. (k-1) . Функция выбирает k наименьшие значения и сортирует их, остальная часть массива имеет оставшиеся значения в неопределенном порядке.

#include <stdio.h>
#include <stdint.h>

static void exchange_bytes(uint8_t *ac, uint8_t *bc, size_t size) {
    while (size-- > 0) { uint8_t t = *ac; *ac++ = *bc; *bc++ = t; }
}

/* select and sort the k smallest elements from an array */
void qsort_select(void *base, size_t nmemb, size_t size,
                  int (*compar)(const void *a, const void *b), size_t k)
{
    struct { uint8_t *base, *last; } stack[64], *sp = stack;
    uint8_t *lb, *ub, *p, *i, *j, *top;

    if (nmemb < 2 || size <= 0)
        return;

    top = (uint8_t *)base + (k < nmemb ? k : nmemb) * size;
    sp->base = (uint8_t *)base;
    sp->last = (uint8_t *)base + (nmemb - 1) * size;
    sp++;
    while (sp > stack) {
        --sp;
        lb = sp->base;
        ub = sp->last;
        while (lb < ub && lb < top) {
            /* select middle element as pivot and exchange with 1st element */
            size_t offset = (ub - lb) >> 1;
            p = lb + offset - offset % size;
            exchange_bytes(lb, p, size);

            /* partition into two segments */
            for (i = lb + size, j = ub;; i += size, j -= size) {
                while (i < j && compar(lb, i) > 0)
                    i += size;
                while (j >= i && compar(j, lb) > 0)
                    j -= size;
                if (i >= j)
                    break;
                exchange_bytes(i, j, size);
            }
            /* move pivot where it belongs */
            exchange_bytes(lb, j, size);

            /* keep processing smallest segment, and stack largest */
            if (j - lb <= ub - j) {
                sp->base = j + size;
                sp->last = ub;
                sp++;
                ub = j - size;
            } else {
                sp->base = lb;
                sp->last = j - size;
                sp++;
                lb = j + size;
            }
        }
    }
}

int int_cmp(const void *a, const void *b) {
    int aa = *(const int *)a;
    int bb = *(const int *)b;
    return (aa > bb) - (aa < bb);
}

#define ARRAY_SIZE  50000

int array[ARRAY_SIZE];

int main(void) {
    int i;
    for (i = 0; i < ARRAY_SIZE; i++) {
        array[i] = ARRAY_SIZE - i;
    }
    qsort_select(array, ARRAY_SIZE, sizeof(*array), int_cmp, 50);
    for (i = 0; i < 50; i++) {
        printf("%d%c", array[i], i + 1 == 50 ? '\n' : ',');
    }
    return 0;
}
...