From e5e887248dcf83add8947744e3dd746a2935bbff Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 16 Oct 2024 12:13:30 +0100 Subject: [PATCH 01/27] Refactored GeneralState to use NetworkState --- .../general_state/tensor_network_state.py | 445 +++--------------- 1 file changed, 66 insertions(+), 379 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index bedfc0f2..6bfaa6f8 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -33,11 +33,12 @@ try: import cuquantum as cq # type: ignore from cuquantum import cutensornet as cutn # type: ignore + from cuquantum.cutensornet.experimental import NetworkState, NetworkOperator # type: ignore except ImportError: warnings.warn("local settings failed to import cuquantum", ImportWarning) -class GeneralState: +class GeneralState: # TODO: Write it as a context manager so that I can call free() """Wrapper of cuTensorNet object for exact simulations via path optimisation.""" def __init__( @@ -65,13 +66,16 @@ def __init__( loglevel: Internal logger output level. """ self._logger = set_logger("GeneralState", loglevel) - self._lib = libhandle + self._lib = libhandle # TODO: Looks like I could omit the libhandle, should I? libhandle.print_device_properties(self._logger) + options = {"handle": self._lib} + # TODO: Consider supporting scratch_fraction of some form of memory limit # Remove end-of-circuit measurements and keep track of them separately # It also resolves implicit swaps + # TODO: Is there any point in keeping the circuit around? self._circuit, self._measurements = _remove_meas_and_implicit_swaps(circuit) - # Identify each qubit with the index of the bond that represents it in the + # Identify each qubit with the index of the bond that represents it in the # tensor network stored in this GeneralState. Qubits are sorted in increasing # lexicographical order, which is the TKET standard. self._qubit_idx_map = {q: i for i, q in enumerate(sorted(self._circuit.qubits))} @@ -80,35 +84,29 @@ def __init__( dim = 2 # We are always dealing with qubits, not qudits qubits_dims = (dim,) * num_qubits # qubit size self._logger.debug(f"Converting a quantum circuit with {num_qubits} qubits.") - data_type = cq.cudaDataType.CUDA_C_64F # for now let that be hard-coded + data_type = "complex128" # for now let that be hard-coded - self._state = cutn.create_state( - self._lib.handle, cutn.StatePurity.PURE, num_qubits, qubits_dims, data_type - ) - self._gate_tensors = [] + self._state = NetworkState(qubits_dims, dtype=data_type, options=options) + + self._gate_tensors = [] # TODO: Do I still need to keep these myself? - # Append all gates to the TN + # Append all gates to the NetworkState for com in self._circuit.get_commands(): try: gate_unitary = com.op.get_unitary() except: raise ValueError( "All commands in the circuit must be unitary gates. The circuit " - f"contains {com}; no unitary matrix could be retrived for it." + f"contains {com}; no unitary matrix could be retrived from it." ) self._gate_tensors.append(_formatted_tensor(gate_unitary, com.op.n_qubits)) gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) - cutn.state_apply_tensor_operator( - handle=self._lib.handle, - tensor_network_state=self._state, - num_state_modes=com.op.n_qubits, - state_modes=gate_qubit_indices, - tensor_data=self._gate_tensors[-1].data.ptr, - tensor_mode_strides=0, - immutable=1, - adjoint=0, - unitary=1, + tensor_id = self._state.apply_tensor_operator( + modes=gate_qubit_indices, + operand=self._gate_tensors[-1], + immutable=True, # TODO: Change for parameterised gates + unitary=True, ) def get_statevector( @@ -133,112 +131,23 @@ def get_statevector( host device (CPU). Arrays are returned in a 1D shape. """ - #################################### - # Configure the TN for contraction # - #################################### - if attributes is None: - attributes = dict() - attribute_pairs = [ - (getattr(cutn.StateAttribute, k), v) for k, v in attributes.items() - ] - - for attr, val in attribute_pairs: - attr_dtype = cutn.state_get_attribute_dtype(attr) - attr_arr = np.asarray(val, dtype=attr_dtype) - cutn.state_configure( - self._lib.handle, - self._state, - attr, - attr_arr.ctypes.data, - attr_arr.dtype.itemsize, - ) - - try: - ###################################### - # Allocate workspace for contraction # - ###################################### - stream = cp.cuda.Stream() - free_mem = self._lib.dev.mem_info[0] - scratch_size = int(scratch_fraction * free_mem) - scratch_space = cp.cuda.alloc(scratch_size) - self._logger.debug( - f"Allocated {scratch_size} bytes of scratch memory on GPU" - ) - work_desc = cutn.create_workspace_descriptor(self._lib.handle) - - cutn.state_prepare( - self._lib.handle, - self._state, - scratch_size, - work_desc, - stream.ptr, - ) - workspace_size_d = cutn.workspace_get_memory_size( - self._lib.handle, - work_desc, - cutn.WorksizePref.RECOMMENDED, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - ) - - if workspace_size_d <= scratch_size: - cutn.workspace_set_memory( - self._lib.handle, - work_desc, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - scratch_space.ptr, - workspace_size_d, - ) - self._logger.debug( - f"Set {workspace_size_d} bytes of workspace memory out of the" - f" allocated scratch space." - ) - - else: - raise MemoryError( - f"Insufficient workspace size on the GPU device {self._lib.dev.id}" - ) - - ################### - # Contract the TN # - ################### - state_vector = cp.empty( - (2,) * self._circuit.n_qubits, dtype="complex128", order="F" - ) - cutn.state_compute( - self._lib.handle, - self._state, - work_desc, - (state_vector.data.ptr,), - stream.ptr, - ) - stream.synchronize() - sv = state_vector.flatten() - if on_host: - sv = cp.asnumpy(sv) - # Apply the phase from the circuit - sv *= np.exp(1j * np.pi * self._circuit.phase) - return sv - - finally: - cutn.destroy_workspace_descriptor(work_desc) # type: ignore - del scratch_space + self._logger.debug("(Statevector) contracting the TN") + state_vector = self._state.compute_state_vector() + sv = state_vector.flatten() # Convert to 1D + if on_host: + sv = cp.asnumpy(sv) + # Apply the phase from the circuit + sv *= np.exp(1j * np.pi * self._circuit.phase) + return sv def expectation_value( self, operator: QubitPauliOperator, - attributes: Optional[dict] = None, - scratch_fraction: float = 0.75, - ) -> complex: + ) -> float: """Calculates the expectation value of the given operator. Args: operator: The operator whose expectation value is to be measured. - attributes: Optional. A dict of cuTensorNet `ExpectationAttribute` keys - and their values. - scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space. Raises: ValueError: If the operator acts on qubits not present in the circuit. @@ -247,179 +156,39 @@ def expectation_value( The expectation value. """ - ############################################ - # Generate the cuTensorNet operator object # - ############################################ - pauli_tensors = { - "X": _formatted_tensor(np.asarray([[0, 1], [1, 0]]), 1), - "Y": _formatted_tensor(np.asarray([[0, -1j], [1j, 0]]), 1), - "Z": _formatted_tensor(np.asarray([[1, 0], [0, -1]]), 1), - "I": _formatted_tensor(np.asarray([[1, 0], [0, 1]]), 1), - } - num_qubits = self._circuit.n_qubits - qubits_dims = (2,) * num_qubits - data_type = cq.cudaDataType.CUDA_C_64F - - tn_operator = cutn.create_network_operator( - self._lib.handle, num_qubits, qubits_dims, data_type - ) + self._logger.debug("(Expectation value) converting operator to NetworkOperator") - self._logger.debug("Adding operator terms:") - for pauli_string, coeff in operator._dict.items(): - if isinstance(coeff, Expr): - numeric_coeff = complex(coeff.evalf()) # type: ignore - else: - numeric_coeff = complex(coeff) # type: ignore - self._logger.debug(f" {numeric_coeff}, {pauli_string}") + paulis = ["I", "X", "Y", "Z"] + pauli_strs = dict() + for pstr, coeff in operator._dict.items(): # Raise an error if the operator acts on qubits that are not in the circuit - if any(q not in self._circuit.qubits for q in pauli_string.map.keys()): + if any(q not in self._circuit.qubits for q in pstr.map.keys()): raise ValueError( - f"The operator is acting on qubits {pauli_string.map.keys()}, " + f"The operator is acting on qubits {pstr.map.keys()}, " "but some of these are not present in the circuit, whose set of " f"qubits is: {self._circuit.qubits}." ) - # Obtain the tensors corresponding to this operator - qubit_pauli_map = { - q: pauli_tensors[pauli.name] for q, pauli in pauli_string.map.items() - } - - num_pauli = len(qubit_pauli_map) - num_modes = (1,) * num_pauli - state_modes = tuple( - (self._qubit_idx_map[qb],) for qb in qubit_pauli_map.keys() - ) - gate_data = tuple(tensor.data.ptr for tensor in qubit_pauli_map.values()) - - cutn.network_operator_append_product( - handle=self._lib.handle, - tensor_network_operator=tn_operator, - coefficient=numeric_coeff, - num_tensors=num_pauli, - num_state_modes=num_modes, - state_modes=state_modes, - tensor_mode_strides=0, - tensor_data=gate_data, - ) - - ###################################################### - # Configure the cuTensorNet expectation value object # - ###################################################### - expectation = cutn.create_expectation( - self._lib.handle, self._state, tn_operator - ) - - if attributes is None: - attributes = dict() - attribute_pairs = [ - (getattr(cutn.ExpectationAttribute, k), v) for k, v in attributes.items() - ] - - for attr, val in attribute_pairs: - attr_dtype = cutn.expectation_get_attribute_dtype(attr) - attr_arr = np.asarray(val, dtype=attr_dtype) - cutn.expectation_configure( - self._lib.handle, - expectation, - attr, - attr_arr.ctypes.data, - attr_arr.dtype.itemsize, - ) - - try: - ###################################### - # Allocate workspace for contraction # - ###################################### - stream = cp.cuda.Stream() - free_mem = self._lib.dev.mem_info[0] - scratch_size = int(scratch_fraction * free_mem) - scratch_space = cp.cuda.alloc(scratch_size) - - self._logger.debug( - f"Allocated {scratch_size} bytes of scratch memory on GPU" - ) - work_desc = cutn.create_workspace_descriptor(self._lib.handle) - cutn.expectation_prepare( - self._lib.handle, - expectation, - scratch_size, - work_desc, - stream.ptr, - ) - workspace_size_d = cutn.workspace_get_memory_size( - self._lib.handle, - work_desc, - cutn.WorksizePref.RECOMMENDED, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - ) - - if workspace_size_d <= scratch_size: - cutn.workspace_set_memory( - self._lib.handle, - work_desc, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - scratch_space.ptr, - workspace_size_d, - ) - self._logger.debug( - f"Set {workspace_size_d} bytes of workspace memory out of the" - f" allocated scratch space." - ) - else: - raise MemoryError( - f"Insufficient workspace size on the GPU device {self._lib.dev.id}" - ) - - ################################# - # Compute the expectation value # - ################################# - expectation_value = np.empty(1, dtype="complex128") - state_norm = np.empty(1, dtype="complex128") - cutn.expectation_compute( - self._lib.handle, - expectation, - work_desc, - expectation_value.ctypes.data, - state_norm.ctypes.data, - stream.ptr, - ) - stream.synchronize() + pauli_list = [pstr[q] for q in self._qubit_idx_map.keys()] + this_pauli_string = "".join(map(lambda x: paulis[x], pauli_list)) + pauli_strs[this_pauli_string] = complex(coeff) - # Note: we can also return `state_norm.item()`, but this should be 1 since - # we are always running unitary circuits - assert np.isclose(state_norm.item(), 1.0) + tn_operator = NetworkOperator.from_pauli_strs(pauli_strs, dtype="complex128") - return expectation_value.item() # type: ignore + self._logger.debug("(Expectation value) contracting the TN") + return self._state.compute_expectation(tn_operator).real - finally: - ##################################################### - # Destroy the Operator and ExpectationValue objects # - ##################################################### - cutn.destroy_workspace_descriptor(work_desc) # type: ignore - cutn.destroy_expectation(expectation) - cutn.destroy_network_operator(tn_operator) - del scratch_space - - def sample( + def sample( # TODO: Support seeds (and test) self, n_shots: int, - attributes: Optional[dict] = None, - scratch_fraction: float = 0.75, ) -> BackendResult: """Obtains samples from the measurements at the end of the circuit. Args: n_shots: The number of samples to obtain. - attributes: Optional. A dict of cuTensorNet `SamplerAttribute` keys and - their values. - scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space. Raises: ValueError: If the circuit contains no measurements. - MemoryError: If there is insufficient workspace on GPU. Returns: A pytket ``BackendResult`` with the data from the shots. """ @@ -435,113 +204,33 @@ def sample( qbit_list, cbit_list = zip(*self._measurements.items()) measured_modes = tuple(self._qubit_idx_map[qb] for qb in qbit_list) - ############################################ - # Configure the cuTensorNet sampler object # - ############################################ - - sampler = cutn.create_sampler( - handle=self._lib.handle, - tensor_network_state=self._state, - num_modes_to_sample=num_measurements, - modes_to_sample=measured_modes, + self._logger.debug("(Sampling) contracting the TN") + samples = self._state.compute_sampling( + nshots=n_shots, + modes=measured_modes, ) - if attributes is None: - attributes = dict() - attribute_pairs = [ - (getattr(cutn.SamplerAttribute, k), v) for k, v in attributes.items() - ] - - for attr, val in attribute_pairs: - attr_dtype = cutn.sampler_get_attribute_dtype(attr) - attr_arr = np.asarray(val, dtype=attr_dtype) - cutn.sampler_configure( - self._lib.handle, - sampler, - attr, - attr_arr.ctypes.data, - attr_arr.dtype.itemsize, - ) + # Convert the data in `samples` to an `OutcomeArray` using `from_readouts` + # which expects a 2D array `samples[QubitId][SampleId]` of 0s and 1s. + self._logger.debug("(Sampling) converting samples to pytket Backend") + readouts = np.empty(shape=(num_measurements, n_shots), dtype=int) + sample_id = 0 + for bitstring, count in samples.items: + outcome = [int(b) for b in bitstring] - try: - ###################################### - # Allocate workspace for contraction # - ###################################### - stream = cp.cuda.Stream() - free_mem = self._lib.dev.mem_info[0] - scratch_size = int(scratch_fraction * free_mem) - scratch_space = cp.cuda.alloc(scratch_size) - - self._logger.debug( - f"Allocated {scratch_size} bytes of scratch memory on GPU" - ) - work_desc = cutn.create_workspace_descriptor(self._lib.handle) - cutn.sampler_prepare( - self._lib.handle, - sampler, - scratch_size, - work_desc, - stream.ptr, - ) - workspace_size_d = cutn.workspace_get_memory_size( - self._lib.handle, - work_desc, - cutn.WorksizePref.RECOMMENDED, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - ) + for _ in range(count): + readouts[sample_id] = outcome # TODO: test endian-ness + sample_id += 1 - if workspace_size_d <= scratch_size: - cutn.workspace_set_memory( - self._lib.handle, - work_desc, - cutn.Memspace.DEVICE, - cutn.WorkspaceKind.SCRATCH, - scratch_space.ptr, - workspace_size_d, - ) - self._logger.debug( - f"Set {workspace_size_d} bytes of workspace memory out of the" - f" allocated scratch space." - ) - else: - raise MemoryError( - f"Insufficient workspace size on the GPU device {self._lib.dev.id}" - ) + shots = OutcomeArray.from_readouts(readouts) - ########################### - # Sample from the circuit # - ########################### - samples = np.empty((num_measurements, n_shots), dtype="int64", order="F") - cutn.sampler_sample( - self._lib.handle, - sampler, - n_shots, - work_desc, - samples.ctypes.data, - stream.ptr, - ) - stream.synchronize() - - # Convert the data in `samples` to an `OutcomeArray` - # `samples` is a 2D numpy array `samples[SampleId][QubitId]`, which is - # the transpose of what `OutcomeArray.from_readouts` expects - shots = OutcomeArray.from_readouts(samples.T) - # We need to specify which bits correspond to which columns in the shots - # table. Since cuTensorNet promises that the ordering of outcomes is - # determined by the ordering we provided as `measured_modes`, which in - # turn corresponds to the ordering of qubits in `qbit_list`, the fact that - # `cbit_list` has the appropriate order in relation to `self._measurements` - # determines this defines the ordering of classical bits we intend. - return BackendResult(c_bits=cbit_list, shots=shots) - - finally: - ############################## - # Destroy the Sampler object # - ############################## - cutn.destroy_workspace_descriptor(work_desc) # type: ignore - cutn.destroy_sampler(sampler) - del scratch_space + # We need to specify which bits correspond to which columns in the shots + # table. Since cuTensorNet promises that the ordering of outcomes is + # determined by the ordering we provided as `measured_modes`, which in + # turn corresponds to the ordering of qubits in `qbit_list`, the fact that + # `cbit_list` has the appropriate order in relation to `self._measurements` + # implies this defines the ordering of classical bits we intend. + return BackendResult(c_bits=cbit_list, shots=shots) def destroy(self) -> None: """Destroy the tensor network and free up GPU memory. @@ -551,19 +240,17 @@ def destroy(self) -> None: `GeneralState` object. GPU memory deallocation is not guaranteed otherwise. """ - cutn.destroy_state(self._state) + self._logger.debug("Freeing memory of NetworkState") + self._state.free() def _formatted_tensor(matrix: NDArray, n_qubits: int) -> cp.ndarray: """Convert a matrix to the tensor format accepted by NVIDIA's API.""" - # Transpose is needed because of the way cuTN stores tensors. - # See https://github.com/NVIDIA/cuQuantum/discussions/124 - # #discussioncomment-8683146 for details. - cupy_matrix = cp.asarray(matrix).T.astype(dtype="complex128", order="F") + cupy_matrix = cp.asarray(matrix).astype(dtype="complex128") # We also need to reshape since a matrix only has 2 bonds, but for an # n-qubit gate we want 2^n bonds for input and another 2^n for output - return cupy_matrix.reshape([2] * (2 * n_qubits), order="F") + return cupy_matrix.reshape([2] * (2 * n_qubits)) def _remove_meas_and_implicit_swaps(circ: Circuit) -> Tuple[Circuit, Dict[Qubit, Bit]]: From 03a0dc325cd4cdddabea6bae8bb922437cd530d6 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 16 Oct 2024 15:38:14 +0100 Subject: [PATCH 02/27] Completed refactor of GeneralState to use NetworkState --- .../general_state/tensor_network_state.py | 44 +++++------ tests/test_general_state.py | 74 +++++++++---------- 2 files changed, 57 insertions(+), 61 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 6bfaa6f8..ffb65ef5 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -37,14 +37,13 @@ except ImportError: warnings.warn("local settings failed to import cuquantum", ImportWarning) - +# TODO: Add the options as argument to be passed to NetworkState class GeneralState: # TODO: Write it as a context manager so that I can call free() """Wrapper of cuTensorNet object for exact simulations via path optimisation.""" def __init__( self, circuit: Circuit, - libhandle: CuTensorNetHandle, loglevel: int = logging.WARNING, ) -> None: """Constructs a tensor network for the output state of a pytket circuit. @@ -52,23 +51,14 @@ def __init__( The qubits are assumed to be initialised in the ``|0>`` state. The resulting object stores the *uncontracted* tensor network. - Note: - A ``libhandle`` is created via a ``with CuTensorNetHandle() as libhandle:`` - statement. The device where the ``GeneralState`` is stored will match the - one specified by the library handle. - Note: The ``circuit`` must not contain any ``CircBox`` or non-unitary command. Args: circuit: A pytket circuit to be converted to a tensor network. - libhandle: An instance of a ``CuTensorNetHandle``. loglevel: Internal logger output level. """ self._logger = set_logger("GeneralState", loglevel) - self._lib = libhandle # TODO: Looks like I could omit the libhandle, should I? - libhandle.print_device_properties(self._logger) - options = {"handle": self._lib} # TODO: Consider supporting scratch_fraction of some form of memory limit # Remove end-of-circuit measurements and keep track of them separately @@ -86,12 +76,13 @@ def __init__( self._logger.debug(f"Converting a quantum circuit with {num_qubits} qubits.") data_type = "complex128" # for now let that be hard-coded - self._state = NetworkState(qubits_dims, dtype=data_type, options=options) + self._state = NetworkState(qubits_dims, dtype=data_type) self._gate_tensors = [] # TODO: Do I still need to keep these myself? + commands = self._circuit.get_commands() # Append all gates to the NetworkState - for com in self._circuit.get_commands(): + for com in commands: try: gate_unitary = com.op.get_unitary() except: @@ -103,12 +94,23 @@ def __init__( gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) tensor_id = self._state.apply_tensor_operator( - modes=gate_qubit_indices, - operand=self._gate_tensors[-1], + gate_qubit_indices, + self._gate_tensors[-1], immutable=True, # TODO: Change for parameterised gates unitary=True, ) + # If the circuit has no gates, apply one identity gate so that CuTensorNet does not panic + # due to no tensor operator in the NetworkState + if len(commands) == 0: + identity_tensor = _formatted_tensor(np.identity(2, dtype="complex128"), 1) + tensor_id = self._state.apply_tensor_operator( + (0, ), + identity_tensor, + immutable=True, + unitary=True, + ) + def get_statevector( self, attributes: Optional[dict] = None, @@ -143,7 +145,7 @@ def get_statevector( def expectation_value( self, operator: QubitPauliOperator, - ) -> float: + ) -> complex: """Calculates the expectation value of the given operator. Args: @@ -174,10 +176,10 @@ def expectation_value( this_pauli_string = "".join(map(lambda x: paulis[x], pauli_list)) pauli_strs[this_pauli_string] = complex(coeff) - tn_operator = NetworkOperator.from_pauli_strs(pauli_strs, dtype="complex128") + tn_operator = NetworkOperator.from_pauli_strings(pauli_strs, dtype="complex128") self._logger.debug("(Expectation value) contracting the TN") - return self._state.compute_expectation(tn_operator).real + return self._state.compute_expectation(tn_operator) def sample( # TODO: Support seeds (and test) self, @@ -211,11 +213,11 @@ def sample( # TODO: Support seeds (and test) ) # Convert the data in `samples` to an `OutcomeArray` using `from_readouts` - # which expects a 2D array `samples[QubitId][SampleId]` of 0s and 1s. + # which expects a 2D array `samples[SampleId][QubitId]` of 0s and 1s. self._logger.debug("(Sampling) converting samples to pytket Backend") - readouts = np.empty(shape=(num_measurements, n_shots), dtype=int) + readouts = np.empty(shape=(n_shots, num_measurements), dtype=int) sample_id = 0 - for bitstring, count in samples.items: + for bitstring, count in samples.items(): outcome = [int(b) for b in bitstring] for _ in range(count): diff --git a/tests/test_general_state.py b/tests/test_general_state.py index f3aadf5b..ff96c2dd 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -7,7 +7,6 @@ from pytket.pauli import QubitPauliString, Pauli from pytket.utils.operators import QubitPauliOperator from pytket.extensions.cutensornet.general_state import GeneralState -from pytket.extensions.cutensornet.structured_state import CuTensorNetHandle @pytest.mark.parametrize( @@ -39,24 +38,23 @@ ], ) def test_get_statevec(circuit: Circuit) -> None: - with CuTensorNetHandle() as libhandle: - state = GeneralState(circuit, libhandle) - sv = state.get_statevector() + state = GeneralState(circuit) + sv = state.get_statevector() - sv_pytket = circuit.get_statevector() - assert np.allclose(sv, sv_pytket, atol=1e-10) + sv_pytket = circuit.get_statevector() + assert np.allclose(sv, sv_pytket, atol=1e-10) - op = QubitPauliOperator( - { - QubitPauliString({q: Pauli.I for q in circuit.qubits}): 1.0, - } - ) + op = QubitPauliOperator( + { + QubitPauliString({q: Pauli.I for q in circuit.qubits}): 1.0, + } + ) - # Calculate the inner product as the expectation value - # of the identity operator: = - state = GeneralState(circuit, libhandle) - ovl = state.expectation_value(op) - assert ovl == pytest.approx(1.0) + # Calculate the inner product as the expectation value + # of the identity operator: = + state = GeneralState(circuit) + ovl = state.expectation_value(op) + assert ovl == pytest.approx(1.0) state.destroy() @@ -83,9 +81,8 @@ def test_sv_toffoli_box_with_implicit_swaps() -> None: Transform.OptimiseCliffords().apply(ket_circ) # Convert and contract - with CuTensorNetHandle() as libhandle: - state = GeneralState(ket_circ, libhandle) - ket_net_vector = state.get_statevector() + state = GeneralState(ket_circ) + ket_net_vector = state.get_statevector() state.destroy() # Compare to pytket statevector @@ -119,24 +116,23 @@ def to_bool_tuple(n_qubits: int, x: int) -> tuple: CnXPairwiseDecomposition().apply(ket_circ) Transform.OptimiseCliffords().apply(ket_circ) - with CuTensorNetHandle() as libhandle: - state = GeneralState(ket_circ, libhandle) - ket_net_vector = state.get_statevector() + state = GeneralState(ket_circ) + ket_net_vector = state.get_statevector() - ket_pytket_vector = ket_circ.get_statevector() - assert np.allclose(ket_net_vector, ket_pytket_vector) + ket_pytket_vector = ket_circ.get_statevector() + assert np.allclose(ket_net_vector, ket_pytket_vector) - # Calculate the inner product as the expectation value - # of the identity operator: = - op = QubitPauliOperator( - { - QubitPauliString({q: Pauli.I for q in ket_circ.qubits}): 1.0, - } - ) + # Calculate the inner product as the expectation value + # of the identity operator: = + op = QubitPauliOperator( + { + QubitPauliString({q: Pauli.I for q in ket_circ.qubits}): 1.0, + } + ) - state = GeneralState(ket_circ, libhandle) - ovl = state.expectation_value(op) - assert ovl == pytest.approx(1.0) + state = GeneralState(ket_circ) + ovl = state.expectation_value(op) + assert ovl == pytest.approx(1.0) state.destroy() @@ -197,9 +193,8 @@ def test_expectation_value(circuit: Circuit, observable: QubitPauliOperator) -> exp_val_tket = observable.state_expectation(circuit.get_statevector()) # Calculate using GeneralState - with CuTensorNetHandle() as libhandle: - state = GeneralState(circuit, libhandle) - exp_val = state.expectation_value(observable) + state = GeneralState(circuit) + exp_val = state.expectation_value(observable) assert np.isclose(exp_val, exp_val_tket) state.destroy() @@ -253,9 +248,8 @@ def test_sampler(circuit: Circuit, measure_all: bool) -> None: circuit.Measure(q, Bit(i)) # Sample using our library - with CuTensorNetHandle() as libhandle: - state = GeneralState(circuit, libhandle) - results = state.sample(n_shots) + state = GeneralState(circuit) + results = state.sample(n_shots) # Verify distribution matches theoretical probabilities for bit_tuple, count in results.get_counts().items(): From 432e7c4d6025a0890eab1755a17a78a59570032c Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 16 Oct 2024 15:48:05 +0100 Subject: [PATCH 03/27] Updated backend --- .../backends/cutensornet_backend.py | 35 +++++++++---------- 1 file changed, 16 insertions(+), 19 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index b09d33c0..5d01afdb 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -24,7 +24,6 @@ from pytket.backends.backend import KwargTypes, Backend, BackendResult from pytket.backends.backendinfo import BackendInfo from pytket.backends.resulthandle import _ResultIdTuple -from pytket.extensions.cutensornet.general import CuTensorNetHandle from pytket.extensions.cutensornet.general_state import ( GeneralState, ) @@ -227,16 +226,15 @@ def process_circuits( if valid_check: self._check_all_circuits(circuit_list) handle_list = [] - with CuTensorNetHandle() as libhandle: - for circuit in circuit_list: - tn = GeneralState(circuit, libhandle) - sv = tn.get_statevector(attributes, scratch_fraction) - res_qubits = [qb for qb in sorted(circuit.qubits)] - handle = ResultHandle(str(uuid4())) - self._cache[handle] = { - "result": BackendResult(q_bits=res_qubits, state=sv) - } - handle_list.append(handle) + for circuit in circuit_list: + tn = GeneralState(circuit) + sv = tn.get_statevector(attributes, scratch_fraction) + res_qubits = [qb for qb in sorted(circuit.qubits)] + handle = ResultHandle(str(uuid4())) + self._cache[handle] = { + "result": BackendResult(q_bits=res_qubits, state=sv) + } + handle_list.append(handle) return handle_list @@ -321,14 +319,13 @@ def process_circuits( if valid_check: self._check_all_circuits(circuit_list) handle_list = [] - with CuTensorNetHandle() as libhandle: - for circuit, circ_shots in zip(circuit_list, all_shots): - tn = GeneralState(circuit, libhandle) - handle = ResultHandle(str(uuid4())) - self._cache[handle] = { - "result": tn.sample(circ_shots, attributes, scratch_fraction) - } - handle_list.append(handle) + for circuit, circ_shots in zip(circuit_list, all_shots): + tn = GeneralState(circuit) + handle = ResultHandle(str(uuid4())) + self._cache[handle] = { + "result": tn.sample(circ_shots) # TODO: recover use of (attributes, scratch_fraction) + } + handle_list.append(handle) return handle_list From d205538d467b991757f5b79899be8b8c22fd7943 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 16 Oct 2024 15:48:41 +0100 Subject: [PATCH 04/27] Cleaning up a bit --- .../cutensornet/general_state/tensor_network_state.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index ffb65ef5..a3b69146 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -78,7 +78,6 @@ def __init__( self._state = NetworkState(qubits_dims, dtype=data_type) - self._gate_tensors = [] # TODO: Do I still need to keep these myself? commands = self._circuit.get_commands() # Append all gates to the NetworkState @@ -90,12 +89,11 @@ def __init__( "All commands in the circuit must be unitary gates. The circuit " f"contains {com}; no unitary matrix could be retrived from it." ) - self._gate_tensors.append(_formatted_tensor(gate_unitary, com.op.n_qubits)) gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) tensor_id = self._state.apply_tensor_operator( gate_qubit_indices, - self._gate_tensors[-1], + _formatted_tensor(gate_unitary, com.op.n_qubits), immutable=True, # TODO: Change for parameterised gates unitary=True, ) @@ -221,7 +219,7 @@ def sample( # TODO: Support seeds (and test) outcome = [int(b) for b in bitstring] for _ in range(count): - readouts[sample_id] = outcome # TODO: test endian-ness + readouts[sample_id] = outcome sample_id += 1 shots = OutcomeArray.from_readouts(readouts) From 535b4b238eb173eb423fb1640c193a2a3044ee1c Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 16 Oct 2024 16:52:02 +0100 Subject: [PATCH 05/27] Supporting parameters on GeneralState and Backend --- .../backends/cutensornet_backend.py | 43 ++++++++----------- .../general_state/tensor_network_state.py | 23 +++++----- 2 files changed, 31 insertions(+), 35 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index 5d01afdb..510af630 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -45,11 +45,6 @@ CustomPass, ) -try: - from cuquantum.cutensornet import StateAttribute, SamplerAttribute # type: ignore -except ImportError: - warnings.warn("local settings failed to import cuquantum", ImportWarning) - from .._metadata import __extension_version__, __extension_name__ @@ -202,33 +197,33 @@ def process_circuits( corresponding get_ method. Note: - Any element from the ``StateAttribute`` enum (see NVIDIA's CuTensorNet + Any element from the ``TNConfig`` enum (see NVIDIA's CuTensorNet API) can be provided as arguments to this method. For instance: - ``process_circuits(..., CONFIG_NUM_HYPER_SAMPLES=100)``. + ``process_circuits(..., tn_config={"num_hyper_samples": 100})``. Args: circuits: List of circuits to be submitted. n_shots: Number of shots in case of shot-based calculation. This should be ``None``, since this backend does not support shots. valid_check: Whether to check for circuit correctness. + tnconfig: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space. Defaults to `0.75`. + scratch space; value between 0 and 1. Defaults to ``0.8``. Returns: Results handle objects. """ - scratch_fraction = float(kwargs.get("scratch_fraction", 0.75)) # type: ignore - attributes = { - k: v for k, v in kwargs.items() if k in StateAttribute._member_names_ - } + scratch_fraction = float(kwargs.get("scratch_fraction", 0.8)) # type: ignore + tnconfig = kwargs.get("tnconfig", dict()) # type: ignore circuit_list = list(circuits) if valid_check: self._check_all_circuits(circuit_list) handle_list = [] for circuit in circuit_list: - tn = GeneralState(circuit) - sv = tn.get_statevector(attributes, scratch_fraction) + tn = GeneralState(circuit, attributes=tnconfig, scratch_fraction=scratch_fraction) + sv = tn.get_statevector() res_qubits = [qb for qb in sorted(circuit.qubits)] handle = ResultHandle(str(uuid4())) self._cache[handle] = { @@ -278,9 +273,9 @@ def process_circuits( corresponding get_ method. Note: - Any element from the ``SamplerAttribute`` enum (see NVIDIA's CuTensorNet + Any element from the ``TNConfig`` enum (see NVIDIA's CuTensorNet API) can be provided as arguments to this method. For instance: - ``process_circuits(..., CONFIG_NUM_HYPER_SAMPLES=100)``. + ``process_circuits(..., tn_config={"num_hyper_samples": 100})``. Args: circuits: List of circuits to be submitted. @@ -288,21 +283,21 @@ def process_circuits( Optionally, this can be a list of shots specifying the number of shots for each circuit separately. valid_check: Whether to check for circuit correctness. + tnconfig: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space. Defaults to `0.75`. + scratch space; value between 0 and 1. Defaults to ``0.8``. Returns: Results handle objects. """ - scratch_fraction = float(kwargs.get("scratch_fraction", 0.75)) # type: ignore - attributes = { - k: v for k, v in kwargs.items() if k in SamplerAttribute._member_names_ - } + scratch_fraction = float(kwargs.get("scratch_fraction", 0.8)) # type: ignore + tnconfig = kwargs.get("tnconfig", dict()) # type: ignore if "seed" in kwargs and kwargs["seed"] is not None: # Current CuTensorNet does not support seeds for Sampler. I created # a feature request in their repository. - raise NotImplementedError( + raise NotImplementedError( # TODO: Support seeds! "The backend does not currently support user-defined seeds." ) @@ -320,10 +315,10 @@ def process_circuits( self._check_all_circuits(circuit_list) handle_list = [] for circuit, circ_shots in zip(circuit_list, all_shots): - tn = GeneralState(circuit) + tn = GeneralState(circuit, attributes=tnconfig, scratch_fraction=scratch_fraction) handle = ResultHandle(str(uuid4())) self._cache[handle] = { - "result": tn.sample(circ_shots) # TODO: recover use of (attributes, scratch_fraction) + "result": tn.sample(circ_shots) } handle_list.append(handle) return handle_list diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index a3b69146..caedb824 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -37,13 +37,14 @@ except ImportError: warnings.warn("local settings failed to import cuquantum", ImportWarning) -# TODO: Add the options as argument to be passed to NetworkState class GeneralState: # TODO: Write it as a context manager so that I can call free() """Wrapper of cuTensorNet object for exact simulations via path optimisation.""" def __init__( self, circuit: Circuit, + attributes: Optional[dict] = None, + scratch_fraction: float = 0.8, loglevel: int = logging.WARNING, ) -> None: """Constructs a tensor network for the output state of a pytket circuit. @@ -56,10 +57,13 @@ def __init__( Args: circuit: A pytket circuit to be converted to a tensor network. + attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. + scratch_fraction: Optional. Fraction of free memory on GPU to allocate as + scratch space; value between 0 and 1. Defaults to ``0.8``. loglevel: Internal logger output level. """ self._logger = set_logger("GeneralState", loglevel) - # TODO: Consider supporting scratch_fraction of some form of memory limit # Remove end-of-circuit measurements and keep track of them separately # It also resolves implicit swaps @@ -76,7 +80,12 @@ def __init__( self._logger.debug(f"Converting a quantum circuit with {num_qubits} qubits.") data_type = "complex128" # for now let that be hard-coded - self._state = NetworkState(qubits_dims, dtype=data_type) + self._state = NetworkState( + qubits_dims, + dtype=data_type, + config=attributes, + options={"memory_limit": f"{int(scratch_fraction*100)}%"}, + ) commands = self._circuit.get_commands() @@ -111,21 +120,13 @@ def __init__( def get_statevector( self, - attributes: Optional[dict] = None, - scratch_fraction: float = 0.75, on_host: bool = True, ) -> Union[cp.ndarray, np.ndarray]: """Contracts the circuit and returns the final statevector. Args: - attributes: Optional. A dict of cuTensorNet `StateAttribute` keys and - their values. - scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space. on_host: Optional. If ``True``, converts cupy ``ndarray`` to numpy ``ndarray``, copying it to host device (CPU). - Raises: - MemoryError: If there is insufficient workspace on GPU. Returns: Either a ``cupy.ndarray`` on a GPU, or a ``numpy.ndarray`` on a host device (CPU). Arrays are returned in a 1D shape. From 11c01d9872437178e8222bd51d65bfc29acd01c2 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Mon, 21 Oct 2024 13:39:31 +0100 Subject: [PATCH 06/27] Adding get_amplitude --- .../general_state/tensor_network_state.py | 31 ++++++++++++++++--- tests/test_general_state.py | 4 +++ 2 files changed, 31 insertions(+), 4 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index caedb824..dfb59967 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -37,6 +37,7 @@ except ImportError: warnings.warn("local settings failed to import cuquantum", ImportWarning) + class GeneralState: # TODO: Write it as a context manager so that I can call free() """Wrapper of cuTensorNet object for exact simulations via path optimisation.""" @@ -112,7 +113,7 @@ def __init__( if len(commands) == 0: identity_tensor = _formatted_tensor(np.identity(2, dtype="complex128"), 1) tensor_id = self._state.apply_tensor_operator( - (0, ), + (0,), identity_tensor, immutable=True, unitary=True, @@ -131,7 +132,6 @@ def get_statevector( Either a ``cupy.ndarray`` on a GPU, or a ``numpy.ndarray`` on a host device (CPU). Arrays are returned in a 1D shape. """ - self._logger.debug("(Statevector) contracting the TN") state_vector = self._state.compute_state_vector() sv = state_vector.flatten() # Convert to 1D @@ -141,6 +141,30 @@ def get_statevector( sv *= np.exp(1j * np.pi * self._circuit.phase) return sv + def get_amplitude( + self, + state: int, + ) -> complex: + """Returns the amplitude of the chosen computational state. + + Notes: + The result is equivalent to ``state.get_statevector[b]``, but this method + is faster when querying a single amplitude (or just a few). + + Args: + state: The integer whose bitstring describes the computational state. + The qubits in the bitstring are in increasing lexicographic order. + + Returns: + The amplitude of the computational state in ``self``. + """ + self._logger.debug("(Amplitude) contracting the TN") + bitstring = format(state, "b").zfill(self._circuit.n_qubits) + amplitude = self._state.compute_amplitude(bitstring) + # Apply the phase from the circuit + amplitude *= np.exp(1j * np.pi * self._circuit.phase) + return complex(amplitude) + def expectation_value( self, operator: QubitPauliOperator, @@ -156,7 +180,6 @@ def expectation_value( Returns: The expectation value. """ - self._logger.debug("(Expectation value) converting operator to NetworkOperator") paulis = ["I", "X", "Y", "Z"] @@ -178,7 +201,7 @@ def expectation_value( tn_operator = NetworkOperator.from_pauli_strings(pauli_strs, dtype="complex128") self._logger.debug("(Expectation value) contracting the TN") - return self._state.compute_expectation(tn_operator) + return complex(self._state.compute_expectation(tn_operator)) def sample( # TODO: Support seeds (and test) self, diff --git a/tests/test_general_state.py b/tests/test_general_state.py index ff96c2dd..9217fc6e 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -56,6 +56,10 @@ def test_get_statevec(circuit: Circuit) -> None: ovl = state.expectation_value(op) assert ovl == pytest.approx(1.0) + # Check that all amplitudes agree + for i in range(len(sv)): + assert sv[i] == state.get_amplitude(i) + state.destroy() From 5856f5387ace001aafabc9a6aee8e747438376b6 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Mon, 21 Oct 2024 13:39:45 +0100 Subject: [PATCH 07/27] Linter --- .../cutensornet/backends/cutensornet_backend.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index 510af630..fe24f881 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -222,13 +222,13 @@ def process_circuits( self._check_all_circuits(circuit_list) handle_list = [] for circuit in circuit_list: - tn = GeneralState(circuit, attributes=tnconfig, scratch_fraction=scratch_fraction) + tn = GeneralState( + circuit, attributes=tnconfig, scratch_fraction=scratch_fraction + ) sv = tn.get_statevector() res_qubits = [qb for qb in sorted(circuit.qubits)] handle = ResultHandle(str(uuid4())) - self._cache[handle] = { - "result": BackendResult(q_bits=res_qubits, state=sv) - } + self._cache[handle] = {"result": BackendResult(q_bits=res_qubits, state=sv)} handle_list.append(handle) return handle_list @@ -315,11 +315,11 @@ def process_circuits( self._check_all_circuits(circuit_list) handle_list = [] for circuit, circ_shots in zip(circuit_list, all_shots): - tn = GeneralState(circuit, attributes=tnconfig, scratch_fraction=scratch_fraction) + tn = GeneralState( + circuit, attributes=tnconfig, scratch_fraction=scratch_fraction + ) handle = ResultHandle(str(uuid4())) - self._cache[handle] = { - "result": tn.sample(circ_shots) - } + self._cache[handle] = {"result": tn.sample(circ_shots)} handle_list.append(handle) return handle_list From 3f94bfa5f1f98b466f37ecc11eba64d968522af8 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Mon, 21 Oct 2024 16:48:01 +0000 Subject: [PATCH 08/27] Fixed test --- tests/test_general_state.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_general_state.py b/tests/test_general_state.py index 9217fc6e..7becf16f 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -58,7 +58,7 @@ def test_get_statevec(circuit: Circuit) -> None: # Check that all amplitudes agree for i in range(len(sv)): - assert sv[i] == state.get_amplitude(i) + assert np.isclose(sv[i], state.get_amplitude(i)) state.destroy() From c30efa83237a94085f7df0c0b007070a2f55ca20 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Tue, 22 Oct 2024 12:57:47 +0000 Subject: [PATCH 09/27] Implemented GeneralBraOpKet to be used to compute arbitrary --- .../cutensornet/general_state/__init__.py | 2 +- .../general_state/tensor_network_state.py | 250 +++++++++++++++++- tests/test_general_state.py | 17 +- 3 files changed, 255 insertions(+), 14 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/__init__.py b/pytket/extensions/cutensornet/general_state/__init__.py index 04c3357b..bb5dadea 100644 --- a/pytket/extensions/cutensornet/general_state/__init__.py +++ b/pytket/extensions/cutensornet/general_state/__init__.py @@ -25,4 +25,4 @@ get_circuit_overlap, ) -from .tensor_network_state import GeneralState +from .tensor_network_state import GeneralState, GeneralBraOpKet diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index dfb59967..37d9a183 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -24,8 +24,8 @@ import numpy as np from sympy import Expr # type: ignore from numpy.typing import NDArray -from pytket.circuit import Circuit, Qubit, Bit, OpType -from pytket.extensions.cutensornet.general import CuTensorNetHandle, set_logger +from pytket.circuit import Circuit, Qubit, Bit, OpType, Op +from pytket.extensions.cutensornet.general import set_logger from pytket.utils import OutcomeArray from pytket.utils.operators import QubitPauliOperator from pytket.backends.backendresult import BackendResult @@ -39,7 +39,7 @@ class GeneralState: # TODO: Write it as a context manager so that I can call free() - """Wrapper of cuTensorNet object for exact simulations via path optimisation.""" + """Wrapper of cuTensorNet's NetworkState for exact simulation of states.""" def __init__( self, @@ -78,14 +78,16 @@ def __init__( num_qubits = self._circuit.n_qubits dim = 2 # We are always dealing with qubits, not qudits qubits_dims = (dim,) * num_qubits # qubit size - self._logger.debug(f"Converting a quantum circuit with {num_qubits} qubits.") data_type = "complex128" # for now let that be hard-coded self._state = NetworkState( qubits_dims, dtype=data_type, config=attributes, - options={"memory_limit": f"{int(scratch_fraction*100)}%"}, + options={ + "memory_limit": f"{int(scratch_fraction*100)}%", + "logger": self._logger, + }, ) commands = self._circuit.get_commands() @@ -108,13 +110,12 @@ def __init__( unitary=True, ) - # If the circuit has no gates, apply one identity gate so that CuTensorNet does not panic - # due to no tensor operator in the NetworkState + # If the circuit has no gates, apply one identity gate so that CuTensorNet does + # not panic due to no tensor operator in the NetworkState if len(commands) == 0: - identity_tensor = _formatted_tensor(np.identity(2, dtype="complex128"), 1) tensor_id = self._state.apply_tensor_operator( (0,), - identity_tensor, + _formatted_tensor(np.identity(2), 1), immutable=True, unitary=True, ) @@ -264,14 +265,241 @@ def destroy(self) -> None: `GeneralState` object. GPU memory deallocation is not guaranteed otherwise. """ - self._logger.debug("Freeing memory of NetworkState") + self._logger.debug("Freeing memory of GeneralState") self._state.free() +class GeneralBraOpKet: # TODO: Write it as a context manager + """Wrapper of cuTensorNet's NetworkState for exact simulation of ````.""" + + def __init__( + self, + bra: Circuit, + ket: Circuit, + attributes: Optional[dict] = None, + scratch_fraction: float = 0.8, + loglevel: int = logging.WARNING, + ) -> None: + """Constructs a tensor network for ````. + + The qubits in ``ket`` and ``bra`` are assumed to be initialised in the ``|0>`` + state. The resulting object stores the *uncontracted* tensor network. + + Note: + The ``circuit`` must not contain any ``CircBox`` or non-unitary command. + The operator is provided when ``contract`` is called. + + Args: + bra: A pytket circuit describing the |bra> state. + ket: A pytket circuit describing the |ket> state. + attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. + scratch_fraction: Optional. Fraction of free memory on GPU to allocate as + scratch space; value between 0 and 1. Defaults to ``0.8``. + loglevel: Internal logger output level. + + Raises: + ValueError: If the circuits for ``ket`` or ``bra`` contain measurements. + ValueError: If the set of qubits of ``ket`` and ``bra`` do not match. + """ + self._logger = set_logger("GeneralBraOpKet", loglevel) + + # Check that the circuits have the same qubits + if set(ket.qubits) != set(bra.qubits): + raise ValueError( + "The circuits given to GeneralBraOpKet must act on the same qubits." + ) + # Remove end-of-circuit measurements and keep track of them separately + # It also resolves implicit swaps + ket, meas = _remove_meas_and_implicit_swaps(ket) + if len(meas) != 0: + raise ValueError( + "The circuits given to a GeneralBraOpKet cannot have measurements." + ) + bra, meas = _remove_meas_and_implicit_swaps(bra) + if len(meas) != 0: + raise ValueError( + "The circuits given to a GeneralBraOpKet cannot have measurements." + ) + # Identify each qubit with the index of the bond that represents it in the + # tensor network stored in this GeneralState. Qubits are sorted in increasing + # lexicographical order, which is the TKET standard. + self.qubit_idx_map = {q: i for i, q in enumerate(sorted(ket.qubits))} + self.n_qubits = ket.n_qubits + + dim = 2 # We are always dealing with qubits, not qudits + qubits_dims = (dim,) * self.n_qubits # qubit size + data_type = "complex128" # for now let that be hard-coded + + self.tn = NetworkState( + qubits_dims, + dtype=data_type, + config=attributes, + options={ + "memory_limit": f"{int(scratch_fraction*100)}%", + "logger": self._logger, + }, + ) + + # Apply all commands from the ket circuit + self._logger.debug("Converting the ket circuit to a NetworkState") + commands = ket.get_commands() + for com in commands: + try: + gate_unitary = com.op.get_unitary() + except: + raise ValueError( + "All commands in the circuit must be unitary gates. The ket circuit" + f" contains {com}; no unitary matrix could be retrived from it." + ) + gate_qubit_indices = tuple(self.qubit_idx_map[qb] for qb in com.qubits) + + tensor_id = self.tn.apply_tensor_operator( + gate_qubit_indices, + _formatted_tensor(gate_unitary, com.op.n_qubits), + immutable=True, # TODO: Change for parameterised gates + unitary=True, + ) + # If the circuit has no gates, apply one identity gate so that CuTensorNet does + # not panic due to no tensor operator in the NetworkState + if len(commands) == 0: + tensor_id = self.tn.apply_tensor_operator( + (0,), + _formatted_tensor(np.identity(2), 1), + immutable=True, + unitary=True, + ) + self.ket_phase = ket.phase + + # Create a placeholder Pauli identity operator, to be replaced when calling + # contract. + self._logger.debug("Creating a placeholder operator and appending it") + self._pauli_op_ids = [] + for mode in range(self.n_qubits): + tensor_id = self.tn.apply_tensor_operator( + (mode,), + _formatted_tensor(np.identity(2), 1), + immutable=False, + unitary=True, + ) + self._pauli_op_ids.append(tensor_id) + + # Apply all commands from the adjoint of the bra circuit + self._logger.debug("Applying the dagger of the bra circuit to the NetworkState") + commands = bra.dagger().get_commands() + for com in commands: + try: + gate_unitary = com.op.get_unitary() + except: + raise ValueError( + "All commands in the circuit must be unitary gates. The bra circuit" + f" contains {com}; no unitary matrix could be retrived from it." + ) + gate_qubit_indices = tuple(self.qubit_idx_map[qb] for qb in com.qubits) + + tensor_id = self.tn.apply_tensor_operator( + gate_qubit_indices, + _formatted_tensor(gate_unitary, com.op.n_qubits), + immutable=True, # TODO: Change for parameterised gates + unitary=True, + ) + # If the circuit has no gates, apply one identity gate so that CuTensorNet does + # not panic due to no tensor operator in the NetworkState + if len(commands) == 0: + tensor_id = self.tn.apply_tensor_operator( + (0,), + _formatted_tensor(np.identity(2), 1), + immutable=True, + unitary=True, + ) + self.bra_phase = bra.phase + + # TODO: Add the list of parameters here + def contract(self, operator: Optional[QubitPauliOperator] = None) -> complex: + """Contract the tensor network to obtain the value of ````. + + Args: + operator: A pytket ``QubitPauliOperator`` describing the operator. If not + given, then the identity operator is used, so it computes inner product. + + Returns: + The value of ```` + + Raises: + ValueError: If ``operator`` acts on qubits that are not in the circuits. + """ + paulis = ["I", "X", "Y", "Z"] + pauli_matrix = { + "I": np.identity(2), + "X": Op.create(OpType.X).get_unitary(), + "Y": Op.create(OpType.Y).get_unitary(), + "Z": Op.create(OpType.Z).get_unitary(), + } + pauli_strs: dict[str, complex] = dict() + + # Some care has to be taken when handling QubitPauliOperators, since identity + # Paulis may be omitted from the dictionary. + if operator is None: + pauli_strs = {"".join("I" for _ in range(self.n_qubits)): complex(1.0)} + else: + for tk_pstr, coeff in operator._dict.items(): + # Raise an error if the operator acts on qubits missing from the circuit + if any(q not in self.qubit_idx_map.keys() for q in tk_pstr.map.keys()): + raise ValueError( + f"The operator is acting on qubits {tk_pstr.map.keys()}, some " + "of these are missing from the set of qubits present in the " + f"circuits: {self.qubit_idx_map.keys()}." + ) + pauli_list = [tk_pstr[q] for q in self.qubit_idx_map.keys()] + this_pauli_string = "".join(map(lambda x: paulis[x], pauli_list)) + pauli_strs[this_pauli_string] = complex(coeff) + + # Calculate the value by iterating over all components of the QubitPauliOperator + value = 0.0 + zero_bitstring = "".join("0" for _ in range(self.n_qubits)) + for pstr, coeff in pauli_strs.items(): + # Update the NetworkState with this Pauli + self._logger.debug(f"Updating the tensors of the Pauli operator {pstr}") + for mode, pauli in enumerate(pstr): + self.tn.update_tensor_operator( + self._pauli_op_ids[mode], + _formatted_tensor(pauli_matrix[pauli], 1), + unitary=True, + ) + + if isinstance(coeff, Expr): + numeric_coeff = complex(coeff.evalf()) # type: ignore + else: + numeric_coeff = complex(coeff) + + # Compute the amplitude of the |0> state. Since NetworkState holds the + # circuit bra.dagger()*operator*ket|0>, the value of the amplitude at <0| + # will be . + self._logger.debug(f"Computing the contribution of Pauli operator {pstr}") + value += numeric_coeff * self.tn.compute_amplitude(zero_bitstring) + + # Apply the phases from the circuits + value *= np.exp(1j * np.pi * self.ket_phase) + value *= np.exp(-1j * np.pi * self.bra_phase) + + return complex(value) + + def destroy(self) -> None: + """Destroy the tensor network and free up GPU memory. + + Note: + Users are required to call `destroy()` when done using a + `GeneralState` object. GPU memory deallocation is not + guaranteed otherwise. + """ + self._logger.debug("Freeing memory of GeneralBraOpKet") + self.tn.free() + + def _formatted_tensor(matrix: NDArray, n_qubits: int) -> cp.ndarray: """Convert a matrix to the tensor format accepted by NVIDIA's API.""" - cupy_matrix = cp.asarray(matrix).astype(dtype="complex128") + cupy_matrix = cp.asarray(matrix, order="C").astype(dtype="complex128") # We also need to reshape since a matrix only has 2 bonds, but for an # n-qubit gate we want 2^n bonds for input and another 2^n for output return cupy_matrix.reshape([2] * (2 * n_qubits)) diff --git a/tests/test_general_state.py b/tests/test_general_state.py index 7becf16f..ecbb23c8 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -6,7 +6,7 @@ from pytket.transform import Transform from pytket.pauli import QubitPauliString, Pauli from pytket.utils.operators import QubitPauliOperator -from pytket.extensions.cutensornet.general_state import GeneralState +from pytket.extensions.cutensornet.general_state import GeneralState, GeneralBraOpKet @pytest.mark.parametrize( @@ -37,7 +37,7 @@ pytest.lazy_fixture("q8_x0h2v5z6"), # type: ignore ], ) -def test_get_statevec(circuit: Circuit) -> None: +def test_basic_circs_state(circuit: Circuit) -> None: state = GeneralState(circuit) sv = state.get_statevector() @@ -60,6 +60,12 @@ def test_get_statevec(circuit: Circuit) -> None: for i in range(len(sv)): assert np.isclose(sv[i], state.get_amplitude(i)) + # Calculate the inner product again, using GeneralBraOpKet + braket = GeneralBraOpKet(circuit, circuit) + ovl = braket.contract() + assert ovl == pytest.approx(1.0) + + braket.destroy() state.destroy() @@ -203,6 +209,13 @@ def test_expectation_value(circuit: Circuit, observable: QubitPauliOperator) -> assert np.isclose(exp_val, exp_val_tket) state.destroy() + # Calculate using GeneralBraOpKet + braket = GeneralBraOpKet(circuit, circuit) + exp_val = braket.contract(observable) + + assert np.isclose(exp_val, exp_val_tket) + braket.destroy() + @pytest.mark.parametrize( "circuit", From 4fe07bf48945ed1d43abb194100da74d67a8f78e Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Tue, 22 Oct 2024 15:43:07 +0000 Subject: [PATCH 10/27] Addded support forparameterised circuits --- .../general_state/tensor_network_state.py | 234 +++++++++++++----- tests/conftest.py | 28 +++ tests/test_general_state.py | 48 ++++ 3 files changed, 250 insertions(+), 60 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 37d9a183..da07b69b 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -22,7 +22,7 @@ except ImportError: warnings.warn("local settings failed to import cupy", ImportWarning) import numpy as np -from sympy import Expr # type: ignore +from sympy import Expr, Symbol # type: ignore from numpy.typing import NDArray from pytket.circuit import Circuit, Qubit, Bit, OpType, Op from pytket.extensions.cutensornet.general import set_logger @@ -68,19 +68,19 @@ def __init__( # Remove end-of-circuit measurements and keep track of them separately # It also resolves implicit swaps - # TODO: Is there any point in keeping the circuit around? - self._circuit, self._measurements = _remove_meas_and_implicit_swaps(circuit) + circuit, self._measurements = _remove_meas_and_implicit_swaps(circuit) # Identify each qubit with the index of the bond that represents it in the # tensor network stored in this GeneralState. Qubits are sorted in increasing # lexicographical order, which is the TKET standard. - self._qubit_idx_map = {q: i for i, q in enumerate(sorted(self._circuit.qubits))} + self._qubit_idx_map = {q: i for i, q in enumerate(sorted(circuit.qubits))} - num_qubits = self._circuit.n_qubits + self._phase = circuit.phase + self.n_qubits = circuit.n_qubits dim = 2 # We are always dealing with qubits, not qudits - qubits_dims = (dim,) * num_qubits # qubit size + qubits_dims = (dim,) * len(circuit.qubits) # qubit size data_type = "complex128" # for now let that be hard-coded - self._state = NetworkState( + self.tn_state = NetworkState( qubits_dims, dtype=data_type, config=attributes, @@ -90,30 +90,40 @@ def __init__( }, ) - commands = self._circuit.get_commands() - + # Maintain a dict of tensor_id->Op for symbolic Ops to be update when the user + # calls any evaluation function with the symbols specified. + self._symbolic_ops: dict[int, Op] = dict() # Append all gates to the NetworkState + commands = circuit.get_commands() for com in commands: - try: - gate_unitary = com.op.get_unitary() - except: - raise ValueError( - "All commands in the circuit must be unitary gates. The circuit " - f"contains {com}; no unitary matrix could be retrived from it." - ) - gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) + is_fixed = len(com.op.free_symbols()) == 0 + if is_fixed: + try: + gate_unitary = com.op.get_unitary() + except: + raise ValueError( + "All commands in the circuit must be unitary gates. The circuit" + f" contains {com}; no unitary matrix could be retrived from it." + ) + else: + # Dummy unitary to be updated later with user specified paramaters + gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") - tensor_id = self._state.apply_tensor_operator( + gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) + tensor_id = self.tn_state.apply_tensor_operator( gate_qubit_indices, _formatted_tensor(gate_unitary, com.op.n_qubits), - immutable=True, # TODO: Change for parameterised gates + immutable=is_fixed, unitary=True, ) + if not is_fixed: + self._symbolic_ops[tensor_id] = com.op + # If the circuit has no gates, apply one identity gate so that CuTensorNet does # not panic due to no tensor operator in the NetworkState if len(commands) == 0: - tensor_id = self._state.apply_tensor_operator( + tensor_id = self.tn_state.apply_tensor_operator( (0,), _formatted_tensor(np.identity(2), 1), immutable=True, @@ -122,29 +132,40 @@ def __init__( def get_statevector( self, + symbol_map: Optional[dict[Symbol, float]] = None, on_host: bool = True, ) -> Union[cp.ndarray, np.ndarray]: """Contracts the circuit and returns the final statevector. Args: + symbol_map: A dictionary where each element of ``sef.free_symbols`` is + assigned a real number. on_host: Optional. If ``True``, converts cupy ``ndarray`` to numpy - ``ndarray``, copying it to host device (CPU). + ``ndarray``, copying it to host device (CPU). Defaults to ``True``. Returns: Either a ``cupy.ndarray`` on a GPU, or a ``numpy.ndarray`` on a host device (CPU). Arrays are returned in a 1D shape. + + Raises: + ValueError: If not every free symbol in the circuits is assigned a + value in ``symbol_map``. """ + if symbol_map is not None: + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + self._logger.debug("(Statevector) contracting the TN") - state_vector = self._state.compute_state_vector() + state_vector = self.tn_state.compute_state_vector() sv = state_vector.flatten() # Convert to 1D if on_host: sv = cp.asnumpy(sv) # Apply the phase from the circuit - sv *= np.exp(1j * np.pi * self._circuit.phase) + sv *= np.exp(1j * np.pi * self._phase) return sv def get_amplitude( self, state: int, + symbol_map: Optional[dict[Symbol, float]] = None, ) -> complex: """Returns the amplitude of the chosen computational state. @@ -155,32 +176,49 @@ def get_amplitude( Args: state: The integer whose bitstring describes the computational state. The qubits in the bitstring are in increasing lexicographic order. + symbol_map: A dictionary where each element of ``sef.free_symbols`` is + assigned a real number. Returns: The amplitude of the computational state in ``self``. + + Raises: + ValueError: If not every free symbol in the circuits is assigned a + value in ``symbol_map``. """ + if symbol_map is not None: + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + self._logger.debug("(Amplitude) contracting the TN") - bitstring = format(state, "b").zfill(self._circuit.n_qubits) - amplitude = self._state.compute_amplitude(bitstring) + bitstring = format(state, "b").zfill(self.n_qubits) + amplitude = self.tn_state.compute_amplitude(bitstring) # Apply the phase from the circuit - amplitude *= np.exp(1j * np.pi * self._circuit.phase) + amplitude *= np.exp(1j * np.pi * self._phase) return complex(amplitude) def expectation_value( self, operator: QubitPauliOperator, + symbol_map: Optional[dict[Symbol, float]] = None, ) -> complex: """Calculates the expectation value of the given operator. Args: operator: The operator whose expectation value is to be measured. - - Raises: - ValueError: If the operator acts on qubits not present in the circuit. + symbol_map: A dictionary where each element of ``sef.free_symbols`` is + assigned a real number. Returns: The expectation value. + + Raises: + ValueError: If the operator acts on qubits not present in the circuit. + ValueError: If not every free symbol in the circuits is assigned a + value in ``symbol_map``. """ + if symbol_map is not None: + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + self._logger.debug("(Expectation value) converting operator to NetworkOperator") paulis = ["I", "X", "Y", "Z"] @@ -188,11 +226,11 @@ def expectation_value( for pstr, coeff in operator._dict.items(): # Raise an error if the operator acts on qubits that are not in the circuit - if any(q not in self._circuit.qubits for q in pstr.map.keys()): + if any(q not in self._qubit_idx_map.keys() for q in pstr.map.keys()): raise ValueError( f"The operator is acting on qubits {pstr.map.keys()}, " "but some of these are not present in the circuit, whose set of " - f"qubits is: {self._circuit.qubits}." + f"qubits is: {self._qubit_idx_map.keys()}." ) pauli_list = [pstr[q] for q in self._qubit_idx_map.keys()] @@ -202,21 +240,30 @@ def expectation_value( tn_operator = NetworkOperator.from_pauli_strings(pauli_strs, dtype="complex128") self._logger.debug("(Expectation value) contracting the TN") - return complex(self._state.compute_expectation(tn_operator)) + return complex(self.tn_state.compute_expectation(tn_operator)) def sample( # TODO: Support seeds (and test) self, n_shots: int, + symbol_map: Optional[dict[Symbol, float]] = None, ) -> BackendResult: """Obtains samples from the measurements at the end of the circuit. Args: n_shots: The number of samples to obtain. - Raises: - ValueError: If the circuit contains no measurements. + symbol_map: A dictionary where each element of ``sef.free_symbols`` is + assigned a real number. + Returns: A pytket ``BackendResult`` with the data from the shots. + + Raises: + ValueError: If the circuit contains no measurements. + ValueError: If not every free symbol in the circuits is assigned a + value in ``symbol_map``. """ + if symbol_map is not None: + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) num_measurements = len(self._measurements) if num_measurements == 0: @@ -230,7 +277,7 @@ def sample( # TODO: Support seeds (and test) measured_modes = tuple(self._qubit_idx_map[qb] for qb in qbit_list) self._logger.debug("(Sampling) contracting the TN") - samples = self._state.compute_sampling( + samples = self.tn_state.compute_sampling( nshots=n_shots, modes=measured_modes, ) @@ -266,7 +313,7 @@ def destroy(self) -> None: guaranteed otherwise. """ self._logger.debug("Freeing memory of GeneralState") - self._state.free() + self.tn_state.free() class GeneralBraOpKet: # TODO: Write it as a context manager @@ -324,7 +371,7 @@ def __init__( # Identify each qubit with the index of the bond that represents it in the # tensor network stored in this GeneralState. Qubits are sorted in increasing # lexicographical order, which is the TKET standard. - self.qubit_idx_map = {q: i for i, q in enumerate(sorted(ket.qubits))} + self._qubit_idx_map = {q: i for i, q in enumerate(sorted(ket.qubits))} self.n_qubits = ket.n_qubits dim = 2 # We are always dealing with qubits, not qudits @@ -341,25 +388,37 @@ def __init__( }, ) + # Maintain a dict of tensor_id->Op for symbolic Ops to be update when the user + # calls any evaluation function with the symbols specified. + self._symbolic_ops: dict[int, Op] = dict() # Apply all commands from the ket circuit self._logger.debug("Converting the ket circuit to a NetworkState") commands = ket.get_commands() for com in commands: - try: - gate_unitary = com.op.get_unitary() - except: - raise ValueError( - "All commands in the circuit must be unitary gates. The ket circuit" - f" contains {com}; no unitary matrix could be retrived from it." - ) - gate_qubit_indices = tuple(self.qubit_idx_map[qb] for qb in com.qubits) + is_fixed = len(com.op.free_symbols()) == 0 + if is_fixed: + try: + gate_unitary = com.op.get_unitary() + except: + raise ValueError( + "All commands in the circuit must be unitary gates. The circuit" + f" contains {com}; no unitary matrix could be retrived from it." + ) + else: + # Dummy unitary to be updated later with user specified paramaters + gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") + gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) tensor_id = self.tn.apply_tensor_operator( gate_qubit_indices, _formatted_tensor(gate_unitary, com.op.n_qubits), - immutable=True, # TODO: Change for parameterised gates + immutable=is_fixed, unitary=True, ) + + if not is_fixed: + self._symbolic_ops[tensor_id] = com.op + # If the circuit has no gates, apply one identity gate so that CuTensorNet does # not panic due to no tensor operator in the NetworkState if len(commands) == 0: @@ -388,21 +447,30 @@ def __init__( self._logger.debug("Applying the dagger of the bra circuit to the NetworkState") commands = bra.dagger().get_commands() for com in commands: - try: - gate_unitary = com.op.get_unitary() - except: - raise ValueError( - "All commands in the circuit must be unitary gates. The bra circuit" - f" contains {com}; no unitary matrix could be retrived from it." - ) - gate_qubit_indices = tuple(self.qubit_idx_map[qb] for qb in com.qubits) + is_fixed = len(com.op.free_symbols()) == 0 + if is_fixed: + try: + gate_unitary = com.op.get_unitary() + except: + raise ValueError( + "All commands in the circuit must be unitary gates. The circuit" + f" contains {com}; no unitary matrix could be retrived from it." + ) + else: + # Dummy unitary to be updated later with user specified paramaters + gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") + gate_qubit_indices = tuple(self._qubit_idx_map[qb] for qb in com.qubits) tensor_id = self.tn.apply_tensor_operator( gate_qubit_indices, _formatted_tensor(gate_unitary, com.op.n_qubits), - immutable=True, # TODO: Change for parameterised gates + immutable=is_fixed, unitary=True, ) + + if not is_fixed: + self._symbolic_ops[tensor_id] = com.op + # If the circuit has no gates, apply one identity gate so that CuTensorNet does # not panic due to no tensor operator in the NetworkState if len(commands) == 0: @@ -414,13 +482,18 @@ def __init__( ) self.bra_phase = bra.phase - # TODO: Add the list of parameters here - def contract(self, operator: Optional[QubitPauliOperator] = None) -> complex: + def contract( + self, + operator: Optional[QubitPauliOperator] = None, + symbol_map: Optional[dict[Symbol, float]] = None, + ) -> complex: """Contract the tensor network to obtain the value of ````. Args: operator: A pytket ``QubitPauliOperator`` describing the operator. If not given, then the identity operator is used, so it computes inner product. + symbol_map: A dictionary where each element of ``sef.free_symbols`` is + assigned a real number. Returns: The value of ```` @@ -428,6 +501,9 @@ def contract(self, operator: Optional[QubitPauliOperator] = None) -> complex: Raises: ValueError: If ``operator`` acts on qubits that are not in the circuits. """ + if symbol_map is not None: + _update_tensors(self.tn, self._symbolic_ops, symbol_map) + paulis = ["I", "X", "Y", "Z"] pauli_matrix = { "I": np.identity(2), @@ -444,13 +520,13 @@ def contract(self, operator: Optional[QubitPauliOperator] = None) -> complex: else: for tk_pstr, coeff in operator._dict.items(): # Raise an error if the operator acts on qubits missing from the circuit - if any(q not in self.qubit_idx_map.keys() for q in tk_pstr.map.keys()): + if any(q not in self._qubit_idx_map.keys() for q in tk_pstr.map.keys()): raise ValueError( f"The operator is acting on qubits {tk_pstr.map.keys()}, some " "of these are missing from the set of qubits present in the " - f"circuits: {self.qubit_idx_map.keys()}." + f"circuits: {self._qubit_idx_map.keys()}." ) - pauli_list = [tk_pstr[q] for q in self.qubit_idx_map.keys()] + pauli_list = [tk_pstr[q] for q in self._qubit_idx_map.keys()] this_pauli_string = "".join(map(lambda x: paulis[x], pauli_list)) pauli_strs[this_pauli_string] = complex(coeff) @@ -538,3 +614,41 @@ def _remove_meas_and_implicit_swaps(circ: Circuit) -> Tuple[Circuit, Dict[Qubit, pure_circ.add_phase(circ.phase) return pure_circ, measure_map # type: ignore + + +def _update_tensors( + tn: NetworkState, + symbolic_ops: dict[int, Op], + symbol_map: dict[Symbol, float], +) -> None: + """Updates the tensors with the specified values for symbols. + + Args: + tn: The NetworkState that we intend to update. + symbolic_ops: A dictionary mapping ``tensor_id`` to the parameterised Op. + symbol_map: A dictionary mapping symbols to real values. + + Raises: + ValueError: If not every free symbol in the circuits is assigned a + value in ``symbol_map``. + """ + for tensor_id, op in symbolic_ops.items(): + subs_params = op.params.copy() + if any(symb not in symbol_map.keys() for symb in op.free_symbols()): + raise ValueError( + f"Missing values for some of the free symbols {op.free_symbols()}. " + f"Symbols given: {symbol_map}." + ) + + for symb in op.free_symbols(): + subs_params = [ + p.subs(symb, symbol_map[symb]) if isinstance(p, Expr) else p + for p in subs_params + ] + gate_unitary = Op.create(op.type, subs_params).get_unitary() + + tn.update_tensor_operator( + tensor_id, + _formatted_tensor(gate_unitary, op.n_qubits), + unitary=True, + ) diff --git a/tests/conftest.py b/tests/conftest.py index 5297f495..650fc8be 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,5 +1,6 @@ import pytest import numpy as np +from sympy import Symbol from scipy.stats import unitary_group # type: ignore from pytket.circuit import Circuit, OpType, Unitary2qBox, ToffoliBox, Qubit from pytket.passes import DecomposeBoxes, CnXPairwiseDecomposition @@ -180,6 +181,20 @@ def q4_lcu1() -> Circuit: return circuit +@pytest.fixture +def q4_lcu1_parameterised() -> Circuit: + a, b, c = Symbol("a"), Symbol("b"), Symbol("c") + circuit = Circuit(4) + circuit.Ry(a, 3).Ry(0.27, 2).CX(2, 3).Ry(b, 2).Ry(0.21, 3) + circuit.Ry(0.12, 0).Ry(a, 1) + circuit.add_gate(OpType.CnX, [0, 1, 2]).add_gate(OpType.CnX, [0, 1, 3]) + circuit.X(0).X(1).add_gate(OpType.CnY, [0, 1, 2]).add_gate(OpType.CnY, [0, 1, 3]).X( + 0 + ).X(1) + circuit.Ry(-b, 0).Ry(-c, 1) + return circuit + + @pytest.fixture def q4_multicontrols() -> Circuit: circ = Circuit(4) @@ -229,6 +244,19 @@ def q5_h0s1rz2ry3tk4tk13() -> Circuit: return circuit +@pytest.fixture +def q5_h0s1rz2ry3tk4tk13_parameterised() -> Circuit: + a, b, c = Symbol("a"), Symbol("b"), Symbol("c") + circuit = Circuit(5) + circuit.H(0) + circuit.S(1) + circuit.Rz(a * c, 2) + circuit.Ry(b + a, 3) + circuit.TK1(a, b, c, 4) + circuit.TK2(a - b, c - a, (a + b) * c, 1, 3) + return circuit + + @pytest.fixture def q8_x0h2v5z6() -> Circuit: circuit = Circuit(8) diff --git a/tests/test_general_state.py b/tests/test_general_state.py index ecbb23c8..b963ec8a 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -1,6 +1,7 @@ import random import numpy as np import pytest +from sympy import Symbol from pytket.circuit import Circuit, ToffoliBox, Qubit, Bit from pytket.passes import DecomposeBoxes, CnXPairwiseDecomposition from pytket.transform import Transform @@ -288,3 +289,50 @@ def test_sampler(circuit: Circuit, measure_all: bool) -> None: assert np.isclose(count / n_shots, prob, atol=0.01) state.destroy() + + +@pytest.mark.parametrize( + "circuit", + [ + pytest.lazy_fixture("q4_lcu1_parameterised"), # type: ignore + pytest.lazy_fixture("q5_h0s1rz2ry3tk4tk13_parameterised"), # type: ignore + ], +) +@pytest.mark.parametrize( + "symbol_map", + [ + {Symbol("a"): 0.3, Symbol("b"): 0.42, Symbol("c"): -0.13}, + {Symbol("a"): 5.3, Symbol("b"): 1.42, Symbol("c"): -0.07, Symbol("d"): 0.53}, + ], +) +def test_parameterised(circuit: Circuit, symbol_map: dict[Symbol, float]) -> None: + state = GeneralState(circuit) + sv = state.get_statevector(symbol_map) + + circuit.symbol_substitution(symbol_map) + sv_pytket = circuit.get_statevector() + assert np.allclose(sv, sv_pytket, atol=1e-10) + + op = QubitPauliOperator( + { + QubitPauliString({q: Pauli.I for q in circuit.qubits}): 1.0, + } + ) + + # Calculate the inner product as the expectation value + # of the identity operator: = + state = GeneralState(circuit) + ovl = state.expectation_value(op) + assert ovl == pytest.approx(1.0) + + # Check that all amplitudes agree + for i in range(len(sv)): + assert np.isclose(sv[i], state.get_amplitude(i)) + + # Calculate the inner product again, using GeneralBraOpKet + braket = GeneralBraOpKet(circuit, circuit) + ovl = braket.contract() + assert ovl == pytest.approx(1.0) + + braket.destroy() + state.destroy() From e26cb324af251397f323f864ab95453303797ea3 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Tue, 22 Oct 2024 15:53:16 +0000 Subject: [PATCH 11/27] Pleasing linter --- pytket/extensions/cutensornet/backends/cutensornet_backend.py | 1 - .../cutensornet/general_state/tensor_network_state.py | 2 -- 2 files changed, 3 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index fe24f881..d6c4b692 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -15,7 +15,6 @@ """Methods to allow tket circuits to be run on the cuTensorNet simulator.""" from abc import abstractmethod -import warnings from typing import List, Union, Optional, Sequence from uuid import uuid4 diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index da07b69b..e8f8135d 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -31,8 +31,6 @@ from pytket.backends.backendresult import BackendResult try: - import cuquantum as cq # type: ignore - from cuquantum import cutensornet as cutn # type: ignore from cuquantum.cutensornet.experimental import NetworkState, NetworkOperator # type: ignore except ImportError: warnings.warn("local settings failed to import cuquantum", ImportWarning) From 7ad5fc2a80b0dc7a33662ac9d82a890da59618ee Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 09:54:41 +0000 Subject: [PATCH 12/27] Added support for RNG seeds --- .../cutensornet/backends/cutensornet_backend.py | 12 ++++-------- .../general_state/tensor_network_state.py | 8 +++++++- tests/test_cutensornet_backend.py | 16 ++++++++++++++++ 3 files changed, 27 insertions(+), 9 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index d6c4b692..94ee1b96 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -281,6 +281,8 @@ def process_circuits( n_shots: Number of shots in case of shot-based calculation. Optionally, this can be a list of shots specifying the number of shots for each circuit separately. + seed: An optional RNG seed. Different calls to ``process_circuits`` with the + same seed will generate the same list of shot outcomes. valid_check: Whether to check for circuit correctness. tnconfig: Optional. A dict of cuTensorNet ``TNConfig`` keys and their values. @@ -292,13 +294,7 @@ def process_circuits( """ scratch_fraction = float(kwargs.get("scratch_fraction", 0.8)) # type: ignore tnconfig = kwargs.get("tnconfig", dict()) # type: ignore - - if "seed" in kwargs and kwargs["seed"] is not None: - # Current CuTensorNet does not support seeds for Sampler. I created - # a feature request in their repository. - raise NotImplementedError( # TODO: Support seeds! - "The backend does not currently support user-defined seeds." - ) + seed = kwargs.get("seed", None) if n_shots is None: raise ValueError( @@ -318,7 +314,7 @@ def process_circuits( circuit, attributes=tnconfig, scratch_fraction=scratch_fraction ) handle = ResultHandle(str(uuid4())) - self._cache[handle] = {"result": tn.sample(circ_shots)} + self._cache[handle] = {"result": tn.sample(circ_shots, seed=seed)} handle_list.append(handle) return handle_list diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index e8f8135d..61f79375 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -240,10 +240,11 @@ def expectation_value( self._logger.debug("(Expectation value) contracting the TN") return complex(self.tn_state.compute_expectation(tn_operator)) - def sample( # TODO: Support seeds (and test) + def sample( self, n_shots: int, symbol_map: Optional[dict[Symbol, float]] = None, + seed: Optional[int] = None, ) -> BackendResult: """Obtains samples from the measurements at the end of the circuit. @@ -251,6 +252,8 @@ def sample( # TODO: Support seeds (and test) n_shots: The number of samples to obtain. symbol_map: A dictionary where each element of ``sef.free_symbols`` is assigned a real number. + seed: An optional RNG seed. Different calls to ``sample`` with the same + seed will generate the same list of shot outcomes. Returns: A pytket ``BackendResult`` with the data from the shots. @@ -275,9 +278,12 @@ def sample( # TODO: Support seeds (and test) measured_modes = tuple(self._qubit_idx_map[qb] for qb in qbit_list) self._logger.debug("(Sampling) contracting the TN") + if seed is not None: + seed = abs(seed) # Must be a positive integer samples = self.tn_state.compute_sampling( nshots=n_shots, modes=measured_modes, + seed=seed, ) # Convert the data in `samples` to an `OutcomeArray` using `from_readouts` diff --git a/tests/test_cutensornet_backend.py b/tests/test_cutensornet_backend.py index 25c9d31c..80816b22 100644 --- a/tests/test_cutensornet_backend.py +++ b/tests/test_cutensornet_backend.py @@ -34,6 +34,22 @@ def test_sampler_bell() -> None: assert np.isclose(res.get_counts()[(1, 1)] / n_shots, 0.5, atol=0.01) +def test_sampler_bell_seed() -> None: + n_shots = 1000 + c = Circuit(2, 2) + c.H(0) + c.CX(0, 1) + c.measure_all() + b = CuTensorNetShotsBackend() + c = b.get_compiled_circuit(c) + res1 = b.run_circuit(c, n_shots=n_shots, seed=1234) + res2 = b.run_circuit(c, n_shots=n_shots, seed=1234) + res3 = b.run_circuit(c, n_shots=n_shots, seed=4321) + print(type(res1.get_shots())) + assert np.all(res1.get_shots() == res2.get_shots()) + assert np.any(res1.get_shots() != res3.get_shots()) + + def test_config_options() -> None: c = Circuit(2, 2) c.H(0) From 7240e551c697d331651aa27b0302b292aaf80941 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 09:59:29 +0000 Subject: [PATCH 13/27] Fixed test with old configuration option --- tests/test_cutensornet_backend.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_cutensornet_backend.py b/tests/test_cutensornet_backend.py index 80816b22..e58cf098 100644 --- a/tests/test_cutensornet_backend.py +++ b/tests/test_cutensornet_backend.py @@ -62,7 +62,7 @@ def test_config_options() -> None: h = b1.process_circuit( c, scratch_fraction=0.3, - CONFIG_NUM_HYPER_SAMPLES=100, + tn_config={"num_hyper_samples": 100}, ) assert np.allclose( b1.get_result(h).get_state(), np.asarray([1, 0, 0, 1]) * 1 / np.sqrt(2) @@ -76,7 +76,7 @@ def test_config_options() -> None: c, n_shots=n_shots, scratch_fraction=0.3, - CONFIG_NUM_HYPER_SAMPLES=100, + tn_config={"num_hyper_samples": 100}, ) assert res.get_shots().shape == (n_shots, 2) assert np.isclose(res.get_counts()[(0, 0)] / n_shots, 0.5, atol=0.01) From 13f74d83d24ad66632889ec17f49efd26388e16b Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 11:08:24 +0000 Subject: [PATCH 14/27] GeneralState and GeneralBraOpKet as context managers --- .../backends/cutensornet_backend.py | 14 ++-- .../general_state/tensor_network_state.py | 40 +++++++-- tests/test_general_state.py | 82 ++++++++----------- 3 files changed, 72 insertions(+), 64 deletions(-) diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index 94ee1b96..4d8a0386 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -221,10 +221,10 @@ def process_circuits( self._check_all_circuits(circuit_list) handle_list = [] for circuit in circuit_list: - tn = GeneralState( + with GeneralState( circuit, attributes=tnconfig, scratch_fraction=scratch_fraction - ) - sv = tn.get_statevector() + ) as tn: + sv = tn.get_statevector() res_qubits = [qb for qb in sorted(circuit.qubits)] handle = ResultHandle(str(uuid4())) self._cache[handle] = {"result": BackendResult(q_bits=res_qubits, state=sv)} @@ -310,11 +310,11 @@ def process_circuits( self._check_all_circuits(circuit_list) handle_list = [] for circuit, circ_shots in zip(circuit_list, all_shots): - tn = GeneralState( - circuit, attributes=tnconfig, scratch_fraction=scratch_fraction - ) handle = ResultHandle(str(uuid4())) - self._cache[handle] = {"result": tn.sample(circ_shots, seed=seed)} + with GeneralState( + circuit, attributes=tnconfig, scratch_fraction=scratch_fraction + ) as tn: + self._cache[handle] = {"result": tn.sample(circ_shots, seed=seed)} handle_list.append(handle) return handle_list diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 61f79375..8123b81e 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -14,7 +14,7 @@ from __future__ import annotations import logging -from typing import Union, Optional, Tuple, Dict +from typing import Union, Optional, Tuple, Dict, Any import warnings try: @@ -37,7 +37,12 @@ class GeneralState: # TODO: Write it as a context manager so that I can call free() - """Wrapper of cuTensorNet's NetworkState for exact simulation of states.""" + """Wrapper of cuTensorNet's NetworkState for exact simulation of states. + + Note: + Preferably used as ``with GeneralState(...) as state:`` so that GPU memory is + automatically released after execution. + """ def __init__( self, @@ -312,16 +317,27 @@ def destroy(self) -> None: """Destroy the tensor network and free up GPU memory. Note: - Users are required to call `destroy()` when done using a - `GeneralState` object. GPU memory deallocation is not - guaranteed otherwise. + The preferred approach is to use a context manager as in + ``with GeneralState(...) as state:``. Otherwise, the user must release + memory explicitly by calling ``destroy()``. """ self._logger.debug("Freeing memory of GeneralState") self.tn_state.free() + def __enter__(self) -> GeneralState: + return self + + def __exit__(self, exc_type: Any, exc_value: Any, exc_tb: Any) -> None: + self.destroy() + class GeneralBraOpKet: # TODO: Write it as a context manager - """Wrapper of cuTensorNet's NetworkState for exact simulation of ````.""" + """Wrapper of cuTensorNet's NetworkState for exact simulation of ````. + + Note: + Preferably used as ``with GeneralBraOpKet(...) as braket:`` so that GPU memory + is automatically released after execution. + """ def __init__( self, @@ -568,13 +584,19 @@ def destroy(self) -> None: """Destroy the tensor network and free up GPU memory. Note: - Users are required to call `destroy()` when done using a - `GeneralState` object. GPU memory deallocation is not - guaranteed otherwise. + The preferred approach is to use a context manager as in + ``with GeneralBraOpKet(...) as braket:``. Otherwise, the user must release + memory explicitly by calling ``destroy()``. """ self._logger.debug("Freeing memory of GeneralBraOpKet") self.tn.free() + def __enter__(self) -> GeneralBraOpKet: + return self + + def __exit__(self, exc_type: Any, exc_value: Any, exc_tb: Any) -> None: + self.destroy() + def _formatted_tensor(matrix: NDArray, n_qubits: int) -> cp.ndarray: """Convert a matrix to the tensor format accepted by NVIDIA's API.""" diff --git a/tests/test_general_state.py b/tests/test_general_state.py index b963ec8a..77ea5998 100644 --- a/tests/test_general_state.py +++ b/tests/test_general_state.py @@ -39,11 +39,7 @@ ], ) def test_basic_circs_state(circuit: Circuit) -> None: - state = GeneralState(circuit) - sv = state.get_statevector() - sv_pytket = circuit.get_statevector() - assert np.allclose(sv, sv_pytket, atol=1e-10) op = QubitPauliOperator( { @@ -51,23 +47,23 @@ def test_basic_circs_state(circuit: Circuit) -> None: } ) - # Calculate the inner product as the expectation value - # of the identity operator: = - state = GeneralState(circuit) - ovl = state.expectation_value(op) - assert ovl == pytest.approx(1.0) + with GeneralState(circuit) as state: + sv = state.get_statevector() + assert np.allclose(sv, sv_pytket, atol=1e-10) - # Check that all amplitudes agree - for i in range(len(sv)): - assert np.isclose(sv[i], state.get_amplitude(i)) + # Calculate the inner product as the expectation value + # of the identity operator: = + ovl = state.expectation_value(op) + assert ovl == pytest.approx(1.0) - # Calculate the inner product again, using GeneralBraOpKet - braket = GeneralBraOpKet(circuit, circuit) - ovl = braket.contract() - assert ovl == pytest.approx(1.0) + # Check that all amplitudes agree + for i in range(len(sv)): + assert np.isclose(sv[i], state.get_amplitude(i)) - braket.destroy() - state.destroy() + # Calculate the inner product again, using GeneralBraOpKet + with GeneralBraOpKet(circuit, circuit) as braket: + ovl = braket.contract() + assert ovl == pytest.approx(1.0) def test_sv_toffoli_box_with_implicit_swaps() -> None: @@ -92,9 +88,8 @@ def test_sv_toffoli_box_with_implicit_swaps() -> None: Transform.OptimiseCliffords().apply(ket_circ) # Convert and contract - state = GeneralState(ket_circ) - ket_net_vector = state.get_statevector() - state.destroy() + with GeneralState(ket_circ) as state: + ket_net_vector = state.get_statevector() # Compare to pytket statevector ket_pytket_vector = ket_circ.get_statevector() @@ -127,8 +122,8 @@ def to_bool_tuple(n_qubits: int, x: int) -> tuple: CnXPairwiseDecomposition().apply(ket_circ) Transform.OptimiseCliffords().apply(ket_circ) - state = GeneralState(ket_circ) - ket_net_vector = state.get_statevector() + with GeneralState(ket_circ) as state: + ket_net_vector = state.get_statevector() ket_pytket_vector = ket_circ.get_statevector() assert np.allclose(ket_net_vector, ket_pytket_vector) @@ -141,12 +136,10 @@ def to_bool_tuple(n_qubits: int, x: int) -> tuple: } ) - state = GeneralState(ket_circ) - ovl = state.expectation_value(op) + with GeneralState(ket_circ) as state: + ovl = state.expectation_value(op) assert ovl == pytest.approx(1.0) - state.destroy() - @pytest.mark.parametrize( "circuit", @@ -204,18 +197,16 @@ def test_expectation_value(circuit: Circuit, observable: QubitPauliOperator) -> exp_val_tket = observable.state_expectation(circuit.get_statevector()) # Calculate using GeneralState - state = GeneralState(circuit) - exp_val = state.expectation_value(observable) + with GeneralState(circuit) as state: + exp_val = state.expectation_value(observable) assert np.isclose(exp_val, exp_val_tket) - state.destroy() # Calculate using GeneralBraOpKet - braket = GeneralBraOpKet(circuit, circuit) - exp_val = braket.contract(observable) + with GeneralBraOpKet(circuit, circuit) as braket: + exp_val = braket.contract(observable) assert np.isclose(exp_val, exp_val_tket) - braket.destroy() @pytest.mark.parametrize( @@ -266,8 +257,8 @@ def test_sampler(circuit: Circuit, measure_all: bool) -> None: circuit.Measure(q, Bit(i)) # Sample using our library - state = GeneralState(circuit) - results = state.sample(n_shots) + with GeneralState(circuit) as state: + results = state.sample(n_shots) # Verify distribution matches theoretical probabilities for bit_tuple, count in results.get_counts().items(): @@ -288,8 +279,6 @@ def test_sampler(circuit: Circuit, measure_all: bool) -> None: assert np.isclose(count / n_shots, prob, atol=0.01) - state.destroy() - @pytest.mark.parametrize( "circuit", @@ -321,18 +310,15 @@ def test_parameterised(circuit: Circuit, symbol_map: dict[Symbol, float]) -> Non # Calculate the inner product as the expectation value # of the identity operator: = - state = GeneralState(circuit) - ovl = state.expectation_value(op) - assert ovl == pytest.approx(1.0) + with GeneralState(circuit) as state: + ovl = state.expectation_value(op) + assert ovl == pytest.approx(1.0) - # Check that all amplitudes agree - for i in range(len(sv)): - assert np.isclose(sv[i], state.get_amplitude(i)) + # Check that all amplitudes agree + for i in range(len(sv)): + assert np.isclose(sv[i], state.get_amplitude(i)) # Calculate the inner product again, using GeneralBraOpKet - braket = GeneralBraOpKet(circuit, circuit) - ovl = braket.contract() + with GeneralBraOpKet(circuit, circuit) as braket: + ovl = braket.contract() assert ovl == pytest.approx(1.0) - - braket.destroy() - state.destroy() From e760ba1e383c77c7d705c6d70bb0de52fbc09766 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 11:09:28 +0000 Subject: [PATCH 15/27] Removed solved TODO comments --- .../cutensornet/general_state/tensor_network_state.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 8123b81e..65ef3c27 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -36,7 +36,7 @@ warnings.warn("local settings failed to import cuquantum", ImportWarning) -class GeneralState: # TODO: Write it as a context manager so that I can call free() +class GeneralState: """Wrapper of cuTensorNet's NetworkState for exact simulation of states. Note: @@ -331,7 +331,7 @@ def __exit__(self, exc_type: Any, exc_value: Any, exc_tb: Any) -> None: self.destroy() -class GeneralBraOpKet: # TODO: Write it as a context manager +class GeneralBraOpKet: """Wrapper of cuTensorNet's NetworkState for exact simulation of ````. Note: From 35fee383c9b7fae778244e2b3e4da5d72bd6fbd0 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 14:10:38 +0000 Subject: [PATCH 16/27] Updated docs --- docs/api.rst | 5 ----- docs/changelog.rst | 8 ++++++++ docs/modules/general_state.rst | 23 +++++------------------ docs/modules/structured_state.rst | 4 ++++ 4 files changed, 17 insertions(+), 23 deletions(-) diff --git a/docs/api.rst b/docs/api.rst index 1e91e2c8..48ed0e2c 100644 --- a/docs/api.rst +++ b/docs/api.rst @@ -1,11 +1,6 @@ API documentation ----------------- -.. autoclass:: pytket.extensions.cutensornet.CuTensorNetHandle - - .. automethod:: destroy - - .. toctree:: modules/general_state.rst modules/structured_state.rst diff --git a/docs/changelog.rst b/docs/changelog.rst index a1fc3c73..ac0fd654 100644 --- a/docs/changelog.rst +++ b/docs/changelog.rst @@ -3,6 +3,14 @@ Changelog ~~~~~~~~~ +Unreleased +---------- + +* New API: ``GeneralBraOpKet`` for exact calculation of arbitrary ```` values. +* New feature: ``GeneralState`` has a new method ``get_amplitude`` to obtain the amplitude of computational basis states. +* New feature: ``GeneralState`` and ``CuTensorNetShotsBackend`` now support RNG seeds for sampling. +* Deprecated ``TensorNetwork`` object. It is still available for the sake backwards compatibility, but it has been removed from doc pages. + 0.9.0 (October 2024) --------------------- diff --git a/docs/modules/general_state.rst b/docs/modules/general_state.rst index 0c552812..bf557a0b 100644 --- a/docs/modules/general_state.rst +++ b/docs/modules/general_state.rst @@ -4,30 +4,17 @@ General state (exact) simulation .. automodule:: pytket.extensions.cutensornet.general_state .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralState() - .. automethod:: __init__ .. automethod:: get_statevector + .. automethod:: get_amplitude .. automethod:: expectation_value .. automethod:: sample .. automethod:: destroy -cuQuantum `contract` API interface -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. autoclass:: pytket.extensions.cutensornet.general_state.TensorNetwork - -.. autoclass:: pytket.extensions.cutensornet.general_state.PauliOperatorTensorNetwork - -.. autoclass:: pytket.extensions.cutensornet.general_state.ExpectationValueTensorNetwork - -.. autofunction:: pytket.extensions.cutensornet.general_state.tk_to_tensor_network - -.. autofunction:: pytket.extensions.cutensornet.general_state.measure_qubits_state - -.. autofunction:: pytket.extensions.cutensornet.general_state.get_operator_expectation_value - -.. autofunction:: pytket.extensions.cutensornet.general_state.get_circuit_overlap - +.. autoclass:: pytket.extensions.cutensornet.general_state.GeneralBraOpKet() + .. automethod:: __init__ + .. automethod:: contract + .. automethod:: destroy Pytket backend ~~~~~~~~~~~~~~ diff --git a/docs/modules/structured_state.rst b/docs/modules/structured_state.rst index 94cb320f..ce1ea6d8 100644 --- a/docs/modules/structured_state.rst +++ b/docs/modules/structured_state.rst @@ -3,6 +3,10 @@ Structured state evolution .. automodule:: pytket.extensions.cutensornet.structured_state +.. autoclass:: pytket.extensions.cutensornet.CuTensorNetHandle + + .. automethod:: destroy + Simulation ~~~~~~~~~~ From ea12cbb24bfd727722151cf140845a3c31760828 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 14:25:25 +0000 Subject: [PATCH 17/27] Fixed problems in docs --- docs/changelog.rst | 5 +++-- docs/modules/general_state.rst | 2 ++ docs/modules/structured_state.rst | 3 +++ 3 files changed, 8 insertions(+), 2 deletions(-) diff --git a/docs/changelog.rst b/docs/changelog.rst index ac0fd654..d2c71682 100644 --- a/docs/changelog.rst +++ b/docs/changelog.rst @@ -6,10 +6,11 @@ Changelog Unreleased ---------- -* New API: ``GeneralBraOpKet`` for exact calculation of arbitrary ```` values. +* New API: ``GeneralBraOpKet`` for exact calculation of arbitrary ```` values. Can be used to calculate inner products, expectation values and arbitrary matrix elements. +* New feature: both ``GeneralState`` and ``GeneralBraOpKet`` admit circuits with parameterised gates. * New feature: ``GeneralState`` has a new method ``get_amplitude`` to obtain the amplitude of computational basis states. * New feature: ``GeneralState`` and ``CuTensorNetShotsBackend`` now support RNG seeds for sampling. -* Deprecated ``TensorNetwork`` object. It is still available for the sake backwards compatibility, but it has been removed from doc pages. +* Deprecated ``TensorNetwork`` object. It is still available for the sake of backwards compatibility, but it has been removed from doc pages. 0.9.0 (October 2024) --------------------- diff --git a/docs/modules/general_state.rst b/docs/modules/general_state.rst index bf557a0b..3d996e20 100644 --- a/docs/modules/general_state.rst +++ b/docs/modules/general_state.rst @@ -4,6 +4,7 @@ General state (exact) simulation .. automodule:: pytket.extensions.cutensornet.general_state .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralState() + .. automethod:: __init__ .. automethod:: get_statevector .. automethod:: get_amplitude @@ -12,6 +13,7 @@ General state (exact) simulation .. automethod:: destroy .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralBraOpKet() + .. automethod:: __init__ .. automethod:: contract .. automethod:: destroy diff --git a/docs/modules/structured_state.rst b/docs/modules/structured_state.rst index ce1ea6d8..aeafb5df 100644 --- a/docs/modules/structured_state.rst +++ b/docs/modules/structured_state.rst @@ -3,6 +3,9 @@ Structured state evolution .. automodule:: pytket.extensions.cutensornet.structured_state +Library handle +~~~~~~~~~~~~~~ + .. autoclass:: pytket.extensions.cutensornet.CuTensorNetHandle .. automethod:: destroy From beab8612f3db225dfefa1679d9b2d9bd71670a49 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 14:38:50 +0000 Subject: [PATCH 18/27] More fixes to docs --- docs/modules/general_state.rst | 2 - .../general_state/tensor_network_state.py | 71 +++++++++---------- 2 files changed, 32 insertions(+), 41 deletions(-) diff --git a/docs/modules/general_state.rst b/docs/modules/general_state.rst index 3d996e20..b0e2620f 100644 --- a/docs/modules/general_state.rst +++ b/docs/modules/general_state.rst @@ -5,7 +5,6 @@ General state (exact) simulation .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralState() - .. automethod:: __init__ .. automethod:: get_statevector .. automethod:: get_amplitude .. automethod:: expectation_value @@ -14,7 +13,6 @@ General state (exact) simulation .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralBraOpKet() - .. automethod:: __init__ .. automethod:: contract .. automethod:: destroy diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 65ef3c27..0eaaf730 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -39,9 +39,22 @@ class GeneralState: """Wrapper of cuTensorNet's NetworkState for exact simulation of states. + Constructs a tensor network for the output state of a pytket circuit. + The qubits are assumed to be initialised in the ``|0>`` state. + The object stores the *uncontracted* tensor network. + Note: Preferably used as ``with GeneralState(...) as state:`` so that GPU memory is automatically released after execution. + The ``circuit`` must not contain any ``CircBox`` or non-unitary command. + + Args: + circuit: A pytket circuit to be converted to a tensor network. + attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. + scratch_fraction: Optional. Fraction of free memory on GPU to allocate as + scratch space; value between 0 and 1. Defaults to ``0.8``. + loglevel: Internal logger output level. """ def __init__( @@ -51,22 +64,6 @@ def __init__( scratch_fraction: float = 0.8, loglevel: int = logging.WARNING, ) -> None: - """Constructs a tensor network for the output state of a pytket circuit. - - The qubits are assumed to be initialised in the ``|0>`` state. - The resulting object stores the *uncontracted* tensor network. - - Note: - The ``circuit`` must not contain any ``CircBox`` or non-unitary command. - - Args: - circuit: A pytket circuit to be converted to a tensor network. - attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and - their values. - scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space; value between 0 and 1. Defaults to ``0.8``. - loglevel: Internal logger output level. - """ self._logger = set_logger("GeneralState", loglevel) # Remove end-of-circuit measurements and keep track of them separately @@ -332,11 +329,29 @@ def __exit__(self, exc_type: Any, exc_value: Any, exc_tb: Any) -> None: class GeneralBraOpKet: - """Wrapper of cuTensorNet's NetworkState for exact simulation of ````. + """Constructs a tensor network for ````. + + The qubits in ``ket`` and ``bra`` are assumed to be initialised in the ``|0>`` + state. The object stores the *uncontracted* tensor network. Note: Preferably used as ``with GeneralBraOpKet(...) as braket:`` so that GPU memory is automatically released after execution. + The ``circuit`` must not contain any ``CircBox`` or non-unitary command. + The operator is provided when ``contract`` is called. + + Args: + bra: A pytket circuit describing the |bra> state. + ket: A pytket circuit describing the |ket> state. + attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and + their values. + scratch_fraction: Optional. Fraction of free memory on GPU to allocate as + scratch space; value between 0 and 1. Defaults to ``0.8``. + loglevel: Internal logger output level. + + Raises: + ValueError: If the circuits for ``ket`` or ``bra`` contain measurements. + ValueError: If the set of qubits of ``ket`` and ``bra`` do not match. """ def __init__( @@ -347,28 +362,6 @@ def __init__( scratch_fraction: float = 0.8, loglevel: int = logging.WARNING, ) -> None: - """Constructs a tensor network for ````. - - The qubits in ``ket`` and ``bra`` are assumed to be initialised in the ``|0>`` - state. The resulting object stores the *uncontracted* tensor network. - - Note: - The ``circuit`` must not contain any ``CircBox`` or non-unitary command. - The operator is provided when ``contract`` is called. - - Args: - bra: A pytket circuit describing the |bra> state. - ket: A pytket circuit describing the |ket> state. - attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and - their values. - scratch_fraction: Optional. Fraction of free memory on GPU to allocate as - scratch space; value between 0 and 1. Defaults to ``0.8``. - loglevel: Internal logger output level. - - Raises: - ValueError: If the circuits for ``ket`` or ``bra`` contain measurements. - ValueError: If the set of qubits of ``ket`` and ``bra`` do not match. - """ self._logger = set_logger("GeneralBraOpKet", loglevel) # Check that the circuits have the same qubits From 1c677921fee186e9a08b052923374446a75f58c8 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Wed, 23 Oct 2024 16:05:31 +0100 Subject: [PATCH 19/27] Finished polishing docs --- docs/modules/general_state.rst | 4 +- .../backends/cutensornet_backend.py | 2 +- .../general_state/tensor_network_state.py | 61 ++++++++++--------- 3 files changed, 35 insertions(+), 32 deletions(-) diff --git a/docs/modules/general_state.rst b/docs/modules/general_state.rst index b0e2620f..d9c219f5 100644 --- a/docs/modules/general_state.rst +++ b/docs/modules/general_state.rst @@ -5,10 +5,10 @@ General state (exact) simulation .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralState() - .. automethod:: get_statevector + .. automethod:: sample .. automethod:: get_amplitude + .. automethod:: get_statevector .. automethod:: expectation_value - .. automethod:: sample .. automethod:: destroy .. autoclass:: pytket.extensions.cutensornet.general_state.GeneralBraOpKet() diff --git a/pytket/extensions/cutensornet/backends/cutensornet_backend.py b/pytket/extensions/cutensornet/backends/cutensornet_backend.py index 4d8a0386..397f51c3 100644 --- a/pytket/extensions/cutensornet/backends/cutensornet_backend.py +++ b/pytket/extensions/cutensornet/backends/cutensornet_backend.py @@ -281,9 +281,9 @@ def process_circuits( n_shots: Number of shots in case of shot-based calculation. Optionally, this can be a list of shots specifying the number of shots for each circuit separately. + valid_check: Whether to check for circuit correctness. seed: An optional RNG seed. Different calls to ``process_circuits`` with the same seed will generate the same list of shot outcomes. - valid_check: Whether to check for circuit correctness. tnconfig: Optional. A dict of cuTensorNet ``TNConfig`` keys and their values. scratch_fraction: Optional. Fraction of free memory on GPU to allocate as diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 0eaaf730..f5462301 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -46,15 +46,17 @@ class GeneralState: Note: Preferably used as ``with GeneralState(...) as state:`` so that GPU memory is automatically released after execution. + The ``circuit`` must not contain any ``CircBox`` or non-unitary command. Args: - circuit: A pytket circuit to be converted to a tensor network. + circuit: A pytket circuit to be converted into a tensor network. attributes: Optional. A dict of cuTensorNet ``TNConfig`` keys and their values. scratch_fraction: Optional. Fraction of free memory on GPU to allocate as scratch space; value between 0 and 1. Defaults to ``0.8``. - loglevel: Internal logger output level. + loglevel: Internal logger output level. Use 30 for warnings only, 20 for + verbose and 10 for debug mode. """ def __init__( @@ -138,8 +140,8 @@ def get_statevector( """Contracts the circuit and returns the final statevector. Args: - symbol_map: A dictionary where each element of ``sef.free_symbols`` is - assigned a real number. + symbol_map: A dictionary where each element of the pytket circuit's + ``.free_symbols()`` is assigned a real number. on_host: Optional. If ``True``, converts cupy ``ndarray`` to numpy ``ndarray``, copying it to host device (CPU). Defaults to ``True``. Returns: @@ -147,7 +149,7 @@ def get_statevector( host device (CPU). Arrays are returned in a 1D shape. Raises: - ValueError: If not every free symbol in the circuits is assigned a + ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ if symbol_map is not None: @@ -169,21 +171,21 @@ def get_amplitude( ) -> complex: """Returns the amplitude of the chosen computational state. - Notes: + Note: The result is equivalent to ``state.get_statevector[b]``, but this method is faster when querying a single amplitude (or just a few). Args: state: The integer whose bitstring describes the computational state. The qubits in the bitstring are in increasing lexicographic order. - symbol_map: A dictionary where each element of ``sef.free_symbols`` is - assigned a real number. + symbol_map: A dictionary where each element of the pytket circuit's + ``.free_symbols()`` is assigned a real number. Returns: The amplitude of the computational state in ``self``. Raises: - ValueError: If not every free symbol in the circuits is assigned a + ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ if symbol_map is not None: @@ -204,16 +206,16 @@ def expectation_value( """Calculates the expectation value of the given operator. Args: - operator: The operator whose expectation value is to be measured. - symbol_map: A dictionary where each element of ``sef.free_symbols`` is - assigned a real number. + operator: The operator whose expectation value is to be calculated. + symbol_map: A dictionary where each element of the pytket circuit's + ``.free_symbols()`` is assigned a real number. Returns: The expectation value. Raises: ValueError: If the operator acts on qubits not present in the circuit. - ValueError: If not every free symbol in the circuits is assigned a + ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ if symbol_map is not None: @@ -252,8 +254,8 @@ def sample( Args: n_shots: The number of samples to obtain. - symbol_map: A dictionary where each element of ``sef.free_symbols`` is - assigned a real number. + symbol_map: A dictionary where each element of the pytket circuit's + ``.free_symbols()`` is assigned a real number. seed: An optional RNG seed. Different calls to ``sample`` with the same seed will generate the same list of shot outcomes. @@ -262,7 +264,7 @@ def sample( Raises: ValueError: If the circuit contains no measurements. - ValueError: If not every free symbol in the circuits is assigned a + ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ if symbol_map is not None: @@ -313,10 +315,9 @@ def sample( def destroy(self) -> None: """Destroy the tensor network and free up GPU memory. - Note: - The preferred approach is to use a context manager as in - ``with GeneralState(...) as state:``. Otherwise, the user must release - memory explicitly by calling ``destroy()``. + The preferred approach is to use a context manager as in + ``with GeneralState(...) as state:``. Otherwise, the user must release + memory explicitly by calling ``destroy()``. """ self._logger.debug("Freeing memory of GeneralState") self.tn_state.free() @@ -337,7 +338,9 @@ class GeneralBraOpKet: Note: Preferably used as ``with GeneralBraOpKet(...) as braket:`` so that GPU memory is automatically released after execution. - The ``circuit`` must not contain any ``CircBox`` or non-unitary command. + + The circuits must not contain any ``CircBox`` or non-unitary command. + The operator is provided when ``contract`` is called. Args: @@ -347,7 +350,8 @@ class GeneralBraOpKet: their values. scratch_fraction: Optional. Fraction of free memory on GPU to allocate as scratch space; value between 0 and 1. Defaults to ``0.8``. - loglevel: Internal logger output level. + loglevel: Internal logger output level. Use 30 for warnings only, 20 for + verbose and 10 for debug mode. Raises: ValueError: If the circuits for ``ket`` or ``bra`` contain measurements. @@ -505,11 +509,11 @@ def contract( Args: operator: A pytket ``QubitPauliOperator`` describing the operator. If not given, then the identity operator is used, so it computes inner product. - symbol_map: A dictionary where each element of ``sef.free_symbols`` is - assigned a real number. + symbol_map: A dictionary where each element of both pytket circuits' + ``.free_symbols()`` is assigned a real number. Returns: - The value of ```` + The value of ````. Raises: ValueError: If ``operator`` acts on qubits that are not in the circuits. @@ -576,10 +580,9 @@ def contract( def destroy(self) -> None: """Destroy the tensor network and free up GPU memory. - Note: - The preferred approach is to use a context manager as in - ``with GeneralBraOpKet(...) as braket:``. Otherwise, the user must release - memory explicitly by calling ``destroy()``. + The preferred approach is to use a context manager as in + ``with GeneralBraOpKet(...) as braket:``. Otherwise, the user must release + memory explicitly by calling ``destroy()``. """ self._logger.debug("Freeing memory of GeneralBraOpKet") self.tn.free() From 19344aac64b55c9038184d8a217251e9999009fc Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Thu, 24 Oct 2024 11:56:59 +0000 Subject: [PATCH 20/27] Fixed a bug with parameterised circuits where the user does not provide a symbol_map --- .../general_state/tensor_network_state.py | 20 +++++++++---------- 1 file changed, 9 insertions(+), 11 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index f5462301..cf234aee 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -152,8 +152,7 @@ def get_statevector( ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ - if symbol_map is not None: - _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) self._logger.debug("(Statevector) contracting the TN") state_vector = self.tn_state.compute_state_vector() @@ -188,8 +187,7 @@ def get_amplitude( ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ - if symbol_map is not None: - _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) self._logger.debug("(Amplitude) contracting the TN") bitstring = format(state, "b").zfill(self.n_qubits) @@ -218,8 +216,7 @@ def expectation_value( ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ - if symbol_map is not None: - _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) self._logger.debug("(Expectation value) converting operator to NetworkOperator") @@ -267,8 +264,7 @@ def sample( ValueError: If not every free symbol in the circuit is assigned a value in ``symbol_map``. """ - if symbol_map is not None: - _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) + _update_tensors(self.tn_state, self._symbolic_ops, symbol_map) num_measurements = len(self._measurements) if num_measurements == 0: @@ -518,8 +514,7 @@ def contract( Raises: ValueError: If ``operator`` acts on qubits that are not in the circuits. """ - if symbol_map is not None: - _update_tensors(self.tn, self._symbolic_ops, symbol_map) + _update_tensors(self.tn, self._symbolic_ops, symbol_map) paulis = ["I", "X", "Y", "Z"] pauli_matrix = { @@ -641,7 +636,7 @@ def _remove_meas_and_implicit_swaps(circ: Circuit) -> Tuple[Circuit, Dict[Qubit, def _update_tensors( tn: NetworkState, symbolic_ops: dict[int, Op], - symbol_map: dict[Symbol, float], + symbol_map: Optional[dict[Symbol, float]], ) -> None: """Updates the tensors with the specified values for symbols. @@ -654,6 +649,9 @@ def _update_tensors( ValueError: If not every free symbol in the circuits is assigned a value in ``symbol_map``. """ + if symbol_map is None: + symbol_map = dict() + for tensor_id, op in symbolic_ops.items(): subs_params = op.params.copy() if any(symb not in symbol_map.keys() for symb in op.free_symbols()): From a522831f66a4f0c135fca9c7a4e64923f2d95c8d Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Thu, 24 Oct 2024 12:46:24 +0000 Subject: [PATCH 21/27] Added a notebook tutorial --- examples/general_state_tutorial.ipynb | 1 + examples/python/general_state_tutorial.py | 178 ++++++++++++++++++++++ 2 files changed, 179 insertions(+) create mode 100644 examples/general_state_tutorial.ipynb create mode 100644 examples/python/general_state_tutorial.py diff --git a/examples/general_state_tutorial.ipynb b/examples/general_state_tutorial.ipynb new file mode 100644 index 00000000..d8b559e4 --- /dev/null +++ b/examples/general_state_tutorial.ipynb @@ -0,0 +1 @@ +{"cells": [{"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["import numpy as np\n", "from sympy import Symbol\n", "from scipy.stats import unitary_group # type: ignore\n", "from pytket.circuit import Circuit, OpType, Unitary2qBox, Qubit, Bit\n", "from pytket.passes import DecomposeBoxes\n", "from pytket.utils import QubitPauliOperator\n", "from pytket._tket.pauli import Pauli, QubitPauliString\n", "from pytket.circuit.display import render_circuit_jupyter"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["from pytket.extensions.cutensornet.general_state import (\n", " GeneralState,\n", " GeneralBraOpKet,\n", ")\n", "from pytket.extensions.cutensornet.backends import CuTensorNetShotsBackend"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Introduction
\n", "This notebook is a guide on how to use the features provided in the `general_state` submodule of pytket-cutensornet. This submodule is a thin wrapper of CuTensorNet's `NetworkState`, allowing users to convert pytket circuits into tensor networks and use CuTensorNet's contraction path optimisation algorithm.
\n", "All simulations realised with this submodule are *exact*. Once the pytket circuit has been converted to a tensor network, the computation has two steps:
\n", " 1. *Contraction path optimisation*. Attempts to find an order of contracting pairs of tensors in which the the total number of FLOPs is minimised. No operation on the tensor network occurs at this point. Runs on CPU.
\n", " 2. *Tensor network contraction*. Uses the ordering of contractions found in the previous step evaluate the tensor network. Runs on GPU.
\n", "
\n", "**Reference**: The original contraction path optimisation algorithm that NVIDIA implemented on CuTensorNet: https://arxiv.org/abs/2002.01935"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralState`
\n", "The class `GeneralState` is used to convert a circuit into a tensor network and query information from the final state. Let's walk through a simple example."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["my_circ = Circuit(5)\n", "my_circ.CX(3, 4)\n", "my_circ.H(2)\n", "my_circ.CZ(0, 1)\n", "my_circ.ZZPhase(0.1, 4, 3)\n", "my_circ.TK2(0.3, 0.5, 0.7, 2, 1)\n", "my_circ.Ry(0.2, 0)\n", "my_circ.measure_all()"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["render_circuit_jupyter(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The first step is to convert our pytket circuit into a tensor network. This is straightforward:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state = GeneralState(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The variable `tn_state` now holds a tensor network representation of `my_circ`.
\n", "**Note**: Circuits must not have mid-circuit measurements or classical logic. The measurements at the end of the circuit are stripped and only considered when calling `tn_state.sample(n_shots)`.
\n", "We can now query information from the state. For instance, let's calculate the probability of in the qubits 0 and 3 agreeing in their outcome."]}, {"cell_type": "markdown", "metadata": {}, "source": ["First, let's generate `|x>` computational basis states where `q[0]` and `q[3]` agree on their values. We can do this with some bitwise operators and list comprehension.
\n", "**Note**: Remember that pytket uses \"increasing lexicographic order\" (ILO) for qubits, so `q[0]` is the most significant bit."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["selected_states = [\n", " x for x in range(2**my_circ.n_qubits) if ( # Iterate over all possible states\n", " x & int(\"10000\", 2) == 0 and x & int(\"00010\", 2) == 0 # both qubits are 0 or...\n", " or x & int(\"10000\", 2) != 0 and x & int(\"00010\", 2) != 0 # both qubits are 1\n", " )\n", "]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can now query the amplitude of all of these states and calculate the probability by summing their squared absolute values."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["amplitudes = []\n", "for x in selected_states:\n", " amplitudes.append(\n", " tn_state.get_amplitude(x)\n", " )\n", "probability = sum(abs(a)**2 for a in amplitudes)\n", "print(f\"Probability: {probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Of course, calculating probabilities by considering the amplitudes of all relevant states is not efficient in general, since we may need to calculate a number of amplitudes that scales exponentially with the number of qubits. An alternative is to use expectation values. In particular, all of the states in `selected_states` are +1 eigenvectors of the `ZIIZI` observable and, hence, we can calculate the probability `p` by solving the equation ` = (+1)p + (-1)(1-p)` using the fact that `ZIIZI` only has +1 and -1 eigenvalues."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_ZIIZI = QubitPauliString(my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I])\n", "observable = QubitPauliOperator({string_ZIIZI: 1.0})\n", "expectation_val = tn_state.expectation_value(observable).real\n", "exp_probability = (expectation_val + 1) / 2\n", "assert np.isclose(probability, exp_probability, atol=0.0001)\n", "print(f\"Probability: {exp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Alternatively, we can estimate the probability by sampling."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 100000\n", "outcomes = tn_state.sample(n_shots)\n", "hit_count = 0\n", "for bit_tuple, count in outcomes.get_counts().items():\n", " if bit_tuple[0] == bit_tuple[3]:\n", " hit_count += count\n", "samp_probability = hit_count / n_shots\n", "assert np.isclose(probability, samp_probability, atol=0.01)\n", "print(f\"Probability: {samp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["When we finish doing computations with the `tn_state` we must destroy it to free GPU memory."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state.destroy()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["To avoid forgetting this final step, we recommend users call `GeneralState` (and `GeneralBraOpKet`) as context managers:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralState(my_circ) as my_state:\n", " expectation_val = my_state.expectation_value(observable)\n", "print(expectation_val)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Using this syntax, `my_state` is automatically destroyed when the code exists the `with ...` block."]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Parameterised circuits
\n", "Circuits that only differ on the parameters of their gates have the same tensor network topology and, hence, we may use the same contraction path for all of them."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["a, b, c = Symbol(\"a\"), Symbol(\"b\"), Symbol(\"c\")\n", "param_circ1 = Circuit(5)\n", "param_circ1.Ry(a, 3).Ry(0.27, 4).CX(4, 3).Ry(b, 2).Ry(0.21, 3)\n", "param_circ1.Ry(0.12, 0).Ry(a, 1)\n", "param_circ1.add_gate(OpType.CnX, [0, 1, 4]).add_gate(OpType.CnX, [4, 1, 3])\n", "param_circ1.X(0).X(1).add_gate(OpType.CnY, [0, 1, 2]).add_gate(OpType.CnY, [0, 4, 3]).X(\n", " 0\n", ").X(1)\n", "param_circ1.Ry(-b, 0).Ry(-c, 1)\n", "render_circuit_jupyter(param_circ1)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can pass a parameterised circuit to `GeneralState`. The value of the parameters is provided when calling methods of `GeneralState`. The contraction path is automatically reused on different calls to the same method."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_circs = 5\n", "with GeneralState(param_circ1) as param_state:\n", " for i in range(n_circs):\n", " symbol_map = {s: np.random.random() for s in [a, b, c]}\n", " exp_val = param_state.expectation_value(observable, symbol_map=symbol_map)\n", " print(f\"Expectation value for circuit {i}: {exp_val.real}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralBraOpKet`
\n", "The `GeneralBraOpKet` can be used to calculate any number that can be represented as the result of some `` where `|bra>` and `|ket>` are the final states of pytket circuits, and `op` is a `QubitPauliOperator`. The circuits for `|bra>` and `|ket>` need not be the same."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["x, y, z = Symbol(\"x\"), Symbol(\"y\"), Symbol(\"z\")\n", "param_circ2 = Circuit(5)\n", "param_circ2.H(0)\n", "param_circ2.S(1)\n", "param_circ2.Rz(x * z, 2)\n", "param_circ2.Ry(y + x, 3)\n", "param_circ2.TK1(x, y, z, 4)\n", "param_circ2.TK2(z - y, z - x, (x + y) * z, 1, 3)\n", "symbol_map = {a: 2.1, b: 1.3, c: 0.7, x: 3.0, y: 1.6, z: -8.3}"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can calculate inner products by providing no `op`:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " inner_prod = braket.contract(symbol_map=symbol_map)\n", "with GeneralBraOpKet(bra=param_circ1, ket=param_circ2) as braket:\n", " inner_prod_conj = braket.contract(symbol_map=symbol_map)\n", "assert np.isclose(np.conj(inner_prod), inner_prod_conj)\n", "print(f\" = {inner_prod}\")\n", "print(f\" = {inner_prod_conj}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["And we are not constrained to Hermitian operators:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_XZIXX = QubitPauliString(param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X])\n", "string_IZZYX = QubitPauliString(param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X])\n", "string_ZIZXY = QubitPauliString(param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y])\n", "operator = QubitPauliOperator({string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j+0.3})\n", "with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " value = braket.contract(operator, symbol_map=symbol_map)\n", "print(value)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Backends
\n", "We provide a pytket `Backend` to obtain shots using `GeneralState`."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's consider a more challenging circuit"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["def random_circuit(n_qubits: int, n_layers: int) -> Circuit:\n", " \"\"\"Random quantum volume circuit.\"\"\"\n", " c = Circuit(n_qubits, n_qubits)\n", " for _ in range(n_layers):\n", " qubits = np.random.permutation([i for i in range(n_qubits)])\n", " qubit_pairs = [[qubits[i], qubits[i + 1]] for i in range(0, n_qubits - 1, 2)]\n", " for pair in qubit_pairs:\n", " # Generate random 4x4 unitary matrix.\n", " SU4 = unitary_group.rvs(4) # random unitary in SU4\n", " SU4 = SU4 / (np.linalg.det(SU4) ** 0.25)\n", " SU4 = np.matrix(SU4)\n", " c.add_unitary2qbox(Unitary2qBox(SU4), *pair)\n", " DecomposeBoxes().apply(c)\n", " return c"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's measure only three of the qubits.
\n", "**Note**: The complexity of this simulation increases exponentially with the number of qubits measured. Other factors leading to intractability are circuit depth and qubit connectivity."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 1000\n", "quantum_vol_circ = random_circuit(n_qubits=40, n_layers=5)\n", "quantum_vol_circ.Measure(Qubit(0), Bit(0))\n", "quantum_vol_circ.Measure(Qubit(1), Bit(1))\n", "quantum_vol_circ.Measure(Qubit(2), Bit(2))"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The `CuTensorNetShotsBackend` is used in the same way as any other pytket `Backend`."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["backend = CuTensorNetShotsBackend()\n", "compiled_circ = backend.get_compiled_circuit(quantum_vol_circ)\n", "results = backend.run_circuit(compiled_circ, n_shots=n_shots)\n", "print(results.get_counts())"]}], "metadata": {"kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"}, "language_info": {"codemirror_mode": {"name": "ipython", "version": 3}, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.4"}}, "nbformat": 4, "nbformat_minor": 2} \ No newline at end of file diff --git a/examples/python/general_state_tutorial.py b/examples/python/general_state_tutorial.py new file mode 100644 index 00000000..221f6347 --- /dev/null +++ b/examples/python/general_state_tutorial.py @@ -0,0 +1,178 @@ +import numpy as np +from sympy import Symbol +from scipy.stats import unitary_group # type: ignore +from pytket.circuit import Circuit, OpType, Unitary2qBox, Qubit, Bit +from pytket.passes import DecomposeBoxes +from pytket.utils import QubitPauliOperator +from pytket._tket.pauli import Pauli, QubitPauliString +from pytket.circuit.display import render_circuit_jupyter + +from pytket.extensions.cutensornet.general_state import ( + GeneralState, + GeneralBraOpKet, +) +from pytket.extensions.cutensornet.backends import CuTensorNetShotsBackend + +# # Introduction +# This notebook is a guide on how to use the features provided in the `general_state` submodule of pytket-cutensornet. This submodule is a thin wrapper of CuTensorNet's `NetworkState`, allowing users to convert pytket circuits into tensor networks and use CuTensorNet's contraction path optimisation algorithm. +# All simulations realised with this submodule are *exact*. Once the pytket circuit has been converted to a tensor network, the computation has two steps: +# 1. *Contraction path optimisation*. Attempts to find an order of contracting pairs of tensors in which the the total number of FLOPs is minimised. No operation on the tensor network occurs at this point. Runs on CPU. +# 2. *Tensor network contraction*. Uses the ordering of contractions found in the previous step evaluate the tensor network. Runs on GPU. +# +# **Reference**: The original contraction path optimisation algorithm that NVIDIA implemented on CuTensorNet: https://arxiv.org/abs/2002.01935 + +# # `GeneralState` +# The class `GeneralState` is used to convert a circuit into a tensor network and query information from the final state. Let's walk through a simple example. + +my_circ = Circuit(5) +my_circ.CX(3, 4) +my_circ.H(2) +my_circ.CZ(0, 1) +my_circ.ZZPhase(0.1, 4, 3) +my_circ.TK2(0.3, 0.5, 0.7, 2, 1) +my_circ.Ry(0.2, 0) +my_circ.measure_all() + +render_circuit_jupyter(my_circ) + +# The first step is to convert our pytket circuit into a tensor network. This is straightforward: +tn_state = GeneralState(my_circ) + +# The variable `tn_state` now holds a tensor network representation of `my_circ`. +# **Note**: Circuits must not have mid-circuit measurements or classical logic. The measurements at the end of the circuit are stripped and only considered when calling `tn_state.sample(n_shots)`. +# We can now query information from the state. For instance, let's calculate the probability of in the qubits 0 and 3 agreeing in their outcome. + +# First, let's generate `|x>` computational basis states where `q[0]` and `q[3]` agree on their values. We can do this with some bitwise operators and list comprehension. +# **Note**: Remember that pytket uses "increasing lexicographic order" (ILO) for qubits, so `q[0]` is the most significant bit. +selected_states = [ + x for x in range(2**my_circ.n_qubits) if ( # Iterate over all possible states + x & int("10000", 2) == 0 and x & int("00010", 2) == 0 # both qubits are 0 or... + or x & int("10000", 2) != 0 and x & int("00010", 2) != 0 # both qubits are 1 + ) +] + +# We can now query the amplitude of all of these states and calculate the probability by summing their squared absolute values. +amplitudes = [] +for x in selected_states: + amplitudes.append( + tn_state.get_amplitude(x) + ) +probability = sum(abs(a)**2 for a in amplitudes) +print(f"Probability: {probability}") + +# Of course, calculating probabilities by considering the amplitudes of all relevant states is not efficient in general, since we may need to calculate a number of amplitudes that scales exponentially with the number of qubits. An alternative is to use expectation values. In particular, all of the states in `selected_states` are +1 eigenvectors of the `ZIIZI` observable and, hence, we can calculate the probability `p` by solving the equation ` = (+1)p + (-1)(1-p)` using the fact that `ZIIZI` only has +1 and -1 eigenvalues. +string_ZIIZI = QubitPauliString(my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I]) +observable = QubitPauliOperator({string_ZIIZI: 1.0}) +expectation_val = tn_state.expectation_value(observable).real +exp_probability = (expectation_val + 1) / 2 +assert np.isclose(probability, exp_probability, atol=0.0001) +print(f"Probability: {exp_probability}") + +# Alternatively, we can estimate the probability by sampling. +n_shots = 100000 +outcomes = tn_state.sample(n_shots) +hit_count = 0 +for bit_tuple, count in outcomes.get_counts().items(): + if bit_tuple[0] == bit_tuple[3]: + hit_count += count +samp_probability = hit_count / n_shots +assert np.isclose(probability, samp_probability, atol=0.01) +print(f"Probability: {samp_probability}") + +# When we finish doing computations with the `tn_state` we must destroy it to free GPU memory. +tn_state.destroy() + +# To avoid forgetting this final step, we recommend users call `GeneralState` (and `GeneralBraOpKet`) as context managers: +with GeneralState(my_circ) as my_state: + expectation_val = my_state.expectation_value(observable) +print(expectation_val) + +# Using this syntax, `my_state` is automatically destroyed when the code exists the `with ...` block. + +# # Parameterised circuits +# Circuits that only differ on the parameters of their gates have the same tensor network topology and, hence, we may use the same contraction path for all of them. +a, b, c = Symbol("a"), Symbol("b"), Symbol("c") +param_circ1 = Circuit(5) +param_circ1.Ry(a, 3).Ry(0.27, 4).CX(4, 3).Ry(b, 2).Ry(0.21, 3) +param_circ1.Ry(0.12, 0).Ry(a, 1) +param_circ1.add_gate(OpType.CnX, [0, 1, 4]).add_gate(OpType.CnX, [4, 1, 3]) +param_circ1.X(0).X(1).add_gate(OpType.CnY, [0, 1, 2]).add_gate(OpType.CnY, [0, 4, 3]).X( + 0 +).X(1) +param_circ1.Ry(-b, 0).Ry(-c, 1) +render_circuit_jupyter(param_circ1) + +# We can pass a parameterised circuit to `GeneralState`. The value of the parameters is provided when calling methods of `GeneralState`. The contraction path is automatically reused on different calls to the same method. +n_circs = 5 +with GeneralState(param_circ1) as param_state: + for i in range(n_circs): + symbol_map = {s: np.random.random() for s in [a, b, c]} + exp_val = param_state.expectation_value(observable, symbol_map=symbol_map) + print(f"Expectation value for circuit {i}: {exp_val.real}") + + +# # `GeneralBraOpKet` +# The `GeneralBraOpKet` can be used to calculate any number that can be represented as the result of some `` where `|bra>` and `|ket>` are the final states of pytket circuits, and `op` is a `QubitPauliOperator`. The circuits for `|bra>` and `|ket>` need not be the same. +x, y, z = Symbol("x"), Symbol("y"), Symbol("z") +param_circ2 = Circuit(5) +param_circ2.H(0) +param_circ2.S(1) +param_circ2.Rz(x * z, 2) +param_circ2.Ry(y + x, 3) +param_circ2.TK1(x, y, z, 4) +param_circ2.TK2(z - y, z - x, (x + y) * z, 1, 3) +symbol_map = {a: 2.1, b: 1.3, c: 0.7, x: 3.0, y: 1.6, z: -8.3} + +# We can calculate inner products by providing no `op`: +with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket: + inner_prod = braket.contract(symbol_map=symbol_map) +with GeneralBraOpKet(bra=param_circ1, ket=param_circ2) as braket: + inner_prod_conj = braket.contract(symbol_map=symbol_map) +assert np.isclose(np.conj(inner_prod), inner_prod_conj) +print(f" = {inner_prod}") +print(f" = {inner_prod_conj}") + +# And we are not constrained to Hermitian operators: +string_XZIXX = QubitPauliString(param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X]) +string_IZZYX = QubitPauliString(param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X]) +string_ZIZXY = QubitPauliString(param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y]) +operator = QubitPauliOperator({string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j+0.3}) +with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket: + value = braket.contract(operator, symbol_map=symbol_map) +print(value) + +# # Backends +# We provide a pytket `Backend` to obtain shots using `GeneralState`. + +# Let's consider a more challenging circuit +def random_circuit(n_qubits: int, n_layers: int) -> Circuit: + """Random quantum volume circuit.""" + c = Circuit(n_qubits, n_qubits) + + for _ in range(n_layers): + qubits = np.random.permutation([i for i in range(n_qubits)]) + qubit_pairs = [[qubits[i], qubits[i + 1]] for i in range(0, n_qubits - 1, 2)] + + for pair in qubit_pairs: + # Generate random 4x4 unitary matrix. + SU4 = unitary_group.rvs(4) # random unitary in SU4 + SU4 = SU4 / (np.linalg.det(SU4) ** 0.25) + SU4 = np.matrix(SU4) + c.add_unitary2qbox(Unitary2qBox(SU4), *pair) + + DecomposeBoxes().apply(c) + return c + +# Let's measure only three of the qubits. +# **Note**: The complexity of this simulation increases exponentially with the number of qubits measured. Other factors leading to intractability are circuit depth and qubit connectivity. +n_shots = 1000 +quantum_vol_circ = random_circuit(n_qubits=40, n_layers=5) +quantum_vol_circ.Measure(Qubit(0), Bit(0)) +quantum_vol_circ.Measure(Qubit(1), Bit(1)) +quantum_vol_circ.Measure(Qubit(2), Bit(2)) + +# The `CuTensorNetShotsBackend` is used in the same way as any other pytket `Backend`. +backend = CuTensorNetShotsBackend() +compiled_circ = backend.get_compiled_circuit(quantum_vol_circ) +results = backend.run_circuit(compiled_circ, n_shots=n_shots) +print(results.get_counts()) \ No newline at end of file From 509e12858fa4d19140cc0eecdddc9bf614d41b05 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Thu, 24 Oct 2024 12:48:41 +0000 Subject: [PATCH 22/27] Linter --- examples/general_state_tutorial.ipynb | 2 +- examples/python/general_state_tutorial.py | 40 +++++++++++++++-------- 2 files changed, 28 insertions(+), 14 deletions(-) diff --git a/examples/general_state_tutorial.ipynb b/examples/general_state_tutorial.ipynb index d8b559e4..9a041f99 100644 --- a/examples/general_state_tutorial.ipynb +++ b/examples/general_state_tutorial.ipynb @@ -1 +1 @@ -{"cells": [{"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["import numpy as np\n", "from sympy import Symbol\n", "from scipy.stats import unitary_group # type: ignore\n", "from pytket.circuit import Circuit, OpType, Unitary2qBox, Qubit, Bit\n", "from pytket.passes import DecomposeBoxes\n", "from pytket.utils import QubitPauliOperator\n", "from pytket._tket.pauli import Pauli, QubitPauliString\n", "from pytket.circuit.display import render_circuit_jupyter"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["from pytket.extensions.cutensornet.general_state import (\n", " GeneralState,\n", " GeneralBraOpKet,\n", ")\n", "from pytket.extensions.cutensornet.backends import CuTensorNetShotsBackend"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Introduction
\n", "This notebook is a guide on how to use the features provided in the `general_state` submodule of pytket-cutensornet. This submodule is a thin wrapper of CuTensorNet's `NetworkState`, allowing users to convert pytket circuits into tensor networks and use CuTensorNet's contraction path optimisation algorithm.
\n", "All simulations realised with this submodule are *exact*. Once the pytket circuit has been converted to a tensor network, the computation has two steps:
\n", " 1. *Contraction path optimisation*. Attempts to find an order of contracting pairs of tensors in which the the total number of FLOPs is minimised. No operation on the tensor network occurs at this point. Runs on CPU.
\n", " 2. *Tensor network contraction*. Uses the ordering of contractions found in the previous step evaluate the tensor network. Runs on GPU.
\n", "
\n", "**Reference**: The original contraction path optimisation algorithm that NVIDIA implemented on CuTensorNet: https://arxiv.org/abs/2002.01935"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralState`
\n", "The class `GeneralState` is used to convert a circuit into a tensor network and query information from the final state. Let's walk through a simple example."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["my_circ = Circuit(5)\n", "my_circ.CX(3, 4)\n", "my_circ.H(2)\n", "my_circ.CZ(0, 1)\n", "my_circ.ZZPhase(0.1, 4, 3)\n", "my_circ.TK2(0.3, 0.5, 0.7, 2, 1)\n", "my_circ.Ry(0.2, 0)\n", "my_circ.measure_all()"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["render_circuit_jupyter(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The first step is to convert our pytket circuit into a tensor network. This is straightforward:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state = GeneralState(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The variable `tn_state` now holds a tensor network representation of `my_circ`.
\n", "**Note**: Circuits must not have mid-circuit measurements or classical logic. The measurements at the end of the circuit are stripped and only considered when calling `tn_state.sample(n_shots)`.
\n", "We can now query information from the state. For instance, let's calculate the probability of in the qubits 0 and 3 agreeing in their outcome."]}, {"cell_type": "markdown", "metadata": {}, "source": ["First, let's generate `|x>` computational basis states where `q[0]` and `q[3]` agree on their values. We can do this with some bitwise operators and list comprehension.
\n", "**Note**: Remember that pytket uses \"increasing lexicographic order\" (ILO) for qubits, so `q[0]` is the most significant bit."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["selected_states = [\n", " x for x in range(2**my_circ.n_qubits) if ( # Iterate over all possible states\n", " x & int(\"10000\", 2) == 0 and x & int(\"00010\", 2) == 0 # both qubits are 0 or...\n", " or x & int(\"10000\", 2) != 0 and x & int(\"00010\", 2) != 0 # both qubits are 1\n", " )\n", "]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can now query the amplitude of all of these states and calculate the probability by summing their squared absolute values."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["amplitudes = []\n", "for x in selected_states:\n", " amplitudes.append(\n", " tn_state.get_amplitude(x)\n", " )\n", "probability = sum(abs(a)**2 for a in amplitudes)\n", "print(f\"Probability: {probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Of course, calculating probabilities by considering the amplitudes of all relevant states is not efficient in general, since we may need to calculate a number of amplitudes that scales exponentially with the number of qubits. An alternative is to use expectation values. In particular, all of the states in `selected_states` are +1 eigenvectors of the `ZIIZI` observable and, hence, we can calculate the probability `p` by solving the equation ` = (+1)p + (-1)(1-p)` using the fact that `ZIIZI` only has +1 and -1 eigenvalues."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_ZIIZI = QubitPauliString(my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I])\n", "observable = QubitPauliOperator({string_ZIIZI: 1.0})\n", "expectation_val = tn_state.expectation_value(observable).real\n", "exp_probability = (expectation_val + 1) / 2\n", "assert np.isclose(probability, exp_probability, atol=0.0001)\n", "print(f\"Probability: {exp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Alternatively, we can estimate the probability by sampling."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 100000\n", "outcomes = tn_state.sample(n_shots)\n", "hit_count = 0\n", "for bit_tuple, count in outcomes.get_counts().items():\n", " if bit_tuple[0] == bit_tuple[3]:\n", " hit_count += count\n", "samp_probability = hit_count / n_shots\n", "assert np.isclose(probability, samp_probability, atol=0.01)\n", "print(f\"Probability: {samp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["When we finish doing computations with the `tn_state` we must destroy it to free GPU memory."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state.destroy()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["To avoid forgetting this final step, we recommend users call `GeneralState` (and `GeneralBraOpKet`) as context managers:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralState(my_circ) as my_state:\n", " expectation_val = my_state.expectation_value(observable)\n", "print(expectation_val)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Using this syntax, `my_state` is automatically destroyed when the code exists the `with ...` block."]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Parameterised circuits
\n", "Circuits that only differ on the parameters of their gates have the same tensor network topology and, hence, we may use the same contraction path for all of them."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["a, b, c = Symbol(\"a\"), Symbol(\"b\"), Symbol(\"c\")\n", "param_circ1 = Circuit(5)\n", "param_circ1.Ry(a, 3).Ry(0.27, 4).CX(4, 3).Ry(b, 2).Ry(0.21, 3)\n", "param_circ1.Ry(0.12, 0).Ry(a, 1)\n", "param_circ1.add_gate(OpType.CnX, [0, 1, 4]).add_gate(OpType.CnX, [4, 1, 3])\n", "param_circ1.X(0).X(1).add_gate(OpType.CnY, [0, 1, 2]).add_gate(OpType.CnY, [0, 4, 3]).X(\n", " 0\n", ").X(1)\n", "param_circ1.Ry(-b, 0).Ry(-c, 1)\n", "render_circuit_jupyter(param_circ1)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can pass a parameterised circuit to `GeneralState`. The value of the parameters is provided when calling methods of `GeneralState`. The contraction path is automatically reused on different calls to the same method."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_circs = 5\n", "with GeneralState(param_circ1) as param_state:\n", " for i in range(n_circs):\n", " symbol_map = {s: np.random.random() for s in [a, b, c]}\n", " exp_val = param_state.expectation_value(observable, symbol_map=symbol_map)\n", " print(f\"Expectation value for circuit {i}: {exp_val.real}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralBraOpKet`
\n", "The `GeneralBraOpKet` can be used to calculate any number that can be represented as the result of some `` where `|bra>` and `|ket>` are the final states of pytket circuits, and `op` is a `QubitPauliOperator`. The circuits for `|bra>` and `|ket>` need not be the same."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["x, y, z = Symbol(\"x\"), Symbol(\"y\"), Symbol(\"z\")\n", "param_circ2 = Circuit(5)\n", "param_circ2.H(0)\n", "param_circ2.S(1)\n", "param_circ2.Rz(x * z, 2)\n", "param_circ2.Ry(y + x, 3)\n", "param_circ2.TK1(x, y, z, 4)\n", "param_circ2.TK2(z - y, z - x, (x + y) * z, 1, 3)\n", "symbol_map = {a: 2.1, b: 1.3, c: 0.7, x: 3.0, y: 1.6, z: -8.3}"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can calculate inner products by providing no `op`:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " inner_prod = braket.contract(symbol_map=symbol_map)\n", "with GeneralBraOpKet(bra=param_circ1, ket=param_circ2) as braket:\n", " inner_prod_conj = braket.contract(symbol_map=symbol_map)\n", "assert np.isclose(np.conj(inner_prod), inner_prod_conj)\n", "print(f\" = {inner_prod}\")\n", "print(f\" = {inner_prod_conj}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["And we are not constrained to Hermitian operators:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_XZIXX = QubitPauliString(param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X])\n", "string_IZZYX = QubitPauliString(param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X])\n", "string_ZIZXY = QubitPauliString(param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y])\n", "operator = QubitPauliOperator({string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j+0.3})\n", "with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " value = braket.contract(operator, symbol_map=symbol_map)\n", "print(value)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Backends
\n", "We provide a pytket `Backend` to obtain shots using `GeneralState`."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's consider a more challenging circuit"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["def random_circuit(n_qubits: int, n_layers: int) -> Circuit:\n", " \"\"\"Random quantum volume circuit.\"\"\"\n", " c = Circuit(n_qubits, n_qubits)\n", " for _ in range(n_layers):\n", " qubits = np.random.permutation([i for i in range(n_qubits)])\n", " qubit_pairs = [[qubits[i], qubits[i + 1]] for i in range(0, n_qubits - 1, 2)]\n", " for pair in qubit_pairs:\n", " # Generate random 4x4 unitary matrix.\n", " SU4 = unitary_group.rvs(4) # random unitary in SU4\n", " SU4 = SU4 / (np.linalg.det(SU4) ** 0.25)\n", " SU4 = np.matrix(SU4)\n", " c.add_unitary2qbox(Unitary2qBox(SU4), *pair)\n", " DecomposeBoxes().apply(c)\n", " return c"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's measure only three of the qubits.
\n", "**Note**: The complexity of this simulation increases exponentially with the number of qubits measured. Other factors leading to intractability are circuit depth and qubit connectivity."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 1000\n", "quantum_vol_circ = random_circuit(n_qubits=40, n_layers=5)\n", "quantum_vol_circ.Measure(Qubit(0), Bit(0))\n", "quantum_vol_circ.Measure(Qubit(1), Bit(1))\n", "quantum_vol_circ.Measure(Qubit(2), Bit(2))"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The `CuTensorNetShotsBackend` is used in the same way as any other pytket `Backend`."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["backend = CuTensorNetShotsBackend()\n", "compiled_circ = backend.get_compiled_circuit(quantum_vol_circ)\n", "results = backend.run_circuit(compiled_circ, n_shots=n_shots)\n", "print(results.get_counts())"]}], "metadata": {"kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"}, "language_info": {"codemirror_mode": {"name": "ipython", "version": 3}, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.4"}}, "nbformat": 4, "nbformat_minor": 2} \ No newline at end of file +{"cells": [{"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["import numpy as np\n", "from sympy import Symbol\n", "from scipy.stats import unitary_group # type: ignore\n", "from pytket.circuit import Circuit, OpType, Unitary2qBox, Qubit, Bit\n", "from pytket.passes import DecomposeBoxes\n", "from pytket.utils import QubitPauliOperator\n", "from pytket._tket.pauli import Pauli, QubitPauliString\n", "from pytket.circuit.display import render_circuit_jupyter"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["from pytket.extensions.cutensornet.general_state import (\n", " GeneralState,\n", " GeneralBraOpKet,\n", ")\n", "from pytket.extensions.cutensornet.backends import CuTensorNetShotsBackend"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Introduction
\n", "This notebook is a guide on how to use the features provided in the `general_state` submodule of pytket-cutensornet. This submodule is a thin wrapper of CuTensorNet's `NetworkState`, allowing users to convert pytket circuits into tensor networks and use CuTensorNet's contraction path optimisation algorithm.
\n", "All simulations realised with this submodule are *exact*. Once the pytket circuit has been converted to a tensor network, the computation has two steps:
\n", " 1. *Contraction path optimisation*. Attempts to find an order of contracting pairs of tensors in which the the total number of FLOPs is minimised. No operation on the tensor network occurs at this point. Runs on CPU.
\n", " 2. *Tensor network contraction*. Uses the ordering of contractions found in the previous step evaluate the tensor network. Runs on GPU.
\n", "
\n", "**Reference**: The original contraction path optimisation algorithm that NVIDIA implemented on CuTensorNet: https://arxiv.org/abs/2002.01935"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralState`
\n", "The class `GeneralState` is used to convert a circuit into a tensor network and query information from the final state. Let's walk through a simple example."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["my_circ = Circuit(5)\n", "my_circ.CX(3, 4)\n", "my_circ.H(2)\n", "my_circ.CZ(0, 1)\n", "my_circ.ZZPhase(0.1, 4, 3)\n", "my_circ.TK2(0.3, 0.5, 0.7, 2, 1)\n", "my_circ.Ry(0.2, 0)\n", "my_circ.measure_all()"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["render_circuit_jupyter(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The first step is to convert our pytket circuit into a tensor network. This is straightforward:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state = GeneralState(my_circ)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The variable `tn_state` now holds a tensor network representation of `my_circ`.
\n", "**Note**: Circuits must not have mid-circuit measurements or classical logic. The measurements at the end of the circuit are stripped and only considered when calling `tn_state.sample(n_shots)`.
\n", "We can now query information from the state. For instance, let's calculate the probability of in the qubits 0 and 3 agreeing in their outcome."]}, {"cell_type": "markdown", "metadata": {}, "source": ["First, let's generate `|x>` computational basis states where `q[0]` and `q[3]` agree on their values. We can do this with some bitwise operators and list comprehension.
\n", "**Note**: Remember that pytket uses \"increasing lexicographic order\" (ILO) for qubits, so `q[0]` is the most significant bit."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["selected_states = [\n", " x\n", " for x in range(2**my_circ.n_qubits)\n", " if ( # Iterate over all possible states\n", " x & int(\"10000\", 2) == 0\n", " and x & int(\"00010\", 2) == 0 # both qubits are 0 or...\n", " or x & int(\"10000\", 2) != 0\n", " and x & int(\"00010\", 2) != 0 # both qubits are 1\n", " )\n", "]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can now query the amplitude of all of these states and calculate the probability by summing their squared absolute values."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["amplitudes = []\n", "for x in selected_states:\n", " amplitudes.append(tn_state.get_amplitude(x))\n", "probability = sum(abs(a) ** 2 for a in amplitudes)\n", "print(f\"Probability: {probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Of course, calculating probabilities by considering the amplitudes of all relevant states is not efficient in general, since we may need to calculate a number of amplitudes that scales exponentially with the number of qubits. An alternative is to use expectation values. In particular, all of the states in `selected_states` are +1 eigenvectors of the `ZIIZI` observable and, hence, we can calculate the probability `p` by solving the equation ` = (+1)p + (-1)(1-p)` using the fact that `ZIIZI` only has +1 and -1 eigenvalues."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_ZIIZI = QubitPauliString(\n", " my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I]\n", ")\n", "observable = QubitPauliOperator({string_ZIIZI: 1.0})\n", "expectation_val = tn_state.expectation_value(observable).real\n", "exp_probability = (expectation_val + 1) / 2\n", "assert np.isclose(probability, exp_probability, atol=0.0001)\n", "print(f\"Probability: {exp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Alternatively, we can estimate the probability by sampling."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 100000\n", "outcomes = tn_state.sample(n_shots)\n", "hit_count = 0\n", "for bit_tuple, count in outcomes.get_counts().items():\n", " if bit_tuple[0] == bit_tuple[3]:\n", " hit_count += count\n", "samp_probability = hit_count / n_shots\n", "assert np.isclose(probability, samp_probability, atol=0.01)\n", "print(f\"Probability: {samp_probability}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["When we finish doing computations with the `tn_state` we must destroy it to free GPU memory."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["tn_state.destroy()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["To avoid forgetting this final step, we recommend users call `GeneralState` (and `GeneralBraOpKet`) as context managers:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralState(my_circ) as my_state:\n", " expectation_val = my_state.expectation_value(observable)\n", "print(expectation_val)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Using this syntax, `my_state` is automatically destroyed when the code exists the `with ...` block."]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Parameterised circuits
\n", "Circuits that only differ on the parameters of their gates have the same tensor network topology and, hence, we may use the same contraction path for all of them."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["a, b, c = Symbol(\"a\"), Symbol(\"b\"), Symbol(\"c\")\n", "param_circ1 = Circuit(5)\n", "param_circ1.Ry(a, 3).Ry(0.27, 4).CX(4, 3).Ry(b, 2).Ry(0.21, 3)\n", "param_circ1.Ry(0.12, 0).Ry(a, 1)\n", "param_circ1.add_gate(OpType.CnX, [0, 1, 4]).add_gate(OpType.CnX, [4, 1, 3])\n", "param_circ1.X(0).X(1).add_gate(OpType.CnY, [0, 1, 2]).add_gate(OpType.CnY, [0, 4, 3]).X(\n", " 0\n", ").X(1)\n", "param_circ1.Ry(-b, 0).Ry(-c, 1)\n", "render_circuit_jupyter(param_circ1)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can pass a parameterised circuit to `GeneralState`. The value of the parameters is provided when calling methods of `GeneralState`. The contraction path is automatically reused on different calls to the same method."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_circs = 5\n", "with GeneralState(param_circ1) as param_state:\n", " for i in range(n_circs):\n", " symbol_map = {s: np.random.random() for s in [a, b, c]}\n", " exp_val = param_state.expectation_value(observable, symbol_map=symbol_map)\n", " print(f\"Expectation value for circuit {i}: {exp_val.real}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# `GeneralBraOpKet`
\n", "The `GeneralBraOpKet` can be used to calculate any number that can be represented as the result of some `` where `|bra>` and `|ket>` are the final states of pytket circuits, and `op` is a `QubitPauliOperator`. The circuits for `|bra>` and `|ket>` need not be the same."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["x, y, z = Symbol(\"x\"), Symbol(\"y\"), Symbol(\"z\")\n", "param_circ2 = Circuit(5)\n", "param_circ2.H(0)\n", "param_circ2.S(1)\n", "param_circ2.Rz(x * z, 2)\n", "param_circ2.Ry(y + x, 3)\n", "param_circ2.TK1(x, y, z, 4)\n", "param_circ2.TK2(z - y, z - x, (x + y) * z, 1, 3)\n", "symbol_map = {a: 2.1, b: 1.3, c: 0.7, x: 3.0, y: 1.6, z: -8.3}"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We can calculate inner products by providing no `op`:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " inner_prod = braket.contract(symbol_map=symbol_map)\n", "with GeneralBraOpKet(bra=param_circ1, ket=param_circ2) as braket:\n", " inner_prod_conj = braket.contract(symbol_map=symbol_map)\n", "assert np.isclose(np.conj(inner_prod), inner_prod_conj)\n", "print(f\" = {inner_prod}\")\n", "print(f\" = {inner_prod_conj}\")"]}, {"cell_type": "markdown", "metadata": {}, "source": ["And we are not constrained to Hermitian operators:"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["string_XZIXX = QubitPauliString(\n", " param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X]\n", ")\n", "string_IZZYX = QubitPauliString(\n", " param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X]\n", ")\n", "string_ZIZXY = QubitPauliString(\n", " param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y]\n", ")\n", "operator = QubitPauliOperator(\n", " {string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j + 0.3}\n", ")\n", "with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket:\n", " value = braket.contract(operator, symbol_map=symbol_map)\n", "print(value)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["# Backends
\n", "We provide a pytket `Backend` to obtain shots using `GeneralState`."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's consider a more challenging circuit"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["def random_circuit(n_qubits: int, n_layers: int) -> Circuit:\n", " \"\"\"Random quantum volume circuit.\"\"\"\n", " c = Circuit(n_qubits, n_qubits)\n", " for _ in range(n_layers):\n", " qubits = np.random.permutation([i for i in range(n_qubits)])\n", " qubit_pairs = [[qubits[i], qubits[i + 1]] for i in range(0, n_qubits - 1, 2)]\n", " for pair in qubit_pairs:\n", " # Generate random 4x4 unitary matrix.\n", " SU4 = unitary_group.rvs(4) # random unitary in SU4\n", " SU4 = SU4 / (np.linalg.det(SU4) ** 0.25)\n", " SU4 = np.matrix(SU4)\n", " c.add_unitary2qbox(Unitary2qBox(SU4), *pair)\n", " DecomposeBoxes().apply(c)\n", " return c"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Let's measure only three of the qubits.
\n", "**Note**: The complexity of this simulation increases exponentially with the number of qubits measured. Other factors leading to intractability are circuit depth and qubit connectivity."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["n_shots = 1000\n", "quantum_vol_circ = random_circuit(n_qubits=40, n_layers=5)\n", "quantum_vol_circ.Measure(Qubit(0), Bit(0))\n", "quantum_vol_circ.Measure(Qubit(1), Bit(1))\n", "quantum_vol_circ.Measure(Qubit(2), Bit(2))"]}, {"cell_type": "markdown", "metadata": {}, "source": ["The `CuTensorNetShotsBackend` is used in the same way as any other pytket `Backend`."]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": ["backend = CuTensorNetShotsBackend()\n", "compiled_circ = backend.get_compiled_circuit(quantum_vol_circ)\n", "results = backend.run_circuit(compiled_circ, n_shots=n_shots)\n", "print(results.get_counts())"]}], "metadata": {"kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"}, "language_info": {"codemirror_mode": {"name": "ipython", "version": 3}, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.4"}}, "nbformat": 4, "nbformat_minor": 2} \ No newline at end of file diff --git a/examples/python/general_state_tutorial.py b/examples/python/general_state_tutorial.py index 221f6347..bef9b47d 100644 --- a/examples/python/general_state_tutorial.py +++ b/examples/python/general_state_tutorial.py @@ -45,23 +45,27 @@ # First, let's generate `|x>` computational basis states where `q[0]` and `q[3]` agree on their values. We can do this with some bitwise operators and list comprehension. # **Note**: Remember that pytket uses "increasing lexicographic order" (ILO) for qubits, so `q[0]` is the most significant bit. selected_states = [ - x for x in range(2**my_circ.n_qubits) if ( # Iterate over all possible states - x & int("10000", 2) == 0 and x & int("00010", 2) == 0 # both qubits are 0 or... - or x & int("10000", 2) != 0 and x & int("00010", 2) != 0 # both qubits are 1 + x + for x in range(2**my_circ.n_qubits) + if ( # Iterate over all possible states + x & int("10000", 2) == 0 + and x & int("00010", 2) == 0 # both qubits are 0 or... + or x & int("10000", 2) != 0 + and x & int("00010", 2) != 0 # both qubits are 1 ) ] # We can now query the amplitude of all of these states and calculate the probability by summing their squared absolute values. amplitudes = [] for x in selected_states: - amplitudes.append( - tn_state.get_amplitude(x) - ) -probability = sum(abs(a)**2 for a in amplitudes) + amplitudes.append(tn_state.get_amplitude(x)) +probability = sum(abs(a) ** 2 for a in amplitudes) print(f"Probability: {probability}") # Of course, calculating probabilities by considering the amplitudes of all relevant states is not efficient in general, since we may need to calculate a number of amplitudes that scales exponentially with the number of qubits. An alternative is to use expectation values. In particular, all of the states in `selected_states` are +1 eigenvectors of the `ZIIZI` observable and, hence, we can calculate the probability `p` by solving the equation ` = (+1)p + (-1)(1-p)` using the fact that `ZIIZI` only has +1 and -1 eigenvalues. -string_ZIIZI = QubitPauliString(my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I]) +string_ZIIZI = QubitPauliString( + my_circ.qubits, [Pauli.Z, Pauli.I, Pauli.I, Pauli.Z, Pauli.I] +) observable = QubitPauliOperator({string_ZIIZI: 1.0}) expectation_val = tn_state.expectation_value(observable).real exp_probability = (expectation_val + 1) / 2 @@ -133,10 +137,18 @@ print(f" = {inner_prod_conj}") # And we are not constrained to Hermitian operators: -string_XZIXX = QubitPauliString(param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X]) -string_IZZYX = QubitPauliString(param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X]) -string_ZIZXY = QubitPauliString(param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y]) -operator = QubitPauliOperator({string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j+0.3}) +string_XZIXX = QubitPauliString( + param_circ2.qubits, [Pauli.X, Pauli.Z, Pauli.I, Pauli.X, Pauli.X] +) +string_IZZYX = QubitPauliString( + param_circ2.qubits, [Pauli.I, Pauli.Z, Pauli.Z, Pauli.Y, Pauli.X] +) +string_ZIZXY = QubitPauliString( + param_circ2.qubits, [Pauli.Z, Pauli.I, Pauli.Z, Pauli.X, Pauli.Y] +) +operator = QubitPauliOperator( + {string_XZIXX: -1.38j, string_IZZYX: 2.36, string_ZIZXY: 0.42j + 0.3} +) with GeneralBraOpKet(bra=param_circ2, ket=param_circ1) as braket: value = braket.contract(operator, symbol_map=symbol_map) print(value) @@ -144,6 +156,7 @@ # # Backends # We provide a pytket `Backend` to obtain shots using `GeneralState`. + # Let's consider a more challenging circuit def random_circuit(n_qubits: int, n_layers: int) -> Circuit: """Random quantum volume circuit.""" @@ -163,6 +176,7 @@ def random_circuit(n_qubits: int, n_layers: int) -> Circuit: DecomposeBoxes().apply(c) return c + # Let's measure only three of the qubits. # **Note**: The complexity of this simulation increases exponentially with the number of qubits measured. Other factors leading to intractability are circuit depth and qubit connectivity. n_shots = 1000 @@ -175,4 +189,4 @@ def random_circuit(n_qubits: int, n_layers: int) -> Circuit: backend = CuTensorNetShotsBackend() compiled_circ = backend.get_compiled_circuit(quantum_vol_circ) results = backend.run_circuit(compiled_circ, n_shots=n_shots) -print(results.get_counts()) \ No newline at end of file +print(results.get_counts()) From e2d8d6c4c46b9e6eb9298c0e4d3d74ea21ccd5e5 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Thu, 24 Oct 2024 12:51:07 +0000 Subject: [PATCH 23/27] Added notebooks to CI --- examples/check-examples | 6 +++--- examples/ci-tested-notebooks.txt | 3 ++- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/examples/check-examples b/examples/check-examples index a9e947cc..9a381a64 100755 --- a/examples/check-examples +++ b/examples/check-examples @@ -9,7 +9,7 @@ do p2j -o -t ${name}-gen.ipynb python/${name}.py cmp ${name}.ipynb ${name}-gen.ipynb rm ${name}-gen.ipynb - # TODO, add this when GPU is added to CI - # Run script: - # python python/${name}.py + + # Run script: + python python/${name}.py done diff --git a/examples/ci-tested-notebooks.txt b/examples/ci-tested-notebooks.txt index 1cfc10e5..ee440a89 100644 --- a/examples/ci-tested-notebooks.txt +++ b/examples/ci-tested-notebooks.txt @@ -1,2 +1,3 @@ mps_tutorial -ttn_tutorial \ No newline at end of file +ttn_tutorial +general_state_tutorial \ No newline at end of file From 0bc00c410cce51f44bffc84641c328203dd83e3f Mon Sep 17 00:00:00 2001 From: Pablo Andres-Martinez <104848389+PabloAndresCQ@users.noreply.github.com> Date: Thu, 24 Oct 2024 14:21:36 +0100 Subject: [PATCH 24/27] Apply suggestions from code review Co-authored-by: Oliver Backhouse --- .../general_state/tensor_network_state.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index cf234aee..7ce087d8 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -102,11 +102,11 @@ def __init__( if is_fixed: try: gate_unitary = com.op.get_unitary() - except: + except Exception as e: raise ValueError( "All commands in the circuit must be unitary gates. The circuit" f" contains {com}; no unitary matrix could be retrived from it." - ) + ) from e else: # Dummy unitary to be updated later with user specified paramaters gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") @@ -412,11 +412,11 @@ def __init__( if is_fixed: try: gate_unitary = com.op.get_unitary() - except: + except Exception as e: raise ValueError( "All commands in the circuit must be unitary gates. The circuit" f" contains {com}; no unitary matrix could be retrived from it." - ) + ) from e else: # Dummy unitary to be updated later with user specified paramaters gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") @@ -464,11 +464,11 @@ def __init__( if is_fixed: try: gate_unitary = com.op.get_unitary() - except: + except Exception as e: raise ValueError( "All commands in the circuit must be unitary gates. The circuit" f" contains {com}; no unitary matrix could be retrived from it." - ) + ) from e else: # Dummy unitary to be updated later with user specified paramaters gate_unitary = np.identity(2**com.op.n_qubits, dtype="complex128") @@ -598,7 +598,7 @@ def _formatted_tensor(matrix: NDArray, n_qubits: int) -> cp.ndarray: return cupy_matrix.reshape([2] * (2 * n_qubits)) -def _remove_meas_and_implicit_swaps(circ: Circuit) -> Tuple[Circuit, Dict[Qubit, Bit]]: +def _remove_meas_and_implicit_swaps(circ: Circuit) -> tuple[Circuit, dict[Qubit, Bit]]: """Convert a pytket Circuit to an equivalent circuit with no measurements or implicit swaps. The measurements are returned as a map between qubits and bits. From 61e2fcf7b39f0876b50835d5b6cc940029be3d09 Mon Sep 17 00:00:00 2001 From: PabloAndresCQ Date: Thu, 24 Oct 2024 13:29:07 +0000 Subject: [PATCH 25/27] Linting --- .github/workflows/lint.yml | 2 +- .../cutensornet/general_state/tensor_network_state.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/lint.yml b/.github/workflows/lint.yml index 807d038c..84df82a8 100644 --- a/.github/workflows/lint.yml +++ b/.github/workflows/lint.yml @@ -30,4 +30,4 @@ jobs: black --check . - name: Run pylint run: | - pylint --recursive=y --ignore=ttn_tutorial.py,mps_tutorial.py */ + pylint --recursive=y pytket/extensions/cutensornet/ diff --git a/pytket/extensions/cutensornet/general_state/tensor_network_state.py b/pytket/extensions/cutensornet/general_state/tensor_network_state.py index 7ce087d8..e89917d5 100644 --- a/pytket/extensions/cutensornet/general_state/tensor_network_state.py +++ b/pytket/extensions/cutensornet/general_state/tensor_network_state.py @@ -14,7 +14,7 @@ from __future__ import annotations import logging -from typing import Union, Optional, Tuple, Dict, Any +from typing import Union, Optional, Any import warnings try: From c070d80c34632cc228168dfc0c974c9f52a43792 Mon Sep 17 00:00:00 2001 From: Jake Arkinstall <65358059+jake-arkinstall@users.noreply.github.com> Date: Thu, 24 Oct 2024 15:36:33 +0100 Subject: [PATCH 26/27] Migrate testing of example notebooks to nix and run them in CI on the dedicated GPU box --- .github/workflows/build-with-nix.yml | 3 +++ examples/check-examples | 9 +++++++-- flake.nix | 2 +- nix-support/pytket-cutensornet.nix | 28 ++++++++++++++++++++++++++++ 4 files changed, 39 insertions(+), 3 deletions(-) diff --git a/.github/workflows/build-with-nix.yml b/.github/workflows/build-with-nix.yml index ac6ed09c..cc801114 100644 --- a/.github/workflows/build-with-nix.yml +++ b/.github/workflows/build-with-nix.yml @@ -39,3 +39,6 @@ jobs: - name: Test pytket-cutensornet # impure is necessary due to nixgl usage (system-dependent cuda) run: nix run .#tests --impure --accept-flake-config + - name: Test example notebooks + # impure is necessary due to nixgl usage (system-dependent cuda) + run: nix run .#example-tests --impure --accept-flake-config diff --git a/examples/check-examples b/examples/check-examples index 9a381a64..47a470c1 100755 --- a/examples/check-examples +++ b/examples/check-examples @@ -10,6 +10,11 @@ do cmp ${name}.ipynb ${name}-gen.ipynb rm ${name}-gen.ipynb - # Run script: - python python/${name}.py + # run tests are performed in nix, allowing + # us to manage the testing environment in a + # reproducible way. + # + # See /nix-support/pytket-cutensornet.nix, + # in the derivation called + # run-pytket-cutensornet-examples. done diff --git a/flake.nix b/flake.nix index 04b67de9..72eb5a54 100644 --- a/flake.nix +++ b/flake.nix @@ -43,9 +43,9 @@ in { packages = { default = pkgs.pytket-cutensornet; - cupy = pkgs.cupy'; pytket-cutensornet = pkgs.pytket-cutensornet; tests = pkgs.run-pytket-cutensornet-tests; + example-tests = pkgs.run-pytket-cutensornet-examples; }; devShells = { default = pkgs.mkShell { buildInputs = [ pkgs.pytket-cutensornet ]; }; diff --git a/nix-support/pytket-cutensornet.nix b/nix-support/pytket-cutensornet.nix index 4895c30d..1228a537 100644 --- a/nix-support/pytket-cutensornet.nix +++ b/nix-support/pytket-cutensornet.nix @@ -65,4 +65,32 @@ EOF export LD_LIBRARY_PATH; ${test-env}/bin/pytest -s ${../tests}; ''; + run-pytket-cutensornet-examples = let + example-env = super.python3.withPackages(ps: with ps; [ + self.pytket-cutensornet + matplotlib + numpy + networkx + ipython + nbmake + pytest + ]); + nixgl-bin = self.lib.getExe self.nixgl.auto.nixGLNvidia; + in super.writeShellScriptBin "run-pytket-cutensornet-examples" '' + HOME=$(mktemp -d); + export HOME; + NIXGL_PATH="$(${nixgl-bin} printenv LD_LIBRARY_PATH)"; + WSL_PATH="/usr/lib/wsl/lib"; + LD_LIBRARY_PATH="$NIXGL_PATH:$WSL_PATH:$LD_LIBRARY_PATH"; + + example_dir=${../examples}; + set -e; + for name in `cat ''${example_dir}/ci-tested-notebooks.txt`; + do + ${example-env}/bin/pytest \ + --nbmake \ + -p no:cacheprovider \ + $example_dir/$name.ipynb; + done; + ''; } From ff5a88b66fe752298aefac29ba8624e1442425c7 Mon Sep 17 00:00:00 2001 From: Jake Arkinstall <65358059+jake-arkinstall@users.noreply.github.com> Date: Thu, 24 Oct 2024 15:46:18 +0100 Subject: [PATCH 27/27] export LD_LIBRARY_PATH in examples test runner to expose cuda drivers --- nix-support/pytket-cutensornet.nix | 1 + 1 file changed, 1 insertion(+) diff --git a/nix-support/pytket-cutensornet.nix b/nix-support/pytket-cutensornet.nix index 1228a537..c6398cd4 100644 --- a/nix-support/pytket-cutensornet.nix +++ b/nix-support/pytket-cutensornet.nix @@ -82,6 +82,7 @@ EOF NIXGL_PATH="$(${nixgl-bin} printenv LD_LIBRARY_PATH)"; WSL_PATH="/usr/lib/wsl/lib"; LD_LIBRARY_PATH="$NIXGL_PATH:$WSL_PATH:$LD_LIBRARY_PATH"; + export LD_LIBRARY_PATH; example_dir=${../examples}; set -e;