Arndt - Algorithms for Programmers (523138), страница 23
Текст из файла (страница 23)
PERMUTATIONSelse ++v;bit_subset b(v);do{// --- do cycle: --ulong i = z | b.next(); // start of cycleType t = f[i];// save start valueulong g = gray_code(i); // next in cyclefor (ulong k=cl-1; k!=0; --k){f[i] = f[g];i = g;g = gray_code(i);}f[i] = t;// --- end (do cycle) --}while ( b.current() );}}135How fast is it? We use the convention that the speed of the trivial (and completely cache-friendly,therefore running at memory bandwidth) reverse is 1.0, our hereby declared time unit for comparison.A little benchmark looks like:CLOCK defined as 1000 MHz // AMD Athlon 1000MHz with 100MHz DDR RAMmemsize=32768 kiloByte // permuting that much memory (in chunks of doubles)reverse(fr,n2);dt= 0.0997416rel=1 // set to onerevbin_permute(fr,n2);dt=0.594105rel=5.95644reverse(fr,n2);dt= 0.0997483rel=1.00007gray_permute(fr,n2);dt=0.119014rel=1.19323reverse(fr,n2);dt= 0.0997618rel=1.0002inverse_gray_permute(fr,n2);dt=0.11028rel=1.10566reverse(fr,n2);dt= 0.0997424rel=1.00001We repeatedly timed reverse to get an impression how much we can trust the observed numbers.
Thebandwidth of the reverse is about 320MByte/sec which should be compared to the output of a specialmemory testing program, revealing that it actually runs at about 83% of the bandwidth one can getwithout using streaming instructions:avg:avg:avg:avg:avg:avg:avg:avg:avg:avg:avg:avg:avg:33554432335544323355443233554432335544323355443233554432335544323355443233554432335544323355443233554432[ 0]"memcpy"[ 1]"char *"[ 2]"short *"[ 3]"int *"[ 4]"long *"[ 5]"long * (4x unrolled)"[ 6]"int64 *"[ 7]"double *"[ 8]"double * (4x unrolled)"[ 9]"streaming K7"[10]"streaming K7 prefetch"[11]"streaming K7 clear"[12]"long * clear"305.869154.713187.943300.720300.584306.135305.372388.695374.271902.1711082.8681318.875341.456MB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/sMB/s// <--=While the revbin_permute takes about 6 units (due to its memory access pattern that is very problematicwrt.
cache usage) the gray_permute only uses 1.20 units, the inverse_gray_permute even5 only 1.10!This is pretty amazing for such a nontrivial permutation.The described permutation can be used to significantly speed up fast transforms of lengths a power oftwo, notably the Walsh transform, see chapter 5.It is instructive to study the complementary masks that occur for cycles off different lengths.
The systemin the structure of the cycles in gray_permute(f, 128) seems non-obvious:0:1:((5 The2,4,3) #=27, 5,6) #=4observed difference between the forward- and backward version is in fact systematic.CHAPTER 7. PERMUTATIONS2: ( 8, 15, 10, 12) #=43: ( 9, 14, 11, 13) #=44: ( 16, 31, 21, 25, 17, 30,5: ( 18, 28, 23, 26, 19, 29,6: ( 32, 63, 42, 51, 34, 60,7: ( 33, 62, 43, 50, 35, 61,8: ( 36, 56, 47, 53, 38, 59,9: ( 37, 57, 46, 52, 39, 58,10: ( 64,127, 85,102, 68,120,11: ( 65,126, 84,103, 69,121,12: ( 66,124, 87,101, 70,123,13: ( 67,125, 86,100, 71,122,14: ( 72,112, 95,106, 76,119,15: ( 73,113, 94,107, 77,118,16: ( 74,115, 93,105, 78,116,17: ( 75,114, 92,104, 79,117,126 elements in 18 nontrivialcycle lengths:2 ...82 fixed points: [0.
1]20, 24)22, 27)40, 48)41, 49)45, 54)44, 55)80, 96)81, 97)82, 99)83, 98)90,108)91,109)88,111)89,110)cycles.136#=8#=8#=8#=8#=8#=8#=8#=8#=8#=8#=8#=8#=8#=8However, one can identify the cycle maxima as the set of different possible values that consist of the bitsof an INV-mask combined with all variations of bits of the VAR-mask (filled in as max for n < 128):-------------------------ldm= 1:1 cycles of length= 2 [2..3]max: .......................11 = 3 cl=2.......................11 = INV......................... = VAR-------------------------ldm= 2:1 cycles of length= 4 [4..7]max: ......................111 = 7 cl=4......................111 = INV.........................
= VARldm= 3:2 cycles of length= 4 [8..f]max: .....................1111 = 15 cl=4max: .....................111. = 14 cl=4.....................111. = INV........................1 = VAR-------------------------ldm= 4:2 cycles of length= 8 [10..1f]max: ....................11111 = 31 cl=8max: ....................111.1 = 29 cl=8....................111.1 = INV.......................1. = VARldm= 5:4 cycles of length= 8 [20..3f]max: ...................111.11 = 59 cl=8max: ...................11111. = 62 cl=8max: ...................111111 = 63 cl=8max: ...................111.1. = 58 cl=8...................111.1.
= INV......................1.1 = VARldm= 6:8 cycles of length= 8 [40..7f]max: ..................111.1.1 = 117 cl=8max: ..................111.11. = 118 cl=8max: ..................111.111 = 119 cl=8max: ..................11111.. = 124 cl=8max: ..................11111.1 = 125 cl=8max: ..................111111. = 126 cl=8max: ..................1111111 = 127 cl=8max: ..................111.1.. = 116 cl=8..................111.1.. = INV.....................1.11 = VARldm= 7:16 cycles of length= 8 [80..ff]max: [...].................111.1... = INV....................1.111 = VAR-------------------------ldm= 8:16 cycles of length=16 [100..1ff]max: [...]................111.1...1 = INV...................1.111. = VARldm= 9:32 cycles of length=16 [200..3ff]max: [...]...............111.1...1.
= INV..................1.111.1 = VARCHAPTER 7. PERMUTATIONS7.90:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:137The reversed Gray code permutation[* ][*][*][*][*][*][*][*][ *][*][*][*][*][*][*][*]0:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:[*][*][*][*][*][* ][*][*][*][*][*][*][*][*][*][ *]Figure 7.8: Permutation matrices of the reversed Gray code permutation (left) and its inverse (right).If the length-n array is permuted in the way the upper half of the length-2n array would be permuted bygray_permute() then all cycles are of the same length. The resulting permutation is equivalent to thereversed Gray code permutation:template <typename Type>inline void gray_rev_permute(const Type *f, Type * restrict g, ulong n)// gray_rev_permute() =^=// { reverse(); gray_permute(); }{for (ulong k=0, m=n-1; k<n; ++k, --m) g[gray_code(m)] = f[k];}The routine, its inverse and in-place versions can be found in [FXT: file perm/grayrevpermute.h].All cycles have the same length, gray_rev_permute(f, 64) gives:0: ( 0, 63, 21, 38, 4, 56,1: ( 1, 62, 20, 39, 5, 57,2: ( 2, 60, 23, 37, 6, 59,3: ( 3, 61, 22, 36, 7, 58,4: ( 8, 48, 31, 42, 12, 55,5: ( 9, 49, 30, 43, 13, 54,6: ( 10, 51, 29, 41, 14, 52,7: ( 11, 50, 28, 40, 15, 53,64 elements in8 nontrivialcycle length is == 8No fixed points.16, 32)17, 33)18, 35)19, 34)26, 44)27, 45)24, 47)25, 46)cycles.#=8#=8#=8#=8#=8#=8#=8#=8If 64 is added to the cycle elements then the cycles in the upper half of the array as in ofgray_permute(f, 128) are reproduced (this is by construction).Let G denote the Gray code permutation, Ḡ the reversed Gray code permutation.
Symbolically one canwrite©ªG(n) =. . . , Ḡ(n/8), Ḡ(n/4), Ḡ(n/2)(7.14)©ª−1−1−1−1G (n) =. . . , Ḡ (n/8), Ḡ (n/4), Ḡ (n/2)(7.15)Now let r the reversion and h the permutation that swaps the upper and the lower half of an array. ThenḠ =Ḡ=−1Ḡ G =Grr G−1G−1 Ḡ = r = Xn−1(7.16)(7.17)(7.18)G Ḡ−1Ḡ G−1 = h = Xn/2(7.19)−1=CHAPTER 7. PERMUTATIONS138Throughout it is assumed that the array length n is a power of two.7.100:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:The green code permutation[ *][* ][*][*][*][*][*][*][*][*][*][*][*][*][*][*]0:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:[ *][*][*][*][*][* ][*][*][*][*][*][*][*][*][*][*]Figure 7.9: Permutation matrices of the green code permutation (left) and its inverse (right).The green code permutation is obtained by replacing the calls to gray_code by green_code in the Graycode permutation. An additional step is required: the highest bit must be masked out in order to keepindices in the range.[FXT: green permute in perm/greenpermute.h]:template <typename Type>inline void green_permute(const Type *f, Type * restrict g, ulong n){for (ulong k=0; k<n; ++k){ulong r = green_code(k);r &= (n-1);g[r] = f[k];}}The inverse istemplate <typename Type>inline void inverse_green_permute(const Type *f, Type * restrict g, ulong n){for (ulong k=0; k<n; ++k){ulong r = green_code(k);r &= (n-1);g[k] = f[r];}}7.11The reversed green code permutationThe analogue to gray_rev_permute for the green permutation istemplate <typename Type>inline void green_rev_permute(const Type *f, Type * restrict g, ulong n){for (ulong k=0, m=n-1; k<n; ++k, --m)CHAPTER 7.
PERMUTATIONS0:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:139[* ][ *][*][*][*][*][*][*][*][*][*][*][*][*][*][*]0:1:2:3:4:5:6:7:8:9:10:11:12:13:14:15:[*][*][*][*][*][*][*][*][*][*][* ][*][*][*][*][ *]Figure 7.10: Permutation matrices of the reversed green code permutation (left) and its inverse (right).{}}ulong r = green_code(m);r &= (n-1);g[r] = f[k];The inverse istemplate <typename Type>inline void inverse_green_rev_permute(const Type *f, Type * restrict g, ulong n){for (ulong k=0, m=n-1; k<n; ++k, --m){ulong r = green_code(m);r &= (n-1);g[k] = f[r];}}Let E denote the green code permutation, Ē the reversed green code permutation and R the revbinpermutation. Further assume that the array length n is a power of two.