Fri 13 Dec 2024 04:07:12 PM EST

diff --git a/src/bytecode.c b/src/bytecode.c
index 76355954a7c3bc59bbfa2b888697ef6c26974927..9c95884795361a8df7b43ad024ac15eaf9140bf0 100644 (file)
--- a/src/bytecode.c
+++ b/src/bytecode.c
@@ -1,4 +1,5 @@
 #include "bytecode.h"
+void qansel_normalize(cpx_mtx_t* stateVector, unsigned int qubitCountPow2);
 
 const char* qansel_instruction_to_string(unsigned char instr)
 {
@@ -448,6 +449,22 @@ void qansel_instruction
 unsigned char qansel_rand_b(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned int qubitCountPow2, unsigned char qubit)
 {
        cpx_t n;
+       cpx_mtx_t tmp;
+       if (ctx->noise != 0)
+       {
+               cpx_mtx_init(&tmp, stateVector->rows, stateVector->cols);
+               cpx_mtx_copy(&tmp, stateVector);
+               for (int i = 0; i < tmp.rows; i++)
+               {
+                       for (int j = 0; j < tmp.cols; j++)
+                       {
+                               cpx_mtx_get(stateVector, i, j, &n);
+                               n.real += n.real >= 0 ? ctx->noise : -(ctx->noise);
+                               cpx_mtx_set(stateVector, i, j, &n);
+                       }
+               }
+               qansel_normalize(stateVector, qubitCountPow2);
+       }
        float prob0 = 0;
        for (unsigned int i = 0; i < qubitCountPow2; i++)
        {
@@ -456,15 +473,20 @@ unsigned char qansel_rand_b(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned
                if (bit == 0) prob0 += cpx_magsqr(&n);
        }
        float r;
-       if (ctx->noise > 0)
+       /*if (ctx->noise > 0)
        {
                r = (ctx->hidden_variable_set) ? qansel_rand_h() : qansel_rand_t(ctx);
                if (r < ctx->noise)
                {
                        prob0 = 0.5;    
                }
-       }
+       }*/
        r = (ctx->hidden_variable_set) ? qansel_rand_h() : qansel_rand_t(ctx);
+       if (ctx->noise != 0)
+       {
+               cpx_mtx_copy(stateVector, &tmp);
+               cpx_mtx_free(&tmp);
+       }
        return r < prob0 ? 0 : 1;       
 }
 
@@ -500,40 +522,43 @@ void qansel_sample_stats(unsigned short* stats, int allMode, int shots, int bitC
        }
 }
 
-unsigned char qansel_measure(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char qubit)
+void qansel_normalize(cpx_mtx_t* stateVector, unsigned int qubitCountPow2)
 {
        cpx_t n;
-       unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount);
-       unsigned char newBit = qansel_rand_b(ctx, stateVector, qubitCountPow2, qubit);
        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);
-               }
+               probTot += cpx_magsqr(&n);
                cpx_mtx_set(stateVector, 0, i, &n);
        }
        float multiplier = probTot == 0 ? 1 : 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;
-               }
+               n.real *= multiplier;
+               n.imaginary *= multiplier;
                cpx_mtx_set(stateVector, 0, i, &n);
        }
+}
 
+unsigned char qansel_measure(QAnselContext* ctx, cpx_mtx_t* stateVector, unsigned char qubitCount, unsigned char qubit)
+{
+       cpx_t n;
+       unsigned int qubitCountPow2 = (unsigned int)pow(2, qubitCount);
+       unsigned char newBit = qansel_rand_b(ctx, stateVector, qubitCountPow2, qubit);
+       for (unsigned int i = 0; i < qubitCountPow2; i++)
+       {
+               unsigned char bit = (i >> qubit) & 1;
+               if (bit != newBit)
+               {
+                       cpx_mtx_get(stateVector, 0, i, &n);
+                       n.real = 0;
+                       n.imaginary = 0;
+                       cpx_mtx_set(stateVector, 0, i, &n);
+               }
+       }
+       qansel_normalize(stateVector, qubitCountPow2);
        return newBit;
 }
 
@@ -794,9 +819,10 @@ void qansel_born(QAnselContext* ctx, cpx_mtx_t* stateVector, int PC, int qubitCo
                return;
        }
 
-       cpx_mtx_t oldStateVector;
-       cpx_mtx_init(&oldStateVector, stateVector->rows, stateVector->cols);
-       cpx_mtx_copy(&oldStateVector, stateVector);
+       cpx_mtx_t tmp;
+       cpx_mtx_init(&tmp, stateVector->rows, stateVector->cols);
+       cpx_mtx_copy(&tmp, stateVector);
+
        unsigned short stats[65536];
        for (unsigned int i = 0; i < qubitCountPow2; i++)
        {
@@ -807,14 +833,14 @@ void qansel_born(QAnselContext* ctx, cpx_mtx_t* stateVector, int PC, int qubitCo
                unsigned int bitidx = 0;
                for (unsigned int j = 0; j < qubitCount; j++)
                {
-                       unsigned char newBit = qansel_measure(ctx, stateVector, qubitCount, j);
+                       unsigned char newBit = qansel_measure(ctx, &tmp, qubitCount, j);
                        bitidx |= newBit << j;
                }
                stats[bitidx]++;
-               cpx_mtx_copy(stateVector, &oldStateVector);
+               cpx_mtx_copy(&tmp, stateVector);
        }
        qansel_sample_stats(stats, QANSEL_ALL_QUANTUM, ctx->bsampling_shots, qubitCount, q0);
-       cpx_mtx_free(&oldStateVector);
+       cpx_mtx_free(&tmp);
 }
 
 void qansel_density_or_print(cpx_mtx_t* stateVector, unsigned char* bitVector, unsigned char density, int bitCount, int qubitCount, unsigned char a0)