9 #include "shrinkstar.h"
11 #include "randomstabilizerstate.h"
13 #define ZEROTHRESHOLD (0.00000001)
15 // Note: Make sure you seed the random number generator before calling this!
16 // i.e. srand((unsigned)time(NULL));
17 void randomstabilizerstate(int n, int* k, int** h, int*** G, int*** GBar, int* Q, int** D, int*** J, double** Pd)
19 // vector and matrix pointers should all be passed unallocated! (no freeing of memory performed)
21 float *randomX, *randomXcopy; // d*n matrix that is in array-form for LAPACKE
22 float *randomXsingvals; // singular values of randomX
23 float *U, *VT, *superb;
28 //srand((unsigned)time(NULL));
30 double randPd = (double)rand()/(double)(RAND_MAX);
34 cumprob += Pd[n-1][i]; // cumprob goes up to 1.0 once you count all Pd[n][0:n+1]
35 if(cumprob > randPd) {
42 //printf("d=%d total=%f randPd=%f\n", d, total, randPd);
43 //printf("!!!d=%d\n", d);
45 *k = n; // we will get it to be k=n-d by caling shrinkstar d times below
47 randomX = calloc(n*d,sizeof(float));
48 randomXcopy = calloc(n*d,sizeof(float));
49 randomXsingvals = calloc(d,sizeof(float));
50 superb = calloc((d),sizeof(float));
51 U = calloc((d*d),sizeof(float));
52 VT = calloc((n*n),sizeof(float));
57 while(numsingvals != d) {
60 randomX[i*n+j] = (float)(rand()%2);
61 randomXcopy[i*n+j] = randomX[i*n+j];
62 //printf("%1.0f ", randomX[i*n+j]);
67 info = LAPACKE_sgesvd( LAPACK_ROW_MAJOR, 'N', 'N', d, n, randomXcopy, n, randomXsingvals, U, n, VT, n, superb );
72 if(fabs(randomXsingvals[i]) >= ZEROTHRESHOLD)
76 //printf("Number of singular values: %d\n", numsingvals);
79 *G = calloc(n, sizeof(int*)); *GBar = calloc(n, sizeof(int*));
81 *h = calloc(n, sizeof(int));
84 (*G)[i] = calloc(n, sizeof(int));
85 (*GBar)[i] = calloc(n, sizeof(int));
91 int* xi = calloc(n, sizeof(int));
94 xi[j] = (int)randomX[i*n+j];
95 //printf("%1.0f ", randomX[i*n+j]);
99 shrinkstar(n, k, *h, *G, *GBar, xi, 0);
102 //printf("n-d=%d\n", n-d);
103 //printf("k=%d\n", *k);
104 //printMatrix(*G, n, n);
109 *Q = (rand()%8); // Q \in Z/8Z
111 *D = calloc(*k, sizeof(int));
113 (*D)[i] = (rand()%4)*2; // D_a \in {0, 2, 4, 6}
115 *J = calloc(*k, sizeof(int*));
116 for(i=0; i<*k; i++) {
117 (*J)[i] = calloc(*k, sizeof(int));
118 for(j=(i+1); j<*k; j++)
119 (*J)[i][j] = (rand()%2)*4; // J_{a,b} \in {0, 4} for a!=b
120 (*J)[i][i] = (2*(*D)[i])%8;
122 for(i=0; i<*k; i++) {
123 for(j=(i+1); j<*k; j++) {
124 (*J)[j][i] = (*J)[i][j];
128 free(randomX); free(randomXcopy); free(randomXsingvals); free(superb); free(U); free(VT);