8.5 Selecting the Mth Largest 343 requires looking at all N elements,if only to find those that are still alive,while the bisections are dominated by the N that occur in the first round.Minimizing O(N logM N)+O(N log2 M)thus yields the result M2Vlog2 N (8.5.1) The square root of the logarithm is so slowly varying that secondary considerations of machine timing become important.We use M=64 as a convenient constant value. Two minor additional tricks in the following routine,selip,are(i)augmenting the set of M random values by an M+Ist,the arithmetic mean,and(ii)choosing 8 the M random values"on the fly"in a pass through the data,by a method that makes later values no less likely to be chosen than earlier ones.(The underlying idea is to give element m>M an M/m chance of being brought into the set.You can prove by induction that this yields the desired result.) 3 #include "nrutil.h" #define M 64 RECIPES #define BIG 1.0e30 #define FREEALL free_vector(sel,1,M+2);free_lvector(isel,1,M+2); float selip(unsigned long k,unsigned long n,float arr[]) Returns the kth smallest value in the array arr [1..n].The input array is not altered. computer, Ameri Press. ART void shel1(unsigned long n,f1oata[▣)； unsigned long i,j,j1,jm,ju,kk,mm,nlo,nxtmm,*isel; Programs float ahi,alo,sum,*sel; if (k<1 lI k>n ll n <0)nrerror("bad input to selip"); isel=lvector(1,M+2); sel=vector(1,M+2); kk=k; to dir ahi-BIG; alo =-BIG; for (;;) Main iteration loop,until desired ele- 1920 SCIENTIFIC COMPUTING(ISBN mm=nlo=0; ment is isolated. sum=0.0; nxtmm=M+1; for(1=1;i<=n;1++) Make a pass through the whole array. 10.6211 if (arr[i]>=aloarr[i]<ahi){ Consider only elements in the current brackets. idge.org Numerical Recipes 43108 mm++; if (arr[i]=alo)nlo++; In case of ties for low bracket. Now use statistical procedure for selecting m in-range elements with equal probability,even without knowing in advance how many there are! (outside if (mm <M)sel[mm]-arr[i]; North Software. else if (mm ==nxtmm){ nxtmm=mm+mm/M; sel[1 ((i+mm+kk)%M)]=arr[i];The operation provides a some- what random number. sum +=arr[i]; if (kk <nlo){ Desired element is tied for lower bound; FREEALL return it. return alo; 2 else if (mm <M){ All in-range elements were kept.So re- shell(mm,sel); turn answer by direct method. ahi sel[kk];8.5 Selecting the Mth Largest 343 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). requires looking at all N elements, if only to find those that are still alive, while the bisections are dominated by the N that occur in the first round. Minimizing O(N logM N) + O(N log2 M) thus yields the result M ∼ 2 √log2 N (8.5.1) The square root of the logarithm is so slowly varying that secondary considerations of machine timing become important. We use M = 64 as a convenient constant value. Two minor additional tricks in the following routine, selip, are (i) augmenting the set of M random values by an M + 1st, the arithmetic mean, and (ii) choosing the M random values “on the fly” in a pass through the data, by a method that makes later values no less likely to be chosen than earlier ones. (The underlying idea is to give element m>M an M/m chance of being brought into the set. You can prove by induction that this yields the desired result.) #include "nrutil.h" #define M 64 #define BIG 1.0e30 #define FREEALL free_vector(sel,1,M+2);free_lvector(isel,1,M+2); float selip(unsigned long k, unsigned long n, float arr[]) Returns the kth smallest value in the array arr[1..n]. The input array is not altered. { void shell(unsigned long n, float a[]); unsigned long i,j,jl,jm,ju,kk,mm,nlo,nxtmm,*isel; float ahi,alo,sum,*sel; if (k < 1 || k > n || n <= 0) nrerror("bad input to selip"); isel=lvector(1,M+2); sel=vector(1,M+2); kk=k; ahi=BIG; alo = -BIG; for (;;) { Main iteration loop, until desired ele￾mm=nlo=0; ment is isolated. sum=0.0; nxtmm=M+1; for (i=1;i<=n;i++) { Make a pass through the whole array. if (arr[i] >= alo && arr[i] <= ahi) { Consider only elements in the current brackets. mm++; if (arr[i] == alo) nlo++; In case of ties for low bracket. Now use statistical procedure for selecting m in-range elements with equal probability, even without knowing in advance how many there are! if (mm <= M) sel[mm]=arr[i]; else if (mm == nxtmm) { nxtmm=mm+mm/M; sel[1 + ((i+mm+kk) % M)]=arr[i]; The % operation provides a some- } what random number. sum += arr[i]; } } if (kk <= nlo) { Desired element is tied for lower bound; FREEALL return it. return alo; } else if (mm <= M) { All in-range elements were kept. So re￾shell(mm,sel); turn answer by direct method. ahi = sel[kk];