~ubuntu-branches/ubuntu/karmic/scilab/karmic

« back to all changes in this revision

Viewing changes to routines/graphics/qsort.c

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2002-03-21 16:57:43 UTC
  • Revision ID: james.westby@ubuntu.com-20020321165743-e9mv12c1tb1plztg
Tags: upstream-2.6
ImportĀ upstreamĀ versionĀ 2.6

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*------------------------------------------------------------------------
 
2
 *      $NetBSD: qsort.c,v 1.5 1995/12/28 08:52:36 thorpej Exp $       
 
3
 *-
 
4
 * Copyright (c) 1992, 1993
 
5
 *      The Regents of the University of California.  All rights reserved.
 
6
 *
 
7
 * Redistribution and use in source and binary forms, with or without
 
8
 * modification, are permitted provided that the following conditions
 
9
 * are met:
 
10
 * 1. Redistributions of source code must retain the above copyright
 
11
 *    notice, this list of conditions and the following disclaimer.
 
12
 * 2. Redistributions in binary form must reproduce the above copyright
 
13
 *    notice, this list of conditions and the following disclaimer in the
 
14
 *    documentation and/or other materials provided with the distribution.
 
15
 * 3. All advertising materials mentioning features or use of this software
 
16
 *    must display the following acknowledgement:
 
17
 *      This product includes software developed by the University of
 
18
 *      California, Berkeley and its contributors.
 
19
 * 4. Neither the name of the University nor the names of its contributors
 
20
 *    may be used to endorse or promote products derived from this software
 
21
 *    without specific prior written permission.
 
22
 *
 
23
 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
 
24
 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 
25
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
 
26
 * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
 
27
 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 
28
 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
 
29
 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
 
30
 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
 
31
 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
 
32
 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
 
33
 * SUCH DAMAGE.
 
34
 *------------------------------------------------------------------------
 
35
 *  Modified for Scilab Jean-Philippe Chancelier 
 
36
 *  to keep a permutation index 
 
37
 *------------------------------------------------------------------------
 
38
 */
 
39
#if defined(THINK_C) || defined (__MWERKS__)  
 
40
#include <types.h>
 
41
#else
 
42
#ifndef __ABSC__
 
43
#include <sys/types.h>
 
44
#endif
 
45
#endif
 
46
#include <stdlib.h>
 
47
 
 
48
#define Min(a, b)       ((a) < (b) ? (a) : (b))
 
49
 
 
50
/*
 
51
 * Qsort routine from Bentley & McIlroy's "Engineering a Sort Function".
 
52
 */
 
53
 
 
54
#define swap(a, b) swapcode(a, b, 1 )
 
55
#define swapind(a, b)  if ( flag==1) swapcodeind(a,b,1)
 
56
#define vecswap(a, b, n) if ((n) > 0) swapcode(a, b, n/es)
 
57
#define vecswapind(a, b, n) if ((n) > 0 && flag == 1) swapcodeind(a,b,n/es1) 
 
58
 
 
59
#define med3(res,tabres,a, b, c, xa,xb,xc,cmp) cmp(a, b) < 0 ? \
 
60
         (cmp(b, c) < 0 ? (res=b,tabres=xb) :  \
 
61
          (cmp(a, c) < 0 ? (res=c,tabres=xc) : (res=a,tabres=xa) )) \
 
62
         :(cmp(b, c) > 0 ? (res=b,tabres=xb) : (cmp(a, c) < 0 ? (res=a,tabres=xa) : (res=c,tabres=xc) ))
 
63
 
 
64
/**************************************************
 
65
 * L'index de permutation est gard'e dans tab avec size ds es1
 
66
 * on utilise aussi escode et es1code qui peuvent etre differents de 
 
67
 * es et es1 
 
68
 * par exemple si on veut pour une matrice faire un sort de chacune 
 
69
 * de ces lignes la matrice etant stockee par colonne a[i+n*j] 
 
70
 * on va faire n appels a sciqsort((char *) &a[i],....)
 
71
 * les increments es et es1 doivent etre de la taille n*sizeof(elt-de-a)
 
72
 * par contre swap code ne doit permuter que les elements de la ligne i
 
73
 * dans ce cas on donnera a escode et es1code la valeur N
 
74
 * sinon cas standard on leur donnera la valeur 1;
 
75
 **************************************************/
 
76
 
 
77
void sciqsort(a,tab, flag, n, es, es1, cmp ,swapcode,swapcodeind)
 
78
  char *a;
 
79
  char *tab;
 
80
  int flag;
 
81
  int n, es,es1;
 
82
  int (*cmp)();  int (*swapcode)();  int (*swapcodeind)();
 
