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############################################################################### 33from collections.abc import Iterable 34import warnings 35from copy import deepcopy 36 37import numpy as np 38from scipy.interpolate import CubicSpline 39 40from qutip.qobj import Qobj 41from qutip.qobjevo import QobjEvo 42from qutip.operators import identity 43from qutip.qip.operations.gates import expand_operator, globalphase 44from qutip.tensor import tensor 45from qutip.mesolve import mesolve 46from qutip.mcsolve import mcsolve 47from qutip.qip.circuit import QubitCircuit 48from qutip.qip.noise import ( 49 Noise, RelaxationNoise, DecoherenceNoise, 50 ControlAmpNoise, RandomNoise, process_noise) 51from qutip.qip.pulse import Pulse, Drift, _merge_qobjevo, _fill_coeff 52 53 54__all__ = ['Processor'] 55 56 57class Processor(object): 58 """ 59 A simulator of a quantum device based on the QuTiP solver 60 :func:`qutip.mesolve`. It is defined by the available driving Hamiltonian 61 and the decoherence time for each component systems. The processor can 62 simulate the evolution under the given control pulses. Noisy evolution is 63 supported by :class:`.Noise` and can be added to the processor. 64 65 Parameters 66 ---------- 67 N: int 68 The number of component systems. 69 70 t1: list or float, optional 71 Characterize the decoherence of amplitude damping for 72 each qubit. A list of size `N` or a float for all qubits. 73 74 t2: list of float, optional 75 Characterize the decoherence of dephasing for 76 each qubit. A list of size `N` or a float for all qubits. 77 78 dims: list, optional 79 The dimension of each component system. 80 Default value is a 81 qubit system of ``dim=[2,2,2,...,2]`` 82 83 spline_kind: str, optional 84 Type of the coefficient interpolation. Default is "step_func" 85 Note that they have different requirement for the length of ``coeff``. 86 87 - "step_func": 88 The coefficient will be treated as a step function. E.g. 89 ``tlist=[0,1,2]`` and ``coeff=[3,2]``, means that the coefficient is 90 3 in t=[0,1) and 2 in t=[2,3). It requires 91 ``len(coeff)=len(tlist)-1`` or ``len(coeff)=len(tlist)``, but in the 92 second case the last element of `coeff` has no effect. 93 94 - "cubic": Use cubic interpolation for the coefficient. It requires 95 ``len(coeff)=len(tlist)`` 96 97 Attributes 98 ---------- 99 N: int 100 The number of component systems. 101 102 pulses: list of :class:`.Pulse` 103 A list of control pulses of this device 104 105 t1: float or list 106 Characterize the decoherence of amplitude damping of 107 each qubit. 108 109 t2: float or list 110 Characterize the decoherence of dephasing for 111 each qubit. 112 113 noise: :class:`.Noise`, optional 114 A list of noise objects. They will be processed when creating the 115 noisy :class:`qutip.QobjEvo` from the processor or run the simulation. 116 117 drift: :class:`qutip.qip.pulse.Drift` 118 A `Drift` object representing the drift Hamiltonians. 119 120 dims: list 121 The dimension of each component system. 122 Default value is a 123 qubit system of ``dim=[2,2,2,...,2]`` 124 125 spline_kind: str 126 Type of the coefficient interpolation. 127 See parameters of :class:`.Processor` for details. 128 """ 129 def __init__(self, N, t1=None, t2=None, 130 dims=None, spline_kind="step_func"): 131 self.N = N 132 self.pulses = [] 133 self.t1 = t1 134 self.t2 = t2 135 self.noise = [] 136 self.drift = Drift() 137 if dims is None: 138 self.dims = [2] * N 139 else: 140 self.dims = dims 141 self.pulse_mode = "discrete" 142 self.spline_kind = spline_kind 143 144 @property 145 def num_qubits(self): 146 return self.N 147 148 @num_qubits.setter 149 def num_qubits(self, value): 150 self.N = value 151 152 def add_drift(self, qobj, targets, cyclic_permutation=False): 153 """ 154 Add a drift Hamiltonians. The drift Hamiltonians are intrinsic 155 of the quantum system and cannot be controlled by external field. 156 157 Parameters 158 ---------- 159 qobj: :class:`qutip.Qobj` 160 The drift Hamiltonian. 161 targets: list 162 The indices of the target qubits 163 (or subquantum system of other dimensions). 164 """ 165 if not isinstance(qobj, Qobj): 166 raise TypeError("The drift Hamiltonian must be a qutip.Qobj.") 167 if not qobj.isherm: 168 raise ValueError("The drift Hamiltonian must be Hermitian.") 169 170 num_qubits = len(qobj.dims[0]) 171 if targets is None: 172 targets = list(range(num_qubits)) 173 if not isinstance(targets, list): 174 targets = [targets] 175 if cyclic_permutation: 176 for i in range(self.N): 177 temp_targets = [(t + i) % self.N for t in targets] 178 self.drift.add_drift(qobj, temp_targets) 179 else: 180 self.drift.add_drift(qobj, targets) 181 182 def add_control(self, qobj, targets=None, cyclic_permutation=False, 183 label=None): 184 """ 185 Add a control Hamiltonian to the processor. It creates a new 186 :class:`.Pulse` 187 object for the device that is turned off 188 (``tlist = None``, ``coeff = None``). To activate the pulse, one 189 can set its `tlist` and `coeff`. 190 191 Parameters 192 ---------- 193 qobj: :class:`qutip.Qobj` 194 The Hamiltonian for the control pulse.. 195 196 targets: list, optional 197 The indices of the target qubits 198 (or subquantum system of other dimensions). 199 200 cyclic_permutation: bool, optional 201 If true, the Hamiltonian will be expanded for 202 all cyclic permutation of the target qubits. 203 204 label: str, optional 205 The label (name) of the pulse 206 """ 207 # Check validity of ctrl 208 if not isinstance(qobj, Qobj): 209 raise TypeError("The control Hamiltonian must be a qutip.Qobj.") 210 if not qobj.isherm: 211 raise ValueError("The control Hamiltonian must be Hermitian.") 212 213 num_qubits = len(qobj.dims[0]) 214 if targets is None: 215 targets = list(range(num_qubits)) 216 if not isinstance(targets, list): 217 targets = [targets] 218 if cyclic_permutation: 219 for i in range(self.N): 220 temp_targets = [(t + i) % self.N for t in targets] 221 if label is not None: 222 temp_label = label + "_" + str(temp_targets) 223 temp_label = label 224 self.pulses.append( 225 Pulse(qobj, temp_targets, spline_kind=self.spline_kind, 226 label=temp_label)) 227 else: 228 self.pulses.append( 229 Pulse(qobj, targets, spline_kind=self.spline_kind, label=label) 230 ) 231 232 def find_pulse(self, pulse_name): 233 if isinstance(pulse_name, str): 234 try: 235 return self.pulses[self.pulse_dict[pulse_name]] 236 except (KeyError): 237 raise KeyError( 238 "Pulse name {} undefined. " 239 "Please define it in the attribute " 240 "`pulse_dict`.".format(pulse_name)) 241 elif isinstance(pulse_name, int): 242 return self.pulses[pulse_name] 243 else: 244 raise TypeError( 245 "pulse_name is either a string or an integer, not " 246 "{}".format(type(pulse_name)) 247 ) 248 249 @property 250 def ctrls(self): 251 """ 252 A list of Hamiltonians of all pulses. 253 """ 254 result = [] 255 for pulse in self.pulses: 256 result.append(pulse.get_ideal_qobj(self.dims)) 257 return result 258 259 @property 260 def coeffs(self): 261 """ 262 A list of the coefficients for all control pulses. 263 """ 264 if not self.pulses: 265 return None 266 coeffs_list = [pulse.coeff for pulse in self.pulses] 267 return coeffs_list 268 269 @coeffs.setter 270 def coeffs(self, coeffs_list): 271 for i, coeff in enumerate(coeffs_list): 272 self.pulses[i].coeff = coeff 273 274 @property 275 def pulse_mode(self): 276 if self.spline_kind == "step_func": 277 return "discrete" 278 elif self.spline_kind == "cubic": 279 return "continuous" 280 else: 281 raise ValueError( 282 "Saved spline_kind not understood.") 283 284 @pulse_mode.setter 285 def pulse_mode(self, mode): 286 if mode == "discrete": 287 spline_kind = "step_func" 288 elif mode == "continuous": 289 spline_kind = "cubic" 290 else: 291 raise ValueError( 292 "Pulse mode must be either discrete or continuous.") 293 294 self.spline_kind = spline_kind 295 for pulse in self.pulses: 296 pulse.spline_kind = spline_kind 297 298 def get_full_tlist(self, tol=1.0e-10): 299 """ 300 Return the full tlist of the ideal pulses. 301 If different pulses have different time steps, 302 it will collect all the time steps in a sorted array. 303 304 Returns 305 ------- 306 full_tlist: array-like 1d 307 The full time sequence for the ideal evolution. 308 """ 309 full_tlist = [pulse.tlist 310 for pulse in self.pulses if pulse.tlist is not None] 311 if not full_tlist: 312 return None 313 full_tlist = np.unique(np.sort(np.hstack(full_tlist))) 314 # account for inaccuracy in float-point number 315 full_tlist = np.concatenate( 316 (full_tlist[:1], full_tlist[1:][np.diff(full_tlist) > tol])) 317 return full_tlist 318 319 def get_full_coeffs(self, full_tlist=None): 320 """ 321 Return the full coefficients in a 2d matrix form. 322 Each row corresponds to one pulse. If the `tlist` are 323 different for different pulses, the length of each row 324 will be same as the `full_tlist` (see method 325 `get_full_tlist`). Interpolation is used for 326 adding the missing coefficient according to `spline_kind`. 327 328 Returns 329 ------- 330 coeffs: array-like 2d 331 The coefficients for all ideal pulses. 332 """ 333 # TODO add tests 334 self._is_pulses_valid() 335 if not self.pulses: 336 return np.array((0, 0), dtype=float) 337 if full_tlist is None: 338 full_tlist = self.get_full_tlist() 339 coeffs_list = [] 340 for pulse in self.pulses: 341 if pulse.tlist is None and pulse.coeff is None: 342 coeffs_list.append(np.zeros(len(full_tlist))) 343 continue 344 if not isinstance(pulse.coeff, (bool, np.ndarray)): 345 raise ValueError( 346 "get_full_coeffs only works for " 347 "NumPy array or bool coeff.") 348 if isinstance(pulse.coeff, bool): 349 if pulse.coeff: 350 coeffs_list.append(np.ones(len(full_tlist))) 351 else: 352 coeffs_list.append(np.zeros(len(full_tlist))) 353 continue 354 if self.spline_kind == "step_func": 355 arg = {"_step_func_coeff": True} 356 coeffs_list.append( 357 _fill_coeff(pulse.coeff, pulse.tlist, full_tlist, arg)) 358 elif self.spline_kind == "cubic": 359 coeffs_list.append( 360 _fill_coeff(pulse.coeff, pulse.tlist, full_tlist, {})) 361 else: 362 raise ValueError("Unknown spline kind.") 363 return np.array(coeffs_list) 364 365 def set_all_tlist(self, tlist): 366 """ 367 Set the same `tlist` for all the pulses. 368 369 Parameters 370 ---------- 371 tlist: array-like, optional 372 A list of time at which the time-dependent coefficients are 373 applied. See :class:`.Pulse` for detailed information` 374 """ 375 if isinstance(tlist, list) and len(tlist) == len(self.pulses): 376 for i, pulse in enumerate(self.pulses): 377 pulse.tlist = tlist[i] 378 else: 379 for pulse in self.pulses: 380 pulse.tlist = tlist 381 382 def add_pulse(self, pulse): 383 """ 384 Add a new pulse to the device. 385 386 Parameters 387 ---------- 388 pulse: :class:`.Pulse` 389 `Pulse` object to be added. 390 """ 391 if isinstance(pulse, Pulse): 392 if pulse.spline_kind is None: 393 pulse.spline_kind = self.spline_kind 394 self.pulses.append(pulse) 395 else: 396 raise ValueError("Invalid input, pulse must be a Pulse object") 397 398 def remove_pulse(self, indices=None, label=None): 399 """ 400 Remove the control pulse with given indices. 401 402 Parameters 403 ---------- 404 indices: int or list of int 405 The indices of the control Hamiltonians to be removed. 406 label: str 407 The label of the pulse 408 """ 409 if indices is not None: 410 if not isinstance(indices, Iterable): 411 indices = [indices] 412 indices.sort(reverse=True) 413 for ind in indices: 414 del self.pulses[ind] 415 else: 416 for ind, pulse in enumerate(self.pulses): 417 if pulse.label == label: 418 del self.pulses[ind] 419 420 def _is_pulses_valid(self): 421 """ 422 Check if the pulses are in the correct shape. 423 424 Returns: bool 425 If they are valid or not 426 """ 427 for i, pulse in enumerate(self.pulses): 428 if pulse.coeff is None or isinstance(pulse.coeff, bool): 429 # constant pulse 430 continue 431 if pulse.tlist is None: 432 raise ValueError( 433 "Pulse id={} is invalid. " 434 "Please define a tlist for the pulse.".format(i)) 435 if pulse.tlist is not None and pulse.coeff is None: 436 raise ValueError( 437 "Pulse id={} is invalid. " 438 "Please define a coeff for the pulse.".format(i)) 439 coeff_len = len(pulse.coeff) 440 tlist_len = len(pulse.tlist) 441 if pulse.spline_kind == "step_func": 442 if coeff_len == tlist_len-1 or coeff_len == tlist_len: 443 pass 444 else: 445 raise ValueError( 446 "The length of tlist and coeff of the pulse " 447 "labelled {} is invalid. " 448 "It's either len(tlist)=len(coeff) or " 449 "len(tlist)-1=len(coeff) for coefficients " 450 "as step function".format(i)) 451 else: 452 if coeff_len == tlist_len: 453 pass 454 else: 455 raise ValueError( 456 "The length of tlist and coeff of the pulse " 457 "labelled {} is invalid. " 458 "It should be either len(tlist)=len(coeff)".format(i)) 459 return True 460 461 def add_noise(self, noise): 462 """ 463 Add a noise object to the processor 464 465 Parameters 466 ---------- 467 noise: :class:`.Noise` 468 The noise object defined outside the processor 469 """ 470 if isinstance(noise, Noise): 471 self.noise.append(noise) 472 else: 473 raise TypeError("Input is not a Noise object.") 474 475 def save_coeff(self, file_name, inctime=True): 476 """ 477 Save a file with the control amplitudes in each timeslot. 478 479 Parameters 480 ---------- 481 file_name: string 482 Name of the file. 483 484 inctime: bool, optional 485 True if the time list should be included in the first column. 486 """ 487 self._is_pulses_valid() 488 coeffs = np.array(self.get_full_coeffs()) 489 if inctime: 490 shp = coeffs.T.shape 491 data = np.empty((shp[0], shp[1] + 1), dtype=np.float64) 492 data[:, 0] = self.get_full_tlist() 493 data[:, 1:] = coeffs.T 494 else: 495 data = coeffs.T 496 497 np.savetxt(file_name, data, delimiter='\t', fmt='%1.16f') 498 499 def read_coeff(self, file_name, inctime=True): 500 """ 501 Read the control amplitudes matrix and time list 502 saved in the file by `save_amp`. 503 504 Parameters 505 ---------- 506 file_name: string 507 Name of the file. 508 509 inctime: bool, optional 510 True if the time list in included in the first column. 511 512 Returns 513 ------- 514 tlist: array_like 515 The time list read from the file. 516 517 coeffs: array_like 518 The pulse matrix read from the file. 519 """ 520 data = np.loadtxt(file_name, delimiter='\t') 521 if not inctime: 522 self.coeffs = data.T 523 return self.coeffs 524 else: 525 tlist = data[:, 0] 526 self.set_all_tlist(tlist) 527 self.coeffs = data[:, 1:].T 528 return self.get_full_tlist, self.coeffs 529 530 def get_noisy_pulses(self, device_noise=False, drift=False): 531 """ 532 It takes the pulses defined in the `Processor` and 533 adds noise according to `Processor.noise`. It does not modify the 534 pulses saved in `Processor.pulses` but returns a new list. 535 The length of the new list of noisy pulses might be longer 536 because of drift Hamiltonian and device noise. They will be 537 added to the end of the pulses list. 538 539 Parameters 540 ---------- 541 device_noise: bool, optional 542 If true, include pulse independent noise such as single qubit 543 Relaxation. Default is False. 544 drift: bool, optional 545 If true, include drift Hamiltonians. Default is False. 546 547 Returns 548 ------- 549 noisy_pulses: list of :class:`.Pulse` 550 A list of noisy pulses. 551 """ 552 pulses = deepcopy(self.pulses) 553 noisy_pulses = process_noise( 554 pulses, self.noise, self.dims, t1=self.t1, t2=self.t2, 555 device_noise=device_noise) 556 if drift: 557 noisy_pulses += [self.drift] 558 return noisy_pulses 559 560 def get_qobjevo(self, args=None, noisy=False): 561 """ 562 Create a :class:`qutip.QobjEvo` representation of the evolution. 563 It calls the method `get_noisy_pulses` and create the `QobjEvo` 564 from it. 565 566 Parameters 567 ---------- 568 args: dict, optional 569 Arguments for :class:`qutip.QobjEvo` 570 noisy: bool, optional 571 If noise are included. Default is False. 572 573 Returns 574 ------- 575 qobjevo: :class:`qutip.QobjEvo` 576 The :class:`qutip.QobjEvo` representation of the unitary evolution. 577 c_ops: list of :class:`qutip.QobjEvo` 578 A list of lindblad operators is also returned. if ``noisy==Flase``, 579 it is always an empty list. 580 """ 581 # TODO test it for non array-like coeff 582 # check validity 583 self._is_pulses_valid() 584 585 if args is None: 586 args = {} 587 else: 588 args = args 589 # set step function 590 591 if not noisy: 592 dynamics = self.pulses 593 else: 594 dynamics = self.get_noisy_pulses( 595 device_noise=True, drift=True) 596 597 qu_list = [] 598 c_ops = [] 599 for pulse in dynamics: 600 if noisy: 601 qu, new_c_ops = pulse.get_noisy_qobjevo(dims=self.dims) 602 c_ops += new_c_ops 603 else: 604 qu = pulse.get_ideal_qobjevo(dims=self.dims) 605 qu_list.append(qu) 606 607 final_qu = _merge_qobjevo(qu_list) 608 final_qu.args.update(args) 609 610 # bring all c_ops to the same tlist, won't need it in QuTiP 5 611 full_tlist = self.get_full_tlist() 612 temp = [] 613 for c_op in c_ops: 614 temp.append(_merge_qobjevo([c_op], full_tlist)) 615 c_ops = temp 616 617 if noisy: 618 return final_qu, c_ops 619 else: 620 return final_qu, [] 621 622 def run_analytically(self, init_state=None, qc=None): 623 """ 624 Simulate the state evolution under the given `qutip.QubitCircuit` 625 with matrice exponentiation. It will calculate the propagator 626 with matrix exponentiation and return a list of :class:`qutip.Qobj`. 627 This method won't include noise or collpase. 628 629 Parameters 630 ---------- 631 qc: :class:`.QubitCircuit`, optional 632 Takes the quantum circuit to be implemented. If not given, use 633 the quantum circuit saved in the processor by ``load_circuit``. 634 635 init_state: :class:`qutip.Qobj`, optional 636 The initial state of the qubits in the register. 637 638 Returns 639 ------- 640 evo_result: :class:`qutip.Result` 641 An instance of the class 642 :class:`qutip.Result` will be returned. 643 """ 644 if init_state is not None: 645 U_list = [init_state] 646 else: 647 U_list = [] 648 tlist = self.get_full_tlist() 649 # TODO replace this by get_complete_coeff 650 coeffs = self.get_full_coeffs() 651 for n in range(len(tlist)-1): 652 H = sum([coeffs[m, n] * self.ctrls[m] 653 for m in range(len(self.ctrls))]) 654 dt = tlist[n + 1] - tlist[n] 655 U = (-1j * H * dt).expm() 656 U = self.eliminate_auxillary_modes(U) 657 U_list.append(U) 658 659 try: # correct_global_phase are defined for ModelProcessor 660 if self.correct_global_phase and self.global_phase != 0: 661 U_list.append(globalphase(self.global_phase, N=self.N)) 662 except AttributeError: 663 pass 664 665 return U_list 666 667 def run(self, qc=None): 668 """ 669 Calculate the propagator of the evolution by matrix exponentiation. 670 This method won't include noise or collpase. 671 672 Parameters 673 ---------- 674 qc: :class:`.QubitCircuit`, optional 675 Takes the quantum circuit to be implemented. If not given, use 676 the quantum circuit saved in the processor by `load_circuit`. 677 678 Returns 679 ------- 680 U_list: list 681 The propagator matrix obtained from the physical implementation. 682 """ 683 if qc: 684 self.load_circuit(qc) 685 return self.run_analytically(qc=qc, init_state=None) 686 687 def run_state(self, init_state=None, analytical=False, states=None, 688 noisy=True, solver="mesolve", **kwargs): 689 """ 690 If `analytical` is False, use :func:`qutip.mesolve` to 691 calculate the time of the state evolution 692 and return the result. Other arguments of mesolve can be 693 given as keyword arguments. 694 695 If `analytical` is True, calculate the propagator 696 with matrix exponentiation and return a list of matrices. 697 Noise will be neglected in this option. 698 699 Parameters 700 ---------- 701 init_state: Qobj 702 Initial density matrix or state vector (ket). 703 704 analytical: bool 705 If True, calculate the evolution with matrices exponentiation. 706 707 states: :class:`qutip.Qobj`, optional 708 Old API, same as init_state. 709 710 solver: str 711 "mesolve" or "mcsolve" 712 713 **kwargs 714 Keyword arguments for the qutip solver. 715 716 Returns 717 ------- 718 evo_result: :class:`qutip.Result` 719 If ``analytical`` is False, an instance of the class 720 :class:`qutip.Result` will be returned. 721 722 If ``analytical`` is True, a list of matrices representation 723 is returned. 724 """ 725 if states is not None: 726 warnings.warn( 727 "states will be deprecated and replaced by init_state", 728 DeprecationWarning) 729 if init_state is None and states is None: 730 raise ValueError("Qubit state not defined.") 731 elif init_state is None: 732 # just to keep the old parameters `states`, 733 # it is replaced by init_state 734 init_state = states 735 if analytical: 736 if kwargs or self.noise: 737 raise warnings.warn( 738 "Analytical matrices exponentiation" 739 "does not process noise or" 740 "any keyword arguments.") 741 return self.run_analytically(init_state=init_state) 742 743 # kwargs can not contain H or tlist 744 if "H" in kwargs or "tlist" in kwargs: 745 raise ValueError( 746 "`H` and `tlist` are already specified by the processor " 747 "and can not be given as a keyword argument") 748 749 # construct qobjevo for unitary evolution 750 if "args" in kwargs: 751 noisy_qobjevo, sys_c_ops = self.get_qobjevo( 752 args=kwargs["args"], noisy=noisy) 753 else: 754 noisy_qobjevo, sys_c_ops = self.get_qobjevo(noisy=noisy) 755 756 # add collpase operators into kwargs 757 if "c_ops" in kwargs: 758 if isinstance(kwargs["c_ops"], (Qobj, QobjEvo)): 759 kwargs["c_ops"] += [kwargs["c_ops"]] + sys_c_ops 760 else: 761 kwargs["c_ops"] += sys_c_ops 762 else: 763 kwargs["c_ops"] = sys_c_ops 764 765 # choose solver: 766 if solver == "mesolve": 767 evo_result = mesolve( 768 H=noisy_qobjevo, rho0=init_state, 769 tlist=noisy_qobjevo.tlist, **kwargs) 770 elif solver == "mcsolve": 771 evo_result = mcsolve( 772 H=noisy_qobjevo, psi0=init_state, 773 tlist=noisy_qobjevo.tlist, **kwargs) 774 775 return evo_result 776 777 def load_circuit(self, qc): 778 """ 779 Translate an :class:`.QubitCircuit` to its 780 corresponding Hamiltonians. (Defined in subclasses) 781 """ 782 raise NotImplementedError("Use the function in the sub-class") 783 784 def eliminate_auxillary_modes(self, U): 785 """ 786 Eliminate the auxillary modes like the cavity modes in cqed. 787 (Defined in subclasses) 788 """ 789 return U 790 791 def get_operators_labels(self): 792 """ 793 Get the labels for each Hamiltonian. 794 It is used in the method``plot_pulses``. 795 It is a 2-d nested list, in the plot, 796 a different color will be used for each sublist. 797 """ 798 label_list = [] 799 for pulse in self.pulses: 800 label_list.append(pulse.label) 801 return [label_list] 802 803 def plot_pulses( 804 self, title=None, figsize=(12, 6), dpi=None, 805 show_axis=False, rescale_pulse_coeffs=True, 806 num_steps=1000): 807 """ 808 Plot the ideal pulse coefficients. 809 810 Parameters 811 ---------- 812 title: str, optional 813 Title for the plot. 814 815 figsize: tuple, optional 816 The size of the figure. 817 818 dpi: int, optional 819 The dpi of the figure. 820 821 show_axis: bool, optional 822 If the axis are shown. 823 824 rescale_pulse_coeffs: bool, optional 825 Rescale the hight of each pulses. 826 827 num_steps: int, optional 828 Number of time steps in the plot. 829 830 Returns 831 ------- 832 fig: matplotlib.figure.Figure 833 The `Figure` object for the plot. 834 835 ax: matplotlib.axes._subplots.AxesSubplot 836 The axes for the plot. 837 838 Notes 839 ----- 840 ``plot_pulses`` only works for array_like coefficients 841 """ 842 import matplotlib.pyplot as plt 843 import matplotlib.gridspec as gridspec 844 color_list = plt.rcParams['axes.prop_cycle'].by_key()['color'] 845 846 # create a axis for each pulse 847 fig = plt.figure(figsize=figsize, dpi=dpi) 848 grids = gridspec.GridSpec(len(self.pulses), 1) 849 grids.update(wspace=0., hspace=0.) 850 851 tlist = np.linspace(0., self.get_full_tlist()[-1], num_steps) 852 dt = tlist[1] - tlist[0] 853 854 # make sure coeffs start and end with zero, for ax.fill 855 tlist = np.hstack(([-dt*1.e-20], tlist, [tlist[-1] + dt*1.e-20])) 856 coeffs = [] 857 for pulse in self.pulses: 858 coeffs.append(_pulse_interpolate(pulse, tlist)) 859 860 pulse_ind = 0 861 axis = [] 862 for i, label_group in enumerate(self.get_operators_labels()): 863 for j, label in enumerate(label_group): 864 grid = grids[pulse_ind] 865 ax = plt.subplot(grid) 866 axis.append(ax) 867 ax.fill(tlist, coeffs[pulse_ind], color_list[i], alpha=0.7) 868 ax.plot(tlist, coeffs[pulse_ind], color_list[i]) 869 if rescale_pulse_coeffs: 870 ymax = np.max(np.abs(coeffs[pulse_ind])) * 1.1 871 else: 872 ymax = np.max(np.abs(coeffs)) * 1.1 873 if ymax != 0.: 874 ax.set_ylim((-ymax, ymax)) 875 876 # disable frame and ticks 877 if not show_axis: 878 ax.set_xticks([]) 879 ax.spines['bottom'].set_visible(False) 880 ax.spines['top'].set_visible(False) 881 ax.spines['right'].set_visible(False) 882 ax.spines['left'].set_visible(False) 883 ax.set_yticks([]) 884 ax.set_ylabel(label, rotation=0) 885 pulse_ind += 1 886 if i == 0 and j == 0 and title is not None: 887 ax.set_title(title) 888 fig.tight_layout() 889 return fig, axis 890 891 892def _pulse_interpolate(pulse, tlist): 893 """ 894 A function that calls Scipy interpolation routine. Used for plotting. 895 """ 896 if pulse.tlist is None and pulse.coeff is None: 897 coeff = np.zeros(len(tlist)) 898 return coeff 899 if isinstance(pulse.coeff, bool): 900 if pulse.coeff: 901 coeff = np.ones(len(tlist)) 902 else: 903 coeff = np.zeros(len(tlist)) 904 return coeff 905 coeff = pulse.coeff 906 if len(coeff) == len(pulse.tlist)-1: # for discrete pulse 907 coeff = np.concatenate([coeff, [0]]) 908 909 from scipy import interpolate 910 if pulse.spline_kind == "step_func": 911 kind = "previous" 912 else: 913 kind = "cubic" 914 inter = interpolate.interp1d( 915 pulse.tlist, coeff, kind=kind, 916 bounds_error=False, fill_value=0.0) 917 return inter(tlist) 918