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);
49 *k = n; // we will get it to be k=n-d by calling shrinkstar d times below
51 randomX = calloc(n*d,sizeof(float));
52 randomXcopy = calloc(n*d,sizeof(float));
53 randomXsingvals = calloc(d,sizeof(float));
54 superb = calloc((d),sizeof(float));
55 U = calloc((d*d),sizeof(float));
56 VT = calloc((n*n),sizeof(float));
61 while(numsingvals != d) {
64 randomX[i*n+j] = (float)(rand()%2);
65 randomXcopy[i*n+j] = randomX[i*n+j];
66 //printf("%1.0f ", randomX[i*n+j]);
71 info = LAPACKE_sgesvd( LAPACK_ROW_MAJOR, 'N', 'N', d, n, randomXcopy, n, randomXsingvals, U, n, VT, n, superb );
76 if(fabs(randomXsingvals[i]) >= ZEROTHRESHOLD)
80 //printf("Number of singular values: %d\n", numsingvals);
83 *G = calloc(n, sizeof(int*)); *GBar = calloc(n, sizeof(int*));
84 *h = calloc(n, sizeof(int));
87 (*G)[i] = calloc(n, sizeof(int));
88 (*GBar)[i] = calloc(n, sizeof(int));
94 int* xi = calloc(n, sizeof(int));
97 xi[j] = (int)randomX[i*n+j];
98 //printf("%1.0f ", randomX[i*n+j]);
101 //printVector(xi, n);
102 shrinkstar(n, k, *h, *G, *GBar, xi, 0);
105 //printf("n-d=%d\n", n-d);
106 //printf("k=%d\n", *k);
107 //printMatrix(*G, n, n);
112 *Q = (rand()%8); // Q \in Z/8Z
114 *D = calloc(*k, sizeof(int));
116 (*D)[i] = (rand()%4)*2; // D_a \in {0, 2, 4, 6}
118 *J = calloc(*k, sizeof(int*));
119 for(i=0; i<*k; i++) {
120 (*J)[i] = calloc(*k, sizeof(int));
121 for(j=(i+1); j<*k; j++)
122 (*J)[i][j] = (rand()%2)*4; // J_{a,b} \in {0, 4} for a!=b
123 (*J)[i][i] = (2*(*D)[i])%8;
125 for(i=0; i<*k; i++) {
126 for(j=(i+1); j<*k; j++) {
127 (*J)[j][i] = (*J)[i][j];
131 free(randomX); free(randomXcopy); free(randomXsingvals); free(superb); free(U); free(VT);