diff --git a/qiskit_dynamics/pulse/backend_parser/__init__.py b/qiskit_dynamics/pulse/backend_parser/__init__.py new file mode 100644 index 000000000..d6dde3096 --- /dev/null +++ b/qiskit_dynamics/pulse/backend_parser/__init__.py @@ -0,0 +1,19 @@ +# -*- coding: utf-8 -*- + +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +""" +String model parser functionality. + +This module is meant for internal use and may be changed at any point. +""" diff --git a/qiskit_dynamics/pulse/backend_parser/operator_from_string.py b/qiskit_dynamics/pulse/backend_parser/operator_from_string.py new file mode 100644 index 000000000..fff05af08 --- /dev/null +++ b/qiskit_dynamics/pulse/backend_parser/operator_from_string.py @@ -0,0 +1,142 @@ +# -*- coding: utf-8 -*- + +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +# pylint: disable=invalid-name + +"""Generate operators from string. + +This file is meant for internal use and may be changed at any point. +""" + +from typing import Dict +import numpy as np + +from qiskit import QiskitError +import qiskit.quantum_info as qi + + +def operator_from_string( + op_label: str, subsystem_label: int, subsystem_dims: Dict[int, int] +) -> np.ndarray: + r"""Generates a dense operator acting on a single subsystem, tensoring + identities for remaining subsystems. + + The single system operator is specified via a string in ``op_label``, + the list of subsystems and their corresponding dimensions are specified in the + dictionary ``subsystem_dims``, with system label being the keys specified as ``int``s, + and system dimensions the values also specified as ``int``s, and ``subsystem_label`` + indicates which subsystem the operator specified by ``op_label`` acts on. + + Accepted ``op_labels`` are: + - `'X'`: If the target subsystem is two dimensional, the + Pauli :math:`X` operator, and if greater than two dimensional, returns + :math:`a + a^\dagger`, where :math:`a` and :math:`a^\dagger` are the + annihiliation and creation operators, respectively. + - `'Y'`: If the target subsystem is two dimensional, the + Pauli :math:`Y` operator, and if greater than two dimensional, returns + :math:`-i(a - a^\dagger)`, where :math:`a` and :math:`a^\dagger` are the + annihiliation and creation operators, respectively. + - `'Z'`: If the target subsystem is two dimensional, the + Pauli :math:`Z` operator, and if greater than two dimensional, returns + :math:`I - 2 * N`, where :math:`N` is the number operator. + - `'a'`, `'A'`, or `'Sm'`: If two dimensional, the sigma minus operator, and if greater, + generalizes to the operator. + - `'C'`, or `'Sp'`: If two dimensional, sigma plus operator, and if greater, + generalizes to the creation operator. + - `'N'`, or `'O'`: The number operator. + - `'I'`: The identity operator. + + Note that the ordering of tensor factors is reversed. + + Args: + op_label: The string labelling the single system operator. + subsystem_label: Index of the subsystem to apply the operator. + subsystem_dims: Dictionary of subsystem labels and dimensions. + + Returns: + np.ndarray corresponding to the specified operator. + + Raises: + QiskitError: If op_label is invalid. + """ + + # construct single system operator + op_func = __operdict.get(op_label, None) + if op_func is None: + raise QiskitError(f"String {op_label} does not correspond to a known operator.") + + dim = subsystem_dims[subsystem_label] + out = qi.Operator(op_func(dim), input_dims=[dim], output_dims=[dim]) + + # sort subsystem labels and dimensions according to subsystem label + sorted_subsystem_keys, sorted_subsystem_dims = zip( + *sorted(zip(subsystem_dims.keys(), subsystem_dims.values())) + ) + + # get subsystem location in ordered list + subsystem_location = sorted_subsystem_keys.index(subsystem_label) + + # construct full operator + return qi.ScalarOp(sorted_subsystem_dims).compose(out, [subsystem_location]).data + + +# functions for generating individual operators +def a(dim: int) -> np.ndarray: + """Annihilation operator.""" + return np.diag(np.sqrt(np.arange(1, dim, dtype=complex)), 1) + + +def adag(dim: int) -> np.ndarray: + """Creation operator.""" + return a(dim).conj().transpose() + + +def N(dim: int) -> np.ndarray: + """Number operator.""" + return np.diag(np.arange(dim, dtype=complex)) + + +def X(dim: int) -> np.ndarray: + """Generalized X operator, written in terms of raising and lowering operators.""" + return a(dim) + adag(dim) + + +def Y(dim: int) -> np.ndarray: + """Generalized Y operator, written in terms of raising and lowering operators.""" + return -1j * (a(dim) - adag(dim)) + + +def Z(dim: int) -> np.ndarray: + """Generalized Z operator, written as id - 2 * N.""" + return ident(dim) - 2 * N(dim) + + +def ident(dim: int) -> np.ndarray: + """Identity operator.""" + return np.eye(dim, dtype=complex) + + +# operator names +__operdict = { + "X": X, + "Y": Y, + "Z": Z, + "a": a, + "A": a, + "Sm": a, + "Sp": adag, + "C": adag, + "N": N, + "O": N, + "I": ident, +} diff --git a/qiskit_dynamics/pulse/backend_parser/regex_parser.py b/qiskit_dynamics/pulse/backend_parser/regex_parser.py new file mode 100644 index 000000000..5b0eef254 --- /dev/null +++ b/qiskit_dynamics/pulse/backend_parser/regex_parser.py @@ -0,0 +1,369 @@ +# -*- coding: utf-8 -*- + +# This code is part of Qiskit. +# +# (C) Copyright IBM 2018, 2019, 2020, 2021, 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +# pylint: disable=invalid-name + +"""Legacy code for parsing operator strings. + +This file is meant for internal use and may be changed at any point. +""" + +import re +import copy +from typing import List, Dict, Tuple +from collections import namedtuple, OrderedDict + +import numpy as np + +from .operator_from_string import operator_from_string + + +def _regex_parser( + operator_str: List[str], subsystem_dims: Dict[int, int], subsystem_list: List[int] +) -> List[Tuple[np.array, str]]: + """Function wrapper for regex parsing object. + + Args: + operator_str: List of strings in accepted format as described in + string_model_parser.parse_hamiltonian_dict. + subsystem_dims: Dictionary mapping subsystem labels to dimensions. + subsystem_list: List of subsystems on which the operators are to be constructed. + Returns: + List of tuples containing pairs operators and their string coefficients. + """ + + return _HamiltonianParser(h_str=operator_str, subsystem_dims=subsystem_dims).parse( + subsystem_list + ) + + +class _HamiltonianParser: + """Legacy object for parsing string specifications of Hamiltonians.""" + + Token = namedtuple("Token", ("type", "name")) + + str_elements = OrderedDict( + QubOpr=re.compile(r"(?PO|Sp|Sm|X|Y|Z|I)(?P[0-9]+)"), + PrjOpr=re.compile(r"P(?P[0-9]+),(?P[0-9]+),(?P[0-9]+)"), + CavOpr=re.compile(r"(?PA|C|N)(?P[0-9]+)"), + Func=re.compile(r"(?P[a-z]+)\("), + Ext=re.compile(r"\.(?Pdag)"), + Var=re.compile(r"[a-z]+[0-9]*"), + Num=re.compile(r"[0-9.]+"), + MathOrd0=re.compile(r"[*/]"), + MathOrd1=re.compile(r"[+-]"), + BrkL=re.compile(r"\("), + BrkR=re.compile(r"\)"), + ) + + def __init__(self, h_str, subsystem_dims): + """Create new quantum operator generator + + Parameters: + h_str (list): list of Hamiltonian string + subsystem_dims (dict): dimension of subsystems + """ + self.h_str = h_str + self.subsystem_dims = {int(label): int(dim) for label, dim in subsystem_dims.items()} + self.str2qopr = {} + + def parse(self, qubit_list=None): + """Parse and generate Hamiltonian terms.""" + td_hams = [] + tc_hams = [] + + # expand sum + self._expand_sum() + + # convert to reverse Polish notation + for ham in self.h_str: + if len(re.findall(r"\|\|", ham)) > 1: + raise Exception(f"Multiple time-dependent terms in {ham}") + p_td = re.search(r"(?P[\S]+)\|\|(?P[\S]+)", ham) + + # find time-dependent term + if p_td: + coef, token = self._tokenizer(p_td.group("opr"), qubit_list) + if token is None: + continue + # combine coefficient to time-dependent term + if coef: + td = "*".join([coef, p_td.group("ch")]) + else: + td = p_td.group("ch") + token = self._shunting_yard(token) + _td = self._token2qobj(token), td + + td_hams.append(_td) + else: + coef, token = self._tokenizer(ham, qubit_list) + if token is None: + continue + token = self._shunting_yard(token) + + if (coef == "") or (coef is None): + coef = "1." + + _tc = self._token2qobj(token), coef + + tc_hams.append(_tc) + + return tc_hams + td_hams + + def _expand_sum(self): + """Takes a string-based Hamiltonian list and expands the _SUM action items out.""" + sum_str = re.compile(r"_SUM\[(?P[a-z]),(?P[a-z\d{}+-]+),(?P[a-z\d{}+-]+),") + brk_str = re.compile(r"]") + + ham_list = copy.copy(self.h_str) + ham_out = [] + + while any(ham_list): + ham = ham_list.pop(0) + p_sums = list(sum_str.finditer(ham)) + p_brks = list(brk_str.finditer(ham)) + if len(p_sums) != len(p_brks): + raise Exception(f"Missing correct number of brackets in {ham}") + + # find correct sum-bracket correspondence + if any(p_sums) == 0: + ham_out.append(ham) + else: + itr = p_sums[0].group("itr") + _l = int(p_sums[0].group("l")) + _u = int(p_sums[0].group("u")) + for ii in range(len(p_sums) - 1): + if p_sums[ii + 1].end() > p_brks[ii].start(): + break + else: + ii = len(p_sums) - 1 + + # substitute iterator value + _temp = [] + for kk in range(_l, _u + 1): + trg_s = ham[p_sums[0].end() : p_brks[ii].start()] + # generate replacement pattern + pattern = {} + for p in re.finditer(r"\{(?P[a-z0-9*/+-]+)\}", trg_s): + if p.group() not in pattern: + sub = parse_binop(p.group("op_str"), operands={itr: str(kk)}) + if sub.isdecimal(): + pattern[p.group()] = sub + else: + pattern[p.group()] = f"{{{sub}}}" + for key, val in pattern.items(): + trg_s = trg_s.replace(key, val) + _temp.append( + "".join([ham[: p_sums[0].start()], trg_s, ham[p_brks[ii].end() :]]) + ) + ham_list.extend(_temp) + + self.h_str = ham_out + + return ham_out + + def _tokenizer(self, op_str, qubit_list=None): + """Convert string to token and coefficient + Check if the index is in qubit_list + """ + + # generate token + _op_str = copy.copy(op_str) + token_list = [] + prev = "none" + while any(_op_str): + for key, parser in _HamiltonianParser.str_elements.items(): + p = parser.match(_op_str) + if p: + # find quantum operators + if key in ["QubOpr", "CavOpr"]: + _key = key + _name = p.group() + if p.group() not in self.str2qopr: + idx = int(p.group("idx")) + if qubit_list is not None and idx not in qubit_list: + return 0, None + name = p.group("opr") + opr = operator_from_string(name, idx, self.subsystem_dims) + self.str2qopr[p.group()] = opr + elif key == "PrjOpr": + _key = key + _name = p.group() + if p.group() not in self.str2qopr: + idx = int(p.group("idx")) + name = "P" + opr = operator_from_string(name, idx, self.subsystem_dims) + self.str2qopr[p.group()] = opr + elif key in ["Func", "Ext"]: + _name = p.group("name") + _key = key + elif key == "MathOrd1": + _name = p.group() + if prev not in ["QubOpr", "PrjOpr", "CavOpr", "Var", "Num"]: + _key = "MathUnitary" + else: + _key = key + else: + _name = p.group() + _key = key + token_list.append(_HamiltonianParser.Token(_key, _name)) + _op_str = _op_str[p.end() :] + prev = _key + break + else: + raise Exception(f"Invalid input string {op_str} is found") + + # split coefficient + coef = "" + if any(k.type == "Var" for k in token_list): + for ii, _ in enumerate(token_list): + if token_list[ii].name == "*": + if all(k.type != "Var" for k in token_list[ii + 1 :]): + coef = "".join([k.name for k in token_list[:ii]]) + token_list = token_list[ii + 1 :] + break + else: + raise Exception(f"Invalid order of operators and coefficients in {op_str}") + + return coef, token_list + + def _shunting_yard(self, token_list): + """Reformat token to reverse Polish notation""" + stack = [] + queue = [] + while any(token_list): + token = token_list.pop(0) + if token.type in ["QubOpr", "PrjOpr", "CavOpr", "Num"]: + queue.append(token) + elif token.type in ["Func", "Ext"]: + stack.append(token) + elif token.type in ["MathUnitary", "MathOrd0", "MathOrd1"]: + while stack and math_priority(token, stack[-1]): + queue.append(stack.pop(-1)) + stack.append(token) + elif token.type in ["BrkL"]: + stack.append(token) + elif token.type in ["BrkR"]: + while stack[-1].type not in ["BrkL", "Func"]: + queue.append(stack.pop(-1)) + if not any(stack): + raise Exception("Missing correct number of brackets") + pop = stack.pop(-1) + if pop.type == "Func": + queue.append(pop) + else: + raise Exception(f"Invalid token {token.name} is found") + + while any(stack): + queue.append(stack.pop(-1)) + + return queue + + def _token2qobj(self, tokens): + """Generate Hamiltonian term from tokens.""" + stack = [] + for token in tokens: + if token.type in ["QubOpr", "PrjOpr", "CavOpr"]: + stack.append(self.str2qopr[token.name]) + elif token.type == "Num": + stack.append(float(token.name)) + elif token.type in ["MathUnitary"]: + if token.name == "-": + stack.append(-stack.pop(-1)) + elif token.type in ["MathOrd0", "MathOrd1"]: + op2 = stack.pop(-1) + op1 = stack.pop(-1) + if token.name == "+": + stack.append(op1 + op2) + elif token.name == "-": + stack.append(op1 - op2) + elif token.name == "*": + stack.append(op1 @ op2) + elif token.name == "/": + stack.append(op1 / op2) + elif token.type in ["Func", "Ext"]: + if token.name == "dag": + stack.append(np.conjugate(np.transpose(stack.pop(-1)))) + else: + raise Exception(f"Invalid token {token.name} of type Func, Ext.") + else: + raise Exception(f"Invalid token {token.name} is found.") + + if len(stack) > 1: + raise Exception("Invalid mathematical operation in string.") + + return stack[0] + + +def math_priority(o1, o2): + """Check priority of given math operation""" + rank = {"MathUnitary": 2, "MathOrd0": 1, "MathOrd1": 0} + diff_ops = rank.get(o1.type, -1) - rank.get(o2.type, -1) + + if diff_ops > 0: + return False + else: + return True + + +# pylint: disable=dangerous-default-value +def parse_binop(op_str, operands={}, cast_str=True): + """Calculate binary operation in string format""" + oprs = OrderedDict( + sum=r"(?P[a-zA-Z0-9]+)\+(?P[a-zA-Z0-9]+)", + sub=r"(?P[a-zA-Z0-9]+)\-(?P[a-zA-Z0-9]+)", + mul=r"(?P[a-zA-Z0-9]+)\*(?P[a-zA-Z0-9]+)", + div=r"(?P[a-zA-Z0-9]+)\/(?P[a-zA-Z0-9]+)", + non=r"(?P[a-zA-Z0-9]+)", + ) + + for key, regr in oprs.items(): + p = re.match(regr, op_str) + if p: + val0 = operands.get(p.group("v0"), p.group("v0")) + if key == "non": + # substitution + retv = val0 + else: + val1 = operands.get(p.group("v1"), p.group("v1")) + # binary operation + if key == "sum": + if val0.isdecimal() and val1.isdecimal(): + retv = int(val0) + int(val1) + else: + retv = "+".join([str(val0), str(val1)]) + elif key == "sub": + if val0.isdecimal() and val1.isdecimal(): + retv = int(val0) - int(val1) + else: + retv = "-".join([str(val0), str(val1)]) + elif key == "mul": + if val0.isdecimal() and val1.isdecimal(): + retv = int(val0) * int(val1) + else: + retv = "*".join([str(val0), str(val1)]) + elif key == "div": + if val0.isdecimal() and val1.isdecimal(): + retv = int(val0) / int(val1) + else: + retv = "/".join([str(val0), str(val1)]) + else: + retv = 0 + break + else: + raise Exception(f"Invalid string {op_str}") + + if cast_str: + return str(retv) + else: + return retv diff --git a/qiskit_dynamics/pulse/backend_parser/string_model_parser.py b/qiskit_dynamics/pulse/backend_parser/string_model_parser.py new file mode 100644 index 000000000..17c72a2f6 --- /dev/null +++ b/qiskit_dynamics/pulse/backend_parser/string_model_parser.py @@ -0,0 +1,306 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +# pylint: disable=invalid-name + +""" +Functionality for importing qiskit.pulse model string representation. + +This file is meant for internal use and may be changed at any point. +""" + +from typing import Tuple, List, Optional +from collections import OrderedDict + +# required for calls to exec +# pylint: disable=unused-import +import numpy as np + +from qiskit import QiskitError + +from .regex_parser import _regex_parser + + +# valid channel characters +CHANNEL_CHARS = ["U", "D", "M", "A", "u", "d", "m", "a"] + + +def parse_hamiltonian_dict( + hamiltonian_dict: dict, subsystem_list: Optional[List[int]] = None +) -> Tuple[np.ndarray, np.ndarray, List[str], dict]: + r"""Convert Pulse backend Hamiltonian dictionary into concrete array format + with an ordered list of corresponding channels. + + The Pulse backend Hamiltonian dictionary, ``hamiltonian_dict``, must have the + following keys: + + * ``'h_str'``: List of Hamiltonian terms in string format (see below). + * ``'qub'``: Dictionary giving subsystem dimensions. Keys are subsystem labels, + values are their dimensions. + * ``'vars'``: Dictionary with variables appearing in the terms in the ``h_str`` list, + keys are strings giving the variables, values are the values of the variables. + + The optional argument ``subsystem_list`` specifies a subset of subsystems to keep when parsing. + If ``None``, all subsystems are kept. If ``subsystem_list`` is specified, then terms + including subsystems not in the list will be ignored. + + Entries in the list ``hamiltonian_dict['h_str']`` must be formatted as a product of + constants (either numerical constants or variables in ``hamiltonian_dict['vars'].keys()``) + with operators. Operators are indicated with a capital letters followed by an integer + indicating the subsystem the operator acts on. Accepted operator strings are: + + * ``'X'``: If the target subsystem is two dimensional, the + Pauli :math:`X` operator, and if greater than two dimensional, returns + :math:`a + a^\dagger`, where :math:`a` and :math:`a^\dagger` are the + annihiliation and creation operators, respectively. + * ``'Y'``: If the target subsystem is two dimensional, the + Pauli :math:`Y` operator, and if greater than two dimensional, returns + :math:`-i(a - a^\dagger)`, where :math:`a` and :math:`a^\dagger` are the + annihiliation and creation operators, respectively. + * ``'Z'``: If the target subsystem is two dimensional, the + Pauli :math:`Z` operator, and if greater than two dimensional, returns + :math:`I - 2 * N`, where :math:`N` is the number operator. + * ``'a'``, ``'A'``, or ``'Sm'``: If two dimensional, the sigma minus operator, and if greater, + generalizes to the operator. + * ``'C'``, or ``'Sp'``: If two dimensional, sigma plus operator, and if greater, + generalizes to the creation operator. + * ``'N'``, or ``'O'``: The number operator. + * ``'I'``: The identity operator. + + In addition to the above, a term in ``hamiltonian_dict['h_str']`` can be associated with + a channel by ending it with a string of the form ``'||Sxx'``, where ``S`` is a valid channel + label, and ``'xx'`` is an integer. Accepted channel labels are: + + * ``'D'`` or ``'d'`` for drive channels. + * ``'U'`` or ``'u'`` for control channels. + * ``'M'`` or ``'m'`` for measurement channels. + * ``'A'`` or ``'a'`` for acquire channels. + + Finally, summations of terms of the above form can be indicated in + ``hamiltonian_dict['h_str']`` via strings with syntax ``'_SUM[i, lb, ub, aa||S{i}]'``, + where: + + * ``i`` is the summation variable. + * ``lb`` and ``ub`` are the simuation endpoints (inclusive). + * ``aa`` is a valid operator string, possibly including the string ``{i}`` to indicate + operators acting on subsystem ``i``. + * ``S{i}`` is the specification of a channel indexed by ``i``. + + + For example, the following ``hamiltonian_dict`` specifies a single + transmon with 4 levels: + + .. code-block:: python + + hamiltonian_dict = { + "h_str": ["v*np.pi*O0", "alpha*np.pi*O0*O0", "r*np.pi*X0||D0"], + "qub": {"0": 4}, + "vars": {"v": 2.1, "alpha": -0.33, "r": 0.02}, + } + + The following example specifies a two transmon system, with single system terms specified + using the summation format: + + .. code-block:: python + + hamiltonian_dict = { + "h_str": [ + "_SUM[i,0,1,wq{i}/2*(I{i}-Z{i})]", + "_SUM[i,0,1,delta{i}/2*O{i}*O{i}]", + "_SUM[i,0,1,-delta{i}/2*O{i}]", + "_SUM[i,0,1,omegad{i}*X{i}||D{i}]", + "jq0q1*Sp0*Sm1", + "jq0q1*Sm0*Sp1", + "omegad1*X0||U0", + "omegad0*X1||U1" + ], + "qub": {"0": 4, "1": 4}, + "vars": { + "delta0": -2.111793476400394, + "delta1": -2.0894421352015744, + "jq0q1": 0.010495754104003914, + "omegad0": 0.9715458990879812, + "omegad1": 0.9803812537440838, + "wq0": 32.517894442809514, + "wq1": 33.0948996120196, + }, + } + + Args: + hamiltonian_dict: Pulse backend Hamiltonian dictionary. + subsystem_list: List of subsystems to include in the model. If ``None`` all are kept. + + Returns: + Tuple: Model converted into concrete arrays - the static Hamiltonian, + a list of Hamiltonians corresponding to different channels, a list of + channel labels corresponding to the list of time-dependent Hamiltonians, + and a dictionary with subsystem dimensions whose keys are the subsystem labels. + """ + + # raise errors for invalid hamiltonian_dict + hamiltonian_pre_parse_exceptions(hamiltonian_dict) + + # get variables + variables = OrderedDict() + if "vars" in hamiltonian_dict: + variables = OrderedDict(hamiltonian_dict["vars"]) + + # Get qubit subspace dimensions + if subsystem_list is None: + subsystem_list = [int(qubit) for qubit in hamiltonian_dict["qub"]] + else: + # if user supplied, make a copy and sort it + subsystem_list = subsystem_list.copy() + subsystem_list.sort() + + # force keys in hamiltonian['qub'] to be ints + qub_dict = {int(key): val for key, val in hamiltonian_dict["qub"].items()} + + subsystem_dims = {int(qubit): qub_dict[int(qubit)] for qubit in subsystem_list} + + # Parse the Hamiltonian + system = _regex_parser( + operator_str=hamiltonian_dict["h_str"], + subsystem_dims=subsystem_dims, + subsystem_list=subsystem_list, + ) + + # Extract which channels are associated with which Hamiltonian terms. + # Assumes one channel appearing in each term appearing at the end. + channels = [] + for _, ham_str in system: + chan_idx = None + + for c in CHANNEL_CHARS: + # if c in ham_str, and all characters after are digits, treat + # as channel + if c in ham_str: + if all(a.isdigit() for a in ham_str[ham_str.index(c) + 1 :]): + chan_idx = ham_str.index(c) + break + + if chan_idx is None: + channels.append(None) + else: + channels.append(ham_str[chan_idx:]) + + # evaluate coefficients + local_vars = {chan: 1.0 for chan in set(channels) if chan is not None} + local_vars.update(variables) + + evaluated_ops = [] + for op, coeff in system: + # pylint: disable=exec-used + exec(f"evaluated_coeff = {coeff}", globals(), local_vars) + evaluated_ops.append(local_vars["evaluated_coeff"] * op) + + # merge terms based on channel + static_hamiltonian = None + hamiltonian_operators = [] + reduced_channels = [] + + for channel, op in zip(channels, evaluated_ops): + # if None, add it to the static hamiltonian + if channel is None: + if static_hamiltonian is None: + static_hamiltonian = op + else: + static_hamiltonian += op + else: + channel = channel.lower() + if channel in reduced_channels: + hamiltonian_operators[reduced_channels.index(channel)] += op + else: + hamiltonian_operators.append(op) + reduced_channels.append(channel) + + # sort channels/operators according to channel ordering + if len(reduced_channels) > 0: + reduced_channels, hamiltonian_operators = zip( + *sorted(zip(reduced_channels, hamiltonian_operators)) + ) + + return static_hamiltonian, list(hamiltonian_operators), list(reduced_channels), subsystem_dims + + +def hamiltonian_pre_parse_exceptions(hamiltonian_dict: dict): + """Raises exceptions for improperly formatted or unsupported elements of + hamiltonian dict specification. + + Parameters: + hamiltonian_dict: Dictionary specification of hamiltonian. + Returns: + Raises: + QiskitError: If some part of the Hamiltonian dictionary is unsupported or invalid. + """ + + ham_str = hamiltonian_dict.get("h_str", []) + if ham_str in ([], [""]): + raise QiskitError("Hamiltonian dict requires a non-empty 'h_str' entry.") + + if hamiltonian_dict.get("qub", {}) == {}: + raise QiskitError( + "Hamiltonian dict requires non-empty 'qub' entry with subsystem dimensions." + ) + + if hamiltonian_dict.get("osc", {}) != {}: + raise QiskitError("Oscillator-type systems are not supported.") + + # verify that if terms in h_str have the divider ||, then the channels are in the valid format + for term in hamiltonian_dict["h_str"]: + malformed_text = f"""Term '{term}' does not conform to required string format. + Channels may only be specified in the format + 'aa||Cxx', where 'aa' specifies an operator, + C is a valid channel character, + and 'xx' is a string of digits.""" + + # if two vertical bars used together, check if channels in correct format + if term.count("|") == 2 and term.count("||") == 1: + # get the string reserved for channel + channel_str = term[term.index("||") + 2 :] + + # if channel string is empty + if len(channel_str) == 0: + raise QiskitError(malformed_text) + + # if first entry in channel string isn't a valid channel character + if channel_str[0] not in CHANNEL_CHARS: + raise QiskitError(malformed_text) + + # Verify either that: all remaining characters are digits, or, + # if term starts with _SUM[ and ends with ], all remaining characters + # are either digits, or starts and ends with {} + if term[-1] == "]" and len(term) > 5 and term[:5] == "_SUM[": + # drop the closing ] + channel_str = channel_str[:-1] + + # if channel string doesn't contain anything other than channel character + if len(channel_str) == 1: + raise QiskitError(malformed_text) + + # if starts with opening bracket, verify that it ends with closing bracket + if channel_str[1] == "{": + if not channel_str[-1] == "}": + raise QiskitError(malformed_text) + # otherwise verify that the remainder of terms only contains digits + elif any(not c.isdigit() for c in channel_str[1:]): + raise QiskitError(malformed_text) + else: + # if channel string doesn't contain anything other than channel character + if len(channel_str) == 1: + raise QiskitError(malformed_text) + + if any(not c.isdigit() for c in channel_str[1:]): + raise QiskitError(malformed_text) + + # if bars present but not in correct format, raise error + elif term.count("|") != 0: + raise QiskitError(malformed_text) diff --git a/test/__init__.py b/test/__init__.py index 46eac9899..be625db2a 100644 --- a/test/__init__.py +++ b/test/__init__.py @@ -11,5 +11,5 @@ # that they have been altered from the originals. """ -Qiskit ode tests +Qiskit Dynamics tests """ diff --git a/test/dynamics/pulse/__init__.py b/test/dynamics/pulse/__init__.py new file mode 100644 index 000000000..cadd41274 --- /dev/null +++ b/test/dynamics/pulse/__init__.py @@ -0,0 +1,15 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +""" +Pulse integration module tests. +""" diff --git a/test/dynamics/pulse/test_operator_from_string.py b/test/dynamics/pulse/test_operator_from_string.py new file mode 100644 index 000000000..68f76402d --- /dev/null +++ b/test/dynamics/pulse/test_operator_from_string.py @@ -0,0 +1,103 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +# pylint: disable=invalid-name + +""" +Tests for pulse.operator_from_string. +""" + +import numpy as np + +from qiskit_dynamics.pulse.backend_parser.operator_from_string import operator_from_string +from ..common import QiskitDynamicsTestCase + + +class TestOperatorFromString(QiskitDynamicsTestCase): + """Test cases for operator_from_string.""" + + def test_correct_single_ops_dim2(self): + """Test that single operators give correct outputs.""" + + ident = np.eye(2) + a = np.array([[0.0, 1.0], [0.0, 0.0]]) + adag = a.conj().transpose() + N = np.diag(np.arange(2)) + X = a + adag + Y = -1j * (a - adag) + Z = ident - 2 * N + + self.assertAllClose(operator_from_string("X", 0, {0: 2}), X) + self.assertAllClose(operator_from_string("Y", 0, {0: 2}), Y) + self.assertAllClose(operator_from_string("Z", 0, {0: 2}), Z) + self.assertAllClose(operator_from_string("a", 0, {0: 2}), a) + self.assertAllClose(operator_from_string("A", 0, {0: 2}), a) + self.assertAllClose(operator_from_string("Sm", 0, {0: 2}), a) + self.assertAllClose(operator_from_string("C", 0, {0: 2}), adag) + self.assertAllClose(operator_from_string("Sp", 0, {0: 2}), adag) + self.assertAllClose(operator_from_string("N", 0, {0: 2}), N) + self.assertAllClose(operator_from_string("O", 0, {0: 2}), N) + self.assertAllClose(operator_from_string("I", 0, {0: 2}), np.eye(2)) + + def test_correct_single_ops_dim4(self): + """Test that single operators give correct outputs.""" + + sq2 = np.sqrt(2) + sq3 = np.sqrt(3) + ident = np.eye(4) + a = np.array( + [[0.0, 1.0, 0.0, 0.0], [0.0, 0.0, sq2, 0.0], [0.0, 0.0, 0.0, sq3], [0.0, 0.0, 0.0, 0.0]] + ) + adag = a.conj().transpose() + N = np.diag(np.arange(4)) + X = a + adag + Y = -1j * (a - adag) + Z = ident - 2 * N + + self.assertAllClose(operator_from_string("X", 0, {0: 4}), X) + self.assertAllClose(operator_from_string("Y", 0, {0: 4}), Y) + self.assertAllClose(operator_from_string("Z", 0, {0: 4}), Z) + self.assertAllClose(operator_from_string("a", 0, {0: 4}), a) + self.assertAllClose(operator_from_string("A", 0, {0: 4}), a) + self.assertAllClose(operator_from_string("Sm", 0, {0: 4}), a) + self.assertAllClose(operator_from_string("C", 0, {0: 4}), adag) + self.assertAllClose(operator_from_string("Sp", 0, {0: 4}), adag) + self.assertAllClose(operator_from_string("N", 0, {0: 4}), N) + self.assertAllClose(operator_from_string("O", 0, {0: 4}), N) + self.assertAllClose(operator_from_string("I", 0, {0: 4}), ident) + + def test_ident_before(self): + """Test adding identity before the subsystem in question.""" + + out = operator_from_string("X", 2, {0: 4, 1: 2, 2: 2}) + expected = np.kron(np.array([[0.0, 1.0], [1.0, 0.0]]), np.eye(8)) + self.assertAllClose(out, expected) + + def test_ident_after(self): + """Test adding identity after the subsystem in question.""" + + out = operator_from_string("Z", 0, {0: 2, 1: 2, 2: 4}) + expected = np.kron(np.eye(8), np.array([[1.0, 0.0], [0.0, -1.0]])) + self.assertAllClose(out, expected) + + def test_ident_after_unordered(self): + """Test adding identity after the subsystem in question.""" + + out = operator_from_string("Z", 0, {2: 4, 0: 2, 1: 2}) + expected = np.kron(np.eye(8), np.array([[1.0, 0.0], [0.0, -1.0]])) + self.assertAllClose(out, expected) + + def test_ident_before_and_after(self): + """Test adding identity before and after the subsystem in question.""" + + out = operator_from_string("a", 1, {0: 2, 1: 2, 2: 4}) + expected = np.kron(np.kron(np.eye(4), np.array([[0.0, 1.0], [0.0, 0.0]])), np.eye(2)) + self.assertAllClose(out, expected) diff --git a/test/dynamics/pulse/test_string_model_parser.py b/test/dynamics/pulse/test_string_model_parser.py new file mode 100644 index 000000000..378d1938f --- /dev/null +++ b/test/dynamics/pulse/test_string_model_parser.py @@ -0,0 +1,437 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. +# pylint: disable=invalid-name + +""" +Tests for pulse.string_model_parser for importing qiskit.pulse model string representation. +""" + +import numpy as np + +from qiskit import QiskitError +from qiskit.quantum_info.operators import Operator + +from qiskit_dynamics.pulse.backend_parser.string_model_parser import ( + parse_hamiltonian_dict, + hamiltonian_pre_parse_exceptions, +) + +from qiskit_dynamics.type_utils import to_array + +from ..common import QiskitDynamicsTestCase + + +class TestHamiltonianPreParseExceptions(QiskitDynamicsTestCase): + """Tests for preparse formatting exceptions.""" + + def test_no_h_str(self): + """Test no h_str empty raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({}) + self.assertTrue("requires a non-empty 'h_str'" in str(qe.exception)) + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": []}) + self.assertTrue("requires a non-empty 'h_str'" in str(qe.exception)) + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": [""]}) + self.assertTrue("requires a non-empty 'h_str'" in str(qe.exception)) + + def test_empty_qub(self): + """Test qub dict empty raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0|||D0"]}) + self.assertTrue("requires non-empty 'qub'" in str(qe.exception)) + + def test_too_many_vertical_bars(self): + """Test that too many vertical bars raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0|||D0"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + def test_single_vertical_bar(self): + """Test that only a single vertical bar raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0|D0"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + def test_multiple_channel_error(self): + """Test multiple channels raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0||D0*D1"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + def test_divider_no_channel(self): + """Test that divider with no channel raises error.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0||"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0|"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + def test_non_digit_after_channel(self): + """Test that everything after the D or U is an int.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0||Da"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0||D1a"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + def test_nothing_after_channel(self): + """Test that everything after the D or U is an int.""" + + with self.assertRaises(QiskitError) as qe: + hamiltonian_pre_parse_exceptions({"h_str": ["a * X0||U"], "qub": {0: 2}}) + self.assertTrue("does not conform" in str(qe.exception)) + + +class TestParseHamiltonianDict(QiskitDynamicsTestCase): + """Tests for parse_hamiltonian_dict.""" + + def setUp(self): + """Build operators.""" + self.X = Operator.from_label("X").data + self.Y = Operator.from_label("Y").data + self.Z = Operator.from_label("Z").data + + dim = 4 + self.a = np.array(np.diag(np.sqrt(range(1, dim)), k=1), dtype=complex) + self.adag = self.a.conj().transpose() + self.N = np.array(np.diag(range(dim)), dtype=complex) + + def test_only_static_terms(self): + """Test a basic system.""" + ham_dict = {"h_str": ["v*np.pi*Z0"], "qub": {"0": 2}, "vars": {"v": 2.1}} + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + self.assertAllClose(static_ham, 2.1 * np.pi * self.Z) + self.assertTrue(not ham_ops) + self.assertTrue(not channels) + self.assertTrue(subsystem_dims == {0: 2}) + + def test_simple_single_q_system(self): + """Test a basic system.""" + ham_dict = { + "h_str": ["v*np.pi*Z0", "0.02*np.pi*X0||D0"], + "qub": {"0": 2}, + "vars": {"v": 2.1}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + self.assertAllClose(static_ham, 2.1 * np.pi * self.Z) + self.assertAllClose(to_array(ham_ops), [0.02 * np.pi * self.X]) + self.assertTrue(channels == ["d0"]) + self.assertTrue(subsystem_dims == {0: 2}) + + def test_simple_single_q_system_repeat_entries(self): + """Test merging of terms with same channel or no channel.""" + ham_dict = { + "h_str": ["v*np.pi*Z0", "0.02*np.pi*X0||D0", "v*np.pi*Z0", "0.02*np.pi*X0||D0"], + "qub": {"0": 2}, + "vars": {"v": 2.1}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + self.assertAllClose(static_ham, 2 * 2.1 * np.pi * self.Z) + self.assertAllClose(to_array(ham_ops), [2 * 0.02 * np.pi * self.X]) + self.assertTrue(channels == ["d0"]) + self.assertTrue(subsystem_dims == {0: 2}) + + def test_simple_single_q_system_repeat_entries_different_case(self): + """Test merging of terms with same channel or no channel, + with upper and lower case. + """ + ham_dict = { + "h_str": ["v*np.pi*Z0", "0.02*np.pi*X0||D0", "v*np.pi*Z0", "0.02*np.pi*X0||d0"], + "qub": {"0": 2}, + "vars": {"v": 2.1}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + self.assertAllClose(static_ham, 2 * 2.1 * np.pi * self.Z) + self.assertAllClose(to_array(ham_ops), [2 * 0.02 * np.pi * self.X]) + self.assertTrue(channels == ["d0"]) + self.assertTrue(subsystem_dims == {0: 2}) + + def test_simple_two_q_system(self): + """Test a two qubit system.""" + + ham_dict = { + "h_str": [ + "v0*np.pi*Z0", + "v1*np.pi*Z1", + "j*np.pi*X0*Y1", + "0.03*np.pi*X1||D1", + "0.02*np.pi*X0||D0", + ], + "qub": {"0": 2, "1": 2}, + "vars": {"v0": 2.1, "v1": 2.0, "j": 0.02}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + ident = np.eye(2) + self.assertAllClose( + static_ham, + 2.1 * np.pi * np.kron(ident, self.Z) + + 2.0 * np.pi * np.kron(self.Z, ident) + + 0.02 * np.pi * np.kron(self.Y, self.X), + ) + self.assertAllClose( + to_array(ham_ops), + [0.02 * np.pi * np.kron(ident, self.X), 0.03 * np.pi * np.kron(self.X, ident)], + ) + self.assertTrue(channels == ["d0", "d1"]) + self.assertTrue(subsystem_dims == {0: 2, 1: 2}) + + def test_simple_two_q_system_measurement_channel(self): + """Test a two qubit system with a measurement-labelled channel.""" + + ham_dict = { + "h_str": [ + "v0*np.pi*Z0", + "v1*np.pi*Z1", + "j*np.pi*X0*Y1", + "0.03*np.pi*X1||M1", + "0.02*np.pi*X0||D0", + ], + "qub": {"0": 2, "1": 2}, + "vars": {"v0": 2.1, "v1": 2.0, "j": 0.02}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + ident = np.eye(2) + self.assertAllClose( + static_ham, + 2.1 * np.pi * np.kron(ident, self.Z) + + 2.0 * np.pi * np.kron(self.Z, ident) + + 0.02 * np.pi * np.kron(self.Y, self.X), + ) + self.assertAllClose( + to_array(ham_ops), + [0.02 * np.pi * np.kron(ident, self.X), 0.03 * np.pi * np.kron(self.X, ident)], + ) + self.assertTrue(channels == ["d0", "m1"]) + self.assertTrue(subsystem_dims == {0: 2, 1: 2}) + + def test_single_oscillator_system(self): + """Test single oscillator system.""" + + ham_dict = { + "h_str": ["v*np.pi*O0", "alpha*np.pi*O0*O0", "r*np.pi*X0||D0"], + "qub": {"0": 4}, + "vars": {"v": 2.1, "alpha": -0.33, "r": 0.02}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + self.assertAllClose(static_ham, 2.1 * np.pi * self.N - 0.33 * np.pi * self.N * self.N) + self.assertAllClose(to_array(ham_ops), [0.02 * np.pi * (self.a + self.adag)]) + self.assertTrue(channels == ["d0"]) + self.assertTrue(subsystem_dims == {0: 4}) + + def test_two_oscillator_system(self): + """Test a two qubit system.""" + + ham_dict = { + "h_str": [ + "v0*np.pi*O0", + "alpha0*np.pi*O0*O0", + "v1*np.pi*O1", + "alpha1*np.pi*O1*O1", + "j*np.pi*X0*Y1", + "0.03*np.pi*X1||D1", + "0.02*np.pi*X0||D0", + ], + "qub": {"0": 4, "1": 4}, + "vars": {"v0": 2.1, "v1": 2.0, "alpha0": -0.33, "alpha1": -0.33, "j": 0.02}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + ident = np.eye(4) + + self.assertAllClose( + static_ham, + 2.1 * np.pi * np.kron(ident, self.N) + - 0.33 * np.pi * np.kron(ident, self.N * self.N) + + 2.0 * np.pi * np.kron(self.N, ident) + - 0.33 * np.pi * np.kron(self.N * self.N, ident) + + 0.02 * np.pi * np.kron(-1j * (self.a - self.adag), self.a + self.adag), + ) + self.assertAllClose( + to_array(ham_ops), + [ + 0.02 * np.pi * np.kron(ident, self.a + self.adag), + 0.03 * np.pi * np.kron(self.a + self.adag, ident), + ], + ) + self.assertTrue(channels == ["d0", "d1"]) + self.assertTrue(subsystem_dims == {0: 4, 1: 4}) + + def test_single_q_high_dim(self): + """Test single q system but higher dim.""" + ham_dict = { + "h_str": ["v*np.pi*Z0", "0.02*np.pi*X0||D0"], + "qub": {"0": 4}, + "vars": {"v": 2.1}, + } + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict(ham_dict) + + self.assertAllClose(static_ham, 2.1 * np.pi * (np.eye(4) - 2 * self.N)) + self.assertAllClose(to_array(ham_ops), [0.02 * np.pi * (self.a + self.adag)]) + self.assertTrue(channels == ["d0"]) + self.assertTrue(subsystem_dims == {0: 4}) + + def test_dagger(self): + """Test correct parsing of dagger.""" + ham_dict = { + "h_str": ["v*np.pi*dag(A0)"], + "qub": {"0": 4}, + "vars": {"v": 2.1}, + } + + static_ham, _, _, _ = parse_hamiltonian_dict(ham_dict) + self.assertAllClose(static_ham, 2.1 * np.pi * self.adag) + + def test_5q_hamiltonian_reduced(self): + """Test case for complicated Hamiltonian with reductions.""" + ham_dict = { + "h_str": [ + "_SUM[i,0,4,wq{i}/2*(I{i}-Z{i})]", + "_SUM[i,0,4,delta{i}/2*O{i}*O{i}]", + "_SUM[i,0,4,-delta{i}/2*O{i}]", + "_SUM[i,0,4,omegad{i}*X{i}||D{i}]", + "jq1q2*Sp1*Sm2", + "jq1q2*Sm1*Sp2", + "jq3q4*Sp3*Sm4", + "jq3q4*Sm3*Sp4", + "jq0q1*Sp0*Sm1", + "jq0q1*Sm0*Sp1", + "jq2q3*Sp2*Sm3", + "jq2q3*Sm2*Sp3", + "omegad1*X0||U0", + "omegad0*X1||U1", + "omegad2*X1||U2", + "omegad1*X2||U3", + "omegad3*X2||U4", + "omegad4*X3||U6", + "omegad2*X3||U5", + "omegad3*X4||U7", + ], + "qub": {"0": 4, "1": 4, "2": 4, "3": 4, "4": 4}, + "vars": { + "delta0": -2.111793476400394, + "delta1": -2.0894421352015744, + "delta2": -2.1179183671068604, + "delta3": -2.0410045431261215, + "delta4": -2.1119885565086776, + "jq0q1": 0.010495754104003914, + "jq1q2": 0.010781715511200012, + "jq2q3": 0.008920779377814226, + "jq3q4": 0.008985191651087791, + "omegad0": 0.9715458990879812, + "omegad1": 0.9803812537440838, + "omegad2": 0.9494756077681784, + "omegad3": 0.9763998543087951, + "omegad4": 0.9829308019780478, + "wq0": 32.517894442809514, + "wq1": 33.0948996120196, + "wq2": 31.74518096417169, + "wq3": 30.51062025552735, + "wq4": 32.16082685025662, + }, + } + + ident = np.eye(4, dtype=complex) + X = self.a + self.adag + X0 = np.kron(ident, X) + X1 = np.kron(X, ident) + N0 = np.kron(ident, self.N) + N1 = np.kron(self.N, ident) + + # test case for subsystems [0, 1] + + w0 = ham_dict["vars"]["wq0"] + w1 = ham_dict["vars"]["wq1"] + delta0 = ham_dict["vars"]["delta0"] + delta1 = ham_dict["vars"]["delta1"] + j01 = ham_dict["vars"]["jq0q1"] + omegad0 = ham_dict["vars"]["omegad0"] + omegad1 = ham_dict["vars"]["omegad1"] + omegad2 = ham_dict["vars"]["omegad2"] + static_ham_expected = ( + w0 * N0 + + 0.5 * delta0 * (N0 @ N0 - N0) + + w1 * N1 + + 0.5 * delta1 * (N1 @ N1 - N1) + + j01 * (np.kron(self.a, self.adag) + np.kron(self.adag, self.a)) + ) + ham_ops_expected = np.array( + [omegad0 * X0, omegad1 * X1, omegad1 * X0, omegad0 * X1, omegad2 * X1] + ) + channels_expected = ["d0", "d1", "u0", "u1", "u2"] + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict( + ham_dict, subsystem_list=[0, 1] + ) + self.assertAllClose(static_ham, static_ham_expected) + self.assertAllClose(ham_ops, ham_ops_expected) + self.assertTrue(channels == channels_expected) + self.assertTrue(subsystem_dims == {0: 4, 1: 4}) + + # test case for subsystems [3, 4] + + w3 = ham_dict["vars"]["wq3"] + w4 = ham_dict["vars"]["wq4"] + delta3 = ham_dict["vars"]["delta3"] + delta4 = ham_dict["vars"]["delta4"] + j34 = ham_dict["vars"]["jq3q4"] + omegad3 = ham_dict["vars"]["omegad3"] + omegad4 = ham_dict["vars"]["omegad4"] + omegad2 = ham_dict["vars"]["omegad2"] + static_ham_expected = ( + w3 * N0 + + 0.5 * delta3 * (N0 @ N0 - N0) + + w4 * N1 + + 0.5 * delta4 * (N1 @ N1 - N1) + + j34 * (np.kron(self.a, self.adag) + np.kron(self.adag, self.a)) + ) + ham_ops_expected = np.array( + [omegad3 * X0, omegad4 * X1, omegad2 * X0, omegad4 * X0, omegad3 * X1] + ) + channels_expected = ["d3", "d4", "u5", "u6", "u7"] + + static_ham, ham_ops, channels, subsystem_dims = parse_hamiltonian_dict( + ham_dict, subsystem_list=[3, 4] + ) + self.assertAllClose(static_ham, static_ham_expected) + self.assertAllClose(ham_ops, ham_ops_expected) + self.assertTrue(channels == channels_expected) + self.assertTrue(subsystem_dims == {3: 4, 4: 4})