diff --git a/annotated.html b/annotated.html index 7ee4d61..55f4851 100644 --- a/annotated.html +++ b/annotated.html @@ -78,7 +78,8 @@ - + +
 Ccsr_t
 Cint_pair
 Cparams_t
 CONE_VEC_TDistance of a classical or quantum CSS code
 Cparams_t
diff --git a/classes.html b/classes.html index 5eff584..6580238 100644 --- a/classes.html +++ b/classes.html @@ -74,7 +74,7 @@
Class Index
-
C | I | P
+
C | I | O | P
C
@@ -83,6 +83,9 @@
I
int_pair
+
O
+
ONE_VEC_T
+
P
params_t
diff --git a/dir_68267d1309a1af8e8297ef4c3efbcdba.html b/dir_68267d1309a1af8e8297ef4c3efbcdba.html index cca9389..8e0258a 100644 --- a/dir_68267d1309a1af8e8297ef4c3efbcdba.html +++ b/dir_68267d1309a1af8e8297ef4c3efbcdba.html @@ -81,6 +81,8 @@ + + diff --git a/dist__cc_8c.html b/dist__cc_8c.html new file mode 100644 index 0000000..f6997ea --- /dev/null +++ b/dist__cc_8c.html @@ -0,0 +1,345 @@ + + + + + + + +Dist m4ri: src/dist_cc.c File Reference + + + + + + + + + + + +
+
+

Files

file  dist_cc.c [code]
 
file  dist_m4ri.c [code]
 
file  dist_m4ri.h [code]
+ + + + + +
+
Dist m4ri +  0.0.1.alpha +
+
Computing distance of a classical or quantum CSS code
+
+
+ + + + + + + + +
+
+ + +
+ +
+ + + +
+
+Classes | +Typedefs | +Functions
+
+
dist_cc.c File Reference
+
+
+
#include <inttypes.h>
+#include <strings.h>
+#include <stdlib.h>
+#include <time.h>
+#include <m4ri/m4ri.h>
+#include "mmio.h"
+#include "util_m4ri.h"
+#include "util_io.h"
+
+Include dependency graph for dist_cc.c:
+
+
+ + + + + + + + + + + + +
+
+

Go to the source code of this file.

+ + + + + +

+Classes

struct  ONE_VEC_T
 distance of a classical or quantum CSS code More...
 
+ + + + +

+Typedefs

typedef struct ONE_VEC_T one_vec_t
 distance of a classical or quantum CSS code More...
 
+ + + + + + + + + +

+Functions

void one_vec_print (const one_vec_t *const pvec)
 print entire one_vec_t structure by pointer More...
 
int start_CC_recurs (one_vec_t *err, one_vec_t *urr, one_vec_t *const syn[], const int wmax, const int max_col_wt, const csr_t *const mH, const csr_t *const mHT, const csr_t *const mL, const int debug)
 recursively construct codewords More...
 
int do_CC_dist (const csr_t *const mH, const csr_t *mL, const int wmax, const int start, const int debug)
 
+

Typedef Documentation

+ +

◆ one_vec_t

+ +
+
+ + + + +
typedef struct ONE_VEC_T one_vec_t
+
+ +

distance of a classical or quantum CSS code

+
+

The program implements two methods:

    +
  1. Random information set (random window) algorithm (upper bound).
    + This works with any code (LDPC or not). (2) depth-first codeword enumeration (connected cluster) algorithm (Lower bound or actual distance if a codeword is found.)
    +
  2. +
+

A. Dumer, A. A. Kovalev, and L. P. Pryadko "Distance verification..." in IEEE Trans. Inf. Th., vol. 63, p. 4675 (2017). doi: 10.1109/TIT.2017.2690381

+

author: Leonid Pryadko leoni.nosp@m.d.pr.nosp@m.yadko.nosp@m.@ucr.nosp@m..edu, Weilei Zeng

+ +
+
+

Function Documentation

+ +

◆ do_CC_dist()

+ +
+
+ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
int do_CC_dist (const csr_t *const mH,
const csr_tmL,
const int wmax,
const int start,
const int debug 
)
+
+

rewrite of the cluster method function using only sparse matrices try recursive version first
+

+

prepare the 1st error vector and the syndrome

+

go up

+

verify the vector

+

classical code

+

success

+

codeword weight found

+

not found a codeword up to wmax

+ +

Definition at line 256 of file dist_cc.c.

+ +

References csr_t::cols, csr_max_row_wght(), csr_print(), csr_transpose(), ERROR, ONE_VEC_T::max, csr_t::rows, ONE_VEC_T::vec, and ONE_VEC_T::wei.

+ +
+
+ +

◆ one_vec_print()

+ +
+
+ + + + + + + + +
void one_vec_print (const one_vec_t *const pvec)
+
+ +

print entire one_vec_t structure by pointer

+ +

Definition at line 37 of file dist_cc.c.

+ +

References ONE_VEC_T::vec, and ONE_VEC_T::wei.

+ +

Referenced by start_CC_recurs().

+ +
+
+ +

◆ start_CC_recurs()

+ +
+
+ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
int start_CC_recurs (one_vec_terr,
one_vec_turr,
one_vec_t *const syn[],
const int wmax,
const int max_col_wt,
const csr_t *const mH,
const csr_t *const mHT,
const csr_t *const mL,
const int debug 
)
+
+ +

recursively construct codewords

+
Parameters
+ + + + +
errerror vector with sorted components
urrunsorted vector so far
synarray of syndrome vectors with sorted components (indexed by weight of error)
+
+
+

row with the first non-zero syndrome bit

+

all valid positions should be to the right of here

+

not there

+

go up

+

reachable goal?

+

classical code

+

success, just get out fast

+

nothing found

+ +

Definition at line 183 of file dist_cc.c.

+ +

References csr_t::i, one_vec_print(), p, csr_t::p, ONE_VEC_T::vec, and ONE_VEC_T::wei.

