// // 0.0 ... r ... 1.0 // |<------->|<-------->|<------->|<-------->| // p00 p01 p10 p11 // int Mk( Complex a[]) { // measure overall state double r = rand() / (RAND_MAX + 1.0); // 0.0 <= r < 1.0 double p = 0; // cumulative probability int k = 3; for( int i = 0; i < 4; ++i) { p += abs2( a[i]); if( r < p) { k = i; break; } } printf( "k = %i\n", k); return k; } void M10( Complex a[]) { // measure qubits 1 and 0 int k = Mk( a); a[0]=a[1]=a[2]=a[3]=0; a[k] = 1; }