Merge pull request #122 from MSRudolph/dev

v0.5.2

Author: MSRudolph

Project.toml

@@ -1,7 +1,7 @@
 name = "PauliPropagation"
 uuid = "293282d5-3c99-4fb6-92d0-fd3280a19750"
 authors = ["Manuel S. Rudolph"]
-version = "0.5.1"
+version = "0.5.2"
 
 [deps]
 BitIntegers = "c3b6d118-76ef-56ca-8cc7-ebb389d030a1"

examples/3-utility-example.ipynb

@@ -138,93 +138,55 @@ end
 # XOR between two Pauli different non-identity strings gives the third one. Ignores signs or any coefficient.
 _bitpaulimultiply(pstr1::PauliStringType, pstr2::PauliStringType) = pstr1 ⊻ pstr2
 
-
-
 # Shift to the right and truncate the first encoded Pauli string. Just a utility function.
 _paulishiftright(pstr::PauliStringType) = pstr >> 2
 
+# Computing bit shift index from Pauli site index
+# this is is the amount we need to shift to get to the target Pauli
+_bitshiftfromsiteindex(siteindex::Integer) = 2 * (siteindex - 1)
+
+# a site is represented by two bits
+# using Bits.jl implementation of mask
+_paulimask(::Type{T}, n_sites) where T = mask(T, 2 * n_sites)
 
+_pauliwindowmask(::Type{T}, index1::Integer, index2::Integer) where T = _paulimask(T, index2 - index1 + 1) << _bitshiftfromsiteindex(index1)
 
 # This function extracts the Pauli at position `index` from the integer Pauli string.
 function _getpaulibits(pstr::PauliStringType, index::Integer)
-    # we need to shift the integer by 2 * (index - 1), then the first two bits are target Pauli
-    bitindex = 2 * (index - 1)
-
-    # shift to the right
-    shifted_pstr = (pstr >> bitindex)
-
-    # AND with 3 (00000011) to get the first two bits
-    return shifted_pstr & typeof(pstr)(3)
+    return _getpaulibits(pstr, index, index)
 end
 
+# This function extracts the Pauli from `index1` to `index2`.
+function _getpaulibits(pstr::PauliStringType, index1::Integer, index2::Integer)
+    T = typeof(pstr)
 
-# Gets the bit at index `bitindex` in the integer Pauli string.
-function _getbit(pauli::Integer, bitindex::Integer)
-    # return integer with ...000[bit].
+    bitindex = _bitshiftfromsiteindex(index1)
 
-    # shift by bitindex
-    shifted_pauli = (pauli >> bitindex)
+    # shift to the right
+    shifted_pstr = (pstr >> bitindex)
 
-    # AND with 1 to get first bit
-    return shifted_pauli & typeof(pauli)(1)
+    # creates all 1s mask of length n bits
+    # AND to get the first n bits
+    return shifted_pstr & _paulimask(T, index2 - index1 + 1)
 end
 
 
 # This function sets the Pauli at position `index` in the integer Pauli string to `target_pauli`.
 function _setpaulibits(pstr::PauliStringType, target_pauli::PauliType, index::Integer)
-    # we need to shift the integer by 2 * (index - 1), then the first two bits are target Pauli
-    bitindex = 2 * (index - 1)
-
-    # read bits of the pauli
-    b1 = _getbit(target_pauli, 0)
-    b2 = _getbit(target_pauli, 1)
-
-    # insert them into the pstr
-    pstr = _setbit(pstr, b1, bitindex)
-    pstr = _setbit(pstr, b2, bitindex + 1)
-    return pstr
+    return _setpaulibits(pstr, target_pauli, index, index)
 end
 
 
-# Sets a bit at index `bitindex` in the integer Pauli string to the value of `target_bit`.
-function _setbit(pstr::PauliStringType, target_bit::Integer, bitindex::Integer)
-    # set bit at bitindex to bit
-
-    if target_bit == true  # set to one
-        pstr = _setbittoone(pstr, bitindex)
-    else
-        pstr = _setbittozero(pstr, bitindex)
-    end
-    return pstr
-end
-
-
-# Sets a bit at index `bitindex` in the integer Pauli string to 1.
-function _setbittoone(pstr::Integer, bitindex::Integer)
-    # set bit at bitindex to 1
-
-    # shift ...00100...  to bitindex
-    shifted_onebit = (typeof(pstr)(1) << bitindex)
-
-    # OR with pauli string to make sure that that bit is 1
-    return pstr | shifted_onebit
-end
-
-
-
-# Sets a bit at index `bitindex` in the integer Pauli string to 0.
-function _setbittozero(pstr::Integer, bitindex::Integer)
-    # set bit at bitindex to 0
+# This function sets the Pauli from `index1` to `index2` to `target_pstr`.
+function _setpaulibits(pstr::PauliStringType, target_pstr::PauliStringType, index1::Integer, index2::Integer)
+    T = typeof(pstr)
 
