1 /*
2 * $Id: rank.c 223 2006-08-17 06:19:38Z rgb $
3 *
4 * See copyright in copyright.h and the accompanying file COPYING
5 *
6 */
7
8 /*
9 *========================================================================
10 * This code evaluates the binary rank of a bit matrix on the field[0,1],
11 * for example:
12 *
13 * 0 1 1 0
14 * 1 0 1 0
15 * 1 0 1 1
16 * 0 0 1 1
17 *
18 * would be a 4x4 bitwise matrix. The rank of the matrix is determined
19 * as is the rank of any matrix, by gauss elimination with pivoting, but
20 * in this case the "pivots" are 1 bits (the largest element in any
21 * column or row) and the required transpositions and subtractions are
22 * all done by pointer or bitwise (xor) operations.
23 *
24 *========================================================================
25 */
26
27
28 #include <dieharder/libdieharder.h>
29
30
binary_rank(uint ** mtx,int mrows,int ncols)31 int binary_rank(uint **mtx,int mrows,int ncols)
32 {
33
34 int i,j,k,m,n,s_uint,rank;
35 int col_ind,uint_col_max;
36 uint mask,colchk;
37 uint *rowp;
38
39 /*
40 * mtx is a vector of unsigned integers, e.g.:
41 *
42 * mtx[0] = 0110...110
43 * mtx[1] = 1010...001
44 * ...
45 * mtx[ncols-1] = 0010...100
46 *
47 * We go through this vector a row at a time, searching each
48 * column for a pivot (1). When we find a pivot, we swap rows
49 * and eliminate the column bitwise.
50 */
51
52 /*
53 * size of an unsigned int
54 */
55 s_uint = 8*sizeof(uint);
56 /*
57 * row size in uint columns. Note that we have to remember
58 * which uint ** column we are in and therefore have to
59 * convert bit column into uint column regularly.
60 * Subtract 1, because it is zero-basd.
61 */
62 uint_col_max = (ncols-1)/s_uint;
63
64 if(verbose == D_BRANK || verbose == D_ALL){
65 printf("Starting bitmatrix:\n");
66 for(m=0;m<mrows;m++){
67 printf("# br: ");
68 dumpbits(&mtx[m][0],32);
69 }
70 }
71
72 rank = 0;
73 i = 0;
74 mask = 1;
75 /*
76 * j is the column BIT index, which can be
77 * larger than s_uint (or rmax_bits).
78 */
79 for(j = 0;j < ncols && i < mrows;j++){
80 /*
81 * This is the mxt[i][j] index of the
82 * current bit column.
83 */
84 col_ind = j/s_uint;
85
86 /*
87 * This handles the transition to the next
88 * uint when the bit index j passes the uint
89 * boundary. mask picks out the correct bit
90 * column from right to left (why not from
91 * left to right?).
92 */
93 j%s_uint ? (mask <<=1) : (mask = 1);
94
95 /*
96 * Find a pivot element (a 1) if there is one
97 * in this column (column fixed by mask).
98 */
99 if(verbose == D_BRANK || verbose == D_ALL){
100 printf("Checking column mask ");
101 dumpbits(&mask,32);
102 }
103 for(k=i;k<mrows;k++){
104 colchk = mtx[k][col_ind]&mask;
105 if(verbose == D_BRANK || verbose == D_ALL){
106 printf("row %d = ",k);
107 dumpbits(&colchk,32);
108 }
109 if(colchk) break;
110 }
111
112 /*
113 * OK, k either points to a row with a "1" in
114 * the mask column or it equals mrows. In the latter
115 * case, the entire column from i on down is zero
116 * and we move on without incrementing rank.
117 *
118 * Otherwise, we bring the pivot ROW to the ith position
119 * (which may involve doing nothing). k is set to
120 * point to the first location that could still be
121 * a 1, which is always k+1;
122 */
123 if(k < mrows){
124 if(verbose == D_BRANK || verbose == D_ALL){
125 printf("Swapping %d and %d rows. before bitmatrix:\n",i,k);
126 for(m=0;m<mrows;m++){
127 printf("# br: ");
128 dumpbits(&mtx[m][col_ind],32);
129 }
130 }
131 if(i != k){
132 if(verbose == D_BRANK || verbose == D_ALL){
133 printf("before: mtx[%d] = %p mtx[%d = %p\n",i,(void*) mtx[i],k,(void*) mtx[k]);
134 }
135 rowp = mtx[i]; /* rowp is the ADDRESS of the ith row */
136 mtx[i] = mtx[k];
137 mtx[k] = rowp;
138 if(verbose == D_BRANK || verbose == D_ALL){
139 printf("after mtx[%d] = %p mtx[%d = %p\n",i,(void*) mtx[i],k,(void*) mtx[k]);
140 }
141 }
142 if(verbose == D_BRANK || verbose == D_ALL){
143 printf("Swapped %d and %d rows. after bitmatrix:\n",i,k);
144 for(m=0;m<mrows;m++){
145 printf("# br: ");
146 dumpbits(&mtx[m][col_ind],32);
147 }
148 }
149 k++; /* First row that could have a 1 in this column after i */
150
151 /*
152 * Now we eliminate the rest of the column, by rows, starting
153 * at k.
154 */
155 while(k<mrows){
156 /*
157 * if the row also has a 1 in this column...
158 */
159 if(mtx[k][col_ind] & mask){
160 /*
161 * we use exclusive or to eliminate the
162 * rest of the column.
163 */
164 n = uint_col_max - col_ind;
165 if(verbose == D_BRANK || verbose == D_ALL){
166 printf("eliminating against row %2d = ",i);
167 dumpbits(&mtx[i][col_ind],32);
168 printf("eliminating row %2d, before = ",k);
169 dumpbits(&mtx[k][col_ind],32);
170 }
171 while (n >= 0){
172 if(verbose == D_BRANK || verbose == D_ALL){
173 printf("xoring column = %2d\n",n);
174 }
175 mtx[k][n] ^= mtx[i][n];
176 n--;
177 }
178 if(verbose == D_BRANK || verbose == D_ALL){
179 printf("eliminating row %2d, after = ",k);
180 dumpbits(&mtx[k][col_ind],32);
181 printf("\n");
182 }
183 }
184 k++;
185 }
186 if(verbose == D_BRANK || verbose == D_ALL){
187 printf("Eliminated. New bitmatrix:\n");
188 for(m=0;m<mrows;m++){
189 printf("# br: ");
190 dumpbits(&mtx[m][col_ind],32);
191 }
192 }
193 i++;
194 if(verbose == D_BRANK || verbose == D_ALL){
195 printf("NEW RANK = %d\n",i);
196 }
197 }
198 }
199
200 return(i);
201
202 }
203
204