河村@徳山高専 です。 先日、クイックソートの実行時間の調査をお願い致しましたが、 送られてきた実行時間が小さすぎて比較ができませんでした。 適正な時間測定ができるよう修正したプログラムを以下に再投稿します。 -----------------------ここからが再投稿--------------------------- 1995年に投稿したクイックソート(qs6)に色々なご意見を頂き ありがとうございました。結局 qs6 より高速なものの報告はありませんでした。 あれから2年以上経ちましたが、いくつかの UNIX の qsort には いまだにバグに近い動作をするものがあります。 今回、重複キー(郵便番号など)に対してより高速な qs7 を 開発しましたので報告します。 種々の計算機で速度調査を行いたいのですが、ご協力願えないでしょうか。 結果は、ホームページ(http://www.tokuyama.ac.jp/home/~kawamura/)で 報告したいと思います。私本人(河村)は現在、徳山高専の非常勤講師ですので 電子メールなどの返事が遅れる可能性があります。予めご了承下さい。 下記のプログラムで5つのファイルを作り(bat main.c mm.c qs6.c qs7.c )、 ' % ./bat >& ttt &' して、結果と機種名とOS名(version番号)をお送り下さい。 時間がかかりすぎるときは bat 中の -DLOOP=20 を -DLOOP=10 などにしてください。 ---- bat --------------------------------------------------------------------- #!/bin/csh echo `hostname` cc -O -DDEBUG -DLOOP=20 -c main.c cc -O -DDEBUG -o qs7debug main.o qs7.c time ./qs7debug 1 80000 100 echo " " cc -O -c -DLOOP=20 main.c cc -O -o qsort main.o cc -O -o qs6 main.o qs6.c cc -O -o qs7 main.o qs7.c time ./qsort 1 20000 400 # key=乱数 個数=20000 要素=400byte time ./qs6 1 20000 400 time ./qs7 1 20000 400 echo " " time ./qsort 1 80000 100 # key=乱数 個数=80000 要素=100byte time ./qs6 1 80000 100 time ./qs7 1 80000 100 echo " " time ./qsort 1 400000 20 # key=乱数 個数=400000 要素=20byte time ./qs6 1 400000 20 time ./qs7 1 400000 20 echo " " time ./qsort 1000 400000 20 # key=乱数%1000 個数=同上 要素=同上 time ./qs6 1000 400000 20 time ./qs7 1000 400000 20 echo " " time ./qsort 100 400000 20 # key=乱数%100 time ./qs6 100 400000 20 time ./qs7 100 400000 20 echo " " #time ./qsort 2 400000 20 # key=乱数%2 time ./qs6 2 400000 20 time ./qs7 2 400000 20 echo " " time ./qsort -1 400000 20 # key=昇順 time ./qs6 -1 400000 20 time ./qs7 -1 400000 20 echo " " time ./qsort -2 400000 20 # key=降順 time ./qs6 -2 400000 20 time ./qs7 -2 400000 20 ---- main.c ------------------------------------------------------------------ #include unsigned int cmp_cnt, ass_cnt; die( s ) char *s; { fprintf(stderr, "***** %s *****\n", s); exit(1); } typedef struct { int key; int data; } el_t; int cmpfnc( xp, yp ) el_t *xp, *yp; {cmp_cnt++; return xp->key - yp->key;} #define KEY(i) (((el_t*)(vec+(i)*rec_siz))->key) #define DATA(i) (((el_t*)(vec+(i)*rec_siz))->data) main( argc, argv ) int argc; char **argv; { int arr_max, div_val, rec_siz, i, c; char *vec, *chk; if (argc != 4) die("Usage: a.out div_val arr_max rec_siz"); div_val = atoi(argv[1]); /*テストデータの種類を指定する*/ arr_max = atoi(argv[2]); /*要素の個数*/ rec_siz = atoi(argv[3]); /*要素の大きさ*/ printf("[%s] div_val=%d arr_max=%d rec_siz=%d arr=%dkB c=%d\n", argv[0],div_val, arr_max, rec_siz, arr_max*rec_siz/1000, LOOP); if (rec_siz < sizeof(el_t) || 11000 < rec_siz) die("ill rec_siz"); if (rec_siz % 4) die("recsiz must be multiplier of 4"); if ((vec = (char*)malloc((arr_max+5)*rec_siz)) == NULL) die("large arr_max"); srandom( arr_max ); /*乱数の初期化*/ for (c=0; c= 2 ) for (i = 0; i < arr_max; i++) KEY(i)= random()%div_val; #ifdef DEBUG for (i = 0; i < arr_max; i++) DATA(i) = i; /*下の検査のための準備*/ #endif qsort( (char*)vec, arr_max, rec_siz, cmpfnc ); #ifdef DEBUG /*以下でソートできたことを検査する*/ for (i = 1; i < arr_max; i++) if (KEY(i-1) > KEY(i)) die("not sorted"); chk = (char*)calloc(arr_max,1); for (i = 0; i < arr_max; i++) chk[DATA(i)] = 1; for (i = 0; i < arr_max; i++) if (chk[i] != 1) die("chk err"); #endif } printf("ass=%u cmp=%u OK ", ass_cnt/c, cmp_cnt/c); } ---- mm.c ----------------------------------------------------------------- /*8バイト整数があるときは、さらに高速化できる*/ static int mmkind, mmsize, high, low; extern int ass_cnt; #define A ((int*)a) #define B ((int*)b) #define C ((int*)c) #define D ((int*)d) static void mmprepare( void *base, int size ) { #ifdef DEBUG if (sizeof(int) != 4) die("sizeof(int) != 4"); if (size <= 0) die("mmsize <= 0"); #endif if ( ((int)base & (4-1)) == 0 && (size & (4-1)) == 0 ) if (size >= 16) mmkind = 1; else mmkind = 0; else mmkind = -1; mmsize = size; high = (size & (-16)); low = (size & 0x0C ); } static void mmswap( register char *a, register char *b ) { register int s; if (a == b) return; if (mmkind >= 0) { if (mmkind > 0) { register char *t = a + high; do { s = A[0]; A[0] = B[0]; B[0] = s; s = A[1]; A[1] = B[1]; B[1] = s; s = A[2]; A[2] = B[2]; B[2] = s; s = A[3]; A[3] = B[3]; B[3] = s; a += 16; b += 16; }while (a < t); } if (low != 0) { s = A[0]; A[0] = B[0]; B[0] = s; if (low >= 8) { s = A[1]; A[1] = B[1]; B[1] = s; if (low == 12) {s = A[2]; A[2] = B[2]; B[2] = s;}}} }else{ register char *t = a + mmsize; do {s = *a; *a++ = *b; *b++ = s;} while (a < t); } ass_cnt += 2; } static void mmswapblock( register char *a, register char *b, int size ) { register int s; if (mmkind >= 0) { register char *t = a + (size & (-16)); register int lo = (size & 0x0C); if (size >= 16) { do { s = A[0]; A[0] = B[0]; B[0] = s; s = A[1]; A[1] = B[1]; B[1] = s; s = A[2]; A[2] = B[2]; B[2] = s; s = A[3]; A[3] = B[3]; B[3] = s; a += 16; b += 16; }while (a < t); } if (lo != 0) { s = A[0]; A[0] = B[0]; B[0] = s; if (lo >= 8) { s = A[1]; A[1] = B[1]; B[1] = s; if (lo == 12) {s = A[2]; A[2] = B[2]; B[2] = s;}}} }else{ register char *t = a + size; do {s = *a; *a++ = *b; *b++ = s;} while (a < t); } ass_cnt += (size/mmsize)*2; } static void mmrot3( register char *a, register char *b, register char *c ) { register int s; if (mmkind >= 0) { if (mmkind > 0) { register char *t = a + high; do { s = A[0]; A[0] = B[0]; B[0] = C[0]; C[0] = s; s = A[1]; A[1] = B[1]; B[1] = C[1]; C[1] = s; s = A[2]; A[2] = B[2]; B[2] = C[2]; C[2] = s; s = A[3]; A[3] = B[3]; B[3] = C[3]; C[3] = s; a += 16; b += 16; c += 16; }while (a < t); } if (low != 0) { s = A[0]; A[0] = B[0]; B[0] = C[0]; C[0] = s; if (low >= 8) { s = A[1]; A[1] = B[1]; B[1] = C[1]; C[1] = s; if (low == 12) {s = A[2]; A[2] = B[2]; B[2] = C[2]; C[2] = s;}}} }else{ register char *t = a + mmsize; do {s = *a; *a++ = *b; *b++ = *c; *c++ = s;} while (a < t); } ass_cnt += 3; } ---- qs6.c ----------------------------------------------------------------- /*****************************************************/ /* */ /* qs6 (Quick sort function) */ /* */ /* by 河村 知行 1995.4.21 */ /* 〒745 山口県徳山市久米 徳山高専 情報電子工学科 */ /* kawamura@tokuyama.ac.jp */ /*****************************************************/ #include "mm.c" typedef struct { char *LL, *RR; } stack_node; /*L,l,R,rを積むスタックの構造体*/ #define PUSH(ll,rr) {top->LL = (ll); top->RR = (rr); ++top;} /*L,l,R,rを積む*/ #define POP(ll,rr) {--top; ll = top->LL; rr = top->RR;} /*L,l,R,rを戻す*/ #define med3(a,b,c) ((*cmp)(a,b)<0 ? \ ((*cmp)(b,c)<0 ? b : ((*cmp)(a,c)<0 ? c : a)) : \ ((*cmp)(b,c)>0 ? b : ((*cmp)(a,c)<0 ? a : c)) ) void qsort (base, nel, size, cmp) char* base; int nel; int size; int (*cmp)(); { register char *l, *r, *m; /*l,r:左右の集団の端 m:配列の中央の位置*/ register int t, eq_l, eq_r; /*eq_l:左の集団が全てsに等しいことを示す*/ char *L = base; /*現在分割している区間の左端の要素の先頭 */ char *R = &base[size * (nel - 1)]; /*現在分割している区間の右端の要素の先頭 */ int chklim = 63; /*昇(降)順検査をする要素数の下限*/ stack_node stack[32], *top = stack; /*32ビットマシンでは32で十分*/ if (nel <= 1) return; /*要素数が1以下のときは即終了する*/ mmprepare( base, size ); goto start; nxt: if (stack == top) return; /*スタックが空になったとき終了する*/ POP(L,R); /*スタックに積んである値を取り出す*/ for (;;) { start: if (L + size == R) {if ((*cmp)(L,R) > 0) mmswap(L,R); goto nxt;}/*要素数2*/ l = L; r = R; t = (r - l + size) / size; /*要素数*/ m = l + size * (t >> 1); /*配列の中央を計算*/ if (t >= 60) { register char *m1; register char *m3; if (t >= 200) { t = size*(t>>3); /*8等分したときのバイト数*/ { register char *p1 = l + t; register char *p2 = p1 + t; register char *p3 = p2 + t; m1 = med3( p1, p2, p3 ); p1 = m + t; p2 = p1 + t; p3 = p2 + t; m3 = med3( p1, p2, p3 ); } }else{ t = size*(t>>2); /*4等分したときのバイト数*/ m1 = l + t; m3 = m + t; } m = med3( m1, m, m3 ); } if ((t = (*cmp)(l,m)) < 0) { /*3-5-?*/ if ((t = (*cmp)(m,r)) < 0) { /*3-5-7の場合*/ if (chklim && nel >= chklim) { /*既に昇順か検査する*/ char *p; chklim = 0; for (p=l; p 0) goto fail; goto nxt; } fail: goto loopA; /*3-5-7*/ } if (t > 0) { if ((*cmp)(l,r) <= 0) {mmswap(m,r); goto loopA;} /*3-5-4*/ mmrot3(r,m,l); goto loopA; /*3-5-2*/ } goto loopB; /*3-5-5*/ } if (t > 0) { /*7-5-?*/ if ((t = (*cmp)(m,r)) > 0) { /*7-5-3の場合*/ if (chklim && nel >= chklim) { /*既に降順か検査する*/ char *p; chklim = 0; for (p=l; p 0) {mmswap(l,r); goto loopB;} /*5-5-3*/ /*5-5-5 のときに分割の種類を決める*/ /*5-5-5*/ for (;;) { if ((l += size) == r) goto nxt; /*5-5-5*/ if (l == m) continue; if ((t = (*cmp)(l,m)) > 0) {mmswap(l,r); l = L; goto loopA;} /*575-5*/ if (t < 0) {mmswap(L,l); l = L; goto loopB;} /*535-5*/ } loopA: eq_l = 1; eq_r = 1; /*タイプAの分割*/ /*左 <= 分割値 < 右*/ for (;;) { for (;;) { if ((l += size) == r) {l -= size; if (l != m) mmswap(m,l); l -= size; goto fin;} if (l == m) continue; if ((t = (*cmp)(l,m)) > 0) {eq_r = 0; break;} if (t < 0) eq_l = 0; } for (;;) { if (l == (r -= size)) {l -= size; if (l != m) mmswap(m,l); l -= size; goto fin;} if (r == m) {m = l; break;} if ((t = (*cmp)(r,m)) < 0) {eq_l = 0; break;} if (t == 0) break; } mmswap(l,r); /*左と右を交換する*/ } loopB: eq_l = 1; eq_r = 1; /*タイプBの分割*/ /*左 < 分割値 <= 右*/ for (;;) { for (;;) { if (l == (r -= size)) {r += size; if (r != m) mmswap(r,m); r += size; goto fin;} if (r == m) continue; if ((t = (*cmp)(r,m)) < 0) {eq_l = 0; break;} if (t > 0) eq_r = 0; } for (;;) { if ((l += size) == r) {r += size; if (r != m) mmswap(r,m); r += size; goto fin;} if (l == m) {m = r; break;} if ((t = (*cmp)(l,m)) > 0) {eq_r = 0; break;} if (t == 0) break; } mmswap(l,r); /*左と右を交換する*/ } fin: if (eq_l == 0) /*左側をソートする必要がある*/ if (eq_r == 0) /*両側をソートする必要がある*/ if (l-L < R-r) {PUSH(r,R); R = l;} /*左から先にソートする*/ else {PUSH(L,l); L = r;} /*右から先にソートする*/ else R = l; /*左側のみソートする必要がある*/ else if (eq_r == 0) L = r; /*右側のみソートする必要がある*/ else goto nxt; /*両側ともソートする必要がない*/ } } /* 次のような改良を試してみたが処理速度の向上はほとんど得られなかった。 1。要素数3のときは別処理を行う 2。要素数3のときは loopの直前で goto nxt; する 3。要素数2の区間はスタックに積まないようにする 4。要素数4を検出してloopを回避する 5。要素数4を検出してloopを回避する(4つに場合分けして処理を減らす) 6。中央に集まった連続した分割値を検出して、再ソートの要素数を減らす*/ ---- qs7.c ----------------------------------------------------------------- /*****************************************************/ /* */ /* qs7 (Quick sort function) */ /* */ /* by 河村 知行 1997.11.28 */ /* 〒745 山口県徳山市久米 徳山高専 情報電子工学科 */ /* kawamura@tokuyama.ac.jp */ /*****************************************************/ #include "mm.c" typedef struct { char *LL, *RR; } stack_node; /*L,Rを積むスタックの構造体*/ #define PUSH(ll,rr) {top->LL = (ll); top->RR = (rr); ++top;} /*L,Rを積む*/ #define POP(ll,rr) {--top; ll = top->LL; rr = top->RR;} /*L,Rを戻す*/ #define min(a,b) ((a)<(b)?(a):(b)) #define med3(a,b,c) ((*cmp)(a,b)<=0 ? \ ((*cmp)(b,c)<=0 ? b : ((*cmp)(a,c)<=0 ? c : a)) : \ ((*cmp)(b,c)>=0 ? b : ((*cmp)(a,c)<=0 ? a : c)) ) void qsort (base, nel, size, cmp) char* base; int nel; int size; int (*cmp)(); { register char *L = base; /*分割中の区間の左端の要素の先頭*/ register char *R = &base[size*(nel-1)]; /*分割中の区間の右端の要素の先頭*/ register char *m; /*中央値の要素の先頭*/ register int t; register char *p; /*作業用変数*/ int n; /*分割中の区間の要素数*/ int chklim = 63; /*昇(降)順検査をする要素数の下限*/ stack_node stack[32], *top = stack; /*32ビットマシンでは32で十分*/ if (nel <= 1) return; /*要素数が1以下のときは即終了する*/ mmprepare( base, size ); loop: if (L + size == R) {if ((*cmp)(L,R) > 0) mmswap(L,R); goto nxt;}/*要素数2*/ n = (R - L + size) / size; /*要素数*/ m = L + size * (n >> 1); /*配列の中央を計算*/ if (n <= 4) { if ((t = (*cmp)(L,m)) < 0) { if ((*cmp)(m,R) > 0) { if ((*cmp)(L,R) <= 0) mmswap(m,R); /*3-5-3,4*/ else mmrot3(R,m,L); /*3-5-2*/ } /*3-5-5,7 のときは何もしない*/ }else if (t > 0) { if ((*cmp)(m,R) >= 0) mmswap(L,R); /*7-5-5,3*/ else if ((*cmp)(L,R) <= 0) mmswap(L,m); /*7-5-7,8*/ else mmrot3(L,m,R); /*7-5-6*/ }else{ if ((*cmp)(m,R) > 0) mmswap(L,R); /*5-5-3*/ } /*5-5-5,7 のときは何もしない*/ if (n == 4) { p = L+size; if ((t = (*cmp)(p,m)) < 0) if ((*cmp)(L,p) <= 0) {} /*3-3,4-5-7 のときは何もしない*/ else mmswap(L,p); /*3-2-5-7*/ else if (t > 0) if ((*cmp)(p,R) <= 0) mmswap(p,m); /*3-6,7-5-7*/ else mmrot3(p,m,R); /*3-8-5-7*/ } /*3-5-5-7 のときは何もしない*/ goto nxt; } /* n <= 4 */ if ((n-3) >= 8) { /*(n-3)>=3 であること*/ char *m1,*m2,*m3,*p1,*p2,*p3; int cnt = 0; /*以下で中央値をできるだけ正しい値にする*/ if ((n-3) >= 41) { /*(n-3)>=9 であること*/ t = (n-3)/9*size; /*9等分したときのバイト数*/ p3=(p2=(p1=L+size+size)+t)+t; m1=med3(p1,p2,p3); if (m1==p2) cnt++; p3=(p2=(p1=p3+t )+t)+t; m2=med3(p1,p2,p3); if (m2==p2) cnt++; p3=(p2=(p1=p3+t )+t)+t; m3=med3(p1,p2,p3); if (m3==p2) cnt++; }else{ t = (n-3)/3*size; /*3等分したときのバイト数*/ m3=(m2=(m1=L+size+size)+t)+t; } m = med3( m1, m2, m3 ); if (m == m2) cnt++; if (cnt == 4 && chklim && n >= chklim) { /*4回とも中央値が中央にあるとき*/ chklim = 0; if ((t = (*cmp)(L,R)) < 0) { /*既に昇順か検査する*/ for (p=L; p 0) goto fail; goto nxt; }else if (t > 0) { /*既に降順か検査する*/ for (p=L; p 0) { mmrot3(R,L,m); l = ll = p; r = R; rr = R-size; } else { mmswap(p,m); l = ll = p+size; r = rr = R; } for (;;){ while (ll <= rr && (t = (*cmp)(ll,L)) <= 0) { /*x>mなる要素を探す*/ if (t == 0) {mmswap(l,ll); l+=size;} /*x==mの要素は左へ集める*/ ll+=size; } while (ll <= rr && (t = (*cmp)(L,rr)) <= 0) { /*m>yなる要素を探す*/ if (t == 0) {mmswap(rr,r); r-=size;} /*m==yの要素は右へ集める*/ rr-=size; } if (rr < ll) break; /*繰返しの唯一の終了条件*/ mmswap(ll,rr); /*xとyを交換する*/ ll+=size; rr-=size; } /* x==m y0) mmswapblock(L,ll-t,t); /*x==mを中央へ*/ t = min(r-rr,R-r); if (t>0) mmswapblock(ll,R+size-t,t); /*m==yを中央へ*/ l = L+(ll-l)-size; /*ここでlは y