83
 
 
84
{
 
85
  char *pa, *pb, *pc, *pd, *pl, *pm, *pn;
 
86
  char *taba, *tabb, *tabc, *tabd, *tabl, *tabm, *tabn;
 
87
  int d,dind, r,r1,  swap_cnt;
 
88
 
 
89
loop:   
 
90
  swap_cnt = 0;
 
91
  if (n < 7) {
 
92
    for (pm = a + es, tabm= tab + es1 ; pm < (char *) a + n * es; pm += es, tabm +=es1 )
 
93
      {
 
94
        for (pl = pm, tabl= tabm ; pl > (char *) a && cmp(pl - es, pl) > 0;  pl -= es, tabl -=es1)
 
95
          {
 
96
            swapind(tabl,tabl- es1);
 
97
            swap(pl, pl - es);
 
98
          }
 
99
      }
 
100
    return;
 
101
  }
 
102
  pm = a + (n / 2) * es;
 
103
  tabm = tab + (n / 2)*es1 ;
 
104
  if (n > 7) {
 
105
    pl = a; 
 
106
    tabl = tab;
 
107
    pn = a + (n - 1) * es; 
 
108
    tabn = tab + (n-1) *es1;
 
109
    if (n > 40) {
 
110
      dind= (n/8) *es1;
 
111
      d =   (n/8) *es;
 
112
      med3(pl,tabl,pl, pl + d, pl + 2 * d, tabl, tabl + dind, tabl + 2 * dind, cmp);
 
113
      med3(pm,tabm,pm - d, pm, pm + d, tabm - dind, tabm, tabm + dind, cmp);
 
114
      med3(pn,tabn,pn - 2 * d, pn - d, pn, tabn - 2 * dind, tabn - dind, tabn, cmp);
 
115
    }
 
116
    med3(pm,tabm,pl, pm, pn, tabl, tabm, tabn, cmp);
 
117
  }
 
118
  swapind(tab,tabm);
 
119
  swap(a, pm);
 
120
 
 
121
  pa = pb = a + es;
 
122
  pc = pd = a + (n - 1) * es;
 
123
  
 
124
  taba = tabb = tab + es1;
 
125
  tabc = tabd = tab + (n - 1) * es1;
 
126
 
 
127
  for (;;) {
 
128
    while (pb <= pc && (r = cmp(pb, a)) <= 0) {
 
129
      if (r == 0) {
 
130
        swap_cnt = 1;
 
131
        swapind(taba,tabb);
 
132
        taba +=es1;
 
133
        swap(pa, pb);
 
134
        pa += es;
 
135
      }
 
136
      pb += es;
 
137
      tabb += es1;
 
138
    }
 
139
    while (pb <= pc && (r = cmp(pc, a)) >= 0) {
 
140
      if (r == 0) {
 
141
        swap_cnt = 1;
 
142
        swapind(tabc,tabd);
 
143
        tabd -= es1;
 
144
        swap(pc, pd);
 
145
        pd -= es;
 
146
      }
 
147
      pc -= es;
 
148
      tabc -= es1;
 
149
    }
 
150
    if (pb > pc)
 
151
      break;
 
152
    swapind(tabb,tabc);
 
153
    tabb += es1;
 
154
    tabc -= es1;
 
155
    swap(pb, pc);
 
156
    swap_cnt = 1;
 
157
    pb += es;
 
158
    pc -= es;
 
159
  }
 
160
 
 
161
  if (swap_cnt == 0) {  /* Switch to insertion sort */
 
162
    for (pm = a + es, tabm= tab + es1 ; pm < (char *) a + n * es; pm += es, tabm +=es1)
 
163
      {
 
164
        for (pl = pm, tabl= tabm ; pl > (char *) a && cmp(pl - es, pl) > 0;  pl -= es, tabl -=es1)
 
165
          {
 
166
            swapind(tabl,tabl- es1);
 
167
            swap(pl, pl - es);
 
168
          }
 
169
      }
 
170
    return;
 
171
  }
 
172
 
 
173
  pn = a + n * es;
 
174
  r = Min(pa - (char *)a, pb - pa);
 
175
  vecswap(a, pb - r, r);
 
176
 
 
177
  tabn = tab + n*es1 ;
 
178
  r1 = Min(taba - (char *) tab, tabb - taba);
 
179
  vecswapind(tab, tabb - r1, r1);
 
180
 
 
181
  r = Min(pd - pc, pn - pd - es);
 
182
  vecswap(pb, pn - r, r);
 
183
 
 
184
  r1 = Min(tabd - tabc, tabn - tabd - es1 );
 
185
  vecswapind(tabb, tabn - r1, r1);
 
186
 
 
187
  if ((r = pb - pa) > es )
 
188
    sciqsort(a, tab,flag, r / es, es, es1, cmp,swapcode,swapcodeind);
 
189
  if ((r = pd - pc) > es) { 
 
190
    /* Iterate rather than recurse to save stack space */
 
191
    a = pn - r;
 
192
    tab = tabn - (tabd - tabc);
 
193
    n = r / es;
 
194
    goto loop;
 
195
  }
 
196
  /*            qsort(pn - r, r / es, es, cmp);*/
 
197
}
 
198
 
 
199
 
 
200