+ +
+
+
+ + + + diff --git a/dist__cc_8c__incl.map b/dist__cc_8c__incl.map new file mode 100644 index 0000000..ef490be --- /dev/null +++ b/dist__cc_8c__incl.map @@ -0,0 +1,12 @@ + + + + + + + + + + + + diff --git a/dist__cc_8c__incl.md5 b/dist__cc_8c__incl.md5 new file mode 100644 index 0000000..bd8a849 --- /dev/null +++ b/dist__cc_8c__incl.md5 @@ -0,0 +1 @@ +8f807167e118aa3d987b5b72cb05b6bc \ No newline at end of file diff --git a/dist__cc_8c__incl.png b/dist__cc_8c__incl.png new file mode 100644 index 0000000..5a5b7c9 Binary files /dev/null and b/dist__cc_8c__incl.png differ diff --git a/dist__cc_8c_source.html b/dist__cc_8c_source.html new file mode 100644 index 0000000..1468217 --- /dev/null +++ b/dist__cc_8c_source.html @@ -0,0 +1,410 @@ + + + + + + + +Dist m4ri: src/dist_cc.c Source File + + + + + + + + + + + +
+
+ + + + + + +
+
Dist m4ri +  0.0.1.alpha +
+
Computing distance of a classical or quantum CSS code
+
+
+ + + + + + + + +
+
+ + +
+ +
+ + +
+
+
+
dist_cc.c
+
+
+Go to the documentation of this file.
1 
+
17 // #include <m4ri/config.h>
+
18 #include <inttypes.h>
+
19 #include <strings.h>
+
20 #include <stdlib.h>
+
21 #include <time.h>
+
22 #include <m4ri/m4ri.h>
+
23 
+
24 #include "mmio.h"
+
25 #include "util_m4ri.h"
+
26 #include "util_io.h"
+
27 //#include "dist_m4ri.h"
+
28 // #include "util.h"
+
29 
+
30 typedef struct ONE_VEC_T{
+
31  int wei;
+
32  int max;
+
33  int vec[0];
+ +
35 
+
37 void one_vec_print(const one_vec_t * const pvec){
+
38  printf(" w=%d [ ",pvec->wei);
+
39  for(int i=0; i < pvec->wei; i++)
+
40  printf("%d ",pvec->vec[i]);
+
41  printf("]\n");
+
42 }
+
43 
+
44 
+
45 // v1[:] = v0[:] + mat[row,:]
+
46 static inline int one_csr_row_combine(one_vec_t * const v1, const one_vec_t * const v0,
+
47  const csr_t * const mat, const int row){
+
48 #ifndef NDEBUG
+
49  if ((!v1) || (!v0) || (!mat))
+
50  ERROR("all arguments must be allocated: v1=%p v0=%p mat=%p\n",v1,v0,mat);
+
51  if(v1 == v0)
+
52  ERROR("the two vectors should not be the same !");
+
53  if((row<0) || (row >= mat->rows) ||
+
54  (v1->max < mat->cols) ||
+
55  (v0->max < mat->cols))
+
56  ERROR("this should not happen\n");
+
57 #endif
+
58  int iM, i0=0, i1=0;
+
59  for (iM = mat->p[row]; iM < mat->p[row+1]; iM++){
+
60  const int ic = mat->i[iM];
+
61  while((i0 < v0->wei) && (v0->vec[i0] < ic))
+
62  v1->vec[i1++] = v0->vec[i0++];
+
63  if(i0 >= v0->wei)
+
64  break;
+
65  if(v0->vec[i0]==ic)
+
66  i0++;
+
67  else
+
68  v1->vec[i1++] = ic;
+
69  }
+
70  if(i0 >= v0->wei)
+
71  for ( ; iM < mat->p[row+1]; iM++){
+
72  const int ic = mat->i[iM];
+
73  v1->vec[i1++] = ic;
+
74  }
+
75  else
+
76  while(i0 < v0->wei)
+
77  v1->vec[i1++] = v0->vec[i0++];
+
78 
+
79  v1->wei = i1;
+
80  return i1;
+
81 }
+
82 
+
84 static inline int one_ordered_ins(one_vec_t * const err, const int j){
+
85  int pos=err->wei-1;
+
86  while(j < err->vec[pos]){
+
87  err->vec[pos+1] = err->vec[pos];
+
88  pos--;
+
89  }
+
90 #ifndef NDEBUG
+
91  if (j == err->vec[pos])
+
92  ERROR("Unexpected! vec[%d]=%d is already present!",pos,j);
+
93 #endif
+
94  err->vec[pos+1]=j;
+
95 #ifndef NDEBUG
+
96  if(err->wei >= err->max)
+
97  ERROR("unexpected!");
+
98  for(int i=0; i < err->wei; i++)
+
99  if(err->vec[i] >= err->vec[i+1]){
+
100  printf("check ordering at i=%d! ",i);
+
101  one_vec_print(err);
+
102  ERROR("unexpected");
+
103  }
+
104 #endif
+
105  err->wei ++;
+
106  return pos+1;
+
107 }
+
108 
+
109 
+
113 static inline int one_ordered_search(one_vec_t * const err, const int val){
+
115  int bot=0, top=err->wei , mid=0;
+
116 #ifndef NDEBUG
+
117  if (!top)
+
118  return -1;
+
119 #endif
+
120  while(top - bot > 1){
+
121  mid = (top+bot) >> 1;
+
122 #ifndef NDEBUG
+
123  if (mid>=err->wei)
+
124  ERROR("this should not happen");
+
125 #endif
+
126  if (err->vec[mid] <= val)
+
127  bot = mid;
+
128  else
+
129  top = mid;
+
130  // printf("bot=%d mid=%d top=%d err->vec[mid]=%d val=%d\n",bot,mid,top,err->vec[mid],val);
+
131  }
+
132  if ( err->vec[bot] == val)
+
133  return bot;
+
134  else
+
135  return -1;
+
136 }
+
137 
+
141 static inline int one_ordered_find_del(one_vec_t * const err, const int val){
+
143  int bot=0, top=err->wei, mid=0;
+
144 #ifndef NDEBUG
+
145  if (!top)
+
146  return 0;
+
147 #endif
+
148  while(top - bot > 1){
+
149  mid = (top+bot) >> 1;
+
150  if (err->vec[mid] <= val)
+
151  bot = mid;
+
152  else
+
153  top = mid;
+
154  }
+
155  if ( err->vec[mid] != val)
+
156  return 0;
+
157  // ERROR("unexpected! value val=%d not found",val);
+
158 
+
159  for(int i=mid; i < err->wei; i++)
+
160  err->vec[i] = err->vec[i+1];
+
161  err->wei --;
+
162  return 1;
+
163 }
+
164 
+
165 
+
167 static inline void one_ordered_pos_del(one_vec_t * const err, _maybe_unused const int val, const int pos){
+
168 #ifndef NDEBUG
+
169  if ((pos<0) || (pos >= err->wei) || (err->wei == 0) || (err->vec[pos] != val))
+
170  ERROR("this should not happen!");
+
171 #endif
+
172  err->wei --;
+
173  for(int i=pos; i < err->wei; i++)
+
174  err->vec[i] = err->vec[i+1];
+
175 }
+
176 
+
183 int start_CC_recurs(one_vec_t *err, one_vec_t *urr, one_vec_t * const syn[],
+
184  const int wmax, const int max_col_wt,
+
185  const csr_t * const mH, const csr_t * const mHT, const csr_t * const mL,
+
186  const int debug){
+
187  const int w=err->wei;
+
188  int row = syn[w]->vec[0];
+
189 #ifndef NDEBUG
+
190  if(debug&64){
+
191  printf("starting CC recurs w=%d row=%d:\n urr: ",w,row);
+
192  one_vec_print(urr);
+
193  printf(" err: ");
+
194  one_vec_print(err);
+
195  for(int i=0; i <= w; i++){
+
196  printf("i=%d ",i);
+
197  one_vec_print(syn[i]);
+
198  }
+
199  }
+
200 #endif
+
201  const int col_min=urr->vec[0];
+
202  for(int i1 = mH->p[row]; i1 < mH->p[row+1]; i1++){
+
203  const int col = mH->i[i1];
+
204  if(col > col_min){
+
205  int pos = one_ordered_search(err, col);
+
206  if(pos == -1){
+
207  urr->vec[w] = col;
+
208  urr->wei++;
+
209 #ifndef NDEBUG
+
210  if(debug&64){
+
211  printf(" pos=%d urr: ",pos);
+
212  one_vec_print(urr);
+
213  }
+
214 #endif
+
215  pos = one_ordered_ins(err,col);
+
216  int swei = one_csr_row_combine(syn[w+1],syn[w], mHT, col);
+
217  int result = 0;
+
218  if (err->wei < wmax){
+
219  if (swei){
+
220  if(swei <= (wmax - err->wei)*max_col_wt){
+
221  result = start_CC_recurs(err,urr,syn,wmax,max_col_wt,mH,mHT,mL,debug);
+
222  if(result == 1)
+
223  return 1;
+
224  }
+
225  }
+
226  // swei == 0 means it is a degenerate vector
+
227  }
+
228  else{ // wei == wmax
+
229  if(!swei){
+
230  if((!mL) ||
+
231  (sparse_syndrome_non_zero(mL, err->wei, err->vec))){
+
232  if(debug&32){
+
233  printf("swei=%d *** success ***\n",swei);
+
234  one_vec_print(syn[w+1]);
+
235  }
+
236  return 1;
+
237  }
+
238  }
+
239  }
+
240  urr->wei--;
+
241  one_ordered_pos_del(err,col,pos);
+
242  if(debug&32){
+
243  printf("xerr: ");
+
244  one_vec_print(err);
+
245  }
+
246  }
+
247  }
+
248  }
+
249  if(debug&32)
+
250  printf("exiting CC recurs\n\n");
+
251  return 0;
+
252 }
+
253 
+
256 int do_CC_dist(const csr_t * const mH, const csr_t * mL,
+
257  const int wmax, const int start, const int debug){
+
258 
+
259  const int nchk = mH->rows, nvar = mH->cols;
+
260  if((start<-1) || (start>=nvar))
+
261  ERROR("invalid start=%d for H[%d,%d]\n",start, nchk, nvar);
+
262 
+
263  csr_t * const mHT = csr_transpose(NULL,mH);
+
264  if(debug&32){
+
265  if((mHT->cols<150)||(debug&2048))
+
266  csr_print(mHT,"HT");
+
267  }
+
268  int max_col_W = csr_max_row_wght(mHT);
+
269 
+
270  one_vec_t *err = calloc(1, sizeof(one_vec_t)+sizeof(int)*wmax);
+
271  one_vec_t *urr = calloc(1, sizeof(one_vec_t)+sizeof(int)*wmax);
+
272  one_vec_t **syn = calloc(wmax+1, sizeof(one_vec_t *));
+
273  if((!syn) || (!err) || (!urr))
+
274  ERROR("memory allocation");
+
275  err->max = wmax;
+
276  for(int i=0; i <= wmax; i++){
+
277  syn[i]=calloc(1, sizeof(one_vec_t)+sizeof(int)*mH->rows);
+
278  if(!syn[i]) ERROR("i=%d memory allocation",i);
+
279  syn[i]->max = mH->rows;
+
280  }
+
281  int result = 0;
+
282  for(int w=1; w <= wmax; w++){ /* cluster weight */
+
283  int beg = 0, end = nvar - wmax ;
+
284  if (start >= 0)
+
285  beg = end = start;
+
286  if(debug&2)
+
287  printf("# recursively searching for w=%d codewords wmax=%d beg=%d end=%d\n",w,wmax,beg,end);
+
288  for(int i = beg; i <= end; i++){ /* start column position */
+
290  err->vec[0] = urr->vec[0] = i;
+
291  err->wei = urr->wei = 1;
+
292  int swei = one_csr_row_combine(syn[1], syn[0], mHT, i);
+
293  if (1<w){
+
294  if (swei){
+
295  result = start_CC_recurs(err,urr,syn,w,max_col_W,mH,mHT,mL,debug);
+
296  if(result == 1)
+
297  break;
+
298  }
+
299  }
+
300  else{ // w==1
+
301  if(!swei){
+
302  if((!mL) ||
+
303  (sparse_syndrome_non_zero(mL, err->wei, err->vec))){
+
304  result = 1;
+
305  break;
+
306  }
+
307  }
+
308  }
+
309  err->wei = urr->wei = 0;
+
310  }
+
311  if(result == 1)
+
312  break;
+
313  }
+
314 
+
315  if(result==1){
+
316  result = err->wei;
+
317  if(debug&16){
+
318  printf("# wmax=%d found cw of weight %d: [",wmax,result);
+
319  int max = ((result<50) || (debug&2048)) ? result : 50 ;
+
320  for(int i=0; i< max; i++)
+
321  printf("%d%s",err->vec[i], i+1!=max?" ": (result==max ? "]\n" : "...]\n"));
+
322  }
+
323  }
+
324  else
+
325  result = -wmax;
+
327  for(int i=0; i<= wmax; i++)
+
328  free(syn[i]);
+
329  free(syn);
+
330  free(err);
+
331  free(urr);
+
332  csr_free(mHT);
+
333  return result;
+
334 }
+
335 
+
336 
+
void one_vec_print(const one_vec_t *const pvec)
print entire one_vec_t structure by pointer
Definition: dist_cc.c:37
+
int do_CC_dist(const csr_t *const mH, const csr_t *mL, const int wmax, const int start, const int debug)
Definition: dist_cc.c:256
+
int start_CC_recurs(one_vec_t *err, one_vec_t *urr, one_vec_t *const syn[], const int wmax, const int max_col_wt, const csr_t *const mH, const csr_t *const mHT, const csr_t *const mL, const int debug)
recursively construct codewords
Definition: dist_cc.c:183
+
struct ONE_VEC_T one_vec_t
distance of a classical or quantum CSS code
+
#define _maybe_unused
Definition: mmio.c:17
+ +
distance of a classical or quantum CSS code
Definition: dist_cc.c:30
+
int wei
Definition: dist_cc.c:31
+
int vec[0]
Definition: dist_cc.c:33
+
int max
Definition: dist_cc.c:32
+ +
int * p
Definition: util_m4ri.h:107
+
int * i
Definition: util_m4ri.h:108
+
int rows
Definition: util_m4ri.h:103
+
int cols
Definition: util_m4ri.h:104
+
params_t *const p
Definition: util_io.c:32
+ +
csr_t * csr_free(csr_t *p)
Definition: util_m4ri.c:461
+
void csr_print(const csr_t *const smat, const char str[])
Definition: util_m4ri.c:562
+
csr_t * csr_transpose(csr_t *dst, const csr_t *const p)
Definition: util_m4ri.c:135
+
int csr_max_row_wght(const csr_t *const p)
return max row weight of CSR matrix p TODO: add code for List of Pairs
Definition: util_m4ri.c:117
+ +
#define ERROR(fmt,...)
Definition: util_m4ri.h:13
+
+ + + + diff --git a/dist__m4ri_8c.html b/dist__m4ri_8c.html index 8f56ec4..7842c3c 100644 --- a/dist__m4ri_8c.html +++ b/dist__m4ri_8c.html @@ -75,6 +75,7 @@
+Macros | Functions
dist_m4ri.c File Reference
@@ -110,130 +111,38 @@

Go to the source code of this file.

+ + + +

+Macros

#define NEW   1
 
- - - - - - - - - - - + + + + + +

Functions

int do_dist_clus (const csr_t *const P, const mzd_t *const G, int debug, int wmax, int start, const int rankG)
 lower bound on the minimum distance by cluster enumeration More...
 
mzp_t * do_skip_pivs (const size_t rank, const mzp_t *const pivs)
 prepare an ordered pivot-skip list of length n-rank More...
 
int do_RW_dist (const csr_t *const spaH0, const csr_t *const spaL0, const int steps, const int wmin, const int classical, const int swait, const int debug)
 Random Information Set search for small-E logical operators. More...
 
int do_dist_rnd (csr_t *spaG0, mzd_t *matP0, int debug, int steps, int wmin)
 
mzp_t * do_skip_pivs (mzp_t *ans, const int rank, const mzp_t *const pivs)
 distance of a classical or quantum CSS code More...
 
int do_RW_dist (const csr_t *const spaH0, const csr_t *const spaL0, const int steps, const int wmin, const int classical, const int debug)
 Random Information Set search for small-E logical operators. More...
 
-

Function Documentation

- -

◆ do_dist_clus()

- -
-
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
int do_dist_clus (const csr_t *const P,
const mzd_t *const G,
int debug,
int wmax,
int start,
const int rankG 
)
-
- -

lower bound on the minimum distance by cluster enumeration

-

WARNING: only intended for LDPC codes

- -

Definition at line 141 of file dist_m4ri.c.

- -

References csr_t::cols, csr_transpose(), and csr_t::rows.

- -
-
- -

◆ do_dist_rnd()

+

Macro Definition Documentation

+ +

◆ NEW

- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - +
int do_dist_rnd (csr_tspaG0,
mzd_t * matP0,
int debug,
int steps,
int wmin 
)#define NEW   1
-

Definition at line 336 of file dist_m4ri.c.

-
- -

◆ do_RW_dist()

+

Function Documentation

+ +

◆ do_RW_dist()

@@ -268,12 +177,6 @@

const int  classical, - - - - const int  - swait, - @@ -300,7 +203,7 @@

Returns
minimum weight of a CW found (or -weigt if early termination condition is reached).

whether to verify logical ops as a vector or individually

-

actual vector

+

actual vector in sparse form

  1. Construct random column permutation P
@@ -310,31 +213,36 @@

Definition at line 223 of file dist_m4ri.c.

+

Definition at line 79 of file dist_m4ri.c.

-

References csr_t::cols, ERROR, and mzd_from_csr().

+

References csr_t::cols, ERROR, and mzd_from_csr().

- -

◆ do_skip_pivs()

+ +

◆ do_skip_pivs()

@@ -342,7 +250,13 @@

mzp_t* do_skip_pivs ( - const size_t  + mzp_t *  + ans, + + + + + const int  rank, @@ -359,11 +273,22 @@

+

distance of a classical or quantum CSS code

+
+

The program implements two methods:

    +
  1. Random information set (random window) algorithm (upper bound).
    + This works with any code (LDPC or not). (2) depth-first codeword enumeration (connected cluster) algorithm (Lower bound or actual distance if a codeword is found.)
    +
  2. +
+

A. Dumer, A. A. Kovalev, and L. P. Pryadko "Distance verification..." in IEEE Trans. Inf. Th., vol. 63, p. 4675 (2017). doi: 10.1109/TIT.2017.2690381

+

author: Leonid Pryadko leoni.nosp@m.d.pr.nosp@m.yadko.nosp@m.@ucr.nosp@m..edu, Weilei Zeng

prepare an ordered pivot-skip list of length n-rank

position to insert the next result

in skip_pivs

-

Definition at line 186 of file dist_m4ri.c.

+

Definition at line 30 of file dist_m4ri.c.

+ +

References ERROR.

diff --git a/dist__m4ri_8c_source.html b/dist__m4ri_8c_source.html index bf95aec..5415861 100644 --- a/dist__m4ri_8c_source.html +++ b/dist__m4ri_8c_source.html @@ -79,551 +79,285 @@
Go to the documentation of this file.
1 
-
2 /************************************************************************
-
3  * calculate distance of a quantum LDPC code using
-
4  * (1) random information set (random window) algorithm (upper bound)
-
5  * (2) depth-first codeword enumeration (cluster) algorithm (lower bound)
-
6  * A. Dumer, A. A. Kovalev, and L. P. Pryadko "Distance verification..."
-
7  * in IEEE Trans. Inf. Th., vol. 63, p. 4675 (2017).
-
8  * doi: 10.1109/TIT.2017.2690381
-
9  *
-
10  * currently: CSS only
-
11  *
-
12  * author: Leonid Pryadko <leonid.pryadko@ucr.edu>
-
13  ************************************************************************/
-
14 // #include <m4ri/config.h>
-
15 #include <inttypes.h>
-
16 #include <strings.h>
-
17 #include <stdlib.h>
-
18 #include <time.h>
-
19 #include <m4ri/m4ri.h>
-
20 
-
21 #include "mmio.h"
-
22 #include "util_m4ri.h"
-
23 #include "util_io.h"
-
24 #include "dist_m4ri.h"
-
25 // #include "util.h"
-
26 
-
31 static inline void addto_inline(mzd_t *row, const csr_t *spaQ, const int i){
-
32 #ifndef NDEBUG
-
33  if (i>=spaQ->rows)
-
34  ERROR("addto: attempt to get row=%d of %d",i,spaQ->rows);
-
35  if (row->ncols != spaQ->cols)
-
36  ERROR("addto: column number mismatch");
-
37  // mzd_print(row);
-
38 #endif
-
39  for (int j=spaQ->p[i]; j<spaQ->p[i+1]; j++)
-
40  mzd_flip_bit(row, 0,spaQ->i[j]); /* flip that bit */
-
41 }
-
42 
-
54 static inline int prep_neis(const int z0, int * nei, const mzd_t * v, _maybe_unused const mzd_t * s, const csr_t * P){
-
55  int cnt=0, max=P->p[z0+1];
-
56  for (int j=P->p[z0]; j<max; j++){
-
57  if (!mzd_read_bit(v,0,P->i[j])) /* never flip a bit twice */
-
58  nei[cnt++]=P->i[j];
-
59  }
-
60 #ifndef NDEBUG
-
61  // if(prm.debug&256){
-
62  printf("nei=[");
-
63  for(int i=0; i< cnt; i++)
-
64  printf(" %d ",nei[i]);
-
65  printf("]\n");
-
66  // }
-
67 #endif
-
68  return cnt;
-
69 }
-
70 
-
79 static int start_rec(const int w, const int wmax, mzd_t * v, mzd_t * s,
-
80  const csr_t * const P, const csr_t * const PT, const mzd_t * const G, const int rankG){
-
81  int res=0, all_zero=1;
-
82  mzd_t *v0=NULL;
-
83  word * rawrow = s->rows[0];
-
84  rci_t ns=v->ncols;
-
85 #ifndef NDEBUG
-
86  // if(prm.debug & 512){
-
87  printf("w=%d wgh(Pv)=%d \nv=",w,syndrome_bit_count(v, P));
-
88  mzd_print(v);
-
89  printf("s="); mzd_print(s);
-
90  // }
-
91 #endif
-
92 
-
93  for(rci_t i=-1 ; i+1< ns ; ){
-
94  i=nextelement(rawrow,s->width,i+1);
-
95  if (i<0)
-
96  break; /* no more non-zero syndrome bits */
-
97  else if(i<ns){ /* found a syndrome to fix */
-
98  if(w>=wmax)
-
99  return -1; /* failed to find a cluster */
-
100  all_zero=0;
-
101  int nei[maxrow];
-
102  int neisize=prep_neis(i,nei,v,s,P);
-
103  // cout << " row="<< i<< " nei=" << nei<< endl;
-
104  // if(neisize==0) ERROR("neisize=0"); /* this is actually OK: just continue */
-
105  for(int p=0; p< neisize; p++){
-
106  int j=nei[p];
-
107  mzd_write_bit(v,0,j,1);
-
108  addto_inline(s,PT,j);
-
109  res=start_rec(w+1,wmax,v,s,P,PT,G, rankG);
-
110  if(res==1)
-
111  return 1; /* just get out fast */
-
112  addto_inline(s,PT,j);
-
113  mzd_write_bit(v,0,j,0); /* clean up */
-
114  // printf("s="); mzd_print(s);
-
115  }
-
116  }
-
117  else break;
-
118  }
-
119  // printf("done with w=%d all_zero=%d\n",w,all_zero);
-
120  if(all_zero){ // syndrome vector was zero
-
121  if(G){
-
122  // cout<< "found v="<< v<< endl;
-
123  v0=mzd_copy(v0,v);
-
124  int result=do_reduce(v0,G,rankG);
-
125  if (result==-1)
-
126  return -1; /* go back down, the found row was trivial */
-
127 #ifndef NDEBUG
-
128  if (((int)mzd_weight(v)!=wmax)||(0!=syndrome_bit_count(v, P)))
-
129  ERROR(" start_rec: something is wrong %zu != wmax=%d ",mzd_weight(v),wmax);
-
130 #endif
-
131  }
-
132  return 1; /* success */
-
133  }
-
134  return 0; /* keep going */
-
135 }
-
136 
-
141 int do_dist_clus(const csr_t * const P, const mzd_t * const G, int debug, int wmax, int start, const int rankG){
-
142  // input: wmax=max cluster weight,
-
143  // start=initial position (-1 to scan all),
-
144  // P: LDPC parity check
-
145  // PT: LDPC parity check *** TRANSPOSED ***
-
146  // G: systematic form of the dual generator
-
147  // these must be orthogonal (no check)
-
148 
-
149  int nc=P->cols, ns=P->rows;
-
150  csr_t* PT=csr_transpose(NULL,P);
-
151  // csr_out(spaG0T);
-
152 
-
153  mzd_t *v=mzd_init(1,nc), *s=mzd_init(1,ns);
-
154  for(int w=2;w<=wmax;w++){ // cluster weight loop
-
155  if (debug&8)
-
156  printf( "# starting w=%d\n", w);
-
157  int beg=0, end=nc-w-1;
-
158  if (start>=0)
-
159  beg=end=start;
-
160  for(int i=beg;i<=end;i++){ // starting bit loop
-
161  if ((debug&8)&&(w==wmax))
-
162  printf( "# w=%d start=%d\n", w,i);
-
163  mzd_set_ui(v,0);
-
164  mzd_set_ui(s,0);
-
165  mzd_write_bit(v,0,i,1);
-
166  // printf("v="); mzd_print(v);
-
167  addto_inline(s,PT,i); // s+=col[i];
-
168  // printf("s="); mzd_print(s);
-
169  int done=start_rec(1,w,v,s,P,PT,G,rankG);
-
170  if(done==1){
-
171  // cout << "distance="<< weight(v)<< endl;
-
172  // cout << "prod="<< P.get_H()*v<< endl;
-
173  csr_free(PT);
-
174  mzd_free(v);
-
175  return w;
-
176  }
-
177  }
-
178  }
-
179  csr_free(PT);
-
180  mzd_free(v);
-
181  return -wmax; /* failed up to wmax */
-
182 }
-
183 
-
184 
-
186 mzp_t * do_skip_pivs(const size_t rank, const mzp_t * const pivs){
-
187  const rci_t n=pivs->length;
-
188  rci_t j1=rank;
-
189  mzp_t * ans = mzp_copy(NULL,pivs);
-
190  qsort(ans->values, rank, sizeof(pivs->values[0]), cmp_rci_t);
-
191 
-
192  for(rci_t j=0; j<n; j++){
-
193  if(!bsearch(&j, ans->values, rank, sizeof(ans->values[0]), cmp_rci_t)){
-
194  ans->values[j1++]=j;
-
195  }
-
196  }
-
197  assert(j1==n);
-
198 
-
199  int j=rank;
-
200  for(size_t i=0; j<n; j++)
-
201  ans->values[i++] = ans->values[j];
-
202  ans->length = n-rank;
-
203 
-
204  if(prm.debug & 8){
-
205  printf("skip_pivs of len=%d: ",ans->length);
-
206  for(int i=0; i< ans->length; i++)
-
207  printf(" %d%s",ans->values[i],i+1 == ans->length ?"\n":"");
-
208  printf("pivs of len=%d, rank=%zu: ",pivs->length, rank);
-
209  for(size_t i=0; i< rank; i++)
-
210  printf(" %d%s",pivs->values[i],i+1 == rank ?"\n":"");
-
211  }
-
212  return ans;
-
213 }
-
214 
-
215 
-
223 int do_RW_dist(const csr_t * const spaH0, const csr_t * const spaL0,
-
224  const int steps, const int wmin, const int classical, const int swait, const int debug){
-
226  const int nvar = spaH0->cols;
-
227  if(((!classical)&&(spaL0==NULL)) ||
-
228  ((classical)&&(spaL0!=NULL))){
-
229  printf("L0 %s NULL classical=%d\n",spaL0==NULL ? "=" : "!=", classical);
-
230  ERROR("L0 should be non-NULL present only for classical code!\n");
-
231  }
-
232 
-
233  int minW=nvar+1;
-
234 
-
235  if(debug&16)
-
236  printf("running do_RW_dist() with steps=%d wmin=%d classical=%d swait=%d nvar=%d\n",steps, wmin, classical, swait, nvar);
-
237 
-
238  mzd_t * mH = mzd_from_csr(NULL, spaH0);
-
239  // mzd_t *mLt = NULL, *eemLt = NULL; //, *mL = NULL;
-
240  rci_t *ee = malloc(nvar*sizeof(rci_t));
-
242  if((!mH) || (!ee))
-
243  ERROR("memory allocation failed!\n");
-
244  // if(p->debug & 16) mzd_print(mH);
-
247  mzp_t * perm=mzp_init(nvar);
-
248  mzp_t * pivs=mzp_init(nvar);
-
249  if((!pivs) || (!perm))
-
250  ERROR("memory allocation failed!\n");
-
251 
-
252  int iwait=0, ichanged=0;
-
253  for (int ii=0; ii< steps; ii++){
-
254  pivs=mzp_rand(pivs);
-
255  mzp_set_ui(perm,1);
-
256  perm=perm_p_trans(perm,pivs,0);
-
259  int rank=0;
-
260  for(int i=0; i< nvar; i++){
-
261  int col=perm->values[i];
-
262  int ret=gauss_one(mH, col, rank);
-
263  if(ret)
-
264  pivs->values[rank++]=col;
-
265  }
-
267  mzp_t * skip_pivs = do_skip_pivs(rank, pivs);
-
268 
-
275  int k = nvar - rank;
-
276  for (int ir=0; ir< k; ir++){
-
277  int cnt=0;
-
278  const int col = ee[cnt++] = skip_pivs->values[ir];
-
279  for(int ix=0; ix<rank; ix++){
-
280  if(mzd_read_bit(mH,ix,col))
-
281  ee[cnt++] = pivs->values[ix];
-
282  if (cnt >= minW)
-
283  break;
-
284  }
-
285  if(cnt < minW){
-
287  qsort(ee, cnt, sizeof(rci_t), cmp_rci_t);
-
288 
-
290  int nz;
-
291  if (classical)
-
292  nz=1;
-
293  else
-
294  nz = sparse_syndrome_non_zero(spaL0, cnt, ee);
-
295  if(nz){
-
299  // if (cnt < minW){
-
300  minW=cnt;
-
301  if (minW <= wmin){
-
302  minW = - minW;
-
303  }
-
304  }
-
305  }
-
306  }
-
307  if(debug & 16)
-
308  printf(" round=%d of %d minW=%d ichanged=%d iwait=%d\n",
-
309  ii+1, steps, minW, ichanged, iwait);
-
310 
-
311  mzp_free(skip_pivs);
-
312 
-
313  iwait = ichanged > 0 ? 0 : iwait+1 ;
-
314  ichanged=0;
-
315  if((swait > 0)&&(iwait > swait)){
-
316  if(debug & 16)
-
317  printf(" iwait=%d >swait=%d, terminating after %d steps\n", iwait, swait, ii+1);
-
318  break;
-
319  }
-
320 
-
321  }
-
323  //alldone: /** early termination label */
-
324 
-
326  mzp_free(perm);
-
327  mzp_free(pivs);
-
328  free(ee);
-
329 
-
330  mzd_free(mH);
-
331 
-
332  return minW;
-
333 }
-
334 
-
335 
-
336 int do_dist_rnd(csr_t *spaG0, mzd_t *matP0, int debug,int steps, int wmin){
-
337  rci_t n = spaG0-> cols;
-
338 
-
339  mzd_t *matG0perp=NULL;
-
340  int weimin=n, rt=0;
-
341  csr_t *spaG1=NULL;
-
342  mzp_t *piv1=mzp_init(n), *q1=mzp_init(n);
-
343 
-
344  for (int ii=0; ii < steps ; ii++){ /* random window decoding loop */
-
345  mzp_rand(piv1); /* random permutation in terms of pivots */
-
346  mzp_set_ui(q1,1); q1=perm_p_trans(q1,piv1,0); /* corresponding permutation */
-
347  spaG1=csr_apply_perm(spaG1,spaG0,q1); /* G with permuted cols */
-
348  matG0perp=mzd_generator_from_csr(matG0perp,spaG1);
-
349  mzd_apply_p_right(matG0perp, piv1); // permute cols back to order in G0
-
350  // mzd_info(matG0perp,0);
-
351  // generator matrix dual to G0
-
352  if((debug & 1)&&(ii==0))
-
353  printf("# rankG=%d\n",matG0perp->nrows);
-
354  if (((debug &512)||(debug & 2048)) ){
-
355  if (debug & 2048) printf("current G\n");
-
356  mzd_t *matG0=mzd_from_csr(NULL,spaG0);
-
357  if ((debug & 2048) && (n<=150)){
-
358  mzd_print(matG0);
-
359  printf("current Gperp=\n");
-
360  mzd_print(matG0perp);
-
361  printf(" G*Gperp_T=\n");
-
362  }
-
363  mzd_t *prod1;
-
364  // mzd_info(prod1,0);
-
365  // mzd_info(matG0,0);
-
366  // mzd_info(matG0perp,0);
-
367  mzd_t * matG0ptran = mzd_transpose(NULL,matG0perp);
-
368  prod1=csr_mzd_mul(NULL,spaG0,matG0ptran,1);
-
369  mzd_free(matG0ptran);
-
370  if ((debug & 2048) && (n<=150))
-
371  mzd_print(prod1);
-
372  printf("weigt of G0*G0perp_T=%d\n",(int)mzd_weight(prod1));
-
373  if ((debug & 2048) && (n<=150)){
-
374  printf("current P=\n");
-
375  mzd_print(matP0);
-
376  printf(" G*P_T=\n");
-
377  }
-
378  mzd_free(prod1);
-
379  prod1=mzd_init(matG0->nrows,matP0->nrows);
-
380  prod1=_mzd_mul_naive(prod1,matG0,matP0,1);
-
381  if ((debug & 2048) && (n<=150))
-
382  mzd_print(prod1);
-
383  printf("weigt of G*P^T=%d\n",(int)mzd_weight(prod1));
-
384  mzd_free(matG0);
-
385  mzd_free(prod1);
-
386  }
-
387 
-
388  rt = matG0perp->nrows;
-
389  if (rt==matP0->nrows)
-
390  ERROR("This is an empty code, distance infinite!");
-
391  rci_t wei;
-
392  mzd_t *row, *row0=NULL;
-
393  int width = matG0perp -> width;
-
394  if((debug & 1)&&(ii==0))
-
395  printf("# n=%d width=%d\n",n,width);
-
396  // printf("##### here ########\n");
-
397  for(int i=0;i<rt;i++){
-
398  fflush(stdout);
-
399  row=mzd_init_window(matG0perp,i,0,i+1,n);
-
400  wei=mzd_weight(row);
-
401  if((wei<weimin)){
-
402  // if (debug & 1024)
-
403  row0=mzd_copy(row0,row);
-
404  // else
-
405  // row0=row;
-
406  rci_t j=do_reduce(row,matP0,matP0->nrows);
-
407  if (j!=-1){ /* not an empty row */
-
408  if (debug & 1024){
-
409  int wei0=mzd_weight_naive(row0);
-
410  if (wei!=wei0){
-
411  printf("### naive wei=%d vs std_bitcount %d\n",wei0,wei);
-
412  mzd_print(row0);
-
413  }
-
414  }
-
415  if (debug & 3)
-
416  printf("# round=%d i=%d w=%d s=%d " "s0=%d "
-
417  "min=%d\n",ii,i,wei,
-
418  debug &2 ? syndrome_bit_count(row0,spaG0):0,syndrome_bit_count(row,spaG0),
-
419  weimin);
-
420  weimin=wei;
-
421  }
-
422  /* else
-
423  printf("# wei=%d, syndrome should be zero: s0=%d \n" , wei,
-
424  syndrome_bit_count(row0,spaG0));
-
425  */
-
426  }
-
427  mzd_free(row);
-
428  if(weimin<=wmin){ // no need to continue
-
429  ii=steps;
-
430  weimin=-weimin;
-
431  break;
-
432  }
-
433 
-
434  }
-
435  if (row0!=NULL){
-
436  mzd_free(row0);
-
437  row0=NULL;
-
438  }
-
439  if (weimin<0)
-
440  break;
-
441  }
-
442  if (debug & 2)
-
443  printf("# n=%d k=%d weimin=%d\n",n,rt-matP0->nrows,weimin);
-
444 
-
445  return weimin;
-
446 }
-
447 
-
448 #ifdef STANDALONE
-
449 
-
450 int main(int argc, char **argv){
-
451  params_t * const p = &prm;
-
452  var_init(argc,argv,p);
-
453 
-
454  if (p->method & 2){ /* cluster */
- -
456  if(prm.max_row_wgt_G>maxrow)
-
457  ERROR("main: increase maxrow=%d to %d",maxrow,prm.max_row_wgt_G);
-
458 
-
459  }
-
460  const int n=p->nvar;
-
461 
-
462 
-
463 
-
464  if (prm.method & 1){ /* RW method */
-
465 #if 0
- -
467 #else
- -
469 #endif /* 0 */
-
470  if (prm.debug & 1){
-
471  printf("### RW upper bound on the distance: %d\n",prm.dist_max);
-
472  if(prm.dist_max <0)
-
473  printf("### negative distance due to wmax=%d set (early termination)\n",prm.wmax);
-
474  }
-
475  prm.wmax=minint(prm.wmax, abs(prm.dist_max)-1);
-
476  }
-
477 
-
478  if (prm.method & 2){ /* cluster method */
-
479  /* convert G to standard form */
-
480  mzp_t *piv0=mzp_init(n); // mzp_out(piv0);
-
481  mzd_t *matG0=mzd_from_csr(NULL,p->spaG);
-
482  rci_t rankG=mzd_gauss_naive(matG0,piv0,1);
-
483  if(prm.debug & 1)
-
484  printf("# rankG=%d\n",rankG);
-
485  mzd_apply_p_right_trans(matG0,piv0);
-
486  // matG0->nrows=rankG;
-
487 
-
488  mzp_t *q0=perm_p_trans(NULL,piv0,1); // permutation equiv to piv0
-
489  csr_t *spaH0=csr_apply_perm(NULL,p->spaH,q0); // permuted sparse H
-
490  // csr_t* spaG0=csr_apply_perm(NULL,p->spaG,q0); // permuted sparse G -- not needed here
-
491  if(prm.debug & 2048){
-
492  if ((n<=150))
-
493  {
-
494  printf("matG0=\n");
-
495  mzd_print(matG0);
-
496 
-
497  printf("matP0=\n");
-
498  mzd_t *matP0=mzd_from_csr(NULL,spaH0);
-
499  mzd_print(matP0);
-
500  mzd_free(matP0);
-
501  }
-
502  int wei=product_weight_csr_mzd(spaH0, matG0,1);
-
503  printf("weigt of H0*G0_T=%d\n",wei);
-
504  if(wei>0)
-
505  ERROR("expected zero weight product!");
-
506  }
-
507  mzp_free(piv0);
-
508  mzp_free(q0);
-
509 
-
510  int dmin=do_dist_clus(spaH0,matG0,prm.debug,prm.wmax,prm.start,rankG);
-
511  csr_free(spaH0);
-
512  mzd_free(matG0);
-
513  if (dmin>0){
-
514  if (prm.debug & 1)
-
515  printf("### Cluster (actual min-weight codeword found): dmin=%d\n",dmin);
-
516  prm.dist_min = dmin; /* actual distance found */
-
517  prm.dist_max = dmin;
-
518  }
-
519  else if (dmin<0){
-
520  if (prm.debug & 1)
-
521  printf("### Cluster dmin=%d (no codewords of weight up to %d)\n",dmin,-dmin);
-
522  if (-dmin==abs(prm.dist_max)-1)
-
523  prm.dist_min=abs(prm.dist_max); /* OK */
-
524  else
-
525  prm.dist_min=-dmin;
-
526  }
-
527  else
-
528  ERROR("unexpected dmin=0\n");
-
529 
-
530  if (prm.dist_min==abs(prm.dist_max)){
-
531  if(p->method==3)
-
532  printf("success (two distance bounds coincide) d=%d\n",prm.dist_min);
-
533  else
-
534  printf("success (found min-weight codeword) d=%d\n",prm.dist_min);
-
535  }
-
536  else if (prm.dist_max>prm.dist_min)
-
537  printf("distance in the interval (inclusive) %d to %d\n", prm.dist_min,prm.dist_max);
-
538  else
-
539  printf("cluster algorithm failed to find a codeword up to wmax=%d\n",-dmin);
-
540  }
-
541  else{ /* just RW */
-
542  // if(prm.debug &1)
-
543  printf("RW algorithm upper bound for the distance d=%d\n", prm.dist_max);
-
544 
-
545  // }
-
546  }
-
547  var_kill(p);
-
548 
-
549  return 0;
-
550 }
-
551 
-
552 #endif /* STANDALONE */
-
int do_RW_dist(const csr_t *const spaH0, const csr_t *const spaL0, const int steps, const int wmin, const int classical, const int swait, const int debug)
Random Information Set search for small-E logical operators.
Definition: dist_m4ri.c:223
-
int do_dist_clus(const csr_t *const P, const mzd_t *const G, int debug, int wmax, int start, const int rankG)
lower bound on the minimum distance by cluster enumeration
Definition: dist_m4ri.c:141
-
mzp_t * do_skip_pivs(const size_t rank, const mzp_t *const pivs)
prepare an ordered pivot-skip list of length n-rank
Definition: dist_m4ri.c:186
-
int do_dist_rnd(csr_t *spaG0, mzd_t *matP0, int debug, int steps, int wmin)
Definition: dist_m4ri.c:336
+
17 // #include <m4ri/config.h>
+
18 #include <inttypes.h>
+
19 #include <strings.h>
+
20 #include <stdlib.h>
+
21 #include <time.h>
+
22 #include <m4ri/m4ri.h>
+
23 
+
24 #include "mmio.h"
+
25 #include "util_m4ri.h"
+
26 #include "util_io.h"
+
27 #include "dist_m4ri.h"
+
28 
+
30 mzp_t * do_skip_pivs(mzp_t * ans, const int rank, const mzp_t * const pivs){
+
31  const rci_t n=pivs->length;
+
32  int *arr=calloc(rank,sizeof(int));
+
33  if(!arr)
+
34  ERROR("memory allocation");
+
35  for(int i=0; i< rank; i++)
+
36  arr[i] = pivs->values[i];
+
37  qsort(arr, rank, sizeof(arr[0]), cmp_rci_t);
+
38  if((ans) && (ans->length != n-rank)){
+
39  mzp_free(ans);
+
40  ans=NULL;
+
41  }
+
42  if (!ans)
+
43  ans=mzp_init(n-rank);
+
44 
+
45  rci_t j=0, j1=0;
+
46  for(rci_t i1=0; i1 < rank; i1++, j++){
+
47  const rci_t piv=arr[i1];
+
48  while(j<piv)
+
49  ans->values[j1++] = j++;
+
50  }
+
51  while(j<n)
+
52  ans->values[j1++] = j++;
+
53  assert(j1==n-rank);
+
54 
+
55 #ifndef NDEBUG
+
56  if(prm.debug&1024){
+
57  printf("orig pivs of len=%d, rank=%d: ",pivs->length, rank);
+
58  for(int i=0; i< rank; i++)
+
59  printf(" %d%s",pivs->values[i],i+1 == rank ?"\n":"");
+
60  printf("srtd pivs of len=%d: ", rank);
+
61  for(int i=0; i< rank; i++)
+
62  printf(" %d%s",arr[i],i+1 == rank ?"\n":"");
+
63  printf("skip_pivs of len=%d: ",ans->length);
+
64  for(int i=0; i< ans->length; i++)
+
65  printf(" %d%s",ans->values[i],i+1 == ans->length ?"\n":"");
+
66  }
+
67 #endif
+
68  free(arr);
+
69  return ans;
+
70 }
+
71 
+
79 int do_RW_dist(const csr_t * const spaH0, const csr_t * const spaL0,
+
80  const int steps, const int wmin, const int classical, const int debug){
+
82  const int nvar = spaH0->cols;
+
83  if(((!classical)&&(spaL0==NULL)) ||
+
84  ((classical)&&(spaL0!=NULL))){
+
85  printf("L0 %s NULL classical=%d\n",spaL0==NULL ? "=" : "!=", classical);
+
86  ERROR("L0 should be non-NULL present only for classical code!\n");
+
87  }
+
88 
+
89  int minW=nvar+1;
+
90 
+
91  if(debug&2)
+
92  printf("# running do_RW_dist() with steps=%d wmin=%d classical=%d nvar=%d\n",
+
93  steps, wmin, classical, nvar);
+
94 
+
95  mzd_t * mH = mzd_from_csr(NULL, spaH0);
+
96  mzd_t *mHT = NULL;
+
97  mzp_t * skip_pivs = NULL;
+
99  rci_t *ee = malloc(nvar*sizeof(rci_t));
+
100 
+
101  if((!mH) || (!ee))
+
102  ERROR("memory allocation failed!\n");
+
103 
+
105  mzp_t * perm=mzp_init(nvar);
+
106  mzp_t * pivs=mzp_init(nvar);
+
107  if((!pivs) || (!perm))
+
108  ERROR("memory allocation failed!\n");
+
109 
+
110  for (int ii=0; ii< steps; ii++){
+
111  pivs=mzp_rand(pivs);
+
112  mzp_set_ui(perm,1);
+
113  perm=perm_p_trans(perm,pivs,0);
+
116  int rank=0;
+
117  for(int i=0; i< nvar; i++){
+
118  int col=perm->values[i];
+
119  int ret=gauss_one(mH, col, rank);
+
120  if(ret)
+
121  pivs->values[rank++]=col;
+
122  }
+
123 
+
125  skip_pivs = do_skip_pivs(skip_pivs, rank, pivs);
+
126 #ifndef NEW
+
127 # define NEW 1
+
128 #endif
+
129 #if (NEW==1)
+
131  mHT = mzd_transpose(mHT,mH);
+
132 #endif
+
139  int k = nvar - rank;
+
140  for (int ir=0; ir< k; ir++){
+
141  int cnt=0;
+
142  const int col = ee[cnt++] = skip_pivs->values[ir];
+
143 #if (NEW==0)
+
144  for(int ix=0; ix<rank; ix++){
+
145  if(mzd_read_bit(mH,ix,col))
+
146  ee[cnt++] = pivs->values[ix];
+
147  if (cnt >= minW)
+
148  break;
+
149  }
+
150 #elif (NEW==2)
+
156  rci_t ic=col, ix=0;
+
157  while (ix < rank){
+
158  int res = mzd_find_pivot(mH, ix, col, &ix, &ic);
+
159  if((res)&&(ic==col)){
+
160  ee[cnt++] = pivs->values[ix++];
+
161  // printf("cnt=%d j=%d\n",cnt,ix);
+
162  if (cnt >= minW)
+
163  break;
+
164  }
+
165  else
+
166  break;
+
167  }
+
168 #else
+
169  word * rawrow = mHT->rows[col];
+
170  rci_t j=-1;
+
171  while(cnt < minW){
+
172  j=nextelement(rawrow,mHT->width,j);
+
173  if(j==-1) // empty line after simplification
+
174  break;
+
175  // printf("cnt=%d j=%d\n",cnt,j);
+
176  ee[cnt++] = pivs->values[j++];
+
177  }
+
178 #endif /* NEW */
+
179  if (cnt < minW){
+
181  qsort(ee, cnt, sizeof(rci_t), cmp_rci_t);
+
182 #ifndef NDEBUG
+
184  if(sparse_syndrome_non_zero(spaH0, cnt, ee)){
+
185  printf("# cw of weight %d: [",cnt);
+
186  for(int i=0; i<cnt;i++)
+
187  printf("%d%s",ee[i],i+1==cnt?" ":"]\n");
+
188  ERROR("this should not happen: cw not orthogonal to H");
+
189  }
+
190 #endif /* NDEBUG */
+
191 
+
193  int nz;
+
194  if (classical)
+
195  nz=1;
+
196  else
+
197  nz = sparse_syndrome_non_zero(spaL0, cnt, ee);
+
198  if(nz){
+
201  if(debug&16){
+
202  printf("# step=%d row=%d minW=%d found cw of W=%d: [",ii,ir,minW,cnt);
+
203  int max = ((cnt<25) || (debug&2048)) ? cnt : 25 ;
+
204  for(int i=0; i< max; i++)
+
205  printf("%d%s", ee[i], i+1!=max?" ": (cnt==max ? "]\n" : "...]\n"));
+
206  }
+
207  minW=cnt;
+
208  if (minW <= wmin){
+
209  minW = - minW;
+
210  }
+
211  }
+
212  }
+
213  }
+
214  if(debug&8){
+
215  if(ii%1000==999)
+
216  printf("# round=%d of %d minW=%d\n", ii+1, steps, minW);
+
217  }
+
218 
+
219  }
+
221  //alldone: /** early termination label */
+
222 
+
224  if(skip_pivs)
+
225  mzp_free(skip_pivs);
+
226  mzp_free(perm);
+
227  mzp_free(pivs);
+
228  free(ee);
+
229  if(mHT)
+
230  mzd_free(mHT);
+
231  mzd_free(mH);
+
232 
+
233  return minW;
+
234 }
+
235 
+
236 
+
237 #ifdef STANDALONE
+
238 
+
239 int do_CC_dist(const csr_t * const mH, const csr_t * mL,
+
240  const int wmax, const int start, const int debug);
+
241 
+
242 int main(int argc, char **argv){
+
243  params_t * const p = &prm;
+
244 
+
245  var_init(argc,argv,p);
+
246 
+
247  // const int n=p->nvar;
+
248 
+
249  if (prm.method & 1){ /* RW method */
+
250 
+ +
252 
+
253  if (prm.debug&1){
+
254  printf("### RW upper bound on the distance: %d\n",prm.dist_max);
+
255  if(prm.dist_max <0)
+
256  printf("### negative distance due to wmax=%d set (early termination)\n",prm.wmax);
+
257  }
+
258  prm.wmax=minint(prm.wmax, abs(prm.dist_max)-1);
+
259  }
+
260 
+
261  if (prm.method & 2){ /* cluster method */
+
262  int dmin=do_CC_dist(p->spaH,p->spaL,p->wmax,p->start,p->debug);
+
263 
+
264  if (dmin>0){
+
265  if (prm.debug&1)
+
266  printf("### Cluster (actual min-weight codeword found): dmin=%d\n",dmin);
+
267  prm.dist_min = dmin; /* actual distance found */
+
268  prm.dist_max = dmin;
+
269  }
+
270  else if (dmin<0){
+
271  if (prm.debug&1)
+
272  printf("### Cluster dmin=%d (no codewords of weight up to %d)\n",dmin,-dmin);
+
273  if (-dmin==abs(prm.dist_max)-1)
+
274  prm.dist_min=abs(prm.dist_max); /* OK */
+
275  else
+
276  prm.dist_min=-dmin;
+
277  }
+
278  else
+
279  ERROR("unexpected dmin=0\n");
+
280 
+
281  if (prm.dist_min==abs(prm.dist_max)){
+
282  if(p->method==3)
+
283  printf("success (two distance bounds coincide) d=%d\n",prm.dist_min);
+
284  else
+
285  printf("success (found min-weight codeword) d=%d\n",prm.dist_min);
+
286  }
+
287  else if (prm.dist_max>prm.dist_min)
+
288  printf("distance in the interval (inclusive) %d to %d\n", prm.dist_min,prm.dist_max);
+
289  else
+
290  printf("cluster algorithm failed to find a codeword up to wmax=%d\n",-dmin);
+
291  }
+
292  else{ /* just RW */
+
293  printf("RW algorithm upper bound for the distance d=%d\n", prm.dist_max);
+
294  }
+
295  var_kill(p);
+
296 
+
297  return 0;
+
298 }
+
299 
+
300 #endif /* STANDALONE */
+
int do_CC_dist(const csr_t *const mH, const csr_t *mL, const int wmax, const int start, const int debug)
Definition: dist_cc.c:256
+
mzp_t * do_skip_pivs(mzp_t *ans, const int rank, const mzp_t *const pivs)
distance of a classical or quantum CSS code
Definition: dist_m4ri.c:30
+
int do_RW_dist(const csr_t *const spaH0, const csr_t *const spaL0, const int steps, const int wmin, const int classical, const int debug)
Random Information Set search for small-E logical operators.
Definition: dist_m4ri.c:79
-
#define _maybe_unused
Definition: mmio.c:17
- -
int * p
Definition: util_m4ri.h:110
-
int * i
Definition: util_m4ri.h:111
-
int rows
Definition: util_m4ri.h:106
-
int cols
Definition: util_m4ri.h:107
- -
int dist_max
Definition: util_io.h:33
-
csr_t * spaL
Definition: util_io.h:50
-
csr_t * spaG
Definition: util_io.h:49
-
csr_t * spaH
Definition: util_io.h:48
-
int dist_min
Definition: util_io.h:34
-
int nvar
Definition: util_io.h:40
-
int max_row_wgt_G
Definition: util_io.h:35
-
int method
Definition: util_io.h:27
-
int wmax
Definition: util_io.h:29
-
int wmin
Definition: util_io.h:30
-
int swait
Definition: util_io.h:42
-
int steps
Definition: util_io.h:28
-
int start
int maxrow; /* WARNING: this is defined in dist_m4ri.h as static const int *‍/
Definition: util_io.h:37
-
int debug
Definition: util_io.h:24
-
int classical
Definition: util_io.h:25
-
void var_kill(params_t *const p)
Definition: util_io.c:217
-
params_t *const p
Definition: util_io.c:31
+ +
int cols
Definition: util_m4ri.h:104
+ +
int dist_max
Definition: util_io.h:35
+
csr_t * spaL
Definition: util_io.h:52
+
csr_t * spaH
Definition: util_io.h:50
+
int dist_min
Definition: util_io.h:36
+
int method
Definition: util_io.h:29
+
int wmax
Definition: util_io.h:31
+
int wmin
Definition: util_io.h:32
+
int steps
Definition: util_io.h:30
+
int start
int max_row_wt; /* WARNING: this is defined in util_io.h as static const int *‍/
Definition: util_io.h:40
+
int debug
Definition: util_io.h:26
+
int classical
Definition: util_io.h:27
+
void var_kill(params_t *const p)
Definition: util_io.c:257
+
params_t *const p
Definition: util_io.c:32
params_t prm
Definition: util_io.c:4
-
void var_init(int argc, char **argv, params_t *const p)
Definition: util_io.c:33
+
void var_init(int argc, char **argv, params_t *const p)
Definition: util_io.c:34
-
csr_t * csr_free(csr_t *p)
Definition: util_m4ri.c:422
-
mzd_t * csr_mzd_mul(mzd_t *C, const csr_t *S, const mzd_t *B, int clear)
Definition: util_m4ri.c:226
-
int syndrome_bit_count(const mzd_t *const row, const csr_t *const spaQ)
Definition: util_m4ri.c:296
-
mzp_t * perm_p_trans(mzp_t *q, const mzp_t *p, const rci_t start)
Definition: util_m4ri.c:405
-
csr_t * csr_apply_perm(csr_t *dst, const csr_t *const src, const mzp_t *const perm)
Definition: util_m4ri.c:583
+
mzp_t * perm_p_trans(mzp_t *q, const mzp_t *p, const rci_t start)
Definition: util_m4ri.c:444
mzd_t * mzd_from_csr(mzd_t *dst, const csr_t *p)
Definition: util_m4ri.c:158
-
mzd_t * mzd_generator_from_csr(mzd_t *G, const csr_t *const H)
Definition: util_m4ri.c:185
-
int csr_max_row_wght(const csr_t *p)
Definition: util_m4ri.c:117
-
int do_reduce(mzd_t *row, const mzd_t *matG0, const rci_t rankG0)
Definition: util_m4ri.c:617
-
rci_t mzd_gauss_naive(mzd_t *M, mzp_t *q, int full)
Definition: util_m4ri.c:76
-
csr_t * csr_transpose(csr_t *dst, const csr_t *p)
Definition: util_m4ri.c:135
-
size_t mzd_weight(const mzd_t *A)
Definition: util_m4ri.c:46
-
size_t mzd_weight_naive(const mzd_t *A)
Definition: util_m4ri.c:16
-
size_t product_weight_csr_mzd(const csr_t *A, const mzd_t *B, int transpose)
Definition: util_m4ri.c:319
-
#define ERROR(fmt,...)
Definition: util_m4ri.h:16
+
#define ERROR(fmt,...)
Definition: util_m4ri.h:13