1# This file is part of QuTiP: Quantum Toolbox in Python. 2# 3# Copyright (c) 2011 and later, Paul D. Nation and Robert J. Johansson. 4# All rights reserved. 5# 6# Redistribution and use in source and binary forms, with or without 7# modification, are permitted provided that the following conditions are 8# met: 9# 10# 1. Redistributions of source code must retain the above copyright notice, 11# this list of conditions and the following disclaimer. 12# 13# 2. Redistributions in binary form must reproduce the above copyright 14# notice, this list of conditions and the following disclaimer in the 15# documentation and/or other materials provided with the distribution. 16# 17# 3. Neither the name of the QuTiP: Quantum Toolbox in Python nor the names 18# of its contributors may be used to endorse or promote products derived 19# from this software without specific prior written permission. 20# 21# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 22# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 23# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A 24# PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT 25# HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 26# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT 27# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, 28# DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 29# THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 30# (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 31# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 32############################################################################### 33"""The Quantum Object (Qobj) class, for representing quantum states and 34operators, and related functions. 35""" 36 37__all__ = ['Qobj', 'qobj_list_evaluate', 'ptrace', 'dag', 'isequal', 38 'issuper', 'isoper', 'isoperket', 'isoperbra', 'isket', 'isbra', 39 'isherm', 'shape', 'dims'] 40 41import warnings 42import types 43import numbers 44 45try: 46 import builtins 47except ImportError: 48 import __builtin__ as builtins 49 50# import math functions from numpy.math: required for td string evaluation 51from numpy import (arccos, arccosh, arcsin, arcsinh, arctan, arctan2, arctanh, 52 ceil, copysign, cos, cosh, degrees, e, exp, expm1, fabs, 53 floor, fmod, frexp, hypot, isinf, isnan, ldexp, log, log10, 54 log1p, modf, pi, radians, sin, sinh, sqrt, tan, tanh) 55 56import numpy as np 57import scipy.sparse as sp 58import qutip.settings as settings 59from qutip import __version__ 60from qutip.fastsparse import fast_csr_matrix, fast_identity 61from qutip.cy.ptrace import _ptrace 62from qutip.permute import _permute 63from qutip.sparse import ( 64 sp_eigs, sp_expm, sp_fro_norm, sp_max_norm, sp_one_norm, sp_L2_norm, 65) 66from qutip.dimensions import ( 67 type_from_dims, enumerate_flat, collapse_dims_super, 68) 69from qutip.cy.spmath import ( 70 zcsr_transpose, zcsr_adjoint, zcsr_isherm, zcsr_trace, zcsr_proj, 71 zcsr_inner, 72) 73from qutip.cy.spmatfuncs import zcsr_mat_elem 74from qutip.cy.sparse_utils import cy_tidyup 75import sys 76if sys.version_info.major >= 3: 77 from itertools import zip_longest 78elif sys.version_info.major < 3: 79 from itertools import izip_longest 80 zip_longest = izip_longest 81 82# OPENMP stuff 83from qutip.cy.openmp.utilities import use_openmp 84if settings.has_openmp: 85 from qutip.cy.openmp.omp_sparse_utils import omp_tidyup 86 87 88class Qobj(object): 89 """A class for representing quantum objects, such as quantum operators 90 and states. 91 92 The Qobj class is the QuTiP representation of quantum operators and state 93 vectors. This class also implements math operations +,-,* between Qobj 94 instances (and / by a C-number), as well as a collection of common 95 operator/state operations. The Qobj constructor optionally takes a 96 dimension ``list`` and/or shape ``list`` as arguments. 97 98 Parameters 99 ---------- 100 inpt : array_like 101 Data for vector/matrix representation of the quantum object. 102 dims : list 103 Dimensions of object used for tensor products. 104 shape : list 105 Shape of underlying data structure (matrix shape). 106 copy : bool 107 Flag specifying whether Qobj should get a copy of the 108 input data, or use the original. 109 fast : bool 110 Flag for fast qobj creation when running ode solvers. 111 This parameter is used internally only. 112 113 114 Attributes 115 ---------- 116 data : array_like 117 Sparse matrix characterizing the quantum object. 118 dims : list 119 List of dimensions keeping track of the tensor structure. 120 shape : list 121 Shape of the underlying `data` array. 122 type : str 123 Type of quantum object: 'bra', 'ket', 'oper', 'operator-ket', 124 'operator-bra', or 'super'. 125 superrep : str 126 Representation used if `type` is 'super'. One of 'super' 127 (Liouville form) or 'choi' (Choi matrix with tr = dimension). 128 isherm : bool 129 Indicates if quantum object represents Hermitian operator. 130 isunitary : bool 131 Indictaes if quantum object represents unitary operator. 132 iscp : bool 133 Indicates if the quantum object represents a map, and if that map is 134 completely positive (CP). 135 ishp : bool 136 Indicates if the quantum object represents a map, and if that map is 137 hermicity preserving (HP). 138 istp : bool 139 Indicates if the quantum object represents a map, and if that map is 140 trace preserving (TP). 141 iscptp : bool 142 Indicates if the quantum object represents a map that is completely 143 positive and trace preserving (CPTP). 144 isket : bool 145 Indicates if the quantum object represents a ket. 146 isbra : bool 147 Indicates if the quantum object represents a bra. 148 isoper : bool 149 Indicates if the quantum object represents an operator. 150 issuper : bool 151 Indicates if the quantum object represents a superoperator. 152 isoperket : bool 153 Indicates if the quantum object represents an operator in column vector 154 form. 155 isoperbra : bool 156 Indicates if the quantum object represents an operator in row vector 157 form. 158 159 Methods 160 ------- 161 copy() 162 Create copy of Qobj 163 conj() 164 Conjugate of quantum object. 165 cosm() 166 Cosine of quantum object. 167 dag() 168 Adjoint (dagger) of quantum object. 169 dnorm() 170 Diamond norm of quantum operator. 171 dual_chan() 172 Dual channel of quantum object representing a CP map. 173 eigenenergies(sparse=False, sort='low', eigvals=0, tol=0, maxiter=100000) 174 Returns eigenenergies (eigenvalues) of a quantum object. 175 eigenstates(sparse=False, sort='low', eigvals=0, tol=0, maxiter=100000) 176 Returns eigenenergies and eigenstates of quantum object. 177 expm() 178 Matrix exponential of quantum object. 179 full(order='C') 180 Returns dense array of quantum object `data` attribute. 181 groundstate(sparse=False, tol=0, maxiter=100000) 182 Returns eigenvalue and eigenket for the groundstate of a quantum 183 object. 184 inv() 185 Return a Qobj corresponding to the matrix inverse of the operator. 186 matrix_element(bra, ket) 187 Returns the matrix element of operator between `bra` and `ket` vectors. 188 norm(norm='tr', sparse=False, tol=0, maxiter=100000) 189 Returns norm of a ket or an operator. 190 permute(order) 191 Returns composite qobj with indices reordered. 192 proj() 193 Computes the projector for a ket or bra vector. 194 ptrace(sel) 195 Returns quantum object for selected dimensions after performing 196 partial trace. 197 sinm() 198 Sine of quantum object. 199 sqrtm() 200 Matrix square root of quantum object. 201 tidyup(atol=1e-12) 202 Removes small elements from quantum object. 203 tr() 204 Trace of quantum object. 205 trans() 206 Transpose of quantum object. 207 transform(inpt, inverse=False) 208 Performs a basis transformation defined by `inpt` matrix. 209 trunc_neg(method='clip') 210 Removes negative eigenvalues and returns a new Qobj that is 211 a valid density operator. 212 unit(norm='tr', sparse=False, tol=0, maxiter=100000) 213 Returns normalized quantum object. 214 215 """ 216 __array_priority__ = 100 # sets Qobj priority above numpy arrays 217 # Disable ufuncs from acting directly on Qobj. This is necessary because we 218 # define __array__. 219 __array_ufunc__ = None 220 221 def __init__(self, inpt=None, dims=None, shape=None, 222 type=None, isherm=None, copy=True, 223 fast=False, superrep=None, isunitary=None): 224 """ 225 Qobj constructor. 226 """ 227 self._isherm = isherm 228 self._type = type 229 self.superrep = superrep 230 self._isunitary = isunitary 231 232 if fast == 'mc': 233 # fast Qobj construction for use in mcsolve with ket output 234 self._data = inpt 235 self.dims = dims 236 self._isherm = False 237 return 238 239 if fast == 'mc-dm': 240 # fast Qobj construction for use in mcsolve with dm output 241 self._data = inpt 242 self.dims = dims 243 self._isherm = True 244 return 245 246 if isinstance(inpt, Qobj): 247 # if input is already Qobj then return identical copy 248 249 self._data = fast_csr_matrix( 250 (inpt.data.data, inpt.data.indices, inpt.data.indptr), 251 shape=inpt.shape, copy=copy, 252 ) 253 254 if dims is None: 255 # Dimensions of quantum object used for keeping track of tensor 256 # components 257 self.dims = inpt.dims 258 else: 259 self.dims = dims 260 261 self.superrep = inpt.superrep 262 self._isunitary = inpt._isunitary 263 264 elif inpt is None: 265 # initialize an empty Qobj with correct dimensions and shape 266 267 if dims is not None: 268 N, M = np.prod(dims[0]), np.prod(dims[1]) 269 self.dims = dims 270 271 elif shape is not None: 272 N, M = shape 273 self.dims = [[N], [M]] 274 275 else: 276 N, M = 1, 1 277 self.dims = [[N], [M]] 278 279 self._data = fast_csr_matrix(shape=(N, M)) 280 281 elif isinstance(inpt, list) or isinstance(inpt, tuple): 282 # case where input is a list 283 data = np.array(inpt) 284 if len(data.shape) == 1: 285 # if list has only one dimension (i.e [5,4]) 286 data = data.transpose() 287 288 _tmp = sp.csr_matrix(data, dtype=complex) 289 self._data = fast_csr_matrix( 290 (_tmp.data, _tmp.indices, _tmp.indptr), 291 shape=_tmp.shape, 292 ) 293 if dims is None: 294 self.dims = [[int(data.shape[0])], [int(data.shape[1])]] 295 else: 296 self.dims = dims 297 298 elif isinstance(inpt, np.ndarray) or sp.issparse(inpt): 299 # case where input is array or sparse 300 if inpt.ndim == 1: 301 inpt = inpt[:, np.newaxis] 302 303 do_copy = copy 304 if not isinstance(inpt, fast_csr_matrix): 305 _tmp = sp.csr_matrix(inpt, dtype=complex, copy=do_copy) 306 _tmp.sort_indices() # Make sure indices are sorted. 307 do_copy = 0 308 else: 309 _tmp = inpt 310 self._data = fast_csr_matrix( 311 (_tmp.data, _tmp.indices, _tmp.indptr), 312 shape=_tmp.shape, 313 copy=do_copy, 314 ) 315 316 if dims is None: 317 self.dims = [[int(inpt.shape[0])], [int(inpt.shape[1])]] 318 else: 319 self.dims = dims 320 321 elif isinstance(inpt, (int, float, complex, 322 np.integer, np.floating, np.complexfloating)): 323 # if input is int, float, or complex then convert to array 324 _tmp = sp.csr_matrix([[inpt]], dtype=complex) 325 self._data = fast_csr_matrix( 326 (_tmp.data, _tmp.indices, _tmp.indptr), 327 shape=_tmp.shape, 328 ) 329 if dims is None: 330 self.dims = [[1], [1]] 331 else: 332 self.dims = dims 333 334 else: 335 warnings.warn("Initializing Qobj from unsupported type: %s" % 336 builtins.type(inpt)) 337 inpt = np.array([[0]]) 338 _tmp = sp.csr_matrix(inpt, dtype=complex, copy=copy) 339 self._data = fast_csr_matrix( 340 (_tmp.data, _tmp.indices, _tmp.indptr), 341 shape=_tmp.shape, 342 ) 343 self.dims = [[int(inpt.shape[0])], [int(inpt.shape[1])]] 344 345 if type == 'super': 346 # Type is not super, i.e. dims not explicitly passed, but oper-like 347 # shape. 348 if dims is None and self.shape[0] == self.shape[1]: 349 sub_shape = np.sqrt(self.shape[0]) 350 # check if root of shape is int 351 if (sub_shape % 1) != 0: 352 raise Exception('Invalid shape for a super operator.') 353 else: 354 sub_shape = int(sub_shape) 355 self.dims = [[[sub_shape], [sub_shape]]]*2 356 357 if superrep: 358 self.superrep = superrep 359 else: 360 if self.type == 'super' and self.superrep is None: 361 self.superrep = 'super' 362 363 # clear type cache 364 self._type = None 365 366 def copy(self): 367 """Create identical copy""" 368 return Qobj(inpt=self) 369 370 def get_data(self): 371 return self._data 372 373 # Here we perfrom a check of the csr matrix type during setting of Q.data 374 def set_data(self, data): 375 if not isinstance(data, fast_csr_matrix): 376 raise TypeError('Qobj data must be in fast_csr format.') 377 else: 378 self._data = data 379 data = property(get_data, set_data) 380 381 def __add__(self, other): 382 """ 383 ADDITION with Qobj on LEFT [ ex. Qobj+4 ] 384 """ 385 self._isunitary = None 386 387 if isinstance(other, eseries): 388 return other.__radd__(self) 389 390 if not isinstance(other, Qobj): 391 if isinstance(other, (int, float, complex, np.integer, np.floating, 392 np.complexfloating, np.ndarray, list, tuple)) \ 393 or sp.issparse(other): 394 other = Qobj(other) 395 else: 396 return NotImplemented 397 398 if np.prod(other.shape) == 1 and np.prod(self.shape) != 1: 399 # case for scalar quantum object 400 dat = other.data[0, 0] 401 if dat == 0: 402 return self 403 404 out = Qobj() 405 406 if self.type in ['oper', 'super']: 407 out.data = self.data + dat * fast_identity( 408 self.shape[0]) 409 else: 410 out.data = self.data 411 out.data.data = out.data.data + dat 412 413 out.dims = self.dims 414 415 if settings.auto_tidyup: 416 out.tidyup() 417 418 if isinstance(dat, (int, float)): 419 out._isherm = self._isherm 420 else: 421 # We use _isherm here to prevent recalculating on self and 422 # other, relying on that bool(None) == False. 423 out._isherm = (self._isherm and other._isherm) or out.isherm 424 425 out.superrep = self.superrep 426 427 return out 428 429 elif np.prod(self.shape) == 1 and np.prod(other.shape) != 1: 430 # case for scalar quantum object 431 dat = self.data[0, 0] 432 if dat == 0: 433 return other 434 435 out = Qobj() 436 if other.type in ['oper', 'super']: 437 out.data = dat * fast_identity(other.shape[0]) + other.data 438 else: 439 out.data = other.data 440 out.data.data = out.data.data + dat 441 out.dims = other.dims 442 443 if settings.auto_tidyup: 444 out.tidyup() 445 446 if isinstance(dat, complex): 447 out._isherm = out.isherm 448 else: 449 out._isherm = self._isherm 450 451 out.superrep = self.superrep 452 453 return out 454 455 elif self.dims != other.dims: 456 raise TypeError('Incompatible quantum object dimensions') 457 458 elif self.shape != other.shape: 459 raise TypeError('Matrix shapes do not match') 460 461 else: # case for matching quantum objects 462 out = Qobj() 463 out.data = self.data + other.data 464 out.dims = self.dims 465 if settings.auto_tidyup: 466 out.tidyup() 467 468 if self.type in ['ket', 'bra', 'operator-ket', 'operator-bra']: 469 out._isherm = False 470 elif self._isherm is None or other._isherm is None: 471 out._isherm = out.isherm 472 elif not self._isherm and not other._isherm: 473 out._isherm = out.isherm 474 else: 475 out._isherm = self._isherm and other._isherm 476 477 if self.superrep and other.superrep: 478 if self.superrep != other.superrep: 479 msg = ("Adding superoperators with different " + 480 "representations") 481 warnings.warn(msg) 482 483 out.superrep = self.superrep 484 485 return out 486 487 def __radd__(self, other): 488 """ 489 ADDITION with Qobj on RIGHT [ ex. 4+Qobj ] 490 """ 491 return self + other 492 493 def __sub__(self, other): 494 """ 495 SUBTRACTION with Qobj on LEFT [ ex. Qobj-4 ] 496 """ 497 return self + (-other) 498 499 def __rsub__(self, other): 500 """ 501 SUBTRACTION with Qobj on RIGHT [ ex. 4-Qobj ] 502 """ 503 return (-self) + other 504 505 def __mul__(self, other): 506 """ 507 MULTIPLICATION with Qobj on LEFT [ ex. Qobj*4 ] 508 """ 509 self._isunitary = None 510 511 if isinstance(other, Qobj): 512 if self.dims[1] == other.dims[0]: 513 out = Qobj() 514 out.data = self.data * other.data 515 dims = [self.dims[0], other.dims[1]] 516 out.dims = dims 517 if settings.auto_tidyup: 518 out.tidyup() 519 if (settings.auto_tidyup_dims 520 and not isinstance(dims[0][0], list) 521 and not isinstance(dims[1][0], list)): 522 # If neither left or right is a superoperator, 523 # we should implicitly partial trace over 524 # matching dimensions of 1. 525 # Using izip_longest allows for the left and right dims 526 # to have uneven length (non-square Qobjs). 527 # We use None as padding so that it doesn't match anything, 528 # and will never cause a partial trace on the other side. 529 mask = [l == r == 1 530 for l, r in zip_longest(dims[0], dims[1], 531 fillvalue=None)] 532 # To ensure that there are still any dimensions left, we 533 # use max() to add a dimensions list of [1] if all matching 534 # dims are traced out of that side. 535 out.dims = [max([1], 536 [dim for dim, m in zip(dims[0], mask) 537 if not m]), 538 max([1], 539 [dim for dim, m in zip(dims[1], mask) 540 if not m])] 541 542 else: 543 out.dims = dims 544 545 out._isherm = None 546 547 if self.superrep and other.superrep: 548 if self.superrep != other.superrep: 549 msg = ("Multiplying superoperators with different " + 550 "representations") 551 warnings.warn(msg) 552 553 out.superrep = self.superrep 554 555 return out 556 557 elif np.prod(self.shape) == 1: 558 out = Qobj(other) 559 out.data *= self.data[0, 0] 560 out.superrep = other.superrep 561 return out.tidyup() if settings.auto_tidyup else out 562 563 elif np.prod(other.shape) == 1: 564 out = Qobj(self) 565 out.data *= other.data[0, 0] 566 out.superrep = self.superrep 567 return out.tidyup() if settings.auto_tidyup else out 568 569 elif self.shape[1] == 1 and other.shape[0] == 1: 570 out = Qobj() 571 out.data = self.data * other.data 572 out.dims = [self.dims[0], other.dims[1]] 573 return out.tidyup() if settings.auto_tidyup else out 574 575 else: 576 raise TypeError("Incompatible Qobj shapes") 577 578 elif isinstance(other, np.ndarray): 579 if other.dtype == 'object': 580 out = np.empty(other.shape, dtype=object) 581 for i, item in enumerate(other.flat): 582 out.flat[i] = self * item 583 return out 584 else: 585 return self.data * other 586 587 elif isinstance(other, list): 588 # if other is a list, do element-wise multiplication 589 out = np.empty(len(other), dtype=object) 590 out[:] = [self * item for item in other] 591 return out 592 593 elif isinstance(other, eseries): 594 return other.__rmul__(self) 595 596 elif isinstance(other, numbers.Number): 597 out = Qobj() 598 out.data = self.data * other 599 out.dims = self.dims 600 out.superrep = self.superrep 601 if settings.auto_tidyup: 602 out.tidyup() 603 if isinstance(other, complex): 604 out._isherm = out.isherm 605 else: 606 out._isherm = self._isherm 607 608 return out 609 610 else: 611 return NotImplemented 612 613 def __rmul__(self, other): 614 """ 615 MULTIPLICATION with Qobj on RIGHT [ ex. 4*Qobj ] 616 """ 617 if isinstance(other, np.ndarray): 618 if other.dtype == 'object': 619 out = np.empty(other.shape, dtype=object) 620 for i, item in enumerate(other.flat): 621 out.flat[i] = item * self 622 return out 623 else: 624 return other * self.data 625 626 elif isinstance(other, list): 627 # if other is a list, do element-wise multiplication 628 out = np.empty(len(other), dtype=object) 629 out[:] = [item * self for item in other] 630 return out 631 632 elif isinstance(other, eseries): 633 return other.__mul__(self) 634 635 elif isinstance(other, numbers.Number): 636 out = Qobj() 637 out.data = other * self.data 638 out.dims = self.dims 639 out.superrep = self.superrep 640 if settings.auto_tidyup: 641 out.tidyup() 642 if isinstance(other, complex): 643 out._isherm = out.isherm 644 else: 645 out._isherm = self._isherm 646 647 return out 648 649 else: 650 raise TypeError("Incompatible object for multiplication") 651 652 def __truediv__(self, other): 653 return self.__div__(other) 654 655 def __div__(self, other): 656 """ 657 DIVISION (by numbers only) 658 """ 659 if isinstance(other, Qobj): # if both are quantum objects 660 raise TypeError("Incompatible Qobj shapes " + 661 "[division with Qobj not implemented]") 662 663 if isinstance(other, (int, float, complex, 664 np.integer, np.floating, np.complexfloating)): 665 out = Qobj() 666 out.data = self.data / other 667 out.dims = self.dims 668 if settings.auto_tidyup: 669 out.tidyup() 670 if isinstance(other, complex): 671 out._isherm = out.isherm 672 else: 673 out._isherm = self._isherm 674 675 out.superrep = self.superrep 676 677 return out 678 679 else: 680 raise TypeError("Incompatible object for division") 681 682 def __neg__(self): 683 """ 684 NEGATION operation. 685 """ 686 out = Qobj() 687 out.data = -self.data 688 out.dims = self.dims 689 out.superrep = self.superrep 690 if settings.auto_tidyup: 691 out.tidyup() 692 out._isherm = self._isherm 693 out._isunitary = self._isunitary 694 return out 695 696 def __getitem__(self, ind): 697 """ 698 GET qobj elements. 699 """ 700 out = self.data[ind] 701 if sp.issparse(out): 702 return out.toarray() 703 else: 704 return out 705 706 def __eq__(self, other): 707 """ 708 EQUALITY operator. 709 """ 710 if (isinstance(other, Qobj) and 711 self.dims == other.dims and 712 not np.any(np.abs((self.data - other.data).data) > 713 settings.atol)): 714 return True 715 else: 716 return False 717 718 def __ne__(self, other): 719 """ 720 INEQUALITY operator. 721 """ 722 return not (self == other) 723 724 def __pow__(self, n, m=None): # calculates powers of Qobj 725 """ 726 POWER operation. 727 """ 728 if self.type not in ['oper', 'super']: 729 raise Exception("Raising a qobj to some power works only for " + 730 "operators and super-operators (square matrices).") 731 732 if m is not None: 733 raise NotImplementedError("modulo is not implemented for Qobj") 734 735 try: 736 data = self.data ** n 737 except (TypeError, ValueError): 738 raise ValueError('Invalid choice of exponent.') 739 out = Qobj(data, dims=self.dims) 740 out.superrep = self.superrep 741 return out.tidyup() if settings.auto_tidyup else out 742 743 def __abs__(self): 744 return abs(self.data) 745 746 def __str__(self): 747 s = "" 748 t = self.type 749 shape = self.shape 750 if self.type in ['oper', 'super']: 751 s += ("Quantum object: " + 752 "dims = " + str(self.dims) + 753 ", shape = " + str(shape) + 754 ", type = " + t + 755 ", isherm = " + str(self.isherm) + 756 ( 757 ", superrep = {0.superrep}".format(self) 758 if t == "super" and self.superrep != "super" 759 else "" 760 ) + "\n") 761 else: 762 s += ("Quantum object: " + 763 "dims = " + str(self.dims) + 764 ", shape = " + str(shape) + 765 ", type = " + t + "\n") 766 s += "Qobj data =\n" 767 768 if shape[0] > 10000 or shape[1] > 10000: 769 # if the system is huge, don't attempt to convert to a 770 # dense matrix and then to string, because it is pointless 771 # and is likely going to produce memory errors. Instead print the 772 # sparse data string representation 773 s += str(self.data) 774 775 elif all(np.imag(self.data.data) == 0): 776 s += str(np.real(self.full())) 777 778 else: 779 s += str(self.full()) 780 781 return s 782 783 def __repr__(self): 784 # give complete information on Qobj without print statement in 785 # command-line we cant realistically serialize a Qobj into a string, 786 # so we simply return the informal __str__ representation instead.) 787 return self.__str__() 788 789 def __call__(self, other): 790 """ 791 Acts this Qobj on another Qobj either by left-multiplication, 792 or by vectorization and devectorization, as 793 appropriate. 794 """ 795 if not isinstance(other, Qobj): 796 raise TypeError("Only defined for quantum objects.") 797 798 if self.type == "super": 799 if other.type == "ket": 800 other = qutip.states.ket2dm(other) 801 802 if other.type == "oper": 803 return qutip.superoperator.vector_to_operator( 804 self * qutip.superoperator.operator_to_vector(other) 805 ) 806 else: 807 raise TypeError("Can only act super on oper or ket.") 808 809 elif self.type == "oper": 810 if other.type == "ket": 811 return self * other 812 else: 813 raise TypeError("Can only act oper on ket.") 814 815 def __getstate__(self): 816 # defines what happens when Qobj object gets pickled 817 self.__dict__.update({'qutip_version': __version__[:5]}) 818 return self.__dict__ 819 820 def __setstate__(self, state): 821 # defines what happens when loading a pickled Qobj 822 if 'qutip_version' in state.keys(): 823 del state['qutip_version'] 824 (self.__dict__).update(state) 825 826 def _repr_latex_(self): 827 """ 828 Generate a LaTeX representation of the Qobj instance. Can be used for 829 formatted output in ipython notebook. 830 """ 831 t = self.type 832 shape = self.shape 833 s = r'' 834 if self.type in ['oper', 'super']: 835 s += ("Quantum object: " + 836 "dims = " + str(self.dims) + 837 ", shape = " + str(shape) + 838 ", type = " + t + 839 ", isherm = " + str(self.isherm) + 840 ( 841 ", superrep = {0.superrep}".format(self) 842 if t == "super" and self.superrep != "super" 843 else "" 844 )) 845 else: 846 s += ("Quantum object: " + 847 "dims = " + str(self.dims) + 848 ", shape = " + str(shape) + 849 ", type = " + t) 850 851 M, N = self.data.shape 852 853 s += r'\begin{equation*}\left(\begin{array}{*{11}c}' 854 855 def _format_float(value): 856 if value == 0.0: 857 return "0.0" 858 elif abs(value) > 1000.0 or abs(value) < 0.001: 859 return ("%.3e" % value).replace("e", r"\times10^{") + "}" 860 elif abs(value - int(value)) < 0.001: 861 return "%.1f" % value 862 else: 863 return "%.3f" % value 864 865 def _format_element(m, n, d): 866 s = " & " if n > 0 else "" 867 if type(d) == str: 868 return s + d 869 else: 870 if abs(np.imag(d)) < settings.atol: 871 return s + _format_float(np.real(d)) 872 elif abs(np.real(d)) < settings.atol: 873 return s + _format_float(np.imag(d)) + "j" 874 else: 875 s_re = _format_float(np.real(d)) 876 s_im = _format_float(np.imag(d)) 877 if np.imag(d) > 0.0: 878 return (s + "(" + s_re + "+" + s_im + "j)") 879 else: 880 return (s + "(" + s_re + s_im + "j)") 881 882 if M > 10 and N > 10: 883 # truncated matrix output 884 for m in range(5): 885 for n in range(5): 886 s += _format_element(m, n, self.data[m, n]) 887 s += r' & \cdots' 888 for n in range(N - 5, N): 889 s += _format_element(m, n, self.data[m, n]) 890 s += r'\\' 891 892 for n in range(5): 893 s += _format_element(m, n, r'\vdots') 894 s += r' & \ddots' 895 for n in range(N - 5, N): 896 s += _format_element(m, n, r'\vdots') 897 s += r'\\' 898 899 for m in range(M - 5, M): 900 for n in range(5): 901 s += _format_element(m, n, self.data[m, n]) 902 s += r' & \cdots' 903 for n in range(N - 5, N): 904 s += _format_element(m, n, self.data[m, n]) 905 s += r'\\' 906 907 elif M > 10 and N <= 10: 908 # truncated vertically elongated matrix output 909 for m in range(5): 910 for n in range(N): 911 s += _format_element(m, n, self.data[m, n]) 912 s += r'\\' 913 914 for n in range(N): 915 s += _format_element(m, n, r'\vdots') 916 s += r'\\' 917 918 for m in range(M - 5, M): 919 for n in range(N): 920 s += _format_element(m, n, self.data[m, n]) 921 s += r'\\' 922 923 elif M <= 10 and N > 10: 924 # truncated horizontally elongated matrix output 925 for m in range(M): 926 for n in range(5): 927 s += _format_element(m, n, self.data[m, n]) 928 s += r' & \cdots' 929 for n in range(N - 5, N): 930 s += _format_element(m, n, self.data[m, n]) 931 s += r'\\' 932 933 else: 934 # full output 935 for m in range(M): 936 for n in range(N): 937 s += _format_element(m, n, self.data[m, n]) 938 s += r'\\' 939 940 s += r'\end{array}\right)\end{equation*}' 941 return s 942 943 def dag(self): 944 """Adjoint operator of quantum object. 945 """ 946 out = Qobj() 947 out.data = zcsr_adjoint(self.data) 948 out.dims = [self.dims[1], self.dims[0]] 949 out._isherm = self._isherm 950 out.superrep = self.superrep 951 return out 952 953 def dual_chan(self): 954 """Dual channel of quantum object representing a completely positive 955 map. 956 """ 957 # Uses the technique of Johnston and Kribs (arXiv:1102.0948), which 958 # is only valid for completely positive maps. 959 if not self.iscp: 960 raise ValueError("Dual channels are only implemented for CP maps.") 961 J = sr.to_choi(self) 962 tensor_idxs = enumerate_flat(J.dims) 963 J_dual = tensor.tensor_swap(J, *( 964 list(zip(tensor_idxs[0][1], tensor_idxs[0][0])) + 965 list(zip(tensor_idxs[1][1], tensor_idxs[1][0])) 966 )).trans() 967 J_dual.superrep = 'choi' 968 return J_dual 969 970 def conj(self): 971 """Conjugate operator of quantum object. 972 """ 973 out = Qobj() 974 out.data = self.data.conj() 975 out.dims = [self.dims[0], self.dims[1]] 976 return out 977 978 def norm(self, norm=None, sparse=False, tol=0, maxiter=100000): 979 """Norm of a quantum object. 980 981 Default norm is L2-norm for kets and trace-norm for operators. 982 Other ket and operator norms may be specified using the `norm` and 983 argument. 984 985 Parameters 986 ---------- 987 norm : str 988 Which norm to use for ket/bra vectors: L2 'l2', max norm 'max', 989 or for operators: trace 'tr', Frobius 'fro', one 'one', or max 990 'max'. 991 992 sparse : bool 993 Use sparse eigenvalue solver for trace norm. Other norms are not 994 affected by this parameter. 995 996 tol : float 997 Tolerance for sparse solver (if used) for trace norm. The sparse 998 solver may not converge if the tolerance is set too low. 999 1000 maxiter : int 1001 Maximum number of iterations performed by sparse solver (if used) 1002 for trace norm. 1003 1004 Returns 1005 ------- 1006 norm : float 1007 The requested norm of the operator or state quantum object. 1008 1009 1010 Notes 1011 ----- 1012 The sparse eigensolver is much slower than the dense version. 1013 Use sparse only if memory requirements demand it. 1014 1015 """ 1016 if self.type in ['oper', 'super']: 1017 if norm is None or norm == 'tr': 1018 _op = self.data * zcsr_adjoint(self.data) 1019 vals = sp_eigs(_op, True, vecs=False, 1020 sparse=sparse, tol=tol, maxiter=maxiter) 1021 return np.sum(np.sqrt(np.abs(vals))) 1022 elif norm == 'fro': 1023 return sp_fro_norm(self.data) 1024 elif norm == 'one': 1025 return sp_one_norm(self.data) 1026 elif norm == 'max': 1027 return sp_max_norm(self.data) 1028 else: 1029 raise ValueError( 1030 "For matrices, norm must be 'tr', 'fro', 'one', or 'max'.") 1031 else: 1032 if norm is None or norm == 'l2': 1033 return sp_L2_norm(self.data) 1034 elif norm == 'max': 1035 return sp_max_norm(self.data) 1036 else: 1037 raise ValueError("For vectors, norm must be 'l2', or 'max'.") 1038 1039 def proj(self): 1040 """Form the projector from a given ket or bra vector. 1041 1042 Parameters 1043 ---------- 1044 Q : :class:`qutip.Qobj` 1045 Input bra or ket vector 1046 1047 Returns 1048 ------- 1049 P : :class:`qutip.Qobj` 1050 Projection operator. 1051 """ 1052 if self.isket: 1053 _out = zcsr_proj(self.data, 1) 1054 _dims = [self.dims[0], self.dims[0]] 1055 elif self.isbra: 1056 _out = zcsr_proj(self.data, 0) 1057 _dims = [self.dims[1], self.dims[1]] 1058 else: 1059 raise TypeError('Projector can only be formed from a bra or ket.') 1060 1061 return Qobj(_out, dims=_dims) 1062 1063 def tr(self): 1064 """Trace of a quantum object. 1065 1066 Returns 1067 ------- 1068 trace : float 1069 Returns the trace of the quantum object. 1070 1071 """ 1072 return zcsr_trace(self.data, self.isherm) 1073 1074 def purity(self): 1075 """Calculate purity of a quantum object. 1076 1077 Returns 1078 ------- 1079 state_purity : float 1080 Returns the purity of a quantum object. 1081 For a pure state, the purity is 1. 1082 For a mixed state of dimension `d`, 1/d<=purity<1. 1083 1084 """ 1085 rho = self 1086 if (rho.type == "super"): 1087 raise TypeError('Purity is defined on a density matrix or state.') 1088 1089 if rho.type == "ket" or rho.type == "bra": 1090 state_purity = (rho*rho.dag()).tr() 1091 1092 if rho.type == "oper": 1093 state_purity = (rho*rho).tr() 1094 1095 return state_purity 1096 1097 def full(self, order='C', squeeze=False): 1098 """Dense array from quantum object. 1099 1100 Parameters 1101 ---------- 1102 order : str {'C', 'F'} 1103 Return array in C (default) or Fortran ordering. 1104 squeeze : bool {False, True} 1105 Squeeze output array. 1106 1107 Returns 1108 ------- 1109 data : array 1110 Array of complex data from quantum objects `data` attribute. 1111 """ 1112 if squeeze: 1113 return self.data.toarray(order=order).squeeze() 1114 else: 1115 return self.data.toarray(order=order) 1116 1117 def __array__(self, *arg, **kwarg): 1118 """Numpy array from Qobj 1119 For compatibility with np.array 1120 """ 1121 return self.full() 1122 1123 def diag(self): 1124 """Diagonal elements of quantum object. 1125 1126 Returns 1127 ------- 1128 diags : array 1129 Returns array of ``real`` values if operators is Hermitian, 1130 otherwise ``complex`` values are returned. 1131 1132 """ 1133 out = self.data.diagonal() 1134 if np.any(np.imag(out) > settings.atol) or not self.isherm: 1135 return out 1136 else: 1137 return np.real(out) 1138 1139 def expm(self, method='dense'): 1140 """Matrix exponential of quantum operator. 1141 1142 Input operator must be square. 1143 1144 Parameters 1145 ---------- 1146 method : str {'dense', 'sparse'} 1147 Use set method to use to calculate the matrix exponentiation. The 1148 available choices includes 'dense' and 'sparse'. Since the 1149 exponential of a matrix is nearly always dense, method='dense' 1150 is set as default.s 1151 1152 Returns 1153 ------- 1154 oper : :class:`qutip.Qobj` 1155 Exponentiated quantum operator. 1156 1157 Raises 1158 ------ 1159 TypeError 1160 Quantum operator is not square. 1161 1162 """ 1163 if self.dims[0][0] != self.dims[1][0]: 1164 raise TypeError('Invalid operand for matrix exponential') 1165 1166 if method == 'dense': 1167 F = sp_expm(self.data, sparse=False) 1168 1169 elif method == 'sparse': 1170 F = sp_expm(self.data, sparse=True) 1171 1172 else: 1173 raise ValueError("method must be 'dense' or 'sparse'.") 1174 1175 out = Qobj(F, dims=self.dims) 1176 return out.tidyup() if settings.auto_tidyup else out 1177 1178 def check_herm(self): 1179 """Check if the quantum object is hermitian. 1180 1181 Returns 1182 ------- 1183 isherm : bool 1184 Returns the new value of isherm property. 1185 """ 1186 self._isherm = None 1187 return self.isherm 1188 1189 def sqrtm(self, sparse=False, tol=0, maxiter=100000): 1190 """Sqrt of a quantum operator. 1191 1192 Operator must be square. 1193 1194 Parameters 1195 ---------- 1196 sparse : bool 1197 Use sparse eigenvalue/vector solver. 1198 tol : float 1199 Tolerance used by sparse solver (0 = machine precision). 1200 maxiter : int 1201 Maximum number of iterations used by sparse solver. 1202 1203 Returns 1204 ------- 1205 oper : :class:`qutip.Qobj` 1206 Matrix square root of operator. 1207 1208 Raises 1209 ------ 1210 TypeError 1211 Quantum object is not square. 1212 1213 Notes 1214 ----- 1215 The sparse eigensolver is much slower than the dense version. 1216 Use sparse only if memory requirements demand it. 1217 1218 """ 1219 if self.dims[0][0] == self.dims[1][0]: 1220 evals, evecs = sp_eigs(self.data, self.isherm, sparse=sparse, 1221 tol=tol, maxiter=maxiter) 1222 numevals = len(evals) 1223 dV = sp.spdiags(np.sqrt(evals, dtype=complex), 0, numevals, 1224 numevals, format='csr') 1225 if self.isherm: 1226 spDv = dV.dot(evecs.T.conj().T) 1227 else: 1228 spDv = dV.dot(np.linalg.inv(evecs.T)) 1229 1230 out = Qobj(evecs.T.dot(spDv), dims=self.dims) 1231 return out.tidyup() if settings.auto_tidyup else out 1232 1233 else: 1234 raise TypeError('Invalid operand for matrix square root') 1235 1236 def cosm(self): 1237 """Cosine of a quantum operator. 1238 1239 Operator must be square. 1240 1241 Returns 1242 ------- 1243 oper : :class:`qutip.Qobj` 1244 Matrix cosine of operator. 1245 1246 Raises 1247 ------ 1248 TypeError 1249 Quantum object is not square. 1250 1251 Notes 1252 ----- 1253 Uses the Q.expm() method. 1254 1255 """ 1256 if self.dims[0][0] == self.dims[1][0]: 1257 return 0.5 * ((1j * self).expm() + (-1j * self).expm()) 1258 else: 1259 raise TypeError('Invalid operand for matrix square root') 1260 1261 def sinm(self): 1262 """Sine of a quantum operator. 1263 1264 Operator must be square. 1265 1266 Returns 1267 ------- 1268 oper : :class:`qutip.Qobj` 1269 Matrix sine of operator. 1270 1271 Raises 1272 ------ 1273 TypeError 1274 Quantum object is not square. 1275 1276 Notes 1277 ----- 1278 Uses the Q.expm() method. 1279 1280 """ 1281 if self.dims[0][0] == self.dims[1][0]: 1282 return -0.5j * ((1j * self).expm() - (-1j * self).expm()) 1283 else: 1284 raise TypeError('Invalid operand for matrix square root') 1285 1286 def inv(self, sparse=False): 1287 """Matrix inverse of a quantum operator 1288 1289 Operator must be square. 1290 1291 Returns 1292 ------- 1293 oper : :class:`qutip.Qobj` 1294 Matrix inverse of operator. 1295 1296 Raises 1297 ------ 1298 TypeError 1299 Quantum object is not square. 1300 """ 1301 if self.shape[0] != self.shape[1]: 1302 raise TypeError('Invalid operand for matrix inverse') 1303 if sparse: 1304 inv_mat = sp.linalg.inv(sp.csc_matrix(self.data)) 1305 else: 1306 inv_mat = np.linalg.inv(self.full()) 1307 return Qobj(inv_mat, dims=self.dims[::-1]) 1308 1309 def unit(self, inplace=False, 1310 norm=None, sparse=False, 1311 tol=0, maxiter=100000): 1312 """Operator or state normalized to unity. 1313 1314 Uses norm from Qobj.norm(). 1315 1316 Parameters 1317 ---------- 1318 inplace : bool 1319 Do an in-place normalization 1320 norm : str 1321 Requested norm for states / operators. 1322 sparse : bool 1323 Use sparse eigensolver for trace norm. Does not affect other norms. 1324 tol : float 1325 Tolerance used by sparse eigensolver. 1326 maxiter : int 1327 Number of maximum iterations performed by sparse eigensolver. 1328 1329 Returns 1330 ------- 1331 oper : :class:`qutip.Qobj` 1332 Normalized quantum object if not in-place, 1333 else None. 1334 1335 """ 1336 if inplace: 1337 nrm = self.norm(norm=norm, sparse=sparse, 1338 tol=tol, maxiter=maxiter) 1339 1340 self.data /= nrm 1341 elif not inplace: 1342 out = self / self.norm(norm=norm, sparse=sparse, 1343 tol=tol, maxiter=maxiter) 1344 if settings.auto_tidyup: 1345 return out.tidyup() 1346 else: 1347 return out 1348 else: 1349 raise Exception('inplace kwarg must be bool.') 1350 1351 def ptrace(self, sel, sparse=None): 1352 """Partial trace of the quantum object. 1353 1354 Parameters 1355 ---------- 1356 sel : int/list 1357 An ``int`` or ``list`` of components to keep after partial trace. 1358 The order is unimportant; no transposition will be done and the 1359 spaces will remain in the same order in the output. 1360 1361 Returns 1362 ------- 1363 oper : :class:`qutip.Qobj` 1364 Quantum object representing partial trace with selected components 1365 remaining. 1366 1367 Notes 1368 ----- 1369 This function is identical to the :func:`qutip.qobj.ptrace` function 1370 that has been deprecated. 1371 1372 """ 1373 if sparse is None: 1374 if self.isket: 1375 sparse = False 1376 elif (self.data.nnz / (self.shape[0] * self.shape[1])) >= 0.1: 1377 sparse = False 1378 if sparse: 1379 q = Qobj() 1380 q.data, q.dims, _ = _ptrace(self, sel) 1381 return q.tidyup() if settings.auto_tidyup else q 1382 else: 1383 return _ptrace_dense(self, sel) 1384 1385 def permute(self, order): 1386 """Permutes a composite quantum object. 1387 1388 Parameters 1389 ---------- 1390 order : list/array 1391 List specifying new tensor order. 1392 1393 Returns 1394 ------- 1395 P : :class:`qutip.Qobj` 1396 Permuted quantum object. 1397 1398 """ 1399 q = Qobj() 1400 q.data, q.dims = _permute(self, order) 1401 q.data.sort_indices() 1402 return q 1403 1404 def tidyup(self, atol=settings.auto_tidyup_atol): 1405 """Removes small elements from the quantum object. 1406 1407 Parameters 1408 ---------- 1409 atol : float 1410 Absolute tolerance used by tidyup. Default is set 1411 via qutip global settings parameters. 1412 1413 Returns 1414 ------- 1415 oper : :class:`qutip.Qobj` 1416 Quantum object with small elements removed. 1417 1418 """ 1419 if self.data.nnz: 1420 # This does the tidyup and returns True if 1421 # The sparse data needs to be shortened 1422 if use_openmp() and self.data.nnz > 500: 1423 if omp_tidyup(self.data.data, atol, self.data.nnz, 1424 settings.num_cpus): 1425 self.data.eliminate_zeros() 1426 else: 1427 if cy_tidyup(self.data.data, atol, self.data.nnz): 1428 self.data.eliminate_zeros() 1429 return self 1430 else: 1431 return self 1432 1433 def transform(self, inpt, inverse=False, sparse=True): 1434 """Basis transform defined by input array. 1435 1436 Input array can be a ``matrix`` defining the transformation, 1437 or a ``list`` of kets that defines the new basis. 1438 1439 1440 Parameters 1441 ---------- 1442 inpt : array_like 1443 A ``matrix`` or ``list`` of kets defining the transformation. 1444 inverse : bool 1445 Whether to return inverse transformation. 1446 sparse : bool 1447 Use sparse matrices when possible. Can be slower. 1448 1449 Returns 1450 ------- 1451 oper : :class:`qutip.Qobj` 1452 Operator in new basis. 1453 1454 Notes 1455 ----- 1456 This function is still in development. 1457 1458 1459 """ 1460 if isinstance(inpt, list) or (isinstance(inpt, np.ndarray) and 1461 len(inpt.shape) == 1): 1462 if len(inpt) != max(self.shape): 1463 raise TypeError( 1464 'Invalid size of ket list for basis transformation') 1465 if sparse: 1466 S = sp.hstack([psi.data for psi in inpt], 1467 format='csr', dtype=complex).conj().T 1468 else: 1469 S = np.hstack([psi.full() for psi in inpt], 1470 dtype=complex).conj().T 1471 elif isinstance(inpt, Qobj) and inpt.isoper: 1472 S = inpt.data 1473 elif isinstance(inpt, np.ndarray): 1474 S = inpt.conj() 1475 sparse = False 1476 else: 1477 raise TypeError('Invalid operand for basis transformation') 1478 1479 # transform data 1480 if inverse: 1481 if self.isket: 1482 data = (S.conj().T) * self.data 1483 elif self.isbra: 1484 data = self.data.dot(S) 1485 else: 1486 if sparse: 1487 data = (S.conj().T) * self.data * S 1488 else: 1489 data = (S.conj().T).dot(self.data.dot(S)) 1490 else: 1491 if self.isket: 1492 data = S * self.data 1493 elif self.isbra: 1494 data = self.data.dot(S.conj().T) 1495 else: 1496 if sparse: 1497 data = S * self.data * (S.conj().T) 1498 else: 1499 data = S.dot(self.data.dot(S.conj().T)) 1500 1501 out = Qobj(data, dims=self.dims) 1502 out._isherm = self._isherm 1503 out.superrep = self.superrep 1504 1505 if settings.auto_tidyup: 1506 return out.tidyup() 1507 else: 1508 return out 1509 1510 def trunc_neg(self, method="clip"): 1511 """Truncates negative eigenvalues and renormalizes. 1512 1513 Returns a new Qobj by removing the negative eigenvalues 1514 of this instance, then renormalizing to obtain a valid density 1515 operator. 1516 1517 1518 Parameters 1519 ---------- 1520 method : str 1521 Algorithm to use to remove negative eigenvalues. "clip" 1522 simply discards negative eigenvalues, then renormalizes. 1523 "sgs" uses the SGS algorithm (doi:10/bb76) to find the 1524 positive operator that is nearest in the Shatten 2-norm. 1525 1526 Returns 1527 ------- 1528 oper : :class:`qutip.Qobj` 1529 A valid density operator. 1530 1531 """ 1532 if not self.isherm: 1533 raise ValueError("Must be a Hermitian operator to remove negative " 1534 "eigenvalues.") 1535 1536 if method not in ('clip', 'sgs'): 1537 raise ValueError("Method {} not recognized.".format(method)) 1538 1539 eigvals, eigstates = self.eigenstates() 1540 if all([eigval >= 0 for eigval in eigvals]): 1541 # All positive, so just renormalize. 1542 return self.unit() 1543 idx_nonzero = eigvals != 0 1544 eigvals = eigvals[idx_nonzero] 1545 eigstates = eigstates[idx_nonzero] 1546 1547 if method == 'clip': 1548 eigvals[eigvals < 0] = 0 1549 elif method == 'sgs': 1550 eigvals = eigvals[::-1] 1551 eigstates = eigstates[::-1] 1552 1553 acc = 0.0 1554 dim = self.shape[0] 1555 n_eigs = len(eigvals) 1556 1557 for idx in reversed(range(n_eigs)): 1558 if eigvals[idx] + acc / (idx + 1) >= 0: 1559 break 1560 else: 1561 acc += eigvals[idx] 1562 eigvals[idx] = 0.0 1563 1564 eigvals[:idx+1] += acc / (idx + 1) 1565 1566 return sum([ 1567 val * qutip.states.ket2dm(state) 1568 for val, state in zip(eigvals, eigstates) 1569 ], Qobj(np.zeros(self.shape), dims=self.dims) 1570 ).unit() 1571 1572 def matrix_element(self, bra, ket): 1573 """Calculates a matrix element. 1574 1575 Gives the matrix element for the quantum object sandwiched between a 1576 `bra` and `ket` vector. 1577 1578 Parameters 1579 ----------- 1580 bra : :class:`qutip.Qobj` 1581 Quantum object of type 'bra' or 'ket' 1582 1583 ket : :class:`qutip.Qobj` 1584 Quantum object of type 'ket'. 1585 1586 Returns 1587 ------- 1588 elem : complex 1589 Complex valued matrix element. 1590 1591 Notes 1592 ----- 1593 It is slightly more computationally efficient to use a ket 1594 vector for the 'bra' input. 1595 """ 1596 if not self.isoper: 1597 raise TypeError("Can only get matrix elements for an operator.") 1598 1599 else: 1600 if bra.isbra and ket.isket: 1601 return zcsr_mat_elem(self.data, bra.data, ket.data, 1) 1602 1603 elif bra.isket and ket.isket: 1604 return zcsr_mat_elem(self.data, bra.data, ket.data, 0) 1605 else: 1606 err = "Can only calculate matrix elements for bra" 1607 err += " and ket vectors." 1608 raise TypeError(err) 1609 1610 def overlap(self, other): 1611 """Overlap between two state vectors or two operators. 1612 1613 Gives the overlap (inner product) between the current bra or ket Qobj 1614 and and another bra or ket Qobj. It gives the Hilbert-Schmidt overlap 1615 when one of the Qobj is an operator/density matrix. 1616 1617 Parameters 1618 ----------- 1619 other : :class:`qutip.Qobj` 1620 Quantum object for a state vector of type 'ket', 'bra' or density 1621 matrix. 1622 1623 Returns 1624 ------- 1625 overlap : complex 1626 Complex valued overlap. 1627 1628 Raises 1629 ------ 1630 TypeError 1631 Can only calculate overlap between a bra, ket and density matrix 1632 quantum objects. 1633 1634 Notes 1635 ----- 1636 Since QuTiP mainly deals with ket vectors, the most efficient inner 1637 product call is the ket-ket version that computes the product 1638 <self|other> with both vectors expressed as kets. 1639 """ 1640 1641 if isinstance(other, Qobj): 1642 1643 if self.isbra: 1644 if other.isket: 1645 return zcsr_inner(self.data, other.data, 1) 1646 elif other.isbra: 1647 # Since we deal mainly with ket vectors, the bra-bra combo 1648 # is not common, and not optimized. 1649 return zcsr_inner(self.data, other.dag().data, 1) 1650 elif other.isoper: 1651 return (qutip.states.ket2dm(self).dag() * other).tr() 1652 else: 1653 err = "Can only calculate overlap for state vector Qobjs" 1654 raise TypeError(err) 1655 1656 elif self.isket: 1657 if other.isbra: 1658 return zcsr_inner(other.data, self.data, 1) 1659 elif other.isket: 1660 return zcsr_inner(self.data, other.data, 0) 1661 elif other.isoper: 1662 return (qutip.states.ket2dm(self).dag() * other).tr() 1663 else: 1664 err = "Can only calculate overlap for state vector Qobjs" 1665 raise TypeError(err) 1666 1667 elif self.isoper: 1668 if other.isket or other.isbra: 1669 return (self.dag() * qutip.states.ket2dm(other)).tr() 1670 elif other.isoper: 1671 return (self.dag() * other).tr() 1672 else: 1673 err = "Can only calculate overlap for state vector Qobjs" 1674 raise TypeError(err) 1675 raise TypeError("Can only calculate overlap for state vector Qobjs") 1676 1677 def eigenstates(self, sparse=False, sort='low', eigvals=0, 1678 tol=0, maxiter=100000, phase_fix=None): 1679 """Eigenstates and eigenenergies. 1680 1681 Eigenstates and eigenenergies are defined for operators and 1682 superoperators only. 1683 1684 Parameters 1685 ---------- 1686 sparse : bool 1687 Use sparse Eigensolver 1688 1689 sort : str 1690 Sort eigenvalues (and vectors) 'low' to high, or 'high' to low. 1691 1692 eigvals : int 1693 Number of requested eigenvalues. Default is all eigenvalues. 1694 1695 tol : float 1696 Tolerance used by sparse Eigensolver (0 = machine precision). 1697 The sparse solver may not converge if the tolerance is set too low. 1698 1699 maxiter : int 1700 Maximum number of iterations performed by sparse solver (if used). 1701 1702 phase_fix : int, None 1703 If not None, set the phase of each kets so that ket[phase_fix,0] 1704 is real positive. 1705 1706 Returns 1707 ------- 1708 eigvals : array 1709 Array of eigenvalues for operator. 1710 1711 eigvecs : array 1712 Array of quantum operators representing the oprator eigenkets. 1713 Order of eigenkets is determined by order of eigenvalues. 1714 1715 Notes 1716 ----- 1717 The sparse eigensolver is much slower than the dense version. 1718 Use sparse only if memory requirements demand it. 1719 1720 """ 1721 evals, evecs = sp_eigs(self.data, self.isherm, sparse=sparse, 1722 sort=sort, eigvals=eigvals, tol=tol, 1723 maxiter=maxiter) 1724 if self.type == 'super': 1725 new_dims = [self.dims[0], [1]] 1726 new_type = 'operator-ket' 1727 else: 1728 new_dims = [self.dims[0], [1] * len(self.dims[0])] 1729 new_type = 'ket' 1730 ekets = np.empty((len(evecs),), dtype=object) 1731 ekets[:] = [Qobj(vec, dims=new_dims, type=new_type) for vec in evecs] 1732 norms = np.array([ket.norm() for ket in ekets]) 1733 if phase_fix is None: 1734 phase = np.array([1] * len(ekets)) 1735 else: 1736 phase = np.array([np.abs(ket[phase_fix, 0]) / ket[phase_fix, 0] 1737 if ket[phase_fix, 0] else 1 1738 for ket in ekets]) 1739 return evals, ekets / norms * phase 1740 1741 def eigenenergies(self, sparse=False, sort='low', 1742 eigvals=0, tol=0, maxiter=100000): 1743 """Eigenenergies of a quantum object. 1744 1745 Eigenenergies (eigenvalues) are defined for operators or superoperators 1746 only. 1747 1748 Parameters 1749 ---------- 1750 sparse : bool 1751 Use sparse Eigensolver 1752 sort : str 1753 Sort eigenvalues 'low' to high, or 'high' to low. 1754 eigvals : int 1755 Number of requested eigenvalues. Default is all eigenvalues. 1756 tol : float 1757 Tolerance used by sparse Eigensolver (0=machine precision). 1758 The sparse solver may not converge if the tolerance is set too low. 1759 maxiter : int 1760 Maximum number of iterations performed by sparse solver (if used). 1761 1762 Returns 1763 ------- 1764 eigvals : array 1765 Array of eigenvalues for operator. 1766 1767 Notes 1768 ----- 1769 The sparse eigensolver is much slower than the dense version. 1770 Use sparse only if memory requirements demand it. 1771 1772 """ 1773 return sp_eigs(self.data, self.isherm, vecs=False, sparse=sparse, 1774 sort=sort, eigvals=eigvals, tol=tol, maxiter=maxiter) 1775 1776 def groundstate(self, sparse=False, tol=0, maxiter=100000, safe=True): 1777 """Ground state Eigenvalue and Eigenvector. 1778 1779 Defined for quantum operators or superoperators only. 1780 1781 Parameters 1782 ---------- 1783 sparse : bool 1784 Use sparse Eigensolver 1785 tol : float 1786 Tolerance used by sparse Eigensolver (0 = machine precision). 1787 The sparse solver may not converge if the tolerance is set too low. 1788 maxiter : int 1789 Maximum number of iterations performed by sparse solver (if used). 1790 safe : bool (default=True) 1791 Check for degenerate ground state 1792 1793 Returns 1794 ------- 1795 eigval : float 1796 Eigenvalue for the ground state of quantum operator. 1797 eigvec : :class:`qutip.Qobj` 1798 Eigenket for the ground state of quantum operator. 1799 1800 Notes 1801 ----- 1802 The sparse eigensolver is much slower than the dense version. 1803 Use sparse only if memory requirements demand it. 1804 1805 """ 1806 if safe: 1807 evals = 2 1808 else: 1809 evals = 1 1810 grndval, grndvec = sp_eigs(self.data, self.isherm, sparse=sparse, 1811 eigvals=evals, tol=tol, maxiter=maxiter) 1812 if safe: 1813 tol = 1e-15 if tol == 0 else tol 1814 if (grndval[1]-grndval[0]) <= 10*tol: 1815 print("WARNING: Ground state may be degenerate. " 1816 "Use Q.eigenstates()") 1817 new_dims = [self.dims[0], [1] * len(self.dims[0])] 1818 grndvec = Qobj(grndvec[0], dims=new_dims) 1819 grndvec = grndvec / grndvec.norm() 1820 return grndval[0], grndvec 1821 1822 def trans(self): 1823 """Transposed operator. 1824 1825 Returns 1826 ------- 1827 oper : :class:`qutip.Qobj` 1828 Transpose of input operator. 1829 1830 """ 1831 out = Qobj() 1832 out.data = zcsr_transpose(self.data) 1833 out.dims = [self.dims[1], self.dims[0]] 1834 return out 1835 1836 def extract_states(self, states_inds, normalize=False): 1837 """Qobj with states in state_inds only. 1838 1839 Parameters 1840 ---------- 1841 states_inds : list of integer 1842 The states that should be kept. 1843 1844 normalize : True / False 1845 Weather or not the new Qobj instance should be normalized (default 1846 is False). For Qobjs that represents density matrices or state 1847 vectors normalized should probably be set to True, but for Qobjs 1848 that represents operators in for example an Hamiltonian, normalize 1849 should be False. 1850 1851 Returns 1852 ------- 1853 q : :class:`qutip.Qobj` 1854 A new instance of :class:`qutip.Qobj` that contains only the states 1855 corresponding to the indices in `state_inds`. 1856 1857 Notes 1858 ----- 1859 Experimental. 1860 1861 """ 1862 if self.isoper: 1863 q = Qobj(self.data[states_inds, :][:, states_inds]) 1864 elif self.isket: 1865 q = Qobj(self.data[states_inds, :]) 1866 elif self.isbra: 1867 q = Qobj(self.data[:, states_inds]) 1868 else: 1869 raise TypeError("Can only eliminate states from operators or " + 1870 "state vectors") 1871 1872 return q.unit() if normalize else q 1873 1874 def eliminate_states(self, states_inds, normalize=False): 1875 """Creates a new quantum object with states in state_inds eliminated. 1876 1877 Parameters 1878 ---------- 1879 states_inds : list of integer 1880 The states that should be removed. 1881 1882 normalize : True / False 1883 Weather or not the new Qobj instance should be normalized (default 1884 is False). For Qobjs that represents density matrices or state 1885 vectors normalized should probably be set to True, but for Qobjs 1886 that represents operators in for example an Hamiltonian, normalize 1887 should be False. 1888 1889 Returns 1890 ------- 1891 q : :class:`qutip.Qobj` 1892 A new instance of :class:`qutip.Qobj` that contains only the states 1893 corresponding to indices that are **not** in `state_inds`. 1894 1895 Notes 1896 ----- 1897 Experimental. 1898 1899 """ 1900 keep_indices = np.array([s not in states_inds 1901 for s in range(self.shape[0])]).nonzero()[0] 1902 1903 return self.extract_states(keep_indices, normalize=normalize) 1904 1905 def dnorm(self, B=None): 1906 """Calculates the diamond norm, or the diamond distance to another 1907 operator. 1908 1909 Parameters 1910 ---------- 1911 B : :class:`qutip.Qobj` or None 1912 If B is not None, the diamond distance d(A, B) = dnorm(A - B) 1913 between this operator and B is returned instead of the diamond 1914 norm. 1915 1916 Returns 1917 ------- 1918 d : float 1919 Either the diamond norm of this operator, or the diamond distance 1920 from this operator to B. 1921 1922 """ 1923 return mts.dnorm(self, B) 1924 1925 @property 1926 def ishp(self): 1927 # FIXME: this needs to be cached in the same ways as isherm. 1928 if self.type in ["super", "oper"]: 1929 try: 1930 J = sr.to_choi(self) 1931 return J.isherm 1932 except TypeError: 1933 return False 1934 else: 1935 return False 1936 1937 @property 1938 def iscp(self): 1939 # FIXME: this needs to be cached in the same ways as isherm. 1940 if self.type in ["super", "oper"]: 1941 try: 1942 J = ( 1943 self 1944 # We can test with either Choi or chi, since the basis 1945 # transformation between them is unitary and hence 1946 # preserves the CP and TP conditions. 1947 if self.superrep in ('choi', 'chi') 1948 else sr.to_choi(self) 1949 ) 1950 # If J isn't hermitian, then that could indicate either 1951 # that J is not normal, or is normal, 1952 # but has complex eigenvalues. 1953 # In either case, it makes no sense to then demand that the 1954 # eigenvalues be non-negative. 1955 if not J.isherm: 1956 return False 1957 eigs = J.eigenenergies() 1958 return all(eigs >= -settings.atol) 1959 except TypeError: 1960 return False 1961 else: 1962 return False 1963 1964 @property 1965 def istp(self): 1966 import qutip.superop_reps as sr 1967 if self.type in ["super", "oper"]: 1968 try: 1969 # Normalize to a super of type choi or chi. 1970 # We can test with either Choi or chi, since the basis 1971 # transformation between them is unitary and hence 1972 # preserves the CP and TP conditions. 1973 if self.type == "super" and self.superrep in ('choi', 'chi'): 1974 qobj = self 1975 else: 1976 qobj = sr.to_choi(self) 1977 1978 # Possibly collapse dims. 1979 if any( 1980 len(index) > 1 1981 for super_index in qobj.dims for index in super_index 1982 ): 1983 qobj = Qobj(qobj, dims=collapse_dims_super(qobj.dims)) 1984 else: 1985 qobj = qobj 1986 1987 # We use the condition from John Watrous' lecture notes, 1988 # Tr_1(J(Phi)) = identity_2. 1989 tr_oper = qobj.ptrace([0]) 1990 ident = ops.identity(tr_oper.shape[0]) 1991 return isequal(tr_oper, ident) 1992 except TypeError: 1993 return False 1994 else: 1995 return False 1996 1997 @property 1998 def iscptp(self): 1999 from qutip.superop_reps import to_choi 2000 if self.type == "super" or self.type == "oper": 2001 reps = ('choi', 'chi') 2002 q_oper = to_choi(self) if self.superrep not in reps else self 2003 return q_oper.iscp and q_oper.istp 2004 else: 2005 return False 2006 2007 @property 2008 def isherm(self): 2009 if self._isherm is not None: 2010 # used previously computed value 2011 return self._isherm 2012 2013 self._isherm = bool(zcsr_isherm(self.data)) 2014 2015 return self._isherm 2016 2017 @isherm.setter 2018 def isherm(self, isherm): 2019 self._isherm = isherm 2020 2021 def check_isunitary(self): 2022 """ 2023 Checks whether qobj is a unitary matrix 2024 """ 2025 if self.isoper: 2026 eye_data = fast_identity(self.shape[0]) 2027 return not ( 2028 np.any( 2029 np.abs((self.data*self.dag().data - eye_data).data) 2030 > settings.atol 2031 ) 2032 or 2033 np.any( 2034 np.abs((self.dag().data*self.data - eye_data).data) 2035 > settings.atol 2036 ) 2037 ) 2038 else: 2039 return False 2040 2041 @property 2042 def isunitary(self): 2043 if self._isunitary is not None: 2044 # used previously computed value 2045 return self._isunitary 2046 2047 self._isunitary = self.check_isunitary() 2048 2049 return self._isunitary 2050 2051 @isunitary.setter 2052 def isunitary(self, isunitary): 2053 self._isunitary = isunitary 2054 2055 @property 2056 def type(self): 2057 if not self._type: 2058 self._type = type_from_dims(self.dims) 2059 2060 return self._type 2061 2062 @property 2063 def shape(self): 2064 if self.data.shape == (1, 1): 2065 return tuple([np.prod(self.dims[0]), np.prod(self.dims[1])]) 2066 else: 2067 return tuple(self.data.shape) 2068 2069 @property 2070 def isbra(self): 2071 return self.type == 'bra' 2072 2073 @property 2074 def isket(self): 2075 return self.type == 'ket' 2076 2077 @property 2078 def isoperbra(self): 2079 return self.type == 'operator-bra' 2080 2081 @property 2082 def isoperket(self): 2083 return self.type == 'operator-ket' 2084 2085 @property 2086 def isoper(self): 2087 return self.type == 'oper' 2088 2089 @property 2090 def issuper(self): 2091 return self.type == 'super' 2092 2093 @staticmethod 2094 def evaluate(qobj_list, t, args): 2095 """Evaluate a time-dependent quantum object in list format. For 2096 example, 2097 2098 qobj_list = [H0, [H1, func_t]] 2099 2100 is evaluated to 2101 2102 Qobj(t) = H0 + H1 * func_t(t, args) 2103 2104 and 2105 2106 qobj_list = [H0, [H1, 'sin(w * t)']] 2107 2108 is evaluated to 2109 2110 Qobj(t) = H0 + H1 * sin(args['w'] * t) 2111 2112 Parameters 2113 ---------- 2114 qobj_list : list 2115 A nested list of Qobj instances and corresponding time-dependent 2116 coefficients. 2117 t : float 2118 The time for which to evaluate the time-dependent Qobj instance. 2119 args : dictionary 2120 A dictionary with parameter values required to evaluate the 2121 time-dependent Qobj intance. 2122 2123 Returns 2124 ------- 2125 output : :class:`qutip.Qobj` 2126 A Qobj instance that represents the value of qobj_list at time t. 2127 2128 """ 2129 2130 q_sum = 0 2131 if isinstance(qobj_list, Qobj): 2132 q_sum = qobj_list 2133 elif isinstance(qobj_list, list): 2134 for q in qobj_list: 2135 if isinstance(q, Qobj): 2136 q_sum += q 2137 elif (isinstance(q, list) and len(q) == 2 and 2138 isinstance(q[0], Qobj)): 2139 if isinstance(q[1], types.FunctionType): 2140 q_sum += q[0] * q[1](t, args) 2141 elif isinstance(q[1], str): 2142 args['t'] = t 2143 q_sum += q[0] * float(eval(q[1], globals(), args)) 2144 else: 2145 raise TypeError('Unrecognized format for ' + 2146 'specification of time-dependent Qobj') 2147 else: 2148 raise TypeError('Unrecognized format for specification ' + 2149 'of time-dependent Qobj') 2150 else: 2151 raise TypeError( 2152 'Unrecongized format for specification of time-dependent Qobj') 2153 2154 return q_sum 2155 2156 2157# ----------------------------------------------------------------------------- 2158# This functions evaluates a time-dependent quantum object on the list-string 2159# and list-function formats that are used by the time-dependent solvers. 2160# Although not used directly in by those solvers, it can for test purposes be 2161# conventient to be able to evaluate the expressions passed to the solver for 2162# arbitrary value of time. This function provides this functionality. 2163# 2164def qobj_list_evaluate(qobj_list, t, args): 2165 """ 2166 Deprecated: See Qobj.evaluate 2167 """ 2168 warnings.warn("Deprecated: Use Qobj.evaluate", DeprecationWarning) 2169 return Qobj.evaluate(qobj_list, t, args) 2170 2171 2172# ----------------------------------------------------------------------------- 2173# 2174# A collection of tests used to determine the type of quantum objects, and some 2175# functions for increased compatibility with quantum optics toolbox. 2176# 2177def dag(A): 2178 """Adjont operator (dagger) of a quantum object. 2179 2180 Parameters 2181 ---------- 2182 A : :class:`qutip.Qobj` 2183 Input quantum object. 2184 2185 Returns 2186 ------- 2187 oper : :class:`qutip.Qobj` 2188 Adjoint of input operator 2189 2190 Notes 2191 ----- 2192 This function is for legacy compatibility only. It is recommended to use 2193 the ``dag()`` Qobj method. 2194 2195 """ 2196 if not isinstance(A, Qobj): 2197 raise TypeError("Input is not a quantum object") 2198 2199 return A.dag() 2200 2201 2202def ptrace(Q, sel): 2203 """Partial trace of the Qobj with selected components remaining. 2204 2205 Parameters 2206 ---------- 2207 Q : :class:`qutip.Qobj` 2208 Composite quantum object. 2209 sel : int/list 2210 An ``int`` or ``list`` of components to keep after partial trace. 2211 2212 Returns 2213 ------- 2214 oper : :class:`qutip.Qobj` 2215 Quantum object representing partial trace with selected components 2216 remaining. 2217 2218 Notes 2219 ----- 2220 This function is for legacy compatibility only. It is recommended to use 2221 the ``ptrace()`` Qobj method. 2222 2223 """ 2224 if not isinstance(Q, Qobj): 2225 raise TypeError("Input is not a quantum object") 2226 2227 return Q.ptrace(sel) 2228 2229 2230def _ptrace_dense(Q, sel): 2231 rd = np.asarray(Q.dims[0], dtype=np.int32).ravel() 2232 nd = len(rd) 2233 if isinstance(sel, int): 2234 sel = np.array([sel]) 2235 else: 2236 sel = np.asarray(sel) 2237 sel = list(np.sort(sel)) 2238 for x in sel: 2239 if not 0 <= x < len(rd): 2240 raise IndexError("Invalid selection index in ptrace.") 2241 dkeep = (rd[sel]).tolist() 2242 qtrace = list(set(np.arange(nd)) - set(sel)) 2243 dtrace = (rd[qtrace]).tolist() 2244 if len(dkeep) + len(dtrace) != len(rd): 2245 raise ValueError("Duplicate selection index in ptrace.") 2246 if not dtrace: 2247 # If we are keeping all dimensions, no need to construct an ndarray. 2248 return Q.copy() 2249 rd = list(rd) 2250 if isket(Q): 2251 vmat = (Q.full() 2252 .reshape(rd) 2253 .transpose(sel + qtrace) 2254 .reshape([np.prod(dkeep, dtype=np.int32), 2255 np.prod(dtrace, dtype=np.int32)])) 2256 rhomat = vmat.dot(vmat.conj().T) 2257 else: 2258 rhomat = np.trace(Q.full() 2259 .reshape(rd + rd) 2260 .transpose(qtrace + [nd + q for q in qtrace] + 2261 sel + [nd + q for q in sel]) 2262 .reshape([np.prod(dtrace, dtype=np.int32), 2263 np.prod(dtrace, dtype=np.int32), 2264 np.prod(dkeep, dtype=np.int32), 2265 np.prod(dkeep, dtype=np.int32)])) 2266 return Qobj(rhomat, dims=[dkeep, dkeep]) 2267 2268 2269def dims(inpt): 2270 """Returns the dims attribute of a quantum object. 2271 2272 Parameters 2273 ---------- 2274 inpt : :class:`qutip.Qobj` 2275 Input quantum object. 2276 2277 Returns 2278 ------- 2279 dims : list 2280 A ``list`` of the quantum objects dimensions. 2281 2282 Notes 2283 ----- 2284 This function is for legacy compatibility only. Using the `Qobj.dims` 2285 attribute is recommended. 2286 2287 """ 2288 if isinstance(inpt, Qobj): 2289 return inpt.dims 2290 else: 2291 raise TypeError("Input is not a quantum object") 2292 2293 2294def shape(inpt): 2295 """Returns the shape attribute of a quantum object. 2296 2297 Parameters 2298 ---------- 2299 inpt : :class:`qutip.Qobj` 2300 Input quantum object. 2301 2302 Returns 2303 ------- 2304 shape : list 2305 A ``list`` of the quantum objects shape. 2306 2307 Notes 2308 ----- 2309 This function is for legacy compatibility only. Using the `Qobj.shape` 2310 attribute is recommended. 2311 2312 """ 2313 if isinstance(inpt, Qobj): 2314 return Qobj.shape 2315 else: 2316 return np.shape(inpt) 2317 2318 2319def isket(Q): 2320 """ 2321 Determines if given quantum object is a ket-vector. 2322 2323 Parameters 2324 ---------- 2325 Q : :class:`qutip.Qobj` 2326 Quantum object 2327 2328 Returns 2329 ------- 2330 isket : bool 2331 True if qobj is ket-vector, False otherwise. 2332 2333 Examples 2334 -------- 2335 >>> psi = basis(5,2) 2336 >>> isket(psi) 2337 True 2338 2339 Notes 2340 ----- 2341 This function is for legacy compatibility only. Using the `Qobj.isket` 2342 attribute is recommended. 2343 2344 """ 2345 return True if isinstance(Q, Qobj) and Q.isket else False 2346 2347 2348def isbra(Q): 2349 """Determines if given quantum object is a bra-vector. 2350 2351 Parameters 2352 ---------- 2353 Q : :class:`qutip.Qobj` 2354 Quantum object 2355 2356 Returns 2357 ------- 2358 isbra : bool 2359 True if Qobj is bra-vector, False otherwise. 2360 2361 Examples 2362 -------- 2363 >>> psi = basis(5,2) 2364 >>> isket(psi) 2365 False 2366 2367 Notes 2368 ----- 2369 This function is for legacy compatibility only. Using the `Qobj.isbra` 2370 attribute is recommended. 2371 2372 """ 2373 return True if isinstance(Q, Qobj) and Q.isbra else False 2374 2375 2376def isoperket(Q): 2377 """Determines if given quantum object is an operator in column vector form 2378 (operator-ket). 2379 2380 Parameters 2381 ---------- 2382 Q : :class:`qutip.Qobj` 2383 Quantum object 2384 2385 Returns 2386 ------- 2387 isoperket : bool 2388 True if Qobj is operator-ket, False otherwise. 2389 2390 Notes 2391 ----- 2392 This function is for legacy compatibility only. Using the `Qobj.isoperket` 2393 attribute is recommended. 2394 2395 """ 2396 return True if isinstance(Q, Qobj) and Q.isoperket else False 2397 2398 2399def isoperbra(Q): 2400 """Determines if given quantum object is an operator in row vector form 2401 (operator-bra). 2402 2403 Parameters 2404 ---------- 2405 Q : :class:`qutip.Qobj` 2406 Quantum object 2407 2408 Returns 2409 ------- 2410 isoperbra : bool 2411 True if Qobj is operator-bra, False otherwise. 2412 2413 Notes 2414 ----- 2415 This function is for legacy compatibility only. Using the `Qobj.isoperbra` 2416 attribute is recommended. 2417 2418 """ 2419 return True if isinstance(Q, Qobj) and Q.isoperbra else False 2420 2421 2422def isoper(Q): 2423 """Determines if given quantum object is a operator. 2424 2425 Parameters 2426 ---------- 2427 Q : :class:`qutip.Qobj` 2428 Quantum object 2429 2430 Returns 2431 ------- 2432 isoper : bool 2433 True if Qobj is operator, False otherwise. 2434 2435 Examples 2436 -------- 2437 >>> a = destroy(5) 2438 >>> isoper(a) 2439 True 2440 2441 Notes 2442 ----- 2443 This function is for legacy compatibility only. Using the `Qobj.isoper` 2444 attribute is recommended. 2445 2446 """ 2447 return True if isinstance(Q, Qobj) and Q.isoper else False 2448 2449 2450def issuper(Q): 2451 """Determines if given quantum object is a super-operator. 2452 2453 Parameters 2454 ---------- 2455 Q : :class:`qutip.Qobj` 2456 Quantum object 2457 2458 Returns 2459 ------- 2460 issuper : bool 2461 True if Qobj is superoperator, False otherwise. 2462 2463 Notes 2464 ----- 2465 This function is for legacy compatibility only. Using the `Qobj.issuper` 2466 attribute is recommended. 2467 2468 """ 2469 return True if isinstance(Q, Qobj) and Q.issuper else False 2470 2471 2472def isequal(A, B, tol=None): 2473 """Determines if two qobj objects are equal to within given tolerance. 2474 2475 Parameters 2476 ---------- 2477 A : :class:`qutip.Qobj` 2478 Qobj one 2479 B : :class:`qutip.Qobj` 2480 Qobj two 2481 tol : float 2482 Tolerence for equality to be valid 2483 2484 Returns 2485 ------- 2486 isequal : bool 2487 True if qobjs are equal, False otherwise. 2488 2489 Notes 2490 ----- 2491 This function is for legacy compatibility only. Instead, it is recommended 2492 to use the equality operator of Qobj instances instead: A == B. 2493 2494 """ 2495 if tol is None: 2496 tol = settings.atol 2497 2498 if not isinstance(A, Qobj) or not isinstance(B, Qobj): 2499 return False 2500 2501 if A.dims != B.dims: 2502 return False 2503 2504 Adat = A.data 2505 Bdat = B.data 2506 elems = (Adat - Bdat).data 2507 if np.any(np.abs(elems) > tol): 2508 return False 2509 2510 return True 2511 2512 2513def isherm(Q): 2514 """Determines if given operator is Hermitian. 2515 2516 Parameters 2517 ---------- 2518 Q : :class:`qutip.Qobj` 2519 Quantum object 2520 2521 Returns 2522 ------- 2523 isherm : bool 2524 True if operator is Hermitian, False otherwise. 2525 2526 Examples 2527 -------- 2528 >>> a = destroy(4) 2529 >>> isherm(a) 2530 False 2531 2532 Notes 2533 ----- 2534 This function is for legacy compatibility only. Using the `Qobj.isherm` 2535 attribute is recommended. 2536 2537 """ 2538 return True if isinstance(Q, Qobj) and Q.isherm else False 2539 2540 2541# TRAILING IMPORTS 2542# We do a few imports here to avoid circular dependencies. 2543from qutip.eseries import eseries 2544import qutip.superop_reps as sr 2545import qutip.tensor as tensor 2546import qutip.operators as ops 2547import qutip.metrics as mts 2548import qutip.states 2549import qutip.superoperator 2550