diff --git a/complex.c b/complex.c index 9263d74..ebb6781 100644 --- a/complex.c +++ b/complex.c @@ -35,6 +35,10 @@ int quda_complex_eq(complex_t op1, complex_t op2) { return 0; } +float quda_complex_abs_square(complex_t c) { + return c.real*c.real + c.imag*c.imag; +} + float quda_complex_abs(complex_t c) { float res = c.real*c.real + c.imag*c.imag; return sqrt(res); diff --git a/complex.h b/complex.h index 801f564..0a3bb8a 100644 --- a/complex.h +++ b/complex.h @@ -32,6 +32,13 @@ complex_t quda_complex_copy(complex_t c); /* Test equality of two complex numbers (returns 1 if equal, 0 otherwise) */ int quda_complex_eq(complex_t op1, complex_t op2); +/* Complex absolute square + * Equivalent to a complex number times its conjugate + * For a quantum state amplitude, this value represents the probability + * of the quantum state. + */ +float quda_complex_abs_square(complex_t c); + /* Complex modulus (absolute value) */ float quda_complex_abs(complex_t c); diff --git a/quantum_reg.c b/quantum_reg.c index db864da..c29f19b 100644 --- a/quantum_reg.c +++ b/quantum_reg.c @@ -3,6 +3,7 @@ #include "quantum_reg.h" #include +#include int quda_quantum_reg_init(quantum_reg* qreg, int qubits) { qreg->qubits = qubits; @@ -26,6 +27,125 @@ void quda_quantum_reg_delete(quantum_reg* qreg) { free(qreg->states); } +void quda_quantum_bit_set(int target, quantum_reg* qreg) { + int i; + uint64_t mask = 1 << target; + for(i=0;inum_states;i++) { + qreg->states[i].state = qreg->states[i].state | mask; + } + + quda_quantum_reg_coalesce(qreg); +} + +void quda_quantum_bit_reset(int target, quantum_reg* qreg) { + int i; + uint64_t mask = ~(1 << target); + for(i=0;inum_states;i++) { + qreg->states[i].state = qreg->states[i].state & mask; + } + + quda_quantum_reg_coalesce(qreg); +} + +/* Performs a measurement on the quantum register */ +int quda_quantum_reg_measure(quantum_reg* qreg, uint64_t* retval) { + if(retval == NULL) return -2; + float f = quda_rand_float(); + int i; + for(i=0;inum_states;i++) { + if(!quda_complex_eq(qreg->states[i].amplitude,QUDA_COMPLEX_ZERO)) { + f -= quda_complex_abs_square(qreg->states[i].amplitude); + if(f < 0) { + *retval = qreg->states[i].state; + return 0; + } + } + } + + return -1; +} + +/* Performs a real-world quantum measurement. + * The register collapses to the physical state measured with probability 1. + */ +int quda_quantum_reg_measure_and_collapse(quantum_reg* qreg, uint64_t* retval) { + if(retval == NULL) return -2; + float f = quda_rand_float(); + int i; + for(i=0;inum_states;i++) { + if(!quda_complex_eq(qreg->states[i].amplitude,QUDA_COMPLEX_ZERO)) { + f -= quda_complex_abs_square(qreg->states[i].amplitude); + if(f < 0) { + *retval = qreg->states[i].state; + qreg->states[0].state = qreg->states[i].state; + qreg->states[0].amplitude = QUDA_COMPLEX_ONE; + qreg->num_states = 1; + return 0; + } + } + } + + return -1; +} + +/* Measure 1 bit of a quantum register */ +int quda_quantum_bit_measure(int target, quantum_reg* qreg) { + float p = 0; + float f = quda_rand_float(); + uint64_t mask = 1 << target; + int i; + // Accumulate probability that the bit is in state |1> + for(i = 0;inum_states;i++) { + if(qreg->states[i].state & mask) { + p += quda_complex_abs_square(qreg->states[i].amplitude); + // TODO: Determine overhead of this comparison + if(p > quda_rand_float()) return 1; + } + } + + return 0; + //return quda_rand_float() < p ? 1 : 0; +} + +int quda_quantum_bit_measure_and_collapse(int target, quantum_reg* qreg) { + // Measure bit conventionally + // TODO: Can allow probability measurement from conventional to complete and remove it below + int retval = quda_quantum_bit_measure(target,qreg); + + // Collapse states to those possible + uint64_t mask = 1 << target; + float p = 0; + int i; + for(i=0;inum_states;i++) { + // TODO: Ideally, remove nested conditions + // TODO: Actually prune in this loop instead of just invalidating + if(qreg->states[i].state & mask) { + if(retval) { // this is a valid state, accumulate probability to renormalize + p += quda_complex_abs_square(qreg->states[i].amplitude); + } else { // this is an invalid state -- nullify + qreg->states[i].amplitude = QUDA_COMPLEX_ZERO; + } + } else { + if(!retval) { // valid state, accumulate probability + p += quda_complex_abs_square(qreg->states[i].amplitude); + } else { // invalid state -- nullify + qreg->states[i].amplitude = QUDA_COMPLEX_ZERO; + } + } + } + + // TODO: Remove this call for optimization within the above loop + quda_quantum_reg_prune(qreg); + + // Renormalize + float k = sqrt(1.0f/p); + for(i=0;inum_states;i++) { + qreg->states[i].amplitude = quda_complex_rmul(qreg->states[i].amplitude,k); + } + + return retval; +} + void quda_quantum_reg_prune(quantum_reg* qreg) { int i,end; for(i=0,end=qreg->num_states-1;i < end;i++) { @@ -112,6 +232,10 @@ int quda_quantum_reg_trim(quantum_reg* qreg) { return 0; } +float quda_rand_float() { + return rand()/(float)RAND_MAX; +} + int qstate_compare(const void* qstate1, const void* qstate2) { uint64_t diff = ((quantum_state_t*)qstate1)->state - ((quantum_state_t*)qstate2)->state; if(diff == 0) return 0; diff --git a/quantum_reg.h b/quantum_reg.h index 0cd7f8f..b5b87e2 100644 --- a/quantum_reg.h +++ b/quantum_reg.h @@ -29,14 +29,54 @@ typedef struct quantum_reg { */ int quda_quantum_reg_init(quantum_reg* qreg, int qubits); -/* Sets the register to a single physical state with probability 1. */ -void quda_quantum_reg_set(quantum_reg* qreg, uint64_t state); - /* Frees the given register for deletion. * If the passed qreg was dynamically allocated, it must still be freed separately. */ void quda_quantum_reg_delete(quantum_reg* qreg); +/* Sets the register to a single physical state with probability 1. */ +void quda_quantum_reg_set(quantum_reg* qreg, uint64_t state); + +/* Sets a single bit of a quantum register to 1 with probability 1. + * If the system is in a superposition, it is collapsed into the subset + * of possible states allowed by this value. + */ +void quda_quantum_bit_set(int target, quantum_reg* qreg); + +/* Sets a single bit of a quantum register to 0 with probability 1. + * If the system is in a superposition, it is collapsed into the subset + * of possible states allowed by this value. + */ +void quda_quantum_bit_reset(int target, quantum_reg* qreg); + +/* Performs a measurement on the quantum register and stores the state + * in 'retval' if non-NULL. + * Returns 0 on success, -1 on retval NULL, -2 on normalization error. + */ +// TODO: Attempt correction for minor normalization errors (ie floating point precision errors) +int quda_quantum_reg_measure(quantum_reg* qreg, uint64_t* retval); + +/* Performs a real-world quantum measurement and stores the state in + * 'retval' if non-NULL. + * On success, returns 0 and the register collapses to the physical state + * measured with probability 1. + * On failure, does not collapse. Returns -1 on retval NULL, -2 on + * normalization error. + */ +// TODO: Attempt correction for minor normalization errors (ie floating point precision errors) +int quda_quantum_reg_measure_and_collapse(quantum_reg* qreg, uint64_t* retval); + +/* Measure 1 bit of a quantum register */ +int quda_quantum_bit_measure(int target, quantum_reg* qreg); + +/* Perform a real-world quantum measurement of 1 quantum register bit. + * The system collapses into the subset of possible states allowed by + * the value of the measured bit. + * Simultaneously prunes zero-amplitude states. + * Does not coalesce identical states. + */ +int quda_quantum_bit_measure_and_collapse(int target, quantum_reg* qreg); + /* Removes zero-amplitude states from the register. */ void quda_quantum_reg_prune(quantum_reg* qreg); @@ -59,6 +99,11 @@ void quda_quantum_reg_coalesce(quantum_reg* qreg); */ int quda_quantum_reg_trim(quantum_reg* qreg); +/* Generates a float in the range [0,1) */ +// TODO: Look at performance implications of using 'double' here +float quda_rand_float(); + +/* Comparator for sorting the quantum_state array */ int qstate_compare(const void* qstate1, const void* qstate2); #endif // __QUDA_QUANTUM_REG_H