From 7d0db71088b651f8e518957e0730ee2cefbcc87c Mon Sep 17 00:00:00 2001
From: Fe-r-oz <feroz.ahmad.email@gmail.com>
Date: Mon, 23 Sep 2024 15:36:10 +0500
Subject: [PATCH] fix #113 - better bit-wrangling abstraction

---
 src/pauli_frames.jl | 45 +++++++++++++++++----------------------------
 1 file changed, 17 insertions(+), 28 deletions(-)

diff --git a/src/pauli_frames.jl b/src/pauli_frames.jl
index 45e9763c8..fd11dcdc4 100644
--- a/src/pauli_frames.jl
+++ b/src/pauli_frames.jl
@@ -79,18 +79,22 @@ function apply!(frame::PauliFrame, op::sMRX) # TODO implement a faster direct ve
     apply!(frame, sMRZ(op.qubit, op.bit))
 end
 
+function _get_bitmask_idxs(frame::PauliFrame, i::Int)
+    T = eltype(frame.frame.tab.xzs)
+    lowbit = T(1)
+    ibig = _div(T, i-1) + 1
+    ismall = _mod(T, i-1)
+    ismallm = lowbit << ismall
+    return ibig, ismall, ismallm
+end
+
 function apply!(frame::PauliFrame, op::sMZ) # TODO sMY, and faster sMX
     op.bit == 0 && return frame
     i = op.qubit
-    xzs = frame.frame.tab.xzs
-    T = eltype(xzs)
-    lowbit = T(1)
-    ibig = _div(T,i-1)+1
-    ismall = _mod(T,i-1)
-    ismallm = lowbit<<(ismall)
+    ibig, ismall, ismallm = _get_bitmask_idxs(frame,i)
 
     @inbounds @simd for f in eachindex(frame)
-        should_flip = !iszero(xzs[ibig,f] & ismallm)
+        should_flip = !iszero(frame.frame.tab.xzs[ibig,f] & ismallm)
         frame.measurements[f,op.bit] = should_flip
     end
 
@@ -99,22 +103,17 @@ end
 
 function apply!(frame::PauliFrame, op::sMRZ) # TODO sMRY, and faster sMRX
     i = op.qubit
-    xzs = frame.frame.tab.xzs
-    T = eltype(xzs)
-    lowbit = T(1)
-    ibig = _div(T,i-1)+1
-    ismall = _mod(T,i-1)
-    ismallm = lowbit<<(ismall)
+    ibig, ismall, ismallm = _get_bitmask_idxs(frame,i)
 
     if op.bit != 0
         @inbounds @simd for f in eachindex(frame)
-            should_flip = !iszero(xzs[ibig,f] & ismallm)
+            should_flip = !iszero(frame.frame.tab.xzs[ibig,f] & ismallm)
             frame.measurements[f,op.bit] = should_flip
         end
     end
     @inbounds @simd for f in eachindex(frame)
-        xzs[ibig,f] &= ~ismallm
-        rand(Bool) && (xzs[end÷2+ibig,f] ⊻= ismallm)
+        frame.frame.tab.xzs[ibig,f] &= ~ismallm
+        rand(Bool) && (frame.frame.tab.xzs[end÷2+ibig,f] ⊻= ismallm)
     end
 
     return frame
@@ -122,12 +121,7 @@ end
 
 function applynoise!(frame::PauliFrame,noise::UnbiasedUncorrelatedNoise,i::Int)
     p = noise.p
-    T = eltype(frame.frame.tab.xzs)
-
-    lowbit = T(1)
-    ibig = _div(T,i-1)+1
-    ismall = _mod(T,i-1)
-    ismallm = lowbit<<(ismall)
+    ibig, ismall, ismallm = _get_bitmask_idxs(frame,i)
     p = p/3
 
     @inbounds @simd for f in eachindex(frame)
@@ -145,12 +139,7 @@ function applynoise!(frame::PauliFrame,noise::UnbiasedUncorrelatedNoise,i::Int)
 end
 
 function applynoise!(frame::PauliFrame,noise::PauliNoise,i::Int)
-    T = eltype(frame.frame.tab.xzs)
-
-    lowbit = T(1)
-    ibig = _div(T,i-1)+1
-    ismall = _mod(T,i-1)
-    ismallm = lowbit<<(ismall)
+    ibig, ismall, ismallm = _get_bitmask_idxs(frame,i)
 
     @inbounds @simd for f in eachindex(frame)
         r = rand()