From 8cd5a3a3a99cc6a6ee6af228a717c1c3f2542c73 Mon Sep 17 00:00:00 2001 From: server2 Date: Sun, 8 Dec 2024 16:37:51 -0500 Subject: [PATCH] Sun Dec 8 04:37:51 PM EST 2024 --- src/bytecode.c.measurement.backup | 1751 +++++++++++++++++++++++++++++ 1 file changed, 1751 insertions(+) create mode 100644 src/bytecode.c.measurement.backup diff --git a/src/bytecode.c.measurement.backup b/src/bytecode.c.measurement.backup new file mode 100644 index 0000000..8212add --- /dev/null +++ b/src/bytecode.c.measurement.backup @@ -0,0 +1,1751 @@ +#include "bytecode.h" + +const char* qansel_instruction_to_string(unsigned char instr) +{ + switch (instr) + { + case QANSEL_INSTRUCTION_X: return "x"; + case QANSEL_INSTRUCTION_Y: return "y"; + case QANSEL_INSTRUCTION_Z: return "z"; + case QANSEL_INSTRUCTION_H: return "h"; + case QANSEL_INSTRUCTION_S: return "s"; + case QANSEL_INSTRUCTION_T: return "t"; + case QANSEL_INSTRUCTION_SDG: return "sdg"; + case QANSEL_INSTRUCTION_TDG: return "tdg"; + case QANSEL_INSTRUCTION_RX: return "rx"; + case QANSEL_INSTRUCTION_RY: return "ry"; + case QANSEL_INSTRUCTION_RZ: return "rz"; + case QANSEL_INSTRUCTION_U1: return "u1"; + case QANSEL_INSTRUCTION_U2: return "u2"; + case QANSEL_INSTRUCTION_U3: return "u3"; + case QANSEL_INSTRUCTION_CX: return "cx"; + case QANSEL_INSTRUCTION_SWAP: return "swap"; + case QANSEL_INSTRUCTION_CCX: return "ccx"; + case QANSEL_INSTRUCTION_CSWAP: return "cswap"; + case QANSEL_INSTRUCTION_MEASURE: return "measure"; + case QANSEL_INSTRUCTION_DENSITY: return "density"; + case QANSEL_INSTRUCTION_BORN: return "born"; + case QANSEL_INSTRUCTION_BSAMPLE: return "bsample"; + case QANSEL_INSTRUCTION_IF_E: return "if=="; + case QANSEL_INSTRUCTION_IF_NE: return "if!="; + case QANSEL_INSTRUCTION_IF_G: return "if>"; + case QANSEL_INSTRUCTION_IF_GE: return "if>="; + case QANSEL_INSTRUCTION_IF_L: return "if<"; + case QANSEL_INSTRUCTION_IF_LE: return "if<="; + case QANSEL_INSTRUCTION_PRINT: return "print"; + case QANSEL_INSTRUCTION_SET: return "set"; + case QANSEL_INSTRUCTION_RESET: return "reset"; + case QANSEL_INSTRUCTION_BARRIER: return "barrier"; + case QANSEL_INSTRUCTION_EXIT: return "exit"; + } + return "Unknown"; +} + +float qansel_rand_s(float s) +{ + unsigned int tmp; + memcpy(&tmp, &s, sizeof(unsigned int)); + srand(tmp); +} +float qansel_rand_h() +{ + return ((float)rand()) / ((float)RAND_MAX); +} +float qansel_rand_t(QAnselContext* ctx) +{ + if (ctx->hardware_rng != QANSEL_HARDWARE_NONE) + { + unsigned int num = 0; + for (unsigned char i = 0; i < 4; i++) + { + unsigned char r = 0; + switch (ctx->hardware_rng) + { + case QANSEL_HARDWARE_TRUERNG: r = fgetc(ctx->random_file); break; + case QANSEL_HARDWARE_QUANTIS: + if (ctx->pointer >= QANSEL_QUANTIS_CHUNK_SIZE) + { + qansel_quantis_chunk(&(ctx->chunk)); + ctx->pointer = 0; + } + r = ctx->chunk[ctx->pointer]; + ctx->pointer += 1; + break; + case QANSEL_HARDWARE_RDSEED: r = qansel_hardware_rand(); break; + } + num = (num << 8) | r; + } + return ((float)num) / ((float)UINT32_MAX); + } + else + { + return qansel_rand_h(); + } +} + +void qansel_cnot(cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char bitA, unsigned char bitB) +{ + if (bitA >= qubitCount || bitB >= qubitCount) return; + unsigned int retLen = (unsigned int)pow(2, qubitCount); + cpx_mtx_t ret; + cpx_mtx_init(&ret, 1, retLen); + cpx_t n; + for (unsigned int i = 0; i < retLen; i++) + { + unsigned char bitAVal = (i >> bitA) & 1; + unsigned char bitBVal = (i >> bitB) & 1; + unsigned char bitBNew = bitAVal ? !bitBVal : bitBVal; + unsigned int j = (i & ~(1 << bitB)) | (bitBNew << bitB); + cpx_mtx_get(stateVector, 0, i, &n); + cpx_mtx_set(&ret, 0, j, &n); + } + cpx_mtx_free(stateVector); + stateVector->ptr = ret.ptr; + stateVector->rows = ret.rows; + stateVector->cols = ret.cols; +} + +void qansel_swap(cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char bitA, unsigned char bitB) +{ + if (bitA >= qubitCount || bitB >= qubitCount) return; + unsigned int retLen = (unsigned int)pow(2, qubitCount); + cpx_mtx_t ret; + cpx_mtx_init(&ret, 1, retLen); + cpx_t n; + for (unsigned int i = 0; i < retLen; i++) + { + unsigned char bitAVal = (i >> bitA) & 1; + unsigned char bitBVal = (i >> bitB) & 1; + unsigned char bitANew = bitBVal; + unsigned char bitBNew = bitAVal; + unsigned int j = (i & ~((1 << bitA) | (1 << bitB))) | ((bitANew << bitA) | (bitBNew << bitB)); + cpx_mtx_get(stateVector, 0, i, &n); + cpx_mtx_set(&ret, 0, j, &n); + } + cpx_mtx_free(stateVector); + stateVector->ptr = ret.ptr; + stateVector->rows = ret.rows; + stateVector->cols = ret.cols; +} + +void qansel_fredkin(cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char bitA, unsigned char bitB, unsigned char bitC) +{ + if (bitA >= qubitCount || bitB >= qubitCount) return; + unsigned int retLen = (unsigned int)pow(2, qubitCount); + cpx_mtx_t ret; + cpx_mtx_init(&ret, 1, retLen); + cpx_t n; + for (unsigned int i = 0; i < retLen; i++) + { + unsigned char bitAVal = (i >> bitA) & 1; + unsigned char bitBVal = (i >> bitB) & 1; + unsigned char bitCVal = (i >> bitC) & 1; + unsigned char bitBNew = bitAVal ? bitCVal : bitBVal; + unsigned char bitCNew = bitAVal ? bitBVal : bitCVal; + unsigned int j = (i & ~((1 << bitB) | (1 << bitC))) | ((bitBNew << bitB) | (bitCNew << bitC)); + cpx_mtx_get(stateVector, 0, i, &n); + cpx_mtx_set(&ret, 0, j, &n); + } + cpx_mtx_free(stateVector); + stateVector->ptr = ret.ptr; + stateVector->rows = ret.rows; + stateVector->cols = ret.cols; +} + +void qansel_toffoli(cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char bitA, unsigned char bitB, unsigned char bitC) +{ + if (bitA >= qubitCount || bitB >= qubitCount) return; + unsigned int retLen = (unsigned int)pow(2, qubitCount); + cpx_mtx_t ret; + cpx_mtx_init(&ret, 1, retLen); + cpx_t n; + for (unsigned int i = 0; i < retLen; i++) + { + unsigned char bitAVal = (i >> bitA) & 1; + unsigned char bitBVal = (i >> bitB) & 1; + unsigned char bitCVal = (i >> bitC) & 1; + unsigned char bitCNew = (bitAVal && bitBVal) ? !bitCVal : bitCVal; + unsigned int j = (i & ~(1 << bitC)) | (bitCNew << bitC); + cpx_mtx_get(stateVector, 0, i, &n); + cpx_mtx_set(&ret, 0, j, &n); + } + cpx_mtx_free(stateVector); + stateVector->ptr = ret.ptr; + stateVector->rows = ret.rows; + stateVector->cols = ret.cols; +} + +float* qansel_unitary(float theta, float phi, float lambda) +{ + cpx_mtx_t m; + cpx_t a, b, c, d; + a.real = cos(theta/2.0); + a.imaginary = 0; + b.real = -cos(lambda) * sin(theta/2.0); + b.imaginary = sin(lambda) * sin(theta/2.0); + c.real = cos(phi) * sin(theta/2.0); + c.imaginary = sin(phi) * sin(theta/2.0); + d.real = cos(phi + lambda) * cos(theta/2.0); + d.imaginary = sin(phi + lambda) * cos(theta/2.0); + //remember it must be transposed!!! + cpx_mtx_init(&m, 2, 2); + cpx_mtx_set(&m, 0, 0, &a); + cpx_mtx_set(&m, 0, 1, &c); + cpx_mtx_set(&m, 1, 0, &b); + cpx_mtx_set(&m, 1, 1, &d); + return m.ptr; +} + +void qansel_queue_init +( + cpx_mtx_t* queueVector, + int qubitCount, + int dofree +) +{ + for (int i = 0; i < qubitCount; i++) + { + if (dofree) free(queueVector[i].ptr); + int vectLen = queueVector[i].rows * (queueVector[i].cols * 2) + queueVector[i].cols * 2; + queueVector[i].rows = 2; + queueVector[i].cols = 2; + queueVector[i].ptr = malloc(vectLen); + memcpy(queueVector[i].ptr, Identity, vectLen); + } +} + +void qansel_instruction_queue +( + cpx_mtx_t** queueVector, + int qubitCount, + unsigned char instr, + unsigned char index, + float arg0, + float arg1, + float arg2 +) +{ + cpx_mtx_t gate, tmp; + gate.rows = 2; + gate.cols = 2; + float* gate_ptr; + int needToFreeGate = 0; + switch (instr) + { + case QANSEL_INSTRUCTION_H: gate_ptr = Hadamard; break; + case QANSEL_INSTRUCTION_X: gate_ptr = PauliX; break; + case QANSEL_INSTRUCTION_Y: gate_ptr = PauliY; break; + case QANSEL_INSTRUCTION_Z: gate_ptr = PauliZ; break; + case QANSEL_INSTRUCTION_S: gate_ptr = PhaseS; break; + case QANSEL_INSTRUCTION_T: gate_ptr = PhaseT; break; + case QANSEL_INSTRUCTION_SDG: gate_ptr = PhaseSdg; break; + case QANSEL_INSTRUCTION_TDG: gate_ptr = PhaseTdg; break; + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + gate_ptr = qansel_unitary(arg0, arg1, arg2); + needToFreeGate = 1; + break; + default: gate_ptr = Identity; break; + } + unsigned char qubit = qubitCount - (index) - 1; + gate.rows = 2; + gate.cols = 2; + gate.ptr = gate_ptr; + tmp.rows = 2; + tmp.cols = 2; + tmp.ptr = malloc((tmp.rows * 2) * (tmp.cols * 2) * sizeof(float)); + cpx_mtx_dot(tmp.ptr, (*queueVector)[qubit].ptr, gate_ptr, (*queueVector)[qubit].rows, (*queueVector)[qubit].cols, 2, 2); + free((*queueVector)[qubit].ptr); + (*queueVector)[qubit].rows = tmp.rows; + (*queueVector)[qubit].cols = tmp.cols; + (*queueVector)[qubit].ptr = tmp.ptr; + if (needToFreeGate) free(gate_ptr); +} + +void qansel_instruction +( + QAnselContext* ctx, + cpx_mtx_t* stateVector, + int qubitCount, + unsigned char instr, + unsigned char index, + float arg0, + float arg1, + float arg2, + //pass in this to flush the queue + cpx_mtx_t** queueVector +) +{ + cpx_mtx_t tmp; + cpx_mtx_t gate; + gate.rows = 2; + gate.cols = 2; + float* gate_ptr; + int needToFreeGate = 0; + switch (instr) + { + case QANSEL_INSTRUCTION_H: gate_ptr = Hadamard; break; + case QANSEL_INSTRUCTION_X: gate_ptr = PauliX; break; + case QANSEL_INSTRUCTION_Y: gate_ptr = PauliY; break; + case QANSEL_INSTRUCTION_Z: gate_ptr = PauliZ; break; + case QANSEL_INSTRUCTION_S: gate_ptr = PhaseS; break; + case QANSEL_INSTRUCTION_T: gate_ptr = PhaseT; break; + case QANSEL_INSTRUCTION_SDG: gate_ptr = PhaseSdg; break; + case QANSEL_INSTRUCTION_TDG: gate_ptr = PhaseTdg; break; + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + gate_ptr = qansel_unitary(arg0, arg1, arg2); + needToFreeGate = 1; + break; + default: gate_ptr = Identity; break; + } + + cpx_t n; + cpx_mtx_t filter; + cpx_mtx_init(&filter, 2, 2); + unsigned char qubit = qubitCount - (queueVector == NULL ? index : 0) - 1; + if (qubit == 0 && queueVector == NULL) + { + memcpy(filter.ptr, gate_ptr, 8 * sizeof(float)); + } + else if (queueVector == NULL) + { + memcpy(filter.ptr, Identity, 8 * sizeof(float)); + } + else + { + memcpy(filter.ptr, (*queueVector)[0].ptr, 8 * sizeof(float)); + } + + for (unsigned char i = 1; i < qubitCount; i++) + { + if (index != QANSEL_ALL_QUANTUM && queueVector == NULL) + { + if (qubit == i) + { + gate.ptr = gate_ptr; + } + else + { + gate.ptr = Identity; + } + } + else if (queueVector == NULL) + { + gate.ptr = gate_ptr; + } + else + { + gate.ptr = (*queueVector)[i].ptr; + } + + tmp.rows = filter.rows * gate.rows; + tmp.cols = filter.cols * gate.cols; + tmp.ptr = malloc((size_t)tmp.rows * ((size_t)tmp.cols * (size_t)2) * (size_t)sizeof(float)); + if (tmp.ptr == NULL) + { + fprintf(stderr, "QAnsel: Ran out of memory.\n"); + exit(1); + } + #ifdef SPEED_TEST + printf("(%ix%i);(%ix%i) (knk)\n", tmp.rows, tmp.cols, gate.rows, gate.cols); + unsigned long int us1, us2; + us1 = get_time(); + cpx_mtx_knk_metal(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tMetal: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_knk_metal_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tMetal2x2: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_knk_threads(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tThreads: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_knk_threads_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tThreads2x2: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_knk(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tBare: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_knk_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + us2 = get_time(); + printf("\tBare2x2: %lu\n", us2 - us1); + + //us1 = get_time(); + //cpx_mtx_knk(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + //us2 = get_time(); + //printf("\tTranspose: %lu\n", us2 - us1); + #else + if (((ctx->optimization_level) & QANSEL_MODE_METAL) && !((ctx->optimization_level) & QANSEL_MODE_THREADED) && tmp.cols >= 64) + { + cpx_mtx_knk_metal_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + } + else if (((ctx->optimization_level) & QANSEL_MODE_THREADED) && tmp.cols >= 64) + { + cpx_mtx_knk_threads_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + } + else + { + cpx_mtx_knk_2x2(tmp.ptr, filter.ptr, gate.ptr, filter.rows, filter.cols, gate.rows, gate.cols); + } + #endif + + free(filter.ptr); + filter.ptr = tmp.ptr; + filter.rows = tmp.rows; + filter.cols = tmp.cols; + } + + cpx_mtx_init(&tmp, stateVector->rows, stateVector->cols); + + #ifdef SPEED_TEST + printf("%ix%i (dot)\n", tmp.rows, tmp.cols); + unsigned long int us1, us2; + us1 = get_time(); + cpx_mtx_dot_metal(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + us2 = get_time(); + printf("\tMetal: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_dot_threads(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + us2 = get_time(); + printf("\tThreads: %lu\n", us2 - us1); + us1 = get_time(); + cpx_mtx_dot(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + us2 = get_time(); + printf("\tBare: %lu\n", us2 - us1); + #else + if (((ctx->optimization_level) & QANSEL_MODE_METAL) && tmp.cols >= 64) + { + cpx_mtx_dot_metal(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + } + else if (((ctx->optimization_level) & QANSEL_MODE_THREADED) && tmp.cols >= 64) + { + cpx_mtx_dot_threads(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + } + else + { + cpx_mtx_dot(tmp.ptr, stateVector->ptr, filter.ptr, stateVector->rows, stateVector->cols, filter.rows, filter.cols); + } + #endif + free(stateVector->ptr); + stateVector->ptr = tmp.ptr; + free(filter.ptr); + if (needToFreeGate) free(gate_ptr); +} + +unsigned char qansel_measure(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char qubit) +{ + unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount); + cpx_t n; + float prob0 = 0; + for (unsigned int i = 0; i < qubitCountPow2; i++) + { + unsigned char bit = (i >> qubit) & 1; + cpx_mtx_get(stateVector, 0, i, &n); + if (bit == 0) prob0 += cpx_magsqr(&n); + } + float r = (ctx->hidden_variable_set) ? qansel_rand_h() : qansel_rand_t(ctx); + //add any error skew if desired + if (ctx->noise > 0) + { + float probtot; + float prob1 = 1 - prob0; + prob0 += ctx->noise; + prob1 += ctx->noise; + probtot = prob0 + prob1; + prob0 /= probtot; + prob1 /= probtot; + } + unsigned char newBit = r < prob0 ? 0 : 1; + float probTot = 0; + for (unsigned int i = 0; i < qubitCountPow2; i++) + { + unsigned char bit = (i >> qubit) & 1; + cpx_mtx_get(stateVector, 0, i, &n); + if (bit != newBit) + { + n.real = 0; + n.imaginary = 0; + } + else + { + probTot += cpx_magsqr(&n); + } + cpx_mtx_set(stateVector, 0, i, &n); + } + + float multiplier = sqrt(1 / probTot); + for (unsigned int i = 0; i < qubitCountPow2; i++) + { + unsigned char bit = (i >> qubit) & 1; + cpx_mtx_get(stateVector, 0, i, &n); + if (bit == newBit) + { + n.real *= multiplier; + n.imaginary *= multiplier; + } + cpx_mtx_set(stateVector, 0, i, &n); + } + + return newBit; +} + +int qansel_get_instruction_bitmax(unsigned char* ptr, int offset, int* bitmax, int* qbitmax) +{ + unsigned char a0, a1, a2; + *bitmax = 0; + *qbitmax = 0; + switch (ptr[offset]) + { + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + case QANSEL_INSTRUCTION_DENSITY: + a0 = ptr[offset + 1]; + if (a0 > QANSEL_QBOUND_UPPER) return 0; + *qbitmax = a0 + 1; + return 1; + case QANSEL_INSTRUCTION_BARRIER: + case QANSEL_INSTRUCTION_BORN: + case QANSEL_INSTRUCTION_BSAMPLE: + a0 = ptr[offset + 1]; + if (a0 > QANSEL_QBOUND_UPPER && a0 != QANSEL_ALL_QUANTUM) return 0; + if (a0 != QANSEL_ALL_QUANTUM) *qbitmax = a0 + 1; + return 1; + case QANSEL_INSTRUCTION_IF_E: + case QANSEL_INSTRUCTION_IF_NE: + case QANSEL_INSTRUCTION_IF_G: + case QANSEL_INSTRUCTION_IF_GE: + case QANSEL_INSTRUCTION_IF_L: + case QANSEL_INSTRUCTION_IF_LE: + a0 = ptr[offset + 1]; + if ((a0 > QANSEL_CBOUND_UPPER || a0 < QANSEL_CBOUND_LOWER) && a0 != QANSEL_ALL_CLASSIC) return 0; + if (a0 != QANSEL_ALL_CLASSIC) *bitmax = (a0 - QANSEL_CBOUND_LOWER) + 1; + return 1; + case QANSEL_INSTRUCTION_SET: + case QANSEL_INSTRUCTION_RESET: + case QANSEL_INSTRUCTION_PRINT: + a0 = ptr[offset + 1]; + if (a0 == QANSEL_ALL || a0 == QANSEL_ALL_QUANTUM || a0 == QANSEL_ALL_CLASSIC) return 1; + if (a0 <= QANSEL_QBOUND_UPPER) + { + *qbitmax = a0 + 1; + return 1; + } + else if (a0 >= QANSEL_CBOUND_LOWER && a0 <= QANSEL_CBOUND_UPPER) + { + *bitmax = (a0 - QANSEL_CBOUND_LOWER) + 1; + return 1; + } + return 0; + case QANSEL_INSTRUCTION_CX: + case QANSEL_INSTRUCTION_SWAP: + a0 = ptr[offset + 1]; + a1 = ptr[offset + 2]; + if (a0 > QANSEL_QBOUND_UPPER) return 0; + if (a1 > QANSEL_QBOUND_UPPER) return 0; + *qbitmax = (a0 > a1 ? a0 : a1) + 1; + return 1; + case QANSEL_INSTRUCTION_CCX: + case QANSEL_INSTRUCTION_CSWAP: + a0 = ptr[offset + 1]; + a1 = ptr[offset + 2]; + a2 = ptr[offset + 3]; + if (a0 > QANSEL_QBOUND_UPPER || a1 > QANSEL_QBOUND_UPPER || a2 > QANSEL_QBOUND_UPPER) return 0; + *qbitmax = ((a0 > a1) && (a0 > a2) ? a0 : ((a1 > a0) && (a1 > a2) ? a1 : a2)) + 1; + return 1; + case QANSEL_INSTRUCTION_MEASURE: + a0 = ptr[offset + 1]; + a1 = ptr[offset + 2]; + if (a0 > QANSEL_QBOUND_UPPER) return 0; + if (a1 > QANSEL_CBOUND_UPPER || a1 < QANSEL_CBOUND_LOWER) return 0; + *qbitmax = a0 + 1; + *bitmax = (a1 - QANSEL_CBOUND_LOWER) + 1; + return 1; + case QANSEL_INSTRUCTION_EXIT: + return 1; + } + return 0; +} + +int qansel_get_instruction_size(unsigned char instr) +{ + switch (instr) + { + case QANSEL_INSTRUCTION_X: return 1 + 1; + case QANSEL_INSTRUCTION_Y: return 1 + 1; + case QANSEL_INSTRUCTION_Z: return 1 + 1; + case QANSEL_INSTRUCTION_H: return 1 + 1; + case QANSEL_INSTRUCTION_S: return 1 + 1; + case QANSEL_INSTRUCTION_T: return 1 + 1; + case QANSEL_INSTRUCTION_SDG: return 1 + 1; + case QANSEL_INSTRUCTION_TDG: return 1 + 1; + case QANSEL_INSTRUCTION_RX: return 1 + 1 + sizeof(float); + case QANSEL_INSTRUCTION_RY: return 1 + 1 + sizeof(float); + case QANSEL_INSTRUCTION_RZ: return 1 + 1 + sizeof(float); + case QANSEL_INSTRUCTION_U1: return 1 + 1 + sizeof(float); + case QANSEL_INSTRUCTION_U2: return 1 + 1 + sizeof(float) * 2; + case QANSEL_INSTRUCTION_U3: return 1 + 1 + sizeof(float) * 3; + case QANSEL_INSTRUCTION_CX: return 1 + 2; + case QANSEL_INSTRUCTION_SWAP: return 1 + 2; + case QANSEL_INSTRUCTION_CCX: return 1 + 3; + case QANSEL_INSTRUCTION_CSWAP: return 1 + 3; + case QANSEL_INSTRUCTION_MEASURE: return 1 + 2; + case QANSEL_INSTRUCTION_DENSITY: return 1 + 1; + case QANSEL_INSTRUCTION_BSAMPLE: return 1 + 1; + case QANSEL_INSTRUCTION_BORN: return 1 + 1; + case QANSEL_INSTRUCTION_IF_E: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_IF_NE: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_IF_G: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_IF_GE: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_IF_L: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_IF_LE: return 1 + 1 + sizeof(unsigned short); + case QANSEL_INSTRUCTION_SET: return 1 + 1; + case QANSEL_INSTRUCTION_RESET: return 1 + 1; + case QANSEL_INSTRUCTION_PRINT: return 1 + 1; + case QANSEL_INSTRUCTION_BARRIER: return 1 + 1; + case QANSEL_INSTRUCTION_EXIT: return 1; + } + return 0; +} + +void qansel_get_barrier(QBytecode** qbc, int idx) +{ + switch ((*qbc)[idx].op) + { + case QANSEL_INSTRUCTION_EXIT: + case QANSEL_INSTRUCTION_IF_E: + case QANSEL_INSTRUCTION_IF_NE: + case QANSEL_INSTRUCTION_IF_G: + case QANSEL_INSTRUCTION_IF_GE: + case QANSEL_INSTRUCTION_IF_L: + case QANSEL_INSTRUCTION_IF_LE: + (*qbc)[idx].barrier_width = 0; + return; + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + (*qbc)[idx].barrier_width = 1; + (*qbc)[idx].barrier[0] = (*qbc)[idx].bytes[1]; + return; + case QANSEL_INSTRUCTION_MEASURE: + case QANSEL_INSTRUCTION_DENSITY: + (*qbc)[idx].barrier_width = 2; + (*qbc)[idx].barrier[0] = (*qbc)[idx].bytes[1]; + (*qbc)[idx].barrier[1] = (*qbc)[idx].bytes[2]; + return; + case QANSEL_INSTRUCTION_BSAMPLE: + case QANSEL_INSTRUCTION_BORN: + case QANSEL_INSTRUCTION_BARRIER: + case QANSEL_INSTRUCTION_PRINT: + case QANSEL_INSTRUCTION_SET: + case QANSEL_INSTRUCTION_RESET: + switch ((*qbc)[idx].bytes[1]) + { + case QANSEL_ALL: + (*qbc)[idx].barrier_width = QANSEL_QUBITS_MAX * 2; + for (int i = 0; i < QANSEL_QUBITS_MAX; i++) (*qbc)[idx].barrier[i] = i; + for (int i = 0; i < QANSEL_QUBITS_MAX; i++) (*qbc)[idx].barrier[i] = i + QANSEL_CBOUND_LOWER; + return; + case QANSEL_ALL_QUANTUM: + (*qbc)[idx].barrier_width = QANSEL_QUBITS_MAX; + for (int i = 0; i < QANSEL_QUBITS_MAX; i++) (*qbc)[idx].barrier[i] = i; + return; + case QANSEL_ALL_CLASSIC: + (*qbc)[idx].barrier_width = QANSEL_QUBITS_MAX; + for (int i = 0; i < QANSEL_QUBITS_MAX; i++) (*qbc)[idx].barrier[i] = i + QANSEL_CBOUND_LOWER; + return; + default: + (*qbc)[idx].barrier_width = 1; + (*qbc)[idx].barrier[0] = (*qbc)[idx].bytes[1]; + return; + } + case QANSEL_INSTRUCTION_CX: + case QANSEL_INSTRUCTION_SWAP: + (*qbc)[idx].barrier_width = 2; + (*qbc)[idx].barrier[0] = (*qbc)[idx].bytes[1]; + (*qbc)[idx].barrier[1] = (*qbc)[idx].bytes[2]; + return; + case QANSEL_INSTRUCTION_CCX: + case QANSEL_INSTRUCTION_CSWAP: + (*qbc)[idx].barrier_width = 3; + (*qbc)[idx].barrier[0] = (*qbc)[idx].bytes[1]; + (*qbc)[idx].barrier[1] = (*qbc)[idx].bytes[2]; + (*qbc)[idx].barrier[2] = (*qbc)[idx].bytes[3]; + return; + } + fprintf(stderr, "QAnsel (#%i): Unknown error in barrier analysis.\n", idx); + exit(1); +} + +void qansel_born(QAnselContext* ctx, cpx_mtx_t* stateVector, int PC, int qubitCount, unsigned char q0, unsigned char mode) +{ + unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount); + if (q0 == QANSEL_ALL_QUANTUM) + { + float *psisquared = malloc(sizeof(float) * qubitCountPow2); + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + cpx_t n; + cpx_mtx_get(stateVector, 0, j, &n); + psisquared[j] = cpx_magsqr(&n); + } + //add any error skew if desired + if (ctx->noise > 0) + { + float totalprob = 0; + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + psisquared[j] += ctx->noise; + totalprob += psisquared[j]; + } + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + psisquared[j] /= totalprob; + } + } + if (mode == 0) + { + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + unsigned int tmp = j; + for (unsigned char k = 0; k < qubitCount; k++) + { + putchar('0' + (tmp >> (qubitCount - 1) & 1)); + tmp <<= 1; + } + printf("\t%.2f%%\n", psisquared[j] * 100); + } + free(psisquared); + return; + } + float *psisquared_sorted = malloc(sizeof(float) * qubitCountPow2); + unsigned short *sorting = malloc(sizeof(unsigned short) * qubitCountPow2); + unsigned short *rsorting = malloc(sizeof(unsigned short) * qubitCountPow2); + for (int i = 0; i < qubitCountPow2; i++) psisquared_sorted[i] = 0; + for (int i = 0; i < qubitCountPow2; i++) + { + float max = -1; + int maxi = -1; + for (int j = 0; j < qubitCountPow2; j++) + { + if (psisquared[j] >= max) + { + max = psisquared[j]; + maxi = j; + } + } + psisquared_sorted[i] = max; + sorting[i] = maxi; + rsorting[maxi] = i; + psisquared[maxi] = -1; + } + //for (int i = 0; i < qubitCountPow2; i++) printf("]%i -> %i: %f[\n", rsorting[i], sorting[i], psisquared_sorted[i]); + free(psisquared); + unsigned short *stats = malloc(sizeof(unsigned short) * qubitCountPow2); + for (int i = 0; i < qubitCountPow2; i++) stats[i] = 0; + for (int i = 0; i < (ctx->bsampling_shots); i++) + { + float r = (ctx->hidden_variable_set) ? qansel_rand_h() : qansel_rand_t(ctx); + float j = 0; + //printf("--------------------------------\n"); + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + //printf("%f >= %f\n", psisquared_sorted[j], r); + if (psisquared_sorted[j] >= (r - 0.0000001)) + { + stats[j]++; + //printf("%i _/\n", j); + break; + } + r -= psisquared_sorted[j]; + if (j == qubitCountPow2 - 1) + { + stats[j]++; + //printf("%i _/\n", j); + } + } + } + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + unsigned int tmp = j; + for (unsigned char k = 0; k < qubitCount; k++) + { + putchar('0' + (tmp >> (qubitCount - 1) & 1)); + tmp <<= 1; + } + //printf("%i -> %i\n", j, sorting[j]); + printf("\t%i\t%.2f%%\n", stats[rsorting[j]], ((float)stats[rsorting[j]] / (float)(ctx->bsampling_shots)) * (float)100); + } + free(psisquared_sorted); + free(stats); + free(sorting); + free(rsorting); + } + else if (q0 <= QANSEL_QBOUND_UPPER) + { + float prob = 0; + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + cpx_t n; + cpx_mtx_get(stateVector, 0, j, &n); + if ((j >> q0) & 1) + { + prob += cpx_magsqr(&n); + } + } + if (mode == 0) + { + printf("0\t%.2f%%\n", (1 - prob) * 100.0); + printf("1\t%.2f%%\n", prob * 100.0); + return; + } + + unsigned short stats = 0; + for (int i = 0; i < (ctx->bsampling_shots); i++) + { + float r = (ctx->hidden_variable_set) ? qansel_rand_h() : qansel_rand_t(ctx); + stats += r < prob ? 1 : 0; + } + printf("0\t%i\t%.2f%%\n", (ctx->bsampling_shots) - stats, (((float)(ctx->bsampling_shots) - (float)stats) / (float)(ctx->bsampling_shots)) * (float)100); + printf("1\t%i\t%.2f%%\n", stats, ((float)stats / (float)(ctx->bsampling_shots)) * (float)100); + } +} + +void qansel_density_or_print(cpx_mtx_t* stateVector, unsigned char* bitVector, unsigned char density, int bitCount, int qubitCount, unsigned char a0) +{ + unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount); + if (a0 == QANSEL_ALL_QUANTUM || a0 == QANSEL_ALL_CLASSIC || a0 == QANSEL_ALL) + { + if (a0 == QANSEL_ALL_QUANTUM || a0 == QANSEL_ALL) + { + printf("[ "); cpx_mtx_print(stateVector); printf(" ]\n"); + } + if (a0 == QANSEL_ALL_CLASSIC || a0 == QANSEL_ALL) + { + for (int32_t j = bitCount - 1; j >= 0; j--) + { + putchar('0' + bitVector[j]); + } + putchar('\n'); + } + } + else if (a0 >= QANSEL_CBOUND_LOWER && a0 <= QANSEL_CBOUND_UPPER) + { + putchar('0' + bitVector[a0 - QANSEL_CBOUND_LOWER]); + putchar('\n'); + } + else + { + cpx_mtx_t tmp; + cpx_mtx_init(&tmp, 1, 2); + for (unsigned int j = 0; j < qubitCountPow2; j++) + { + if ((j >> a0) & 1) + { + cpx_t a, b; + cpx_mtx_get(&tmp, 0, 1, &a); + cpx_mtx_get(stateVector, 0, j, &b); + a.real += b.real; + a.imaginary += b.imaginary; + cpx_mtx_set(&tmp, 0, 1, &a); + } + else + { + cpx_t a, b; + cpx_mtx_get(&tmp, 0, 0, &a); + cpx_mtx_get(stateVector, 0, j, &b); + a.real += b.real; + a.imaginary += b.imaginary; + cpx_mtx_set(&tmp, 0, 0, &a); + } + } + float multiplier = 0; + cpx_t n; + cpx_mtx_get(&tmp, 0, 0, &n); + multiplier += cpx_magsqr(&n); + cpx_mtx_get(&tmp, 0, 1, &n); + multiplier += cpx_magsqr(&n); + multiplier = sqrt(1 / multiplier); + n.real *= multiplier; + n.imaginary *= multiplier; + cpx_mtx_set(&tmp, 0, 1, &n); + cpx_mtx_get(&tmp, 0, 0, &n); + n.real *= multiplier; + n.imaginary *= multiplier; + cpx_mtx_set(&tmp, 0, 0, &n); + + if (density) + { + cpx_t a, b, c, d, x, y, z, w; + cpx_mtx_get(&tmp, 0, 0, &a); + cpx_mtx_get(&tmp, 0, 1, &b); + cpx_mtx_get(&tmp, 0, 0, &c); + cpx_mtx_get(&tmp, 0, 1, &d); + c.imaginary *= -1; + d.imaginary *= -1; + cpx_mul(&x, &a, &c); + cpx_mul(&y, &a, &d); + cpx_mul(&z, &b, &c); + cpx_mul(&w, &b, &d); + char* sx = cpx_str(&x); + char* sy = cpx_str(&y); + char* sz = cpx_str(&z); + char* sw = cpx_str(&w); + printf("[ %s, %s ]\n", sx, sy); + printf("[ %s, %s ]\n", sz, sw); + free(sx); + free(sy); + free(sz); + free(sw); + } + else + { + printf("[ "); cpx_mtx_print(&tmp); printf(" ]\n"); + } + cpx_mtx_free(&tmp); + } +} + +float qansel_get_float(unsigned char* program, int offset) +{ + float ret; + memcpy(&ret, program + offset, sizeof(float)); + return ret; +} + +short qansel_get_short(unsigned char* program, int offset) +{ + short ret; + memcpy(&ret, program + offset, sizeof(short)); + return ret; +} + +int qansel_get_int(unsigned char* program, int offset) +{ + int ret; + memcpy(&ret, program + offset, sizeof(int)); + return ret; +} + +void qansel_set(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned char* bitVector, int qubitCount, int bitCount, unsigned char q0, unsigned char value) +{ + unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount); + if (q0 == QANSEL_ALL) + { + cpx_mtx_set2(stateVector, 0, 0, 1, 0); + for (unsigned int j = 1; j < qubitCountPow2; j++) + { + cpx_mtx_set2(stateVector, 0, j, 0, 0); + } + for (unsigned char j = 0; j < bitCount; j++) + { + bitVector[j] = 0; + } + } + else if (q0 == QANSEL_ALL_QUANTUM) + { + cpx_mtx_set2(stateVector, 0, 0, 1, 0); + for (unsigned int j = 1; j < qubitCountPow2; j++) + { + cpx_mtx_set2(stateVector, 0, j, 0, 0); + } + } + else if (q0 == QANSEL_ALL_CLASSIC) + { + for (unsigned char j = 0; j < bitCount; j++) + { + bitVector[j] = 0; + } + } + else if (q0 <= QANSEL_QBOUND_UPPER) + { + unsigned char bit = qansel_measure(ctx, stateVector, qubitCount, q0); + if (bit != value) + { + qansel_instruction(ctx, stateVector, qubitCount, QANSEL_INSTRUCTION_X, q0, 0, 0, 0, NULL); + } + } + else if (q0 >= QANSEL_CBOUND_LOWER && q0 <= QANSEL_CBOUND_UPPER) + { + bitVector[q0 - QANSEL_CBOUND_LOWER] = 0; + } +} + +unsigned char qansel_compare(unsigned char* bitVector, int bitCount, int PC, unsigned char a0, short op) +{ + unsigned char ret = 0; + short val; + if (a0 == QANSEL_ALL_CLASSIC) + { + val = 0; + for (int32_t j = bitCount - 1; j >= 0; j--) + { + val = (val << 1) | bitVector[j]; + } + if (val == op) ret |= QANSEL_FLAGS_EQUAL; + if (val > op) ret |= QANSEL_FLAGS_GREATER; + if (val < op) ret |= QANSEL_FLAGS_LESSER; + } + else if (a0 >= QANSEL_CBOUND_LOWER && a0 <= QANSEL_CBOUND_UPPER) + { + val = bitVector[a0 - QANSEL_CBOUND_LOWER]; + if (val == op) ret |= QANSEL_FLAGS_EQUAL; + if (val > op) ret |= QANSEL_FLAGS_GREATER; + if (val < op) ret |= QANSEL_FLAGS_LESSER; + } + return ret; +} + +//computes program efficiency +// points are awarded for how long strings of single qubit +// instructions are +int qansel_efficiency(QBytecode* program, int programSize, int head) +{ + int score = 0; + int tmpscore = 0; + int idx = head; + + /* + printf("--------------------------------\n"); + for (int i = 0; i < programSize; i++) + { + if (i != -1) + { + printf("%02X: %s\n", i, qansel_instruction_to_string(program[i].op)); + } + else + { + printf("%02X: ERR\n", i); + } + } + printf("================================\n"); + printf("%i\n", head); + for (int i = 0; i < programSize; i++) + { + printf("%i -> %i\n", i, program[i].next); + } + printf("================================\n"); + */ + + for (int i = 0; i < programSize; i++) + { + //printf("%02X: %s (%i)\n", i, qansel_instruction_to_string(program[idx].op), idx); + + switch (program[idx].op) + { + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + switch (tmpscore) + { + case 0: tmpscore = 1; break; + default: tmpscore *= 2; + } + break; + case QANSEL_INSTRUCTION_IF_E: + case QANSEL_INSTRUCTION_IF_NE: + case QANSEL_INSTRUCTION_IF_G: + case QANSEL_INSTRUCTION_IF_GE: + case QANSEL_INSTRUCTION_IF_L: + case QANSEL_INSTRUCTION_IF_LE: + break; + default: + if (tmpscore > 1) score += tmpscore; + tmpscore = 0; + } + int oidx = idx; + idx = program[idx].next; + if (idx == -1) break; + } + return score; +} + +int qansel_crawl(QAnselContext* ctx, unsigned char* program, int programSize, int* qubitCount, int* bitCount) +{ + if (ctx->verbose) printf("Crawling program . . .\n"); + int PC = 0; + *qubitCount = 0; + *bitCount = 0; + while (PC < programSize) + { + int next = qansel_get_instruction_size(program[PC]); + if (next == 0) + { + printf("QAnsel (%04X): Invalid instruction 0x%02x.\n", PC, program[PC]); + return 0; + } + int bitmax, qbitmax; + int success = qansel_get_instruction_bitmax(program, PC, &bitmax, &qbitmax); + if (!success) + { + fprintf(stderr, "QAnsel (%04X): Invalid index.\n", PC); + return 0; + } + if (bitmax > *bitCount) *bitCount = bitmax; + if (qbitmax > *qubitCount) *qubitCount = qbitmax; + PC += next; + } + if (*qubitCount > ctx->qubit_count) + { + fprintf(stderr, "QAnsel: Not enough qubits allocated.\n"); + return 0; + } + if (*bitCount > ctx->bit_count) + { + fprintf(stderr, "QAnsel: Not enough classical bits allocated.\n"); + return 0; + } + *qubitCount = ctx->qubit_count; + *bitCount = ctx->bit_count; + return 1; +} + +void qansel_reorder(QAnselContext* ctx, unsigned char* program, int programSize) +{ + if (ctx->verbose) printf("Reordering . . .\n"); + + //break out program for easier manipulation + QBytecode* ramInstr = malloc(0); + int ramInstrLen = 0; + int copyifop = 0; + unsigned char ifop[16]; + int PC = 0; + + while (PC < programSize) + { + int next = qansel_get_instruction_size(program[PC]); + switch (program[PC]) + { + case QANSEL_INSTRUCTION_IF_E: + case QANSEL_INSTRUCTION_IF_NE: + case QANSEL_INSTRUCTION_IF_G: + case QANSEL_INSTRUCTION_IF_GE: + case QANSEL_INSTRUCTION_IF_L: + case QANSEL_INSTRUCTION_IF_LE: + memcpy(ifop, program + PC, next); + copyifop = 1; + PC += next; + continue; + } + ramInstr = realloc(ramInstr, sizeof(QBytecode) * (ramInstrLen + 1)); + ramInstr[ramInstrLen].size = next; + memcpy(ramInstr[ramInstrLen].bytes, program + PC, next); + ramInstr[ramInstrLen].op = program[PC]; + ramInstr[ramInstrLen].use_ifop = copyifop; + if (copyifop) + { + memcpy(ramInstr[ramInstrLen].ifop, ifop, sizeof(ifop)); + copyifop = 0; + } + ramInstrLen++; + PC += next; + } + + for (int i = 0; i < ramInstrLen; i++) + { + if (i == 0) + { + ramInstr[0].prev = -1; + ramInstr[0].next = ramInstrLen > 1 ? 1 : -1; + } + else if (i == ramInstrLen - 1) + { + ramInstr[ramInstrLen - 1].prev = i - 1; + ramInstr[ramInstrLen - 1].next = -1; + } + else + { + ramInstr[i].prev = i - 1; + ramInstr[i].next = i + 1; + } + //we can only rearrange single-qubit instructions + // so treat other instructions as pre-checked + switch (ramInstr[i].op) + { + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + ramInstr[i].checked = 0; + break; + default: ramInstr[i].checked = 1; + } + } + + for (int i = 0; i < ramInstrLen; i++) + { + qansel_get_barrier(&ramInstr, i); + } + + QBytecode* reordered[2]; + reordered[0] = malloc(ramInstrLen * sizeof(QBytecode)); + reordered[1] = malloc(ramInstrLen * sizeof(QBytecode)); + memcpy(reordered[0], ramInstr, ramInstrLen * sizeof(QBytecode)); + memcpy(reordered[1], ramInstr, ramInstrLen * sizeof(QBytecode)); + + int efficiency = qansel_efficiency(ramInstr, ramInstrLen, 0); + int head[2]; + head[0] = 0; + head[1] = 0; + + for (int reorid = 0; reorid < 2; reorid++) + { + while (1) + { + //find the next unchecked instruction + int found = -1; + int idx = head[reorid]; + for (int i = 0; i < ramInstrLen; i++) + { + if (reordered[reorid][idx].checked == 0) + { + found = idx; + break; + } + idx = reordered[reorid][idx].next; + } + if (found == -1) break; + //check if it can be moved up + reordered[reorid][found].checked = 1; + int blocked = 0; + idx = reordered[reorid][idx].prev; + while (idx != -1) + { + for (int i = 0; i < reordered[reorid][idx].barrier_width; i++) + { + //verify no classical bit barrier + if (reordered[reorid][found].use_ifop) + { + switch (reordered[reorid][found].ifop[1]) + { + case QANSEL_ALL_CLASSIC: + case QANSEL_ALL: + blocked = 1; + break; + default: + if (reordered[reorid][found].ifop[1] == reordered[reorid][idx].barrier[i]) + blocked = 1; + } + if (blocked) break; + } + //verify no quantum bit barrier + if (reordered[reorid][found].bytes[1] == reordered[reorid][idx].barrier[i]) + { + blocked = 1; + break; + } + } + if (blocked) break; + idx = reordered[reorid][idx].prev; + } + if (reordered[reorid][found].prev != idx) + { + QBytecode* tmpbc = malloc(ramInstrLen * sizeof(QBytecode)); + memcpy(tmpbc, reordered[reorid], ramInstrLen * sizeof(QBytecode)); + reordered[reorid][ reordered[reorid][found].prev ].next = reordered[reorid][found].next; + if (reordered[reorid][found].next != -1) + { + reordered[reorid][ reordered[reorid][found].next ].prev = reordered[reorid][found].prev; + } + + if (idx == -1) + { + //idx.next.prev = found + //found.next = head + //found.prev = -1 + //head = found + reordered[reorid][head[reorid]].prev = found; + reordered[reorid][found].next = head[reorid]; + reordered[reorid][found].prev = -1; + head[reorid] = found; + //printf("head: %i\n", found); + } + else + { + /* + [idx].next -> .prev <- o + found + + */ + //idx.next.prev = found + //found.next = idx.next + //found.prev = idx + //idx.next = found + reordered[reorid][ reordered[reorid][idx].next ].prev = found; + reordered[reorid][found].next = reordered[reorid][idx].next; + reordered[reorid][found].prev = idx; + reordered[reorid][idx].next = found; + } + if (reorid == 1) + { + /* + printf("%i\n", head[reorid]); + for (int i = 0; i < ramInstrLen; i++) + { + printf("%i -> %i\n", i, reordered[reorid][i].next); + } + printf("===========================\n"); + */ + int tmpeff = qansel_efficiency(reordered[reorid], ramInstrLen, head[reorid]); + head[reorid] = 0; //correct? + if (tmpeff >= efficiency) + { + efficiency = tmpeff; + } + else + { + memcpy(reordered[reorid], tmpbc, ramInstrLen * sizeof(QBytecode)); + } + } + free(tmpbc); + } + } + } + + efficiency = qansel_efficiency(ramInstr, ramInstrLen, 0); + int best = -1; + int besti = -1; + for (int i = 0; i < 2; i++) + { + int tmpeff = qansel_efficiency(reordered[i], ramInstrLen, head[i]); + if (tmpeff > best) + { + besti = i; + best = tmpeff; + } + } + if (best > efficiency) + { + memcpy(ramInstr, reordered[besti], ramInstrLen * sizeof(QBytecode)); + } + for (int i = 0; i < 2; i++) + { + free(reordered[i]); + } + + if (best > efficiency) + { + int copyloc = 0; + besti = head[besti]; + for (int i = 0; i < ramInstrLen; i++) + { + int next; + if (ramInstr[besti].use_ifop) + { + next = qansel_get_instruction_size(ramInstr[besti].ifop[0]); + memcpy(program + copyloc, ramInstr[besti].ifop, next); + copyloc += next; + } + next = qansel_get_instruction_size(ramInstr[besti].op); + memcpy(program + copyloc, ramInstr[besti].bytes, next); + besti = ramInstr[besti].next; + copyloc += next; + } + } + + if (ctx->verbose) + { + if (best > efficiency) + { + printf("Efficiency score increased from %i to %i.\n", efficiency, best); + } + else + { + printf("Could not increase efficiency.\n"); + } + } + + free(ramInstr); +} + +void qansel_run(QAnselContext* ctx, unsigned char* program, int programSize, int qubitCount, int bitCount, unsigned char* outputBitVector) +{ + int useQueue = ((ctx->optimization_level) & QANSEL_MODE_SMART) ? 1 : 0; + int PC = 0; + unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount); + unsigned char bitVector[bitCount]; + memset(bitVector, 0, bitCount); + cpx_mtx_t stateVector; + cpx_mtx_t* queueVector = NULL; + cpx_mtx_init(&stateVector, 1, qubitCountPow2); + cpx_mtx_set2(&stateVector, 0, 0, 1, 0); + if (useQueue) + { + queueVector = malloc(qubitCount * sizeof(cpx_mtx_t)); + for (int i = 0; i < qubitCount; i++) + { + queueVector[i].ptr = malloc(sizeof(Identity)); + memcpy(queueVector[i].ptr, Identity, sizeof(Identity)); + queueVector[i].rows = 2; + queueVector[i].cols = 2; + } + } + if (ctx->display_delay) { ctx->display_delay = display(&stateVector, qubitCount, ctx->display_delay); } + unsigned char skip = 0, a0 = 0, a1 = 0, a2 = 0; + unsigned char flags = 0; + unsigned short tmp = 0; + unsigned char queueFlushed = 1; + float f0, f1, f2; + + while (PC < programSize) + { + //printf("%i;%i\n", PC, programSize); + int next = qansel_get_instruction_size(program[PC]); + if (skip) + { + skip = 0; + } + else + { + unsigned char instr = program[PC]; + + //flush the queue if any non-single qubit instructions are called + if (useQueue) + { + switch (instr) + { + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + queueFlushed = 0; + break; + case QANSEL_INSTRUCTION_IF_E: + case QANSEL_INSTRUCTION_IF_NE: + case QANSEL_INSTRUCTION_IF_G: + case QANSEL_INSTRUCTION_IF_GE: + case QANSEL_INSTRUCTION_IF_L: + case QANSEL_INSTRUCTION_IF_LE: + break; + default: + if (!queueFlushed) + { + qansel_instruction(ctx, &stateVector, qubitCount, instr, 0, 0, 0, 0, &queueVector); + for (int i = 0; i < qubitCount; i++) + { + memcpy(queueVector[i].ptr, Identity, sizeof(Identity)); + } + queueFlushed = 1; + } + break; + } + } + + if (ctx->verbose) + { + printf("%08X: %s\n", PC, qansel_instruction_to_string(instr)); + } + switch (instr) + { + case QANSEL_INSTRUCTION_X: + case QANSEL_INSTRUCTION_Y: + case QANSEL_INSTRUCTION_Z: + case QANSEL_INSTRUCTION_H: + case QANSEL_INSTRUCTION_S: + case QANSEL_INSTRUCTION_T: + case QANSEL_INSTRUCTION_SDG: + case QANSEL_INSTRUCTION_TDG: + case QANSEL_INSTRUCTION_RX: + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_RZ: + case QANSEL_INSTRUCTION_U1: + case QANSEL_INSTRUCTION_U2: + case QANSEL_INSTRUCTION_U3: + switch (instr) + { + case QANSEL_INSTRUCTION_RX: + f0 = qansel_get_float(program, PC + 2); + f1 = -M_PI / 2; + f2 = M_PI / 2; + break; + case QANSEL_INSTRUCTION_RY: + case QANSEL_INSTRUCTION_U1: + f0 = qansel_get_float(program, PC + 2); + f1 = 0; + f2 = 0; + break; + case QANSEL_INSTRUCTION_RZ: + f0 = 0; + f1 = 0; + f2 = qansel_get_float(program, PC + 2); + break; + case QANSEL_INSTRUCTION_U2: + f0 = qansel_get_float(program, PC + 2); + f1 = qansel_get_float(program, PC + 2 + sizeof(float)); + f2 = 0; + break; + case QANSEL_INSTRUCTION_U3: + f0 = qansel_get_float(program, PC + 2); + f1 = qansel_get_float(program, PC + 2 + sizeof(float)); + f2 = qansel_get_float(program, PC + 2 + sizeof(float) * 2); + break; + default: + f0 = 0; + f1 = 0; + f2 = 0; + break; + } + a0 = program[PC + 1]; + if (useQueue) qansel_instruction_queue(&queueVector, qubitCount, instr, a0, f0, f1, f2); + else qansel_instruction(ctx, &stateVector, qubitCount, instr, a0, f0, f1, f2, NULL); + break; + case QANSEL_INSTRUCTION_CX: + a0 = program[PC + 1]; + a1 = program[PC + 2]; + qansel_cnot(&stateVector, qubitCount, a0, a1); + break; + case QANSEL_INSTRUCTION_SWAP: + a0 = program[PC + 1]; + a1 = program[PC + 2]; + qansel_swap(&stateVector, qubitCount, a0, a1); + break; + case QANSEL_INSTRUCTION_CCX: + a0 = program[PC + 1]; + a1 = program[PC + 2]; + a2 = program[PC + 3]; + qansel_toffoli(&stateVector, qubitCount, a0, a1, a2); + break; + case QANSEL_INSTRUCTION_CSWAP: + a0 = program[PC + 1]; + a1 = program[PC + 2]; + a2 = program[PC + 3]; + qansel_fredkin(&stateVector, qubitCount, a0, a1, a2); + break; + case QANSEL_INSTRUCTION_MEASURE: + a0 = program[PC + 1]; + a1 = program[PC + 2] - QANSEL_CBOUND_LOWER; + bitVector[a1] = qansel_measure(ctx, &stateVector, qubitCount, a0); + break; + case QANSEL_INSTRUCTION_BORN: + a0 = program[PC + 1]; + qansel_born(ctx, &stateVector, PC, qubitCount, a0, 0); + break; + case QANSEL_INSTRUCTION_BSAMPLE: + a0 = program[PC + 1]; + qansel_born(ctx, &stateVector, PC, qubitCount, a0, 1); + break; + case QANSEL_INSTRUCTION_DENSITY: + a0 = program[PC + 1]; + qansel_density_or_print(&stateVector, bitVector, 1, bitCount, qubitCount, a0); + break; + case QANSEL_INSTRUCTION_PRINT: + a0 = program[PC + 1]; + qansel_density_or_print(&stateVector, bitVector, 0, bitCount, qubitCount, a0); + break; + case QANSEL_INSTRUCTION_BARRIER: + a0 = program[PC + 1]; + break; + case QANSEL_INSTRUCTION_SET: + a0 = program[PC + 1]; + qansel_set(ctx, &stateVector, bitVector, qubitCount, bitCount, a0, 1); + break; + case QANSEL_INSTRUCTION_RESET: + a0 = program[PC + 1]; + qansel_set(ctx, &stateVector, bitVector, qubitCount, bitCount, a0, 0); + break; + case QANSEL_INSTRUCTION_IF_E: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if (flags & QANSEL_FLAGS_EQUAL) skip = 0; + break; + case QANSEL_INSTRUCTION_IF_NE: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if (!(flags & QANSEL_FLAGS_EQUAL)) skip = 0; + break; + case QANSEL_INSTRUCTION_IF_G: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if (flags & QANSEL_FLAGS_GREATER) skip = 0; + break; + case QANSEL_INSTRUCTION_IF_GE: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if ((flags & QANSEL_FLAGS_GREATER) && (flags & QANSEL_FLAGS_EQUAL)) skip = 0; + break; + case QANSEL_INSTRUCTION_IF_L: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if (flags & QANSEL_FLAGS_LESSER) skip = 0; + break; + case QANSEL_INSTRUCTION_IF_LE: + a0 = program[PC + 1]; + tmp = qansel_get_short(program, PC + 2); + flags = qansel_compare(bitVector, bitCount, PC, a0, tmp); + skip = 1; if ((flags & QANSEL_FLAGS_LESSER) && (flags & QANSEL_FLAGS_EQUAL)) skip = 0; + break; + + } + } + PC += next; + if (ctx->display_delay) { ctx->display_delay = display(&stateVector, qubitCount, ctx->display_delay); } + } + if (outputBitVector != NULL) + { + for (int i = 0; i < bitCount; i++) + { + outputBitVector[i] = bitVector[i]; + } + } + cpx_mtx_free(&stateVector); + if (useQueue) + { + for (int i = 0; i < qubitCount; i++) + { + free(queueVector[i].ptr); + } + free(queueVector); + } + if (ctx->display_delay) { display(NULL, 0, 0); } + +} + +int qanselExecuteBytecode(unsigned char* buff, int sizeofbuff, QAnselContext* ctx) +{ + if (!qanselContextValidate(ctx)) return 0; + qansel_rand_s(qansel_hardware_getseed()); + + if (ctx != NULL && ctx->hidden_variable_set) + { + unsigned int tmp; + memcpy(&tmp, &(ctx->hidden_variable), sizeof(unsigned int)); + srand(tmp); + } + + unsigned short vals; + float valf; + int pos = 0; + + int qubitCount, bitCount; + if (!qansel_crawl(ctx, buff, sizeofbuff, &qubitCount, &bitCount)) + { + return 0; + } + + if (ctx != NULL && ctx->sampling_shots > 0) + { + unsigned short stats[65536]; + for (unsigned int i = 0; i < (1 << bitCount); i++) + { + stats[i] = 0; + } + unsigned char bitVect[bitCount]; + memset(bitVect, 0, bitCount); + for (int i = 0; i < bitCount; i++) bitVect[i] = 0; + unsigned int shots = ctx->sampling_shots; + for (unsigned int i = 0; i < shots; i++) + { + qansel_run(ctx, buff, sizeofbuff, qubitCount, bitCount, bitVect); + unsigned short stat = 0; + for (signed char j = bitCount - 1; j >= 0; j--) + { + stat = (stat << 1) | bitVect[j]; + } + stats[stat]++; + } + unsigned int count = 0; + for (unsigned int i = 0; i < (1 << bitCount); i++) + { + unsigned int tmp = i; + for (unsigned char j = 0; j < bitCount; j++) + { + unsigned char bit = (tmp >> (bitCount - 1) & 1); + if (j == (bitCount - (ctx->sampling_bit) - 1) && bit) + { + count += stats[i]; + } + if ((ctx->sampling_bit) == QANSEL_ALL) + { + putchar('0' + bit); + } + tmp <<= 1; + } + if ((ctx->sampling_bit) == QANSEL_ALL) + { + printf("\t%i\t%.2f%%\n", stats[i], ((float)stats[i] / (float)shots) * (float)100); + } + } + if ((ctx->sampling_bit) != QANSEL_ALL) + { + float prob = ((float)count / (float)shots) * (float)100; + printf("0\t%i\t%.2f%%\n", shots-count, ((float)100)-prob); + printf("1\t%i\t%.2f%%\n", count, prob); + } + } + else + { + qansel_run(ctx, buff, sizeofbuff, qubitCount, bitCount, NULL); + } + return 1; +} -- 2.39.5