From 6dbc68b0461b8cfa3ac9dd300b9ff7381ee9eed6 Mon Sep 17 00:00:00 2001
From: a <a@gmail.com>
Date: Fri, 13 Dec 2024 16:07:12 -0500
Subject: [PATCH] Fri 13 Dec 2024 04:07:12 PM EST

---
 src/bytecode.c | 80 +++++++++++++++++++++++++++++++++-----------------
 1 file changed, 53 insertions(+), 27 deletions(-)

diff --git a/src/bytecode.c b/src/bytecode.c
index 7635595..9c95884 100644
--- 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)
-- 
2.39.5