-    # flip all bits
-    pstr = ~pstr
+    bitindex = _bitshiftfromsiteindex(index1)
 
-    # set target bit to one
-    pstr = _setbittoone(pstr, bitindex)
+    window_mask = _pauliwindowmask(T, index1, index2)
 
-    # flip all bits back, only the target bit is 0
-    pstr = ~pstr
-    return pstr
+    # set bits to target pstr
+    return (pstr & ~window_mask) | (T(target_pstr) << bitindex)
 end
 
 
@@ -234,10 +196,10 @@ end
 
     # length is the number of bits in the integer
     n_bits = min(bitsize(pstr), 2_048)  # for max 1024 qubits.
-    mask = zero(pstr)
+    mask = zero(T)
     for ii in 0:(n_bits-1)
         if ii % 2 == 0
-            mask = _setbittoone(mask, ii)
+            mask = mask | (T(1) << ii)
         end
     end
     return mask

src/PauliAlgebra/bitoperations.jl

@@ -178,16 +178,32 @@ end
 Gets the Paulis on indices `qinds` of a `pstr` in the integer representation.
 """
 function getpauli(pstr::PauliStringType, qinds)
-    new_pstr = identitylike(pstr)
-    # Get the Paulis on the indices `qinds`
-    for (ii, qind) in enumerate(qinds)
-        pauli = getpauli(pstr, qind)
-        new_pstr = setpauli(new_pstr, pauli, ii)
+    pstr_new = zero(pstr)
+    for (i, qind) in enumerate(qinds)
+        pair = _getpaulibits(pstr, qind) # get two bits for pauli at qind
+        pstr_new |= (pair << (2 * (i - 1))) # append pair using bitwise OR at 2i
     end
-    return new_pstr
+    return pstr_new
+end
+
+
+"""
+    getpauli(pstr::PauliStringType, qind1::Int, qind2::Int)
+
+Gets the Paulis from `qind1` to `qind2` of a `pstr` in the integer representation.
+This function is useful for extracting a continuous sub-PauliString.
+"""
+function getpauli(pstr::PauliStringType, qind1::Int, qind2::Int)
+    # Get the Paulis on the indices from `qind1` to `qind2`
+    if qind1 > qind2
+        throw(ArgumentError("`qind1` should be less than or equal to `qind2`. Got `qind1=$qind1` and `qind2=$qind2`."))
+    end
+
+    return _getpaulibits(pstr, qind1, qind2)
 
 end
 
+
 """
     setpauli(pstr::PauliStringType, target_pauli::PauliType, index::Integer)
 
