diff --git a/pygsti/algorithms/contract.py b/pygsti/algorithms/contract.py index f4bf2022e..e1cb760fc 100644 --- a/pygsti/algorithms/contract.py +++ b/pygsti/algorithms/contract.py @@ -19,6 +19,7 @@ from pygsti import tools as _tools from pygsti.modelmembers import operations as _op from pygsti.modelmembers import povms as _povm +from pygsti.modelmembers import states as _state def contract(model, to_what, dataset=None, maxiter=1000000, tol=0.01, use_direct_cp=True, method="Nelder-Mead", @@ -192,7 +193,7 @@ def _contract_to_cp_direct(model, verbosity, tp_also=False, maxiter=100000, tol= printer.log(("--- Contract to %s (direct) ---" % ("CPTP" if tp_also else "CP")), 1) for (opLabel, gate) in model.operations.items(): - new_op = gate.copy() + new_op = gate.to_dense() if(tp_also): for k in range(new_op.shape[1]): new_op[0, k] = 1.0 if k == 0 else 0.0 @@ -294,7 +295,7 @@ def _contract_to_cp_direct(model, verbosity, tp_also=False, maxiter=100000, tol= #else: print "contract_to_cp_direct success in %d iterations" % it #DEBUG printer.log("Direct CP contraction of %s gate gives frobenius diff of %g" % - (opLabel, _tools.frobeniusdist(mdl.operations[opLabel], gate)), 2) + (opLabel, _tools.frobeniusdist(mdl.operations[opLabel].to_dense(), gate.to_dense())), 2) #mdl.log("Choi-Truncate to %s" % ("CPTP" if tp_also else "CP"), { 'maxiter': maxiter } ) distance = mdl.frobeniusdist(model) @@ -302,8 +303,10 @@ def _contract_to_cp_direct(model, verbosity, tp_also=False, maxiter=100000, tol= if tp_also: # TP also constrains prep vectors op_dim = mdl.dim - for rhoVec in list(mdl.preps.values()): - rhoVec[0, 0] = 1.0 / op_dim**0.25 + for lbl, rhoVec in mdl.preps.items(): + new_rho_vec = rhoVec.to_dense() + new_rho_vec[0] = 1.0 / op_dim**0.25 + mdl.preps[lbl] = new_rho_vec mdl._need_to_rebuild = True return distance, mdl @@ -315,13 +318,18 @@ def _contract_to_tp(model, verbosity): #printer.log('', 2) printer.log("--- Contract to TP ---", 1) mdl = model.copy() - for gate in list(mdl.operations.values()): - gate[0, 0] = 1.0 - for k in range(1, gate.shape[1]): gate[0, k] = 0.0 + for lbl, gate in mdl.operations.items(): + new_gate_mx = gate.to_dense() + new_gate_mx[0, 0] = 1.0 + for k in range(1, new_gate_mx.shape[1]): + new_gate_mx[0, k] = 0.0 + mdl.operations[lbl] = new_gate_mx op_dim = mdl.dim - for rhoVec in list(mdl.preps.values()): - rhoVec[0, 0] = 1.0 / op_dim**0.25 + for lbl, rhoVec in mdl.preps.items(): + new_rho_vec = rhoVec.to_dense() + new_rho_vec[0] = 1.0 / op_dim**0.25 + mdl.preps[lbl] = new_rho_vec mdl._need_to_rebuild = True distance = mdl.frobeniusdist(model) @@ -363,20 +371,20 @@ def _contract_to_valid_spam(model, verbosity=0): # rhoVec must be positive semidefinite and trace = 1 for prepLabel, rhoVec in mdl.preps.items(): - vec = rhoVec.copy() + vec = rhoVec.to_dense() #Ensure trace == 1.0 (maybe later have this be optional) firstElTarget = 1.0 / firstElTrace # TODO: make this function more robust # to multiple rhovecs -- can only use on ratio, # so maybe take average of ideal ratios for each rhoVec # and apply that? The function works fine now for just one rhovec. - if abs(firstElTarget - vec[0, 0]) > TOL: - r = firstElTarget / vec[0, 0] + if abs(firstElTarget - vec[0]) > TOL: + r = firstElTarget / vec[0] vec *= r # multiply rhovec by factor for povmLbl in list(mdl.povms.keys()): scaled_effects = [] for ELabel, EVec in mdl.povms[povmLbl].items(): - scaled_effects.append((ELabel, EVec / r)) + scaled_effects.append((ELabel, EVec.to_dense() / r)) # Note: always creates an unconstrained POVM mdl.povms[povmLbl] = _povm.UnconstrainedPOVM(scaled_effects, mdl.evotype, mdl.state_space) @@ -386,11 +394,11 @@ def _contract_to_valid_spam(model, verbosity=0): lowEval = min([ev.real for ev in _np.linalg.eigvals(mx)]) while(lowEval < -TOL): # only element with trace (even for multiple qubits) -- keep this constant and decrease others - idEl = vec[0, 0] - vec /= 1.00001; vec[0, 0] = idEl + idEl = vec[0] + vec /= 1.00001; vec[0] = idEl lowEval = min([ev.real for ev in _np.linalg.eigvals(_tools.ppvec_to_stdmx(vec))]) - diff += _np.linalg.norm(model.preps[prepLabel] - vec) + diff += _np.linalg.norm(model.preps[prepLabel].to_dense() - vec) mdl.preps[prepLabel] = vec # EVec must have eigenvals between 0 and 1 <==> positive semidefinite and trace <= 1 @@ -399,6 +407,7 @@ def _contract_to_valid_spam(model, verbosity=0): for ELabel, EVec in mdl.povms[povmLbl].items(): #if isinstance(EVec, _povm.ComplementPOVMEffect): # continue #don't contract complement vectors + EVec = EVec.to_dense() evals, evecs = _np.linalg.eig(_tools.ppvec_to_stdmx(EVec)) if(min(evals) < 0.0 or max(evals) > 1.0): if all([ev > 1.0 for ev in evals]): @@ -410,7 +419,7 @@ def _contract_to_valid_spam(model, verbosity=0): if ev > 1.0: evals[k] = 1.0 mx = _np.dot(evecs, _np.dot(_np.diag(evals), _np.linalg.inv(evecs))) vec = _tools.stdmx_to_ppvec(mx) - diff += _np.linalg.norm(model.povms[povmLbl][ELabel] - vec) + diff += _np.linalg.norm(model.povms[povmLbl][ELabel].to_dense() - vec) scaled_effects.append((ELabel, vec)) else: scaled_effects.append((ELabel, EVec)) # no scaling @@ -423,12 +432,12 @@ def _contract_to_valid_spam(model, verbosity=0): printer.log("--- Contract to valid SPAM ---", 1) printer.log(("Sum of norm(deltaE) and norm(deltaRho) = %g" % diff), 1) for (prepLabel, rhoVec) in model.preps.items(): - printer.log(" %s: %s ==> %s " % (prepLabel, str(_np.transpose(rhoVec)), - str(_np.transpose(mdl.preps[prepLabel]))), 2) + printer.log(" %s: %s ==> %s " % (prepLabel, str(_np.transpose(rhoVec.to_dense())), + str(_np.transpose(mdl.preps[prepLabel].to_dense()))), 2) for povmLbl, povm in model.povms.items(): printer.log((" %s (POVM)" % povmLbl), 2) for ELabel, EVec in povm.items(): - printer.log(" %s: %s ==> %s " % (ELabel, str(_np.transpose(EVec)), - str(_np.transpose(mdl.povms[povmLbl][ELabel]))), 2) + printer.log(" %s: %s ==> %s " % (ELabel, str(_np.transpose(EVec.to_dense())), + str(_np.transpose(mdl.povms[povmLbl][ELabel].to_dense()))), 2) return mdl # return contracted model diff --git a/pygsti/algorithms/gaugeopt.py b/pygsti/algorithms/gaugeopt.py index e43198106..5321b2388 100644 --- a/pygsti/algorithms/gaugeopt.py +++ b/pygsti/algorithms/gaugeopt.py @@ -637,13 +637,13 @@ def _objective_fn(gauge_group_el, oob_check): for opLbl in mdl.operations: wt = item_weights.get(opLbl, opWeight) ret += wt * (1.0 - _tools.entanglement_fidelity( - target_model.operations[opLbl], mdl.operations[opLbl]))**2 + target_model.operations[opLbl].to_dense(), mdl.operations[opLbl].to_dense()))**2 elif gates_metric == "tracedist": for opLbl in mdl.operations: wt = item_weights.get(opLbl, opWeight) ret += opWeight * _tools.jtracedist( - target_model.operations[opLbl], mdl.operations[opLbl]) + target_model.operations[opLbl].to_dense(), mdl.operations[opLbl].to_dense()) else: raise ValueError("Invalid gates_metric: %s" % gates_metric) @@ -666,9 +666,9 @@ def _objective_fn(gauge_group_el, oob_check): elif spam_metric == "fidelity": for preplabel, prep in mdl.preps.items(): wt = item_weights.get(preplabel, spamWeight) - rhoMx1 = _tools.vec_to_stdmx(prep, mxBasis) + rhoMx1 = _tools.vec_to_stdmx(prep.to_dense(), mxBasis) rhoMx2 = _tools.vec_to_stdmx( - target_model.preps[preplabel], mxBasis) + target_model.preps[preplabel].to_dense(), mxBasis) ret += wt * (1.0 - _tools.fidelity(rhoMx1, rhoMx2))**2 for povmlabel, povm in mdl.povms.items(): @@ -679,9 +679,9 @@ def _objective_fn(gauge_group_el, oob_check): elif spam_metric == "tracedist": for preplabel, prep in mdl.preps.items(): wt = item_weights.get(preplabel, spamWeight) - rhoMx1 = _tools.vec_to_stdmx(prep, mxBasis) + rhoMx1 = _tools.vec_to_stdmx(prep.to_dense(), mxBasis) rhoMx2 = _tools.vec_to_stdmx( - target_model.preps[preplabel], mxBasis) + target_model.preps[preplabel].to_dense(), mxBasis) ret += wt * _tools.tracedist(rhoMx1, rhoMx2) for povmlabel, povm in mdl.povms.items(): @@ -768,11 +768,11 @@ def _cptp_penalty_jac_fill(cp_penalty_vec_grad_to_fill, mdl_pre, mdl_post, dS = _np.rollaxis(dS, 2) # shape (n, d1, d2) for i, (gl, gate) in enumerate(mdl_post.operations.items()): - pre_op = mdl_pre.operations[gl] + pre_op = mdl_pre.operations[gl].to_dense() #get sgn(chi-matrix) == d(|chi|_Tr)/dchi in std basis # so sgnchi == d(|chi_std|_Tr)/dchi_std - chi = _tools.fast_jamiolkowski_iso_std(gate, op_basis) + chi = _tools.fast_jamiolkowski_iso_std(gate.to_dense(), op_basis) assert(_np.linalg.norm(chi - chi.T.conjugate()) < 1e-4), \ "chi should be Hermitian!" @@ -785,7 +785,7 @@ def _cptp_penalty_jac_fill(cp_penalty_vec_grad_to_fill, mdl_pre, mdl_post, # transform). This shuffle op commutes with the derivative, so that # dchi_std/dp := d(M(G'))/dp = M(d(S_inv*G*S)/dp) = M( d(S_inv)*G*S + S_inv*G*dS ) # = M( (-S_inv*dS*S_inv)*G*S + S_inv*G*dS ) = M( S_inv*(-dS*S_inv*G*S) + G*dS ) - left = -1 * _np.dot(dS, gate) # shape (n,d1,d2) + left = -1 * _np.dot(dS, gate.to_dense()) # shape (n,d1,d2) right = _np.swapaxes(_np.dot(pre_op, dS), 0, 1) # shape (d1, n, d2) -> (n,d1,d2) result = _np.swapaxes(_np.dot(S_inv, left + right), 0, 1) # shape (n, d1, d2) diff --git a/pygsti/extras/rpe/rpeconstruction.py b/pygsti/extras/rpe/rpeconstruction.py index 04f6a9f9e..df65d567c 100644 --- a/pygsti/extras/rpe/rpeconstruction.py +++ b/pygsti/extras/rpe/rpeconstruction.py @@ -95,9 +95,9 @@ def create_parameterized_rpe_model(alpha_true, epsilon_true, aux_rot, spam_depol effect_labels=ELabels, effect_expressions=EExpressions) outputModel.operations[loose_axis_gate_label] = \ - _np.dot(_np.dot(_np.linalg.inv(modelAux1.operations[auxiliary_axis_gate_label]), - outputModel.operations[loose_axis_gate_label]), - modelAux1.operations[auxiliary_axis_gate_label]) + _np.linalg.inv(modelAux1.operations[auxiliary_axis_gate_label].to_dense())@\ + outputModel.operations[loose_axis_gate_label].to_dense()@\ + modelAux1.operations[auxiliary_axis_gate_label].to_dense() outputModel = outputModel.depolarize(op_noise=gate_depol, spam_noise=spam_depol) diff --git a/pygsti/extras/rpe/rpetools.py b/pygsti/extras/rpe/rpetools.py index b7c684d38..9a40b3d16 100644 --- a/pygsti/extras/rpe/rpetools.py +++ b/pygsti/extras/rpe/rpetools.py @@ -257,7 +257,7 @@ def extract_alpha(model, rpeconfig_inst): The value of alpha for the input model. """ op_label = rpeconfig_inst.fixed_axis_gate_label - decomp = _decompose_gate_matrix(model.operations[op_label]) + decomp = _decompose_gate_matrix(model.operations[op_label].to_dense()) alphaVal = decomp['pi rotations'] * _np.pi return alphaVal @@ -284,7 +284,7 @@ def extract_epsilon(model, rpeconfig_inst): The value of epsilon for the input model. """ op_label = rpeconfig_inst.loose_axis_gate_label - decomp = _decompose_gate_matrix(model.operations[op_label]) + decomp = _decompose_gate_matrix(model.operations[op_label].to_dense()) epsilonVal = decomp['pi rotations'] * _np.pi return epsilonVal @@ -313,10 +313,10 @@ def extract_theta(model, rpeconfig_inst): The value of theta for the input model. """ op_label = rpeconfig_inst.loose_axis_gate_label - decomp = _decompose_gate_matrix(model.operations[op_label]) + decomp = _decompose_gate_matrix(model.operations[op_label].to_dense()) target_axis = rpeconfig_inst.loose_axis_target - decomp = _decompose_gate_matrix(model.operations[op_label]) + decomp = _decompose_gate_matrix(model.operations[op_label].to_dense()) thetaVal = _np.real_if_close([_np.arccos( _np.dot(decomp['axis of rotation'], target_axis))])[0] if thetaVal > _np.pi / 2: diff --git a/pygsti/modelmembers/instruments/instrument.py b/pygsti/modelmembers/instruments/instrument.py index 2e0a5a3a0..035d5fde5 100644 --- a/pygsti/modelmembers/instruments/instrument.py +++ b/pygsti/modelmembers/instruments/instrument.py @@ -226,7 +226,7 @@ def num_elements(self): ------- int """ - return sum([g.size for g in self.values()]) + return sum([g._ptr.size for g in self.values()]) @property def num_params(self): diff --git a/pygsti/modelmembers/instruments/tpinstrument.py b/pygsti/modelmembers/instruments/tpinstrument.py index a5662ebd2..0510a7e8e 100644 --- a/pygsti/modelmembers/instruments/tpinstrument.py +++ b/pygsti/modelmembers/instruments/tpinstrument.py @@ -86,11 +86,12 @@ def __init__(self, op_matrices, evotype="default", state_space=None, called_from if op_matrices is not None: if isinstance(op_matrices, dict): - matrix_list = [(k, v) for k, v in op_matrices.items()] # gives definite ordering + matrix_list = [(k, v[()]) if (isinstance(v,_np.ndarray) and v.shape == ()) else (k,v) for k, v in op_matrices.items()] # gives definite ordering elif isinstance(op_matrices, list): - matrix_list = op_matrices # assume it's is already an ordered (key,value) list + matrix_list = [(k, v[()]) if (isinstance(v,_np.ndarray) and v.shape == ()) else (k,v) for k, v in op_matrices] else: raise ValueError("Invalid `op_matrices` arg of type %s" % type(op_matrices)) + matrix_list = [(k,v.to_dense()) if isinstance(v, _mm.ModelMember) else (k,v) for k,v in matrix_list] assert(len(matrix_list) > 0 or state_space is not None), \ "Must specify `state_space` when there are no instrument members!" @@ -99,7 +100,7 @@ def __init__(self, op_matrices, evotype="default", state_space=None, called_from evotype = _Evotype.cast(evotype, state_space=state_space) # Create gate objects that are used to parameterize this instrument - MT_mx = sum([v for k, v in matrix_list]) # sum-of-instrument-members matrix + MT_mx = sum([v for _, v in matrix_list]) # sum-of-instrument-members matrix MT = _op.FullTPOp(MT_mx, None, evotype, state_space) self.param_ops.append(MT) @@ -114,14 +115,6 @@ def __init__(self, op_matrices, evotype="default", state_space=None, called_from items = [(k, _TPInstrumentOp(self.param_ops, i, None)) for i, (k, v) in enumerate(matrix_list)] - #DEBUG - #print("POST INIT PARAM GATES:") - #for i,v in enumerate(self.param_ops): - # print(i,":\n",v) - # - #print("POST nINIT ITEMS:") - #for k,v in items: - # print(k,":\n",v) else: assert(state_space is not None), "`state_space` cannot be `None` when there are no members!" evotype = _Evotype.cast(evotype, state_space=state_space) @@ -273,7 +266,7 @@ def num_elements(self): # same as in Instrument CONSOLIDATE? ------- int """ - return sum([g.size for g in self.values()]) + return sum([g._ptr.size for g in self.values()]) @property def num_params(self): # same as in Instrument CONSOLIDATE? diff --git a/pygsti/modelmembers/instruments/tpinstrumentop.py b/pygsti/modelmembers/instruments/tpinstrumentop.py index 59005d2b4..b72ff1c53 100644 --- a/pygsti/modelmembers/instruments/tpinstrumentop.py +++ b/pygsti/modelmembers/instruments/tpinstrumentop.py @@ -122,13 +122,13 @@ def _construct_matrix(self): nEls = self.num_instrument_elements self._ptr.flags.writeable = True if self.index < nEls - 1: - self._ptr[:, :] = _np.asarray(self.relevant_param_ops[1] # i.e. param_ops[self.index + 1] - + self.relevant_param_ops[0]) # i.e. param_ops[0] + self._ptr[:, :] = _np.asarray(self.relevant_param_ops[1].to_dense() # i.e. param_ops[self.index + 1] + + self.relevant_param_ops[0].to_dense()) # i.e. param_ops[0] else: assert(self.index == nEls - 1), \ "Invalid index %d > %d" % (self.index, nEls - 1) - self._ptr[:, :] = _np.asarray(-sum(self.relevant_param_ops) # all instrument param_ops == relevant - - (nEls - 3) * self.relevant_param_ops[0]) + self._ptr[:, :] = _np.asarray(-sum([op.to_dense() for op in self.relevant_param_ops]) # all instrument param_ops == relevant + - (nEls - 3) * self.relevant_param_ops[0].to_dense()) assert(self._ptr.shape == (self.dim, self.dim)) self._ptr.flags.writeable = False diff --git a/pygsti/modelmembers/operations/__init__.py b/pygsti/modelmembers/operations/__init__.py index 579a85523..c82154433 100644 --- a/pygsti/modelmembers/operations/__init__.py +++ b/pygsti/modelmembers/operations/__init__.py @@ -15,7 +15,7 @@ from .composederrorgen import ComposedErrorgen from .composedop import ComposedOp -from .denseop import DenseOperator, DenseOperatorInterface +from .denseop import DenseOperator from .depolarizeop import DepolarizeOp from .eigpdenseop import EigenvalueParamDenseOp from .embeddederrorgen import EmbeddedErrorgen diff --git a/pygsti/modelmembers/operations/denseop.py b/pygsti/modelmembers/operations/denseop.py index fec40c75e..0015a7c69 100644 --- a/pygsti/modelmembers/operations/denseop.py +++ b/pygsti/modelmembers/operations/denseop.py @@ -126,158 +126,16 @@ def check_deriv_wrt_params(operation, deriv_to_check=None, wrt_filter=None, eps= " norm diff = %g" % _np.linalg.norm(fd_deriv - deriv_to_check)) # pragma: no cover - -class DenseOperatorInterface(object): - """ - Adds a numpy-array-mimicing interface onto an operation object. - """ - - def __init__(self): - pass - - @property - def _ptr(self): - raise NotImplementedError("Derived classes must implement the _ptr property!") - - def _ptr_has_changed(self): - """ Derived classes should override this function to handle rep updates - when the `_ptr` property is changed. """ - pass - - def to_array(self): - """ - Return the array used to identify this operation within its range of possible values. - - For instance, if the operation is a unitary operation, this returns a - unitary matrix regardless of the evolution type. The related :meth:`to_dense` - method, in contrast, returns the dense representation of the operation, which - varies by evolution type. - - Note: for efficiency, this doesn't copy the underlying data, so - the caller should copy this data before modifying it. - - Returns - ------- - numpy.ndarray - """ - return _np.asarray(self._ptr) - # *must* be a numpy array for Cython arg conversion - - def to_sparse(self, on_space='minimal'): - """ - Return the operation as a sparse matrix. - - Parameters - ---------- - on_space : {'minimal', 'Hilbert', 'HilbertSchmidt'} - The space that the returned dense operation acts upon. For unitary matrices and bra/ket vectors, - use `'Hilbert'`. For superoperator matrices and super-bra/super-ket vectors use `'HilbertSchmidt'`. - `'minimal'` means that `'Hilbert'` is used if possible given this operator's evolution type, and - otherwise `'HilbertSchmidt'` is used. - - Returns - ------- - scipy.sparse.csr_matrix - """ - return _sps.csr_matrix(self.to_dense(on_space)) - - def __copy__(self): - # We need to implement __copy__ because we defer all non-existing - # attributes to self.base (a numpy array) which *has* a __copy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - cpy.__dict__.update(self.__dict__) - return cpy - - def __deepcopy__(self, memo): - # We need to implement __deepcopy__ because we defer all non-existing - # attributes to self.base (a numpy array) which *has* a __deepcopy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - memo[id(self)] = cpy - for k, v in self.__dict__.items(): - setattr(cpy, k, _copy.deepcopy(v, memo)) - return cpy - - #Access to underlying ndarray - def __getitem__(self, key): - self.dirty = True - return self._ptr.__getitem__(key) - - def __getslice__(self, i, j): - self.dirty = True - return self.__getitem__(slice(i, j)) # Called for A[:] - - def __setitem__(self, key, val): - self.dirty = True - ret = self._ptr.__setitem__(key, val) - self._ptr_has_changed() - return ret - - def __getattr__(self, attr): - #use __dict__ so no chance for recursive __getattr__ - #ret = getattr(self.__dict__['_rep'].base, attr) - ret = getattr(self._ptr, attr) - self.dirty = True - return ret - - def __str__(self): - s = "%s with shape %s\n" % (self.__class__.__name__, str(self._ptr.shape)) - s += _mt.mx_to_string(self._ptr, width=4, prec=2) - return s - - #Mimic array behavior - def __pos__(self): return self._ptr - def __neg__(self): return -self._ptr - def __abs__(self): return abs(self._ptr) - def __add__(self, x): return self._ptr + x - def __radd__(self, x): return x + self._ptr - def __sub__(self, x): return self._ptr - x - def __rsub__(self, x): return x - self._ptr - def __mul__(self, x): return self._ptr * x - def __rmul__(self, x): return x * self._ptr - def __truediv__(self, x): return self._ptr / x - def __rtruediv__(self, x): return x / self._ptr - def __floordiv__(self, x): return self._ptr // x - def __rfloordiv__(self, x): return x // self._ptr - def __pow__(self, x): return self._ptr ** x - def __eq__(self, x): return _np.array_equal(self._ptr, x) - def __len__(self): return len(self._ptr) - def __int__(self): return int(self._ptr) - def __long__(self): return int(self._ptr) - def __float__(self): return float(self._ptr) - def __complex__(self): return complex(self._ptr) - - -class DenseOperator(DenseOperatorInterface, _KrausOperatorInterface, _LinearOperator): +def __str__(self): + s = "%s with shape %s\n" % (self.__class__.__name__, str(self._ptr.shape)) + s += _mt.mx_to_string(self._ptr, width=4, prec=2) + return s +class DenseOperator(_KrausOperatorInterface, _LinearOperator): """ - TODO: update docstring An operator that behaves like a dense super-operator matrix. This class is the common base class for more specific dense operators. - Parameters - ---------- - mx : numpy.ndarray - The operation as a dense process matrix. - - basis : Basis or {'pp','gm','std'} or None - The basis used to construct the Hilbert-Schmidt space representation - of this state as a super-operator. If None, certain functionality, - such as access to Kraus operators, will be unavailable. - - evotype : Evotype or str - The evolution type. The special value `"default"` is equivalent - to specifying the value of `pygsti.evotypes.Evotype.default_evotype`. - - state_space : StateSpace, optional - The state space for this operation. If `None` a default state space - with the appropriate number of qubits is used. - Attributes ---------- base : numpy.ndarray @@ -310,6 +168,25 @@ def from_kraus_operators(cls, kraus_operators, basis='pp', evotype="default", st return cls(superop, basis, evotype, state_space) def __init__(self, mx, basis, evotype, state_space=None): + """ + Parameters + ---------- + mx : numpy.ndarray + The operation as a dense process matrix. + + basis : Basis or {'pp','gm','std'} or None + The basis used to construct the Hilbert-Schmidt space representation + of this state as a super-operator. If None, certain functionality, + such as access to Kraus operators, will be unavailable. + + evotype : Evotype or str + The evolution type. The special value `"default"` is equivalent + to specifying the value of `pygsti.evotypes.Evotype.default_evotype`. + + state_space : StateSpace, optional + The state space for this operation. If `None` a default state space + with the appropriate number of qubits is used. + """ mx = _LinearOperator.convert_to_matrix(mx) state_space = _statespace.default_space_for_dim(mx.shape[0]) if (state_space is None) \ else _statespace.StateSpace.cast(state_space) @@ -317,7 +194,6 @@ def __init__(self, mx, basis, evotype, state_space=None): self._basis = _Basis.cast(basis, state_space.dim) if (basis is not None) else None # for Hilbert-Schmidt space rep = evotype.create_dense_superop_rep(mx, self._basis, state_space) _LinearOperator.__init__(self, rep, evotype) - DenseOperatorInterface.__init__(self) @property def _ptr(self): @@ -349,6 +225,44 @@ def to_dense(self, on_space='minimal'): """ return self._rep.to_dense(on_space) # both types of possible reps implement 'to_dense' + def to_array(self): + """ + Return the array used to identify this operation within its range of possible values. + + For instance, if the operation is a unitary operation, this returns a + unitary matrix regardless of the evolution type. The related :meth:`to_dense` + method, in contrast, returns the dense representation of the operation, which + varies by evolution type. + + Note: for efficiency, this doesn't copy the underlying data, so + the caller should copy this data before modifying it. + + Returns + ------- + numpy.ndarray + """ + return _np.asarray(self._ptr) + # *must* be a numpy array for Cython arg conversion + + def to_sparse(self, on_space='minimal'): + """ + Return the operation as a sparse matrix. + + Parameters + ---------- + on_space : {'minimal', 'Hilbert', 'HilbertSchmidt'} + The space that the returned dense operation acts upon. For unitary matrices and bra/ket vectors, + use `'Hilbert'`. For superoperator matrices and super-bra/super-ket vectors use `'HilbertSchmidt'`. + `'minimal'` means that `'Hilbert'` is used if possible given this operator's evolution type, and + otherwise `'HilbertSchmidt'` is used. + + Returns + ------- + scipy.sparse.csr_matrix + """ + return _sps.csr_matrix(self.to_dense(on_space)) + + def to_memoized_dict(self, mmg_memo): """Create a serializable dict with references to other objects in the memo. @@ -450,31 +364,17 @@ def set_kraus_operators(self, kraus_operators): superop = _bt.change_basis(std_superop, 'std', self._basis) self.set_dense(superop) # this may fail if derived class doesn't allow it + def __str__(self): + s = "%s with shape %s\n" % (self.__class__.__name__, str(self._ptr.shape)) + s += _mt.mx_to_string(self._ptr, width=4, prec=2) + return s -class DenseUnitaryOperator(DenseOperatorInterface, _KrausOperatorInterface, _LinearOperator): +class DenseUnitaryOperator(_KrausOperatorInterface, _LinearOperator): """ - TODO: update docstring An operator that behaves like a dense (unitary) operator matrix. This class is the common base class for more specific dense operators. - Parameters - ---------- - mx : numpy.ndarray - The operation as a dense process matrix. - - basis : Basis or {'pp','gm','std'}, optional - The basis used to construct the Hilbert-Schmidt space representation - of this state as a super-operator. - - evotype : Evotype or str - The evolution type. The special value `"default"` is equivalent - to specifying the value of `pygsti.evotypes.Evotype.default_evotype`. - - state_space : StateSpace, optional - The state space for this operation. If `None` a default state space - with the appropriate number of qubits is used. - Attributes ---------- base : numpy.ndarray @@ -524,11 +424,29 @@ def quick_init(cls, unitary_mx, superop_mx, basis, evotype, state_space): self._basis = basis _LinearOperator.__init__(self, rep, evotype) - DenseOperatorInterface.__init__(self) return self def __init__(self, mx, basis, evotype, state_space): - """ Initialize a new LinearOperator """ + """ + Initialize a new DenseUnitaryOperator + + Parameters + ---------- + mx : numpy.ndarray + The operation as a dense process matrix. + + basis : Basis or {'pp','gm','std'}, optional + The basis used to construct the Hilbert-Schmidt space representation + of this state as a super-operator. + + evotype : Evotype or str + The evolution type. The special value `"default"` is equivalent + to specifying the value of `pygsti.evotypes.Evotype.default_evotype`. + + state_space : StateSpace, optional + The state space for this operation. If `None` a default state space + with the appropriate number of qubits is used. + """ mx = _LinearOperator.convert_to_matrix(mx) state_space = _statespace.default_space_for_udim(mx.shape[0]) if (state_space is None) \ else _statespace.StateSpace.cast(state_space) @@ -554,7 +472,6 @@ def __init__(self, mx, basis, evotype, state_space): self._basis = basis _LinearOperator.__init__(self, rep, evotype) - DenseOperatorInterface.__init__(self) @property def _ptr(self): @@ -591,6 +508,43 @@ def to_dense(self, on_space='minimal'): return self._unitary return self._rep.to_dense(on_space) # both types of possible reps implement 'to_dense' + def to_array(self): + """ + Return the array used to identify this operation within its range of possible values. + + For instance, if the operation is a unitary operation, this returns a + unitary matrix regardless of the evolution type. The related :meth:`to_dense` + method, in contrast, returns the dense representation of the operation, which + varies by evolution type. + + Note: for efficiency, this doesn't copy the underlying data, so + the caller should copy this data before modifying it. + + Returns + ------- + numpy.ndarray + """ + return _np.asarray(self._ptr) + # *must* be a numpy array for Cython arg conversion + + def to_sparse(self, on_space='minimal'): + """ + Return the operation as a sparse matrix. + + Parameters + ---------- + on_space : {'minimal', 'Hilbert', 'HilbertSchmidt'} + The space that the returned dense operation acts upon. For unitary matrices and bra/ket vectors, + use `'Hilbert'`. For superoperator matrices and super-bra/super-ket vectors use `'HilbertSchmidt'`. + `'minimal'` means that `'Hilbert'` is used if possible given this operator's evolution type, and + otherwise `'HilbertSchmidt'` is used. + + Returns + ------- + scipy.sparse.csr_matrix + """ + return _sps.csr_matrix(self.to_dense(on_space)) + def to_memoized_dict(self, mmg_memo): """Create a serializable dict with references to other objects in the memo. @@ -654,3 +608,7 @@ def set_kraus_operators(self, kraus_operators): assert(len(kraus_operators) == 1), "Length of `kraus_operators` must == 1 for a unitary channel!" self.set_dense(kraus_operators[0]) + def __str__(self): + s = "%s with shape %s\n" % (self.__class__.__name__, str(self._ptr.shape)) + s += _mt.mx_to_string(self._ptr, width=4, prec=2) + return s diff --git a/pygsti/modelmembers/operations/fullarbitraryop.py b/pygsti/modelmembers/operations/fullarbitraryop.py index 42f42b3af..912963441 100644 --- a/pygsti/modelmembers/operations/fullarbitraryop.py +++ b/pygsti/modelmembers/operations/fullarbitraryop.py @@ -82,7 +82,7 @@ def num_params(self): int the number of independent parameters. """ - return self.size + return self._ptr.size def to_vector(self): """ diff --git a/pygsti/modelmembers/operations/linearop.py b/pygsti/modelmembers/operations/linearop.py index 7cd00c073..f54be2d23 100644 --- a/pygsti/modelmembers/operations/linearop.py +++ b/pygsti/modelmembers/operations/linearop.py @@ -671,7 +671,7 @@ def deriv_wrt_params(self, wrt_filter=None): Array of derivatives with shape (dimension^2, num_params) """ if self.num_params == 0: - derivMx = _np.zeros((self.size, 0), 'd') + derivMx = _np.zeros((self._ptr.size, 0), 'd') if wrt_filter is None: return derivMx else: diff --git a/pygsti/modelmembers/povms/conjugatedeffect.py b/pygsti/modelmembers/povms/conjugatedeffect.py index df2745e33..add8ec6ca 100644 --- a/pygsti/modelmembers/povms/conjugatedeffect.py +++ b/pygsti/modelmembers/povms/conjugatedeffect.py @@ -17,117 +17,12 @@ from pygsti.modelmembers import term as _term from pygsti.tools import matrixtools as _mt - -class DenseEffectInterface(object): - """ - Adds a numpy-array-mimicing interface onto a POVM effect object. - """ - # Note: this class may not really be necessary, and maybe methods should just be - # placed within ConjugatedStatePOVMEffect? - - @property - def _ptr(self): - raise NotImplementedError("Derived classes must implement the _ptr property!") - - def _ptr_has_changed(self): - """ Derived classes should override this function to handle rep updates - when the `_ptr` property is changed. """ - pass - - @property - def columnvec(self): - """ - Direct access the the underlying data as column vector, i.e, a (dim,1)-shaped array. - """ - bv = self._ptr.view() - bv.shape = (bv.size, 1) # 'base' is by convention a (N,1)-shaped array - return bv - - def __copy__(self): - # We need to implement __copy__ because we defer all non-existing - # attributes to self.columnvec (a numpy array) which *has* a __copy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - cpy.__dict__.update(self.__dict__) - return cpy - - def __deepcopy__(self, memo): - # We need to implement __deepcopy__ because we defer all non-existing - # attributes to self.columnvec (a numpy array) which *has* a __deepcopy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - memo[id(self)] = cpy - for k, v in self.__dict__.items(): - setattr(cpy, k, _copy.deepcopy(v, memo)) - return cpy - - #Access to underlying array - def __getitem__(self, key): - if not isinstance(key, (int, _np.int64)): # don't set dirty flag if returning a single element - self.dirty = True - return self.columnvec.__getitem__(key) - - def __getslice__(self, i, j): - self.dirty = True - return self.__getitem__(slice(i, j)) # Called for A[:] - - def __setitem__(self, key, val): - self.dirty = True - ret = self.columnvec.__setitem__(key, val) - self._ptr_has_changed() - return ret - - def __getattr__(self, attr): - #use __dict__ so no chance for recursive __getattr__ - if '_rep' in self.__dict__: # sometimes in loading __getattr__ gets called before the instance is loaded - ret = getattr(self.columnvec, attr) - else: - raise AttributeError("No attribute:", attr) - self.dirty = True - return ret - - #Mimic array - def __pos__(self): return self.columnvec - def __neg__(self): return -self.columnvec - def __abs__(self): return abs(self.columnvec) - def __add__(self, x): return self.columnvec + x - def __radd__(self, x): return x + self.columnvec - def __sub__(self, x): return self.columnvec - x - def __rsub__(self, x): return x - self.columnvec - def __mul__(self, x): return self.columnvec * x - def __rmul__(self, x): return x * self.columnvec - def __truediv__(self, x): return self.columnvec / x - def __rtruediv__(self, x): return x / self.columnvec - def __floordiv__(self, x): return self.columnvec // x - def __rfloordiv__(self, x): return x // self.columnvec - def __pow__(self, x): return self.columnvec ** x - def __eq__(self, x): return self.columnvec == x - def __len__(self): return len(self.columnvec) - def __int__(self): return int(self.columnvec) - def __long__(self): return int(self.columnvec) - def __float__(self): return float(self.columnvec) - def __complex__(self): return complex(self.columnvec) - - -class ConjugatedStatePOVMEffect(DenseEffectInterface, _POVMEffect): +class ConjugatedStatePOVMEffect(_POVMEffect): """ - TODO: update docstring - A POVM effect vector that behaves like a numpy array. - This class is the common base class for parameterizations of an effect vector - that have a dense representation and can be accessed like a numpy array. + that have a dense representation. - Parameters - ---------- - vec : numpy.ndarray - The POVM effect vector as a dense numpy array. - evotype : {"statevec", "densitymx"} - The evolution type. Attributes ---------- @@ -140,6 +35,20 @@ class ConjugatedStatePOVMEffect(DenseEffectInterface, _POVMEffect): """ def __init__(self, state, called_from_reduce=False): + """ + Parameters + ---------- + vec : numpy.ndarray + The POVM effect vector as a dense numpy array. + + evotype : {"statevec", "densitymx"} + The evolution type. + + called_from_reduce : bool, optional (default False) + Special flag used when pickling. Users should not need to + interact with this flag directly. + """ + self.state = state evotype = state._evotype rep = evotype.create_conjugatedstate_effect_rep(state._rep) diff --git a/pygsti/modelmembers/states/densestate.py b/pygsti/modelmembers/states/densestate.py index b645f7795..faa4724b4 100644 --- a/pygsti/modelmembers/states/densestate.py +++ b/pygsti/modelmembers/states/densestate.py @@ -23,136 +23,12 @@ from pygsti.tools import optools as _ot -class DenseStateInterface(object): +class DenseState(_State): """ - Adds a numpy-array-mimicing interface onto a state object. - """ - - def __init__(self): - pass - - @property - def _ptr(self): - raise NotImplementedError("Derived classes must implement the _ptr property!") - - def _ptr_has_changed(self): - """ Derived classes should override this function to handle rep updates - when the `_ptr` property is changed. """ - pass - - def to_array(self): - """ - Return the array used to identify this state within its range of possible values. - - For instance, if the state is a pure state, this returns a complex pure-state - vector regardless of the evolution type. The related :meth:`to_dense` - method, in contrast, returns the dense representation of the state, which - varies by evolution type. - - Returns - ------- - numpy.ndarray - """ - return self._ptr # *must* be a numpy array for Cython arg conversion - - @property - def columnvec(self): - """ - Direct access the the underlying data as column vector, i.e, a (dim,1)-shaped array. - """ - bv = self._ptr.view() - bv.shape = (bv.size, 1) # 'base' is by convention a (N,1)-shaped array - return bv - - def __copy__(self): - # We need to implement __copy__ because we defer all non-existing - # attributes to self.columnvec (a numpy array) which *has* a __copy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - cpy.__dict__.update(self.__dict__) - return cpy - - def __deepcopy__(self, memo): - # We need to implement __deepcopy__ because we defer all non-existing - # attributes to self.columnvec (a numpy array) which *has* a __deepcopy__ - # implementation that we don't want to use, as it results in just a - # copy of the numpy array. - cls = self.__class__ - cpy = cls.__new__(cls) - memo[id(self)] = cpy - for k, v in self.__dict__.items(): - setattr(cpy, k, _copy.deepcopy(v, memo)) - return cpy - - #Access to underlying array - def __getitem__(self, key): - self.dirty = True - return self.columnvec.__getitem__(key) - - def __getslice__(self, i, j): - self.dirty = True - return self.__getitem__(slice(i, j)) # Called for A[:] - - def __setitem__(self, key, val): - self.dirty = True - ret = self.columnvec.__setitem__(key, val) - self._ptr_has_changed() - return ret - - def __getattr__(self, attr): - #use __dict__ so no chance for recursive __getattr__ - if '_rep' in self.__dict__: # sometimes in loading __getattr__ gets called before the instance is loaded - ret = getattr(self.columnvec, attr) - else: - raise AttributeError("No attribute:", attr) - self.dirty = True - return ret - - #Mimic array - def __pos__(self): return self.columnvec - def __neg__(self): return -self.columnvec - def __abs__(self): return abs(self.columnvec) - def __add__(self, x): return self.columnvec + x - def __radd__(self, x): return x + self.columnvec - def __sub__(self, x): return self.columnvec - x - def __rsub__(self, x): return x - self.columnvec - def __mul__(self, x): return self.columnvec * x - def __rmul__(self, x): return x * self.columnvec - def __truediv__(self, x): return self.columnvec / x - def __rtruediv__(self, x): return x / self.columnvec - def __floordiv__(self, x): return self.columnvec // x - def __rfloordiv__(self, x): return x // self.columnvec - def __pow__(self, x): return self.columnvec ** x - def __eq__(self, x): return self.columnvec == x - def __len__(self): return len(self.columnvec) - def __int__(self): return int(self.columnvec) - def __long__(self): return int(self.columnvec) - def __float__(self): return float(self.columnvec) - def __complex__(self): return complex(self.columnvec) - - def __str__(self): - s = "%s with dimension %d\n" % (self.__class__.__name__, self.dim) - s += _mt.mx_to_string(self.to_dense(on_space='minimal'), width=4, prec=2) - return s - - -class DenseState(DenseStateInterface, _State): - """ - TODO: update docstring - A state preparation vector that is interfaced/behaves as a dense super-ket (a numpy array). + A state preparation vector that is parameterized as a dense super-ket. This class is the common base class for parameterizations of a state vector - that have a dense representation and can be accessed like a numpy array. - - Parameters - ---------- - vec : numpy.ndarray - The SPAM vector as a dense numpy array. - - evotype : {"statevec", "densitymx"} - The evolution type. + that have a dense representation. Attributes ---------- @@ -165,6 +41,21 @@ class DenseState(DenseStateInterface, _State): """ def __init__(self, vec, basis, evotype, state_space): + """ + Parameters + ---------- + vec : numpy.ndarray + The SPAM vector as a dense numpy array. + + evotype : {"statevec", "densitymx"} + The evolution type. + + state_space : `StateSpace`, optional (default None) + The state space for this state. If `None` a default state space + with the appropriate number of qubits is used. + """ + + vec = _State._to_vector(vec) if state_space is None: state_space = _statespace.default_space_for_dim(vec.shape[0]) @@ -175,8 +66,6 @@ def __init__(self, vec, basis, evotype, state_space): rep = evotype.create_dense_state_rep(vec, self._basis, state_space) _State.__init__(self, rep, evotype) - DenseStateInterface.__init__(self) - #assert(self._ptr.flags['C_CONTIGUOUS'] and self._ptr.flags['OWNDATA']) # not true for TPState @property def _ptr(self): @@ -254,8 +143,12 @@ def _is_similar(self, other, rtol, atol): return self._ptr.shape == other._ptr.shape # similar (up to params) if have same data shape + def __str__(self): + s = "%s with dimension %d\n" % (self.__class__.__name__, self.dim) + s += _mt.mx_to_string(self.to_dense(on_space='minimal'), width=4, prec=2) + return s -class DensePureState(DenseStateInterface, _State): +class DensePureState(_State): """ TODO: docstring - a state that is interfaced as a dense ket """ @@ -286,7 +179,6 @@ def __init__(self, purevec, basis, evotype, state_space): self._purevec = purevec; self._basis = basis _State.__init__(self, rep, evotype) - DenseStateInterface.__init__(self) @property def _ptr(self): @@ -368,3 +260,8 @@ def _is_similar(self, other, rtol, atol): the same local structure, i.e., not considering parameter values or submembers """ return self._ptr.shape == other._ptr.shape # similar (up to params) if have same data shape + + def __str__(self): + s = "%s with dimension %d\n" % (self.__class__.__name__, self.dim) + s += _mt.mx_to_string(self.to_dense(on_space='minimal'), width=4, prec=2) + return s \ No newline at end of file diff --git a/pygsti/modelmembers/states/fullstate.py b/pygsti/modelmembers/states/fullstate.py index feaa06d77..260dd2771 100644 --- a/pygsti/modelmembers/states/fullstate.py +++ b/pygsti/modelmembers/states/fullstate.py @@ -76,7 +76,7 @@ def num_params(self): int the number of independent parameters. """ - return self.size + return self._ptr.size def to_vector(self): """ diff --git a/pygsti/models/explicitmodel.py b/pygsti/models/explicitmodel.py index 01eed989b..5eb06affe 100644 --- a/pygsti/models/explicitmodel.py +++ b/pygsti/models/explicitmodel.py @@ -801,14 +801,15 @@ def _tpdist(self): """ penalty = 0.0 for operationMx in list(self.operations.values()): - penalty += abs(operationMx[0, 0] - 1.0)**2 - for k in range(1, operationMx.shape[1]): - penalty += abs(operationMx[0, k])**2 + mx = operationMx.to_dense() + penalty += abs(mx[0, 0] - 1.0)**2 + for k in range(1, mx.shape[1]): + penalty += abs(mx[0, k])**2 op_dim = self.state_space.dim firstEl = 1.0 / op_dim**0.25 for rhoVec in list(self.preps.values()): - penalty += abs(rhoVec[0, 0] - firstEl)**2 + penalty += abs(rhoVec.to_dense()[0] - firstEl)**2 return _np.sqrt(penalty) @@ -1343,8 +1344,8 @@ def kick(self, absmag=1.0, bias=0, seed=None): kicked_gs = self.copy() rndm = _np.random.RandomState(seed) for opLabel, gate in self.operations.items(): - delta = absmag * 2.0 * (rndm.random_sample(gate.shape) - 0.5) + bias - kicked_gs.operations[opLabel] = _op.FullArbitraryOp(kicked_gs.operations[opLabel] + delta) + delta = absmag * 2.0 * (rndm.random_sample(gate.to_dense().shape) - 0.5) + bias + kicked_gs.operations[opLabel] = _op.FullArbitraryOp(kicked_gs.operations[opLabel].to_dense() + delta) #Note: does not alter intruments! return kicked_gs diff --git a/pygsti/objectivefns/objectivefns.py b/pygsti/objectivefns/objectivefns.py index e97c7829d..3d57040cf 100644 --- a/pygsti/objectivefns/objectivefns.py +++ b/pygsti/objectivefns/objectivefns.py @@ -5771,7 +5771,7 @@ def _cptp_penalty(mdl, prefactor, op_basis): a (real) 1D array of length len(mdl.operations). """ return prefactor * _np.sqrt(_np.array([_tools.tracenorm( - _tools.fast_jamiolkowski_iso_std(gate, op_basis) + _tools.fast_jamiolkowski_iso_std(gate.to_dense(), op_basis) ) for gate in mdl.operations.values()], 'd')) @@ -5840,7 +5840,7 @@ def _cptp_penalty_jac_fill(cp_penalty_vec_grad_to_fill, mdl, prefactor, op_basis #get sgn(chi-matrix) == d(|chi|_Tr)/dchi in std basis # so sgnchi == d(|chi_std|_Tr)/dchi_std - chi = _tools.fast_jamiolkowski_iso_std(gate, op_basis) + chi = _tools.fast_jamiolkowski_iso_std(gate.to_dense(), op_basis) assert(_np.linalg.norm(chi - chi.T.conjugate()) < 1e-4), \ "chi should be Hermitian!" @@ -5915,7 +5915,7 @@ def _spam_penalty_jac_fill(spam_penalty_vec_grad_to_fill, mdl, prefactor, op_bas #get sgn(denMx) == d(|denMx|_Tr)/d(denMx) in std basis # dmDim = denMx.shape[0] - denmx = _tools.vec_to_stdmx(prepvec, op_basis) + denmx = _tools.vec_to_stdmx(prepvec.to_dense(), op_basis) assert(_np.linalg.norm(denmx - denmx.T.conjugate()) < 1e-4), \ "denMx should be Hermitian!" @@ -5950,7 +5950,7 @@ def _spam_penalty_jac_fill(spam_penalty_vec_grad_to_fill, mdl, prefactor, op_bas nparams = effectvec.num_params #get sgn(EMx) == d(|EMx|_Tr)/d(EMx) in std basis - emx = _tools.vec_to_stdmx(effectvec, op_basis) + emx = _tools.vec_to_stdmx(effectvec.to_dense(), op_basis) # dmDim = EMx.shape[0] assert(_np.linalg.norm(emx - emx.T.conjugate()) < 1e-4), \ "EMx should be Hermitian!" diff --git a/pygsti/tools/jamiolkowski.py b/pygsti/tools/jamiolkowski.py index 3cea8291d..cfd60a4ac 100644 --- a/pygsti/tools/jamiolkowski.py +++ b/pygsti/tools/jamiolkowski.py @@ -375,7 +375,7 @@ def magnitudes_of_negative_choi_eigenvalues(model): """ ret = [] for (_, gate) in model.operations.items(): - J = jamiolkowski_iso(gate, model.basis, choi_mx_basis=model.basis.create_simple_equivalent('std')) + J = jamiolkowski_iso(gate.to_dense(on_space='HilbertSchmidt'), model.basis, choi_mx_basis=model.basis.create_simple_equivalent('std')) evals = _np.linalg.eigvals(J) # could use eigvalsh, but wary of this since eigh can be wrong... for ev in evals: ret.append(-ev.real if ev.real < 0 else 0.0) diff --git a/pygsti/tools/optools.py b/pygsti/tools/optools.py index a3354bf6a..241d53e02 100644 --- a/pygsti/tools/optools.py +++ b/pygsti/tools/optools.py @@ -2225,8 +2225,8 @@ def project_model(model, target_model, gsDict[p].set_all_parameterizations("full") NpDict[p] = 0 - errgens = [error_generator(model.operations[gl], - target_model.operations[gl], + errgens = [error_generator(model.operations[gl].to_dense(), + target_model.operations[gl].to_dense(), target_model.basis, gen_type, logG_weight) for gl in opLabels] @@ -2256,7 +2256,7 @@ def project_model(model, target_model, lnd_error_gen = _np.tensordot(HBlk.block_data, HGens, (0, 0)) + \ _np.tensordot(otherBlk.block_data, otherGens, ((0, 1), (0, 1))) - targetOp = target_model.operations[gl] + targetOp = target_model.operations[gl].to_dense() if 'H' in projectiontypes: gsDict['H'].operations[gl] = operation_from_error_generator( diff --git a/scripts/api_names.yaml b/scripts/api_names.yaml index c07036445..be4ff4d5a 100644 --- a/scripts/api_names.yaml +++ b/scripts/api_names.yaml @@ -1556,8 +1556,6 @@ objects: prs_as_polys: prs_as_polynomials test_map: null operation.py: - BasedDenseOperatorInterface: - __name__: null CliffordOp: __name__: null ComposedDenseOp: @@ -1606,11 +1604,6 @@ objects: DenseOperator: __name__: null base: null - DenseOperatorInterface: - __name__: null - deriv_wrt_params: null - todense: to_dense - tosparse: to_sparse DepolarizeOp: __name__: null copy: null diff --git a/test/unit/construction/test_modelconstruction.py b/test/unit/construction/test_modelconstruction.py index cff5fad14..d1435f4ec 100644 --- a/test/unit/construction/test_modelconstruction.py +++ b/test/unit/construction/test_modelconstruction.py @@ -106,8 +106,8 @@ def test_build_explicit_model(self): gateset2b = mc.create_explicit_model_from_expressions([('Q0',)], ['Gi', 'Gx', 'Gy'], ["I(Q0)", "X(pi/2,Q0)", "Y(pi/2,Q0)"], effect_labels=['1', '0']) - self.assertArraysAlmostEqual(model.effects['0'], gateset2b.effects['1']) - self.assertArraysAlmostEqual(model.effects['1'], gateset2b.effects['0']) + self.assertArraysAlmostEqual(model.effects['0'].to_dense(), gateset2b.effects['1'].to_dense()) + self.assertArraysAlmostEqual(model.effects['1'].to_dense(), gateset2b.effects['0'].to_dense()) # This is slightly confusing. Single qubit rotations are always stored in "pp" basis internally # UPDATE: now this isn't even allowed, as the 'densitymx' type represents states as *real* vectors. @@ -667,21 +667,21 @@ def _test_leakA(self): leakA_ans = np.array([[0., 1., 0.], [1., 0., 0.], [0., 0., 1.]], 'd') - self.assertArraysAlmostEqual(self.leakA, leakA_ans) + self.assertArraysAlmostEqual(self.leakA.to_dense(), leakA_ans) def _test_rotXa(self): rotXa_ans = np.array([[1., 0., 0., 0.], [0., 1., 0., 0.], [0., 0., 0, -1.], [0., 0., 1., 0]], 'd') - self.assertArraysAlmostEqual(self.rotXa, rotXa_ans) + self.assertArraysAlmostEqual(self.rotXa.to_dense(), rotXa_ans) def _test_rotX2(self): rotX2_ans = np.array([[1., 0., 0., 0.], [0., 1., 0., 0.], [0., 0., -1., 0.], [0., 0., 0., -1.]], 'd') - self.assertArraysAlmostEqual(self.rotX2, rotX2_ans) + self.assertArraysAlmostEqual(self.rotX2.to_dense(), rotX2_ans) def _test_rotLeak(self): rotLeak_ans = np.array([[0.5, 0., 0., -0.5, 0.70710678], @@ -689,7 +689,7 @@ def _test_rotLeak(self): [0., 0., 0., 0., 0.], [0.5, 0., 0., -0.5, -0.70710678], [0.70710678, 0., 0., 0.70710678, 0.]], 'd') - self.assertArraysAlmostEqual(self.rotLeak, rotLeak_ans) + self.assertArraysAlmostEqual(self.rotLeak.to_dense(), rotLeak_ans) def _test_leakB(self): leakB_ans = np.array([[0.5, 0., 0., -0.5, 0.70710678], @@ -697,7 +697,7 @@ def _test_leakB(self): [0., 0., 0., 0., 0.], [-0.5, 0., 0., 0.5, 0.70710678], [0.70710678, 0., 0., 0.70710678, 0.]], 'd') - self.assertArraysAlmostEqual(self.leakB, leakB_ans) + self.assertArraysAlmostEqual(self.leakB.to_dense(), leakB_ans) def _test_rotXb(self): rotXb_ans = np.array([[1., 0., 0., 0., 0., 0.], @@ -706,7 +706,7 @@ def _test_rotXb(self): [0., 0., 0., -1., 0., 0.], [0., 0., 0., 0., 1., 0.], [0., 0., 0., 0., 0., 1.]], 'd') - self.assertArraysAlmostEqual(self.rotXb, rotXb_ans) + self.assertArraysAlmostEqual(self.rotXb.to_dense(), rotXb_ans) def _test_CnotA(self): CnotA_ans = np.array([[1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], @@ -725,7 +725,7 @@ def _test_CnotA(self): [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0], [0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]) - self.assertArraysAlmostEqual(self.CnotA, CnotA_ans) + self.assertArraysAlmostEqual(self.CnotA.to_dense(), CnotA_ans) def _test_CnotB(self): CnotB_ans = np.array([[1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], @@ -745,7 +745,7 @@ def _test_CnotB(self): [0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0]]) - self.assertArraysAlmostEqual(self.CnotB, CnotB_ans) + self.assertArraysAlmostEqual(self.CnotB.to_dense(), CnotB_ans) def test_raises_on_bad_basis(self): with self.assertRaises(AssertionError): diff --git a/test/unit/construction/test_qutrit.py b/test/unit/construction/test_qutrit.py index 166695be8..2f7200a58 100644 --- a/test/unit/construction/test_qutrit.py +++ b/test/unit/construction/test_qutrit.py @@ -11,7 +11,7 @@ def test_ideal_qutrit(self): def test_noisy_qutrit(self): mdl_sim = qutrit.create_qutrit_model(error_scale=0.1, similarity=True, seed=1234, basis='qt') mdl_ideal = qutrit.create_qutrit_model(error_scale=0.1, similarity=True, seed=1234, basis='qt') - self.assertArraysAlmostEqual(mdl_sim['Gi', 'QT'], mdl_ideal['Gi', 'QT']) + self.assertArraysAlmostEqual(mdl_sim['Gi', 'QT'].to_dense(), mdl_ideal['Gi', 'QT'].to_dense()) #just test building a gate in the qutrit basis # Can't do this b/c need a 'T*' triplet space designator for "triplet space" and it doesn't seem diff --git a/test/unit/extras/interpygate/test_construction.py b/test/unit/extras/interpygate/test_construction.py index 5f6a4e5ac..25e6406d7 100644 --- a/test/unit/extras/interpygate/test_construction.py +++ b/test/unit/extras/interpygate/test_construction.py @@ -119,7 +119,7 @@ def test_create_opfactory(self): interpolator_and_args='linear') op = opfactory_linear.create_op([0,np.pi/4]) op.from_vector([1]) - self.assertArraysAlmostEqual(op, self.static_target) + self.assertArraysAlmostEqual(op.to_dense(), self.static_target) if USE_CSAPS: opfactory_spline = interp.InterpolatedOpFactory.create_by_interpolating_physical_process( @@ -128,7 +128,7 @@ def test_create_opfactory(self): interpolator_and_args='spline') op = opfactory_spline.create_op([0,np.pi/4]) op.from_vector([1]) - self.assertArraysAlmostEqual(op, self.static_target) + self.assertArraysAlmostEqual(op.to_dense(), self.static_target) interpolator_and_args = (_linND, {'rescale': True}) opfactory_custom = opfactory_spline = interp.InterpolatedOpFactory.create_by_interpolating_physical_process( @@ -137,7 +137,7 @@ def test_create_opfactory(self): interpolator_and_args=interpolator_and_args) op = opfactory_custom.create_op([0,np.pi/4]) op.from_vector([1]) - self.assertArraysAlmostEqual(op, self.static_target) + self.assertArraysAlmostEqual(op.to_dense(), self.static_target) diff --git a/test/unit/extras/rb/test_theory.py b/test/unit/extras/rb/test_theory.py index 5c22ac779..90b0b429a 100644 --- a/test/unit/extras/rb/test_theory.py +++ b/test/unit/extras/rb/test_theory.py @@ -55,7 +55,7 @@ def setUpClass(cls): Zrot_channel = ot.unitary_to_pauligate(Zrot_unitary) for key in cls.target_model.operations.keys(): - cls.mdl.operations[key] = np.dot(Zrot_channel, cls.target_model.operations[key]) + cls.mdl.operations[key] = np.dot(Zrot_channel, cls.target_model.operations[key].to_dense()) class RBTheoryWeightedInfidelityTester(GaugeTransformBase, InfidelityBase, BaseCase): @@ -72,17 +72,17 @@ def setUpClass(cls): depmap_X = np.array([[1., 0., 0., 0.], [0., lx, 0., 0.], [0., 0., lx, 0.], [0, 0., 0., lx]]) ly = 1. - depol_strength_Y depmap_Y = np.array([[1., 0., 0., 0.], [0., ly, 0., 0.], [0., 0., ly, 0.], [0, 0., 0., ly]]) - cls.mdl.operations['Gx'] = np.dot(depmap_X, cls.target_model.operations['Gx']) - cls.mdl.operations['Gy'] = np.dot(depmap_Y, cls.target_model.operations['Gy']) + cls.mdl.operations['Gx'] = np.dot(depmap_X, cls.target_model.operations['Gx'].to_dense()) + cls.mdl.operations['Gy'] = np.dot(depmap_Y, cls.target_model.operations['Gy'].to_dense()) Gx_weight = 1 Gy_weight = 2 cls.weights = {'Gx': Gx_weight, 'Gy': Gy_weight} - GxAGI = ot.average_gate_infidelity(cls.mdl.operations['Gx'], cls.target_model.operations['Gx']) - GyAGI = ot.average_gate_infidelity(cls.mdl.operations['Gy'], cls.target_model.operations['Gy']) + GxAGI = ot.average_gate_infidelity(cls.mdl.operations['Gx'].to_dense(), cls.target_model.operations['Gx'].to_dense()) + GyAGI = ot.average_gate_infidelity(cls.mdl.operations['Gy'].to_dense(), cls.target_model.operations['Gy'].to_dense()) cls.expected_AGI = (Gx_weight * GxAGI + Gy_weight * GyAGI) / (Gx_weight + Gy_weight) - GxAEI = ot.entanglement_infidelity(cls.mdl.operations['Gx'], cls.target_model.operations['Gx']) - GyAEI = ot.entanglement_infidelity(cls.mdl.operations['Gy'], cls.target_model.operations['Gy']) + GxAEI = ot.entanglement_infidelity(cls.mdl.operations['Gx'].to_dense(), cls.target_model.operations['Gx'].to_dense()) + GyAEI = ot.entanglement_infidelity(cls.mdl.operations['Gy'].to_dense(), cls.target_model.operations['Gy'].to_dense()) cls.expected_EI = (Gx_weight * GxAEI + Gy_weight * GyAEI) / (Gx_weight + Gy_weight) diff --git a/test/unit/modelmembers/test_modelmembergraph.py b/test/unit/modelmembers/test_modelmembergraph.py index eac99a510..6968832bd 100644 --- a/test/unit/modelmembers/test_modelmembergraph.py +++ b/test/unit/modelmembers/test_modelmembergraph.py @@ -21,14 +21,16 @@ def test_explicit_model_comparisons(self): # Change parameter, similar but not equivalent ex_mdl3 = ex_mdl1.copy() - ex_mdl3.operations['Gxpi2', 0][0, 0] = 0.0 + op_mx = ex_mdl3.operations['Gxpi2', 0].to_dense() + op_mx[0,0]= 0.0 + ex_mdl3.operations['Gxpi2', 0].set_dense(op_mx) ex_mmg3 = ex_mdl3.create_modelmember_graph() self.assertTrue(ex_mmg3.is_similar(ex_mmg1)) self.assertFalse(ex_mmg3.is_equivalent(ex_mmg1)) # Change parameterization, not similar or equivalent ex_mdl4 = ex_mdl1.copy() - ex_mdl4.operations['Gxpi2', 0] = _op.StaticArbitraryOp(ex_mdl4.operations['Gxpi2', 0]) + ex_mdl4.operations['Gxpi2', 0] = _op.StaticArbitraryOp(ex_mdl4.operations['Gxpi2', 0].to_dense()) ex_mmg4 = ex_mdl4.create_modelmember_graph() self.assertFalse(ex_mmg4.is_similar(ex_mmg1)) self.assertFalse(ex_mmg4.is_equivalent(ex_mmg1)) diff --git a/test/unit/modelmembers/test_operation.py b/test/unit/modelmembers/test_operation.py index d3b44deef..f8d927f97 100644 --- a/test/unit/modelmembers/test_operation.py +++ b/test/unit/modelmembers/test_operation.py @@ -139,45 +139,9 @@ def test_set_value_raises_on_bad_size(self): with self.assertRaises((ValueError, AssertionError)): self.gate.set_dense(np.zeros((1, 1), 'd')) # bad size - def test_arithmetic(self): - result = self.gate + self.gate - self.assertEqual(type(result), np.ndarray) - result = self.gate + (-self.gate) - self.assertEqual(type(result), np.ndarray) - result = self.gate - self.gate - self.assertEqual(type(result), np.ndarray) - result = self.gate - abs(self.gate) - self.assertEqual(type(result), np.ndarray) - result = 2 * self.gate - self.assertEqual(type(result), np.ndarray) - result = self.gate * 2 - self.assertEqual(type(result), np.ndarray) - result = 2 / self.gate - self.assertEqual(type(result), np.ndarray) - result = self.gate / 2 - self.assertEqual(type(result), np.ndarray) - result = self.gate // 2 - self.assertEqual(type(result), np.ndarray) - result = self.gate**2 - self.assertEqual(type(result), np.ndarray) - result = self.gate.transpose() - self.assertEqual(type(result), np.ndarray) - - M = np.identity(4, 'd') - - result = self.gate + M - self.assertEqual(type(result), np.ndarray) - result = self.gate - M - self.assertEqual(type(result), np.ndarray) - result = M + self.gate - self.assertEqual(type(result), np.ndarray) - result = M - self.gate - self.assertEqual(type(result), np.ndarray) - - class MutableDenseOpBase(DenseOpBase): def test_set_value(self): - M = np.asarray(self.gate) # gate as a matrix + M = np.asarray(self.gate.to_dense()) # gate as a matrix self.gate.set_dense(M) # TODO assert correctness @@ -185,31 +149,9 @@ def test_transform(self): gate_copy = self.gate.copy() T = FullGaugeGroupElement(np.identity(4, 'd')) gate_copy.transform_inplace(T) - self.assertArraysAlmostEqual(gate_copy, self.gate) + self.assertArraysAlmostEqual(gate_copy.to_dense(), self.gate.to_dense()) # TODO test a non-trivial case - def test_element_accessors(self): - e1 = self.gate[1, 1] - e2 = self.gate[1][1] - self.assertAlmostEqual(e1, e2) - - s1 = self.gate[1, :] - s2 = self.gate[1] - s3 = self.gate[1][:] - a1 = self.gate[:] - self.assertArraysAlmostEqual(s1, s2) - self.assertArraysAlmostEqual(s1, s3) - - s4 = self.gate[2:4, 1] - - self.gate[1, 1] = e1 - self.gate[1, :] = s1 - self.gate[1] = s1 - self.gate[2:4, 1] = s4 - - result = len(self.gate) - # TODO assert correctness - def test_depolarize(self): dp = self.gate.depolarize(0.05) dp = self.gate.depolarize([0.05, 0.10, 0.15]) @@ -221,7 +163,7 @@ def test_rotate(self): class ImmutableDenseOpBase(DenseOpBase): def test_raises_on_set_value(self): - M = np.asarray(self.gate) # gate as a matrix + M = np.asarray(self.gate.to_dense()) # gate as a matrix with self.assertRaises(ValueError): self.gate.set_dense(M) @@ -230,24 +172,6 @@ def test_raises_on_transform(self): with self.assertRaises((ValueError, NotImplementedError)): self.gate.transform_inplace(T) - def test_element_accessors(self): - e1 = self.gate[1, 1] - e2 = self.gate[1][1] - self.assertAlmostEqual(e1, e2) - - s1 = self.gate[1, :] - s2 = self.gate[1] - s3 = self.gate[1][:] - a1 = self.gate[:] - self.assertArraysAlmostEqual(s1, s2) - self.assertArraysAlmostEqual(s1, s3) - - s4 = self.gate[2:4, 1] - - result = len(self.gate) - # TODO assert correctness - - class DenseOpTester(ImmutableDenseOpBase, BaseCase): n_params = 0 @@ -271,7 +195,6 @@ def test_convert_to_matrix_raises_on_bad_input(self): with self.assertRaises(ValueError): op.DenseOperator.convert_to_matrix(bad_mx) - class StaticStdOpTester(BaseCase): def test_statevec(self): std_unitaries = itgs.standard_gatename_unitaries() @@ -394,22 +317,6 @@ def test_convert(self): conv = op.convert(self.gate, "full TP", "gm") # TODO assert correctness - def test_first_row_read_only(self): - # check that first row is read-only - e1 = self.gate[0, 0] - with self.assertRaises(ValueError): - self.gate[0, 0] = e1 - with self.assertRaises(ValueError): - self.gate[0][0] = e1 - with self.assertRaises(ValueError): - self.gate[0, :] = [e1, 0, 0, 0] - with self.assertRaises(ValueError): - self.gate[0][:] = [e1, 0, 0, 0] - with self.assertRaises(ValueError): - self.gate[0, 1:2] = [0] - with self.assertRaises(ValueError): - self.gate[0][1:2] = [0] - class AffineShiftOpTester(DenseOpBase, BaseCase): n_params = 3 @@ -419,40 +326,14 @@ def build_gate(): return op.AffineShiftOp(mat) def test_set_dense(self): - M = np.asarray(self.gate) # gate as a matrix + M = np.asarray(self.gate.to_dense()) # gate as a matrix self.gate.set_dense(M) def test_transform(self): gate_copy = self.gate.copy() T = FullGaugeGroupElement(np.identity(4, 'd')) gate_copy.transform_inplace(T) - self.assertArraysAlmostEqual(gate_copy, self.gate) - - def test_element_accessors(self): - e1 = self.gate[1, 1] - e2 = self.gate[1][1] - self.assertAlmostEqual(e1, e2) - - s1 = self.gate[1, :] - s2 = self.gate[1] - s3 = self.gate[1][:] - a1 = self.gate[:] - self.assertArraysAlmostEqual(s1, s2) - self.assertArraysAlmostEqual(s1, s3) - - s4 = self.gate[2:4, 1] - - def test_set_elements(self): - gate_copy = self.gate.copy() - - #allowed sets: - gate_copy[1,0] = 2 - gate_copy[2,0] = 2 - - #unallowed sets: - with self.assertRaises(ValueError): - gate_copy[1,1] = 2 - + self.assertArraysAlmostEqual(gate_copy.to_dense(), self.gate.to_dense()) class StaticOpTester(ImmutableDenseOpBase, BaseCase): n_params = 0 diff --git a/test/unit/modelmembers/test_spamvec.py b/test/unit/modelmembers/test_spamvec.py index 2b80225f0..4b74c105d 100644 --- a/test/unit/modelmembers/test_spamvec.py +++ b/test/unit/modelmembers/test_spamvec.py @@ -91,61 +91,16 @@ def test_vector_conversion(self): self.vec.from_vector(v) deriv = self.vec.deriv_wrt_params() # TODO assert correctness - - def test_element_accessors(self): - a = self.vec[:] - b = self.vec[0] - #with self.assertRaises(ValueError): - # self.vec.shape = (2,2) #something that would affect the shape?? - - self.vec_as_str = str(self.vec) - a1 = self.vec[:] # invoke getslice method - # TODO assert correctness def test_copy(self): vec_copy = self.vec.copy() self.assertArraysAlmostEqual(vec_copy.to_dense(), self.vec.to_dense()) self.assertEqual(type(vec_copy), type(self.vec)) - def test_arithmetic(self): - result = self.vec + self.vec - self.assertEqual(type(result), np.ndarray) - result = self.vec + (-self.vec) - self.assertEqual(type(result), np.ndarray) - result = self.vec - self.vec - self.assertEqual(type(result), np.ndarray) - result = self.vec - abs(self.vec) - self.assertEqual(type(result), np.ndarray) - result = 2 * self.vec - self.assertEqual(type(result), np.ndarray) - result = self.vec * 2 - self.assertEqual(type(result), np.ndarray) - result = 2 / self.vec - self.assertEqual(type(result), np.ndarray) - result = self.vec / 2 - self.assertEqual(type(result), np.ndarray) - result = self.vec // 2 - self.assertEqual(type(result), np.ndarray) - result = self.vec**2 - self.assertEqual(type(result), np.ndarray) - result = self.vec.transpose() - self.assertEqual(type(result), np.ndarray) - - V = np.ones((4, 1), 'd') - - result = self.vec + V - self.assertEqual(type(result), np.ndarray) - result = self.vec - V - self.assertEqual(type(result), np.ndarray) - result = V + self.vec - self.assertEqual(type(result), np.ndarray) - result = V - self.vec - self.assertEqual(type(result), np.ndarray) - - + class MutableDenseStateBase(DenseStateBase): def test_set_value(self): - v = np.asarray(self.vec) + v = np.asarray(self.vec.to_dense()) self.vec.set_dense(v) # TODO assert correctness @@ -162,7 +117,7 @@ def test_depolarize(self): class ImmutableDenseStateBase(DenseStateBase): def test_raises_on_set_value(self): - v = np.asarray(self.vec) + v = np.asarray(self.vec.to_dense()) with self.assertRaises(ValueError): self.vec.set_dense(v) diff --git a/test/unit/objects/test_explicitmodel.py b/test/unit/objects/test_explicitmodel.py index 3c4d39a90..29d1baeb1 100644 --- a/test/unit/objects/test_explicitmodel.py +++ b/test/unit/objects/test_explicitmodel.py @@ -51,14 +51,14 @@ def test_randomize_with_unitary(self): def test_rotate_1q(self): sslbls = ExplicitStateSpace("Q0") - rotXPi = create_operation("X(pi,Q0)", sslbls, "pp") - rotXPiOv2 = create_operation("X(pi/2,Q0)", sslbls, "pp") - rotYPiOv2 = create_operation("Y(pi/2,Q0)", sslbls, "pp") + rotXPi = create_operation("X(pi,Q0)", sslbls, "pp").to_dense() + rotXPiOv2 = create_operation("X(pi/2,Q0)", sslbls, "pp").to_dense() + rotYPiOv2 = create_operation("Y(pi/2,Q0)", sslbls, "pp").to_dense() gateset_rot = self.model.rotate((np.pi / 2, 0, 0)) # rotate all gates by pi/2 about X axis - self.assertArraysAlmostEqual(gateset_rot['Gi'], rotXPiOv2) - self.assertArraysAlmostEqual(gateset_rot['Gx'], rotXPi) - self.assertArraysAlmostEqual(gateset_rot['Gx'], np.dot(rotXPiOv2, rotXPiOv2)) - self.assertArraysAlmostEqual(gateset_rot['Gy'], np.dot(rotXPiOv2, rotYPiOv2)) + self.assertArraysAlmostEqual(gateset_rot['Gi'].to_dense(), rotXPiOv2) + self.assertArraysAlmostEqual(gateset_rot['Gx'].to_dense(), rotXPi) + self.assertArraysAlmostEqual(gateset_rot['Gx'].to_dense(), np.dot(rotXPiOv2, rotXPiOv2)) + self.assertArraysAlmostEqual(gateset_rot['Gy'].to_dense(), np.dot(rotXPiOv2, rotYPiOv2)) def test_rotate_2q(self): gateset_2q_rot = self.gateset_2q.rotate(rotate=list(np.zeros(15, 'd'))) @@ -81,9 +81,9 @@ def test_depolarize(self): [0, 0, 0.9, 0], [0, -0.9, 0, 0]], 'd') gateset_dep = self.model.depolarize(op_noise=0.1) - self.assertArraysAlmostEqual(gateset_dep['Gi'], Gi_dep) - self.assertArraysAlmostEqual(gateset_dep['Gx'], Gx_dep) - self.assertArraysAlmostEqual(gateset_dep['Gy'], Gy_dep) + self.assertArraysAlmostEqual(gateset_dep['Gi'].to_dense(), Gi_dep) + self.assertArraysAlmostEqual(gateset_dep['Gx'].to_dense(), Gx_dep) + self.assertArraysAlmostEqual(gateset_dep['Gy'].to_dense(), Gy_dep) def test_depolarize_with_spam_noise(self): gateset_spam = self.model.depolarize(spam_noise=0.1) diff --git a/test/unit/objects/test_instrument.py b/test/unit/objects/test_instrument.py index b11850a5c..b730612e8 100644 --- a/test/unit/objects/test_instrument.py +++ b/test/unit/objects/test_instrument.py @@ -51,12 +51,13 @@ def setUp(self): self.n_elements = 32 self.model = std.target_model() - E = self.model.povms['Mdefault']['0'] - Erem = self.model.povms['Mdefault']['1'] + E = self.model.povms['Mdefault']['0'].to_dense() + E = E.reshape((len(E),1)) + Erem = self.model.povms['Mdefault']['1'].to_dense() + Erem = Erem.reshape((len(Erem),1)) self.Gmz_plus = np.dot(E, E.T) self.Gmz_minus = np.dot(Erem, Erem.T) # XXX is this used? - self.povm_ident = self.model.povms['Mdefault']['0'] + self.model.povms['Mdefault']['1'] self.instrument = self.constructor({'plus': self.Gmz_plus, 'minus': self.Gmz_minus}) self.model.instruments['Iz'] = self.instrument super(InstrumentInstanceBase, self).setUp() diff --git a/test/unit/objects/test_model.py b/test/unit/objects/test_model.py index 8f5f0e0c9..691760278 100644 --- a/test/unit/objects/test_model.py +++ b/test/unit/objects/test_model.py @@ -611,12 +611,12 @@ def test_transform(self): # TODO is this needed? for opLabel in cp.operations: - self.assertArraysAlmostEqual(cp[opLabel], Tinv @ self.model[opLabel] @ T) + self.assertArraysAlmostEqual(cp[opLabel].to_dense(), Tinv @ self.model[opLabel].to_dense() @ T) for prepLabel in cp.preps: - self.assertArraysAlmostEqual(cp[prepLabel], Tinv @ self.model[prepLabel]) + self.assertArraysAlmostEqual(cp[prepLabel].to_dense(), Tinv @ self.model[prepLabel].to_dense()) for povmLabel in cp.povms: for effectLabel, eVec in cp.povms[povmLabel].items(): - self.assertArraysAlmostEqual(eVec, np.transpose(T) @ self.model.povms[povmLabel][effectLabel]) + self.assertArraysAlmostEqual(eVec.to_dense(), np.transpose(T) @ self.model.povms[povmLabel][effectLabel].to_dense()) def test_gpindices(self): # Test instrument construction with elements whose gpindices @@ -629,7 +629,9 @@ def test_gpindices(self): v = mdl.to_vector() Np = mdl.num_params gate_with_gpindices = FullArbitraryOp(np.identity(4, 'd')) - gate_with_gpindices[0, :] = v[0:4] + updated_op = gate_with_gpindices.to_dense() + updated_op[0,:] = v[0:4] + gate_with_gpindices.set_dense(updated_op) gate_with_gpindices.set_gpindices(np.concatenate( (np.arange(0, 4), np.arange(Np, Np + 12))), mdl) # manually set gpindices mdl.operations['Gnew2'] = gate_with_gpindices @@ -645,7 +647,7 @@ def test_check_paramvec_raises_on_error(self): self.model._check_paramvec(debug=True) # param vec is now out of sync! def test_probs_warns_on_nan_in_input(self): - self.model['rho0'][:] = np.nan + self.model['rho0'].set_dense(np.array([np.nan, np.nan, np.nan, np.nan])) with self.assertWarns(Warning): self.model.probabilities(self.gatestring1) @@ -658,11 +660,11 @@ def test_set_operation_matrix(self): Gi_test_matrix = np.random.random((4, 4)) Gi_test_matrix[0, :] = [1, 0, 0, 0] # so TP mode works self.model["Gi"] = Gi_test_matrix # set operation matrix - self.assertArraysAlmostEqual(self.model['Gi'], Gi_test_matrix) + self.assertArraysAlmostEqual(self.model['Gi'].to_dense(), Gi_test_matrix) Gi_test_dense_op = FullArbitraryOp(Gi_test_matrix) self.model["Gi"] = Gi_test_dense_op # set gate object - self.assertArraysAlmostEqual(self.model['Gi'], Gi_test_matrix) + self.assertArraysAlmostEqual(self.model['Gi'].to_dense(), Gi_test_matrix) class StaticModelTester(StaticModelBase, StandardMethodBase, BaseCase): @@ -671,7 +673,7 @@ def test_set_operation_matrix(self): Gi_test_matrix = np.random.random((4, 4)) Gi_test_dense_op = FullArbitraryOp(Gi_test_matrix) self.model["Gi"] = Gi_test_dense_op # set gate object - self.assertArraysAlmostEqual(self.model['Gi'], Gi_test_matrix) + self.assertArraysAlmostEqual(self.model['Gi'].to_dense(), Gi_test_matrix) def test_bulk_fill_dprobs_with_high_smallness_threshold(self): self.skipTest("TODO should probably warn user?") diff --git a/test/unit/tools/fixtures.py b/test/unit/tools/fixtures.py index 0d84edc48..06526400d 100644 --- a/test/unit/tools/fixtures.py +++ b/test/unit/tools/fixtures.py @@ -69,6 +69,7 @@ def mdl_lsgst(self): resource_alloc={'mem_limit': self.CM + 1024**3}, verbosity=0 ) + print(models[-1]) return models[-1] diff --git a/test/unit/tools/test_edesigntools.py b/test/unit/tools/test_edesigntools.py index 0a888d87b..29f038c56 100644 --- a/test/unit/tools/test_edesigntools.py +++ b/test/unit/tools/test_edesigntools.py @@ -125,8 +125,10 @@ def setUp(self): #E0 = target_model.effects['0'] #E1 = target_model.effects['1'] # Alternate indexing that uses POVM label explicitly - E0 = self.target_model['Mdefault']['0'] # 'Mdefault' = POVM label, '0' = effect label - E1 = self.target_model['Mdefault']['1'] + E0 = self.target_model['Mdefault']['0'].to_dense() + E0 = E0.reshape(len(E0),1) # 'Mdefault' = POVM label, '0' = effect label + E1 = self.target_model['Mdefault']['1'].to_dense() + E1 = E1.reshape(len(E1),1) Gmz_plus = _np.dot(E0,E0.T) #note effect vectors are stored as column vectors Gmz_minus = _np.dot(E1,E1.T) self.target_model_inst[('Iz',0)] = TPInstrument({'p0': Gmz_plus, 'p1': Gmz_minus}) diff --git a/test/unit/tools/test_jamiolkowski.py b/test/unit/tools/test_jamiolkowski.py index 5eadbef28..48d20742f 100644 --- a/test/unit/tools/test_jamiolkowski.py +++ b/test/unit/tools/test_jamiolkowski.py @@ -31,13 +31,14 @@ def setUp(self): #Build a test gate -- old # X(pi,Qhappy)*LX(pi,0,2) self.testGate = create_operation("LX(pi,0,2)", self.sslbls, self.stdSmall) - self.testGateGM_mx = bt.change_basis(self.testGate, self.stdSmall, self.gmSmall) - self.expTestGate_mx = bt.flexible_change_basis(self.testGate, self.stdSmall, self.std) + self.testGatemx = self.testGate.to_dense() + self.testGateGM_mx = bt.change_basis(self.testGatemx, self.stdSmall, self.gmSmall) + self.expTestGate_mx = bt.flexible_change_basis(self.testGatemx, self.stdSmall, self.std) self.expTestGateGM_mx = bt.change_basis(self.expTestGate_mx, self.std, self.gm) def checkBasis(self, cmb): #Op with Jamio map on gate in std and gm bases - Jmx1 = j.jamiolkowski_iso(self.testGate, op_mx_basis=self.stdSmall, + Jmx1 = j.jamiolkowski_iso(self.testGatemx, op_mx_basis=self.stdSmall, choi_mx_basis=cmb) Jmx2 = j.jamiolkowski_iso(self.testGateGM_mx, op_mx_basis=self.gmSmall, choi_mx_basis=cmb) @@ -64,7 +65,7 @@ def checkBasis(self, cmb): #Reverse transform without specifying stateSpaceDims, then contraction, should yield same result revExpTestGate_mx = j.jamiolkowski_iso_inv(Jmx1, choi_mx_basis=cmb, op_mx_basis=self.std) self.assertArraysAlmostEqual(bt.resize_std_mx(revExpTestGate_mx, 'contract', self.std, self.stdSmall), - self.testGate) + self.testGatemx) def test_std_basis(self): #mx_dim = sum([ int(np.sqrt(d)) for d in ]) diff --git a/test/unit/tools/test_optools.py b/test/unit/tools/test_optools.py index a93c9bc5d..4ba71d241 100644 --- a/test/unit/tools/test_optools.py +++ b/test/unit/tools/test_optools.py @@ -131,7 +131,7 @@ def test_log_diff_model_projection(self): proj2, _ = ot.project_model(pm1, self.target_model, [ptype], gen_type, logG_weight=0) pm2 = proj2[0] for pm1_op, pm2_op in zip(pm1.operations.values(), pm2.operations.values()): - self.assertArraysAlmostEqual(pm1_op, pm2_op) + self.assertArraysAlmostEqual(pm1_op.to_dense(), pm2_op.to_dense()) def test_logTiG_model_projection(self): gen_type = 'logTiG' @@ -141,7 +141,7 @@ def test_logTiG_model_projection(self): proj2, _ = ot.project_model(pm1, self.target_model, [ptype], gen_type, logG_weight=0) pm2 = proj2[0] for pm1_op, pm2_op in zip(pm1.operations.values(), pm2.operations.values()): - self.assertArraysAlmostEqual(pm1_op, pm2_op) + self.assertArraysAlmostEqual(pm1_op.to_dense(), pm2_op.to_dense()) def test_logGTi_model_projection(self): gen_type = 'logGTi' @@ -151,7 +151,7 @@ def test_logGTi_model_projection(self): proj2, _ = ot.project_model(pm1, self.target_model, [ptype], gen_type, logG_weight=0) pm2 = proj2[0] for pm1_op, pm2_op in zip(pm1.operations.values(), pm2.operations.values()): - self.assertArraysAlmostEqual(pm1_op, pm2_op) + self.assertArraysAlmostEqual(pm1_op.to_dense(), pm2_op.to_dense()) def test_raises_on_basis_mismatch(self): with self.assertRaises(ValueError): @@ -171,14 +171,11 @@ def test_std_errgens(self): # Note: bases must have first element == identity for projectionType in projectionTypes: - #REMOVE ot.std_scale_factor(4, projectionType) for basisName in basisNames: - #REMOVE ot.std_error_generators(4, projectionType, basisName) ot.elementary_errorgens_dual(4, projectionType, basisName) def test_std_errgens_raise_on_bad_projection_type(self): with self.assertRaises(AssertionError): - #REMOVE ot.std_error_generators(4, "foobar", 'gm') ot.elementary_errorgens_dual(4, "foobar", 'gm') def test_lind_errgens(self): @@ -291,17 +288,15 @@ def test_err_gen(self): basisNames = ['std', 'gm', 'pp'] # , 'qt'] #dim must == 3 for qt for (lbl, gateTarget), gate in zip(self.target_model.operations.items(), self.mdl_datagen.operations.values()): + gate = gate.to_dense() + gateTarget = gateTarget.to_dense() + errgen = ot.error_generator(gate, gateTarget, self.target_model.basis, 'logG-logT') altErrgen = ot.error_generator(gate, gateTarget, self.target_model.basis, 'logTiG') altErrgen2 = ot.error_generator(gate, gateTarget, self.target_model.basis, 'logGTi') with self.assertRaises(ValueError): ot.error_generator(gate, gateTarget, self.target_model.basis, 'adsf') - #OLD: tested above - #for projectionType in projectionTypes: - # for basisName in basisNames: - # ot.std_errorgen_projections(errgen, projectionType, basisName) - originalGate = ot.operation_from_error_generator(errgen, gateTarget, self.target_model.basis, 'logG-logT') altOriginalGate = ot.operation_from_error_generator(altErrgen, gateTarget, self.target_model.basis, 'logTiG') altOriginalGate2 = ot.operation_from_error_generator(altErrgen, gateTarget, self.target_model.basis, 'logGTi') @@ -313,25 +308,25 @@ def test_err_gen(self): @fake_minimize def test_err_gen_nonunitary(self): - errgen_nonunitary = ot.error_generator(self.mdl_datagen.operations['Gxi'], - self.mdl_datagen.operations['Gxi'], + errgen_nonunitary = ot.error_generator(self.mdl_datagen.operations['Gxi'].to_dense(), + self.mdl_datagen.operations['Gxi'].to_dense(), self.mdl_datagen.basis) # Perfect match, should get all 0s - self.assertArraysAlmostEqual(np.zeros_like(self.mdl_datagen.operations['Gxi']), errgen_nonunitary) + self.assertArraysAlmostEqual(np.zeros_like(self.mdl_datagen.operations['Gxi'].to_dense()), errgen_nonunitary) def test_err_gen_not_near_gate(self): # Both should warn with self.assertWarns(UserWarning): - errgen_notsmall = ot.error_generator(self.mdl_datagen.operations['Gxi'], self.target_model.operations['Gix'], + errgen_notsmall = ot.error_generator(self.mdl_datagen.operations['Gxi'].to_dense(), self.target_model.operations['Gix'].to_dense(), self.target_model.basis, 'logTiG') with self.assertWarns(UserWarning): - errgen_notsmall = ot.error_generator(self.mdl_datagen.operations['Gxi'], self.target_model.operations['Gix'], + errgen_notsmall = ot.error_generator(self.mdl_datagen.operations['Gxi'].to_dense(), self.target_model.operations['Gix'].to_dense(), self.target_model.basis, 'logGTi') def test_err_gen_raises_on_bad_type(self): with self.assertRaises(ValueError): - ot.error_generator(self.mdl_datagen.operations['Gxi'], self.target_model.operations['Gxi'], + ot.error_generator(self.mdl_datagen.operations['Gxi'].to_dense(), self.target_model.operations['Gxi'].to_dense(), self.target_model.basis, 'foobar') def test_err_gen_assert_shape_raises_on_ndims_too_high(self):