@@ -198,6 +214,24 @@ function setpauli(pstr::PauliStringType, target_pauli::PauliType, index::Integer
     return _setpaulibits(pstr, target_pauli, index)
 end
 
+
+"""
+    setpauli(pstr::PauliStringType, target_paulis::PauliStringType, index1::Integer, index2::Integer)
+
+Sets the Paulis from `index1` to `index2` of an integer Pauli string to `target_paulis`.
+"""
+function setpauli(pstr::PauliStringType, target_paulis::PauliStringType, index1::Integer, index2::Integer)
+
+    if index1 > index2
+        throw(ArgumentError("`index1` should be less than or equal to `index2`. Got `index1=$index1` and `index2=$index2`."))
+    end
+
+    # TODO: check that `target_paulis` is of the correct length
+
+    return _setpaulibits(pstr, target_paulis, index1, index2)
+end
+
+
 """
     setpauli(pstr::PauliStringType, target_pauli::Symbol, index::Integer)
 
@@ -291,3 +325,39 @@ function _getprettystr(psum::Dict, nqubits::Int; max_lines=20)
     return str
 
 end
+
+
+# Conversion of PauliString and PauliSum to different coefficient types
+function Base.convert(::Type{PauliString{TT1,CT1}}, pstr::PauliString{TT2,CT2}) where {TT1,TT2,CT1,CT2}
+    if TT1 != TT2
+        throw(ArgumentError("Cannot change term type from $TT2 to $TT1"))
+    end
+    return PauliString(pstr.nqubits, convert(TT1, pstr.term), convert(CT1, pstr.coeff))
+end
+
+function Base.convert(::Type{PauliSum{TT1,CT1}}, psum::PauliSum{TT2,CT2}) where {TT1,TT2,CT1,CT2}
+    if TT1 != TT2
+        throw(ArgumentError("Cannot change term type from $TT2 to $TT1"))
+    end
+    return PauliSum(pstr.nqubits, convert(Dict{TT1,CT1}, psum.terms))
+
+end
+
+# # Examples
+# ```julia
+# convertcoefftype(Float64, PauliString(2, :X, 1, 1+0im))
+# ```
+# """
+function convertcoefftype(::Type{CT1}, pstr::PauliString{TT,CT2}) where {TT,CT1,CT2}
+    return PauliString(pstr.nqubits, pstr.term, convert(CT1, pstr.coeff))
+end
+
+# # Examples
+# ```julia
+# psum = PauliSum(PauliString(2, :X, 1, 1+0im))
+# convertcoefftype(Float64, psum)
+# ```
+# """
+function convertcoefftype(::Type{CT1}, psum::PauliSum{TT,CT2}) where {TT,CT1,CT2}
+    return PauliSum(psum.nqubits, convert(Dict{TT,CT1}, psum.terms))
+end

src/PauliAlgebra/utils.jl

@@ -40,6 +40,7 @@ export
     pauliprod,
     commutes,
     commutator,
+    convertcoefftype,
     getinttype
 
 include("PauliTransferMatrix/PauliTransferMatrix.jl")

src/PauliPropagation.jl

@@ -66,6 +66,10 @@ If |x><x| is a computational basis state, it we compute Tr[psum * |x><x|] = <x|p
 For example, `overlapwithcomputational(psum, [1,2,4])` returns the overlap with `|1101000...>`.
 """
 function overlapwithcomputational(psum::PauliSum, onebitinds)
+    if length(psum) == 0
+        return 0.0
+    end
+
     val = zero(numcoefftype(psum))
     for (pstr, coeff) in psum
         val += tonumber(coeff) * _calcsignwithones(pstr, onebitinds)

src/stateoverlap.jl

@@ -112,6 +112,45 @@ end
 
 end
 
+
+@testset "Get Pauli for a range of indices" begin
+    nq = 12
+    pstr = PauliString(nq, [:X, :Z], [3, 9])
+
+    @test getpauli(pstr.term, range(2, 10)) == getpauli(pstr.term, range(2, 10))
+    @test getpauli(pstr.term, range(2, 10)) == getpauli(pstr.term, 2, 10)
+end
+
+@testset "Set Pauli Continuous" begin
+    nq = 12
+    pstr = PauliString(nq, [:X, :Z], [3, 9])
+    target_pstr = PauliString(nq, [:Y, :X], [4, 5])
+    @test setpauli(pstr.term, target_pstr.term, 2, 10) == setpauli(pstr.term, target_pstr.term, range(2, 10))
+end
+
+
+@testset "Convert PauliString/PauliSum coeff types" begin
+
+    @testset "PauliString" begin
+        pstr = PauliString(2, :X, 1, 1 + 0im)
+        @test convertcoefftype(Float64, pstr) == PauliString(2, :X, 1, 1.)
+
+        pstr = PauliString(2, [:X, :Y], [1, 2], 1)
+        @test convertcoefftype(ComplexF64, pstr) == PauliString(2, [:X, :Y], [1, 2], 1 + 0im)
+    end
+
+    @testset "PauliSum" begin
+        psum = PauliSum(ComplexF64, 2)
+        add!(psum, [:X, :Y], [1, 2], 2 + 0im)
+        add!(psum, :Z, 1, 1 + 0im)
+        expected_psum = PauliSum(2)
+        add!(expected_psum, [:X, :Y], [1, 2], 2)
+        add!(expected_psum, :Z, 1, 1)
+        @test convertcoefftype(Float64, psum) == expected_psum
+    end
+
+end
+
 @testset "Set Pauli for `PauliString` type" begin
 
     # nqubits = 4