Quick Start

In this chapter, the basic behaviour of QuBricks is demonstrated in the context of some simple problems. For specific documentation on methods and their usage, please refer to the API documentation is subsequent chapters.

In the following, we assume that the following has been run in your Python environment.

>>> from qubricks import *
>>> from qubricks.wall import *
>>> import numpy as np

Note

QuBricks current only works in Python 2, since it depends on python-parameters which in turn is compatible only with Python 2.

Getting Started

In this section, we motivate the use of QuBricks in simple quantum systems. We will make use of the submodule qubricks.wall, which brings together many of the “bricks” that make up QuBricks into usable standalone tools.

Consider a single isolated electron in a magnetic field. Suppose that this magnetic field was composed of stable magnetic field along the Z axis \(B_z\), and a noisy magnetic field along the X axis \(B_x(t)\). The Hamiltonian describing the mechanics is then:

\[\begin{split}\begin{aligned} H & = & \mu_B \left(\begin{array}{cc} B_z & B_x(t)\\ B_x(t) & -B_z \end{array}\right),\end{aligned}\end{split}\]

where \(B_z\) is the magnetic field along \(z\) and \(B_x(t) = \bar{B}_{x} + \tilde{B}_x(t)\) is a noise field along \(x\) centred on some nominal value \(\bar{B}_{x}\) with a time-dependent noisy component \(\tilde{B_x}(t)\).

Let us assume that the noise in \(B_x(t)\) is white, so that:

\[\begin{split}\left< \tilde{B}_x(t) \tilde{B}_x(t^\prime) \right> = D\delta(t-t^\prime),\end{split}\]

We can model such white noise using a Lindblad superoperator.

This is a simple system which can be analytically solved. Evolution under this Hamiltonian will lead to the electron’s spin gyrating around an axis between Z and X (i.e. at an angle \(\theta = \tan^{-1}(B_z/J_x)\) from the x-axis) at a frequency of \(2\sqrt{B_z^2 + B_x^2}\). The effect of high frequency noise in \(B_x\) is to progressively increase the mixedness in the Z quadrature until such time as measurements of Z are unbiased. For example, when \(B_z=0\), the return probability for an initially up state is given by: \(p = \frac{1}{2}(1+\cos{2B_{x}t})\). Since \(\left< \tilde{B_{x}}^2 \right>_{\textrm{avg}} = D/t\), we find by taylor expanding that: \(\left<p\right> = 1 - Dt\). A more careful analysis would have found that:

\[\begin{split}\left<p\right> = \frac{1}{2}(1+exp(-2Dt)) .\end{split}\]

It is possible to find approximations for a general \(\theta\), but we leave that as an exercise. Alternatively, you can take advantage of QuBricks to simulate these results for us.

Figure 1: Dynamics of a single electron spin under a magnetic field aligned 45 degrees in the XZ plane.

For example, suppose we wanted to examine the case where \(B_z = B_x\), as shown in figure 1. This can be simulated using the following QuBricks code:

q = SimpleQuantumSystem(
  hamiltonian={
    'B':[[1,0],[0,-1]],
    'J':[[0,1],[1,0]]
  },
  parameters={
    'B':(40,'neV'),
    'J':(40,'neV'),
    'T_2': (400,'ns'),
    'D': (lambda T_2,c_hbar: 0.5*c_hbar**2/T_2, '{mu}J^2*ns')
  },
  measurements={
    'E': ExpectationMeasurement( [[1,0],[0,-1]], [[0,1],[1,0]], [[0,-1j],[1j,0]] ),
  },
  derivative_ops={
    'J_noise': LindbladStateOperator(coefficient='D',operator=[[0,1],[1,0]])
  }
)

ts = np.linspace(0,1e-6,1000)
r = q.measure.E.integrate(psi_0s=[ [1,0] ],times=ts, operators=['evolution','J_noise'])

We could then plot the results using matplotlib:

import matplotlib.pyplot as plt
from mplstyles import SampleStyle

style = SampleStyle()

with style:
	# If B=0, plot theoretical exponential decay curve
	if q.p.B == 0:
	    p_D = lambda t,D,c_hbar: np.exp(-2*D/c_hbar**2*t)
	    plt.plot(ts*1e9,q.p.range(p_D,t=ts),linestyle='--')

	plt.plot(r['time'][0]*1e9,r['expectation'][0,:,0],label="$\\left<Z\\right>$")
	plt.plot(r['time'][0]*1e9,r['expectation'][0,:,1],label="$\\left<X\\right>$")
	plt.plot(r['time'][0]*1e9,r['expectation'][0,:,2],label="$\\left<Y\\right>$")

	# Formatting options
	plt.grid()
	plt.legend(loc=0)
	plt.hlines([np.exp(-1),-np.exp(-1)],*plt.xlim())
	plt.xlabel('$\mathrm{Time}\,(ns)$')
	plt.ylabel("$E_Z$")

	style.savefig('results.pdf', polish=False, )

This would result in the following plot:

The above code takes advantage of several attributes and methods of QuantumSystem instances which may not be entirely clear. At this point, you can look them up in the API reference in subsequent chapters.

Advanced Usage

For more fine-grained control, one can subclass QuantumSystem, Measurement, StateOperator and Basis as necessary. For more information about which methods are available on these objects, refer to the API documentation below. The templates for subclassing are shown below.

QuantumSystem

from qubricks import QuantumSystem

class CustomSystem(QuantumSystem):
	'''
	Refer to the API documentation for `QuantumSystem` for more information.
	'''

	def init(self, **kwargs):
		'''
		This method can be used by subclasses to initialise the state
		of the `QuantumSystem` instance. Any excess kwargs beyond `parameters`
		passed to q `QuantumSystem` instance will be passed to this method.
		'''
		pass

	def init_parameters(self):
		'''
		This method can be used by subclasses to add any additional
		parameters required to describe the quantum system. If this 
		method returns a dictionary, then it is used to update the
		parameters stored in the `Parameters` instance. This would
		be equivalent to:
		
		>>> system.p << system.init_parameters()
		
		Parameters may, of course, be set directly by this method, using 
		(for example):
		
		>>> self.p.x = 1
		'''
		pass
	
	def init_hamiltonian(self):
		'''
		This method can be used by subclasses to initialise the Hamiltonian
		to be used to describe the Quantum System. The Hamiltonian can
		either be set directly in this method, using:
		
		>>> self.hamiltonian = <Operator or OperatorSet>
		
		Alternatively, if this method returns an Operator or OperatorSet
		object, then it will be set as the Hamiltonian for this QuantumSystem
		instance.
		'''
		pass

	def init_bases(self):
		'''
		This method can be used by subclasses to initialise the bases to be
		used by this instance of `QuantumSystem`. Bases can be added directly
		using the `add_basis` method; or, if this method returns a dictionary 
		of `Basis` objects (indexed by string names), then they will be added
		as bases of this system.
		'''
		pass

	def init_states(self):
		'''
		This method can be used by subclasses to initialise named states, ensembles,
		and subspaces. This can be done directly, using the corresponding 
		`add_state` and `add_subspace` methods. If a dictionary is returned, then
		it is assumed to be a dictionary of states indexed by names, which are then
		added to the system using `add_state`. Note that this assumes that the 
		states are represented in the standard basis of this `QuantumSystem` object.
		'''
		pass

	def init_measurements(self):
		'''
		This method can be used by subclasses to initialise the measurements that
		will be used with this `QuantumSystem` instance. This can be done directly, 
		using `add_measurement`; or, if this method returns a dictionary of `Measurement`
		objects indexed by string names, then they will be added as potential measurements
		of this quantum system.
		'''
		pass
	
	def init_derivative_ops(self):
		'''
		This method can be used by subclasses to initialise the `StateOperator` objects
		use to describe the time derivative of the evolution of the quantum system
		described by this object. Derivative operators may be added directly using
		`add_derivative_op`, or, if a dictionary of `StateOperator` objects is returned
		indexed with string names, then they are added as derivative operators of this object.
		If the operators depend on the hamiltonian or other properties of the quantum system,
		then the operators should be implemented in `get_derivative_ops` instead.
		
		This method should also initialise the default_derivative_ops property.
		'''
		pass
	
	def get_derivative_ops(self, components=None):
		'''
		This method can be used by subclasses to specify the `StateOperator` objects
		use to describe the time derivative of the evolution of the quantum system
		described by this object. These operators are added just before integration
		the operators described in `init_derivative_ops` and the 'evolution' operator
		describing Schroedinger evolution. Any properties of this
		`QuantumSystem` instance should not change before integration.
		
		:param components: The components activated in the Hamiltonian for this integration.
		:type components: iterable
		'''
		pass
	

Measurement

from qubricks import Measurement

class CustomMeasurement(Measurement):
	'''
	Refer to the API documentation `Measurement` for more details.
	'''

	def init(self, *args, **kwargs):
		'''
		This method should initialise the Measurement instance in whatever way
		is necessary to prepare the instance for use. Note that any arguments 
		passed to the Measurement constructor will also be passed to this method.
		There is no restriction on the method signature for the init method.
		'''
		raise NotImplementedError("Measurement.init has not been implemented.")
	
	def measure(self, data=None, params={}, int_kwargs={}, **kwargs):
		'''
		This method should return the value of a measurement as a numpy array with
		data type and shape as specified in `result_type` and `result_shape` respectively.

		.. note:: It is possible to return types other than numpy array and still
			be compatible with iteration (see MeasurementWrapper) provided you overload
			the `iterate_results_init` and `iterate_results_add` methods.

		Implementations of `measure` will typically be provided by integration data
		by a `MeasurementWrapper` instance (which will be a structured numpy array
		as returned by `Integrator.integrate) as the value for the `data` keyword.
		A consistent set of values for `times` and `initial` will also be passed
		as keywords inside `int_kwargs`.

		.. note:: If an implementation of `measure` omits the `data` keyword, QuBricks
			assumes that all integration required by the `measure` operator will be
			performed internally. It can use the reference to a QuantumSystem
			instance at `Measurement.system` for this purpose. If the `data` keyword
			is present (for testing/etc), but pre-computed integration data is undesired,
			override the `is_independent` method to return `True`. If external data
			is *required*, then simply remove the default value of `data`.

		Apart from the required keywords of `data` and `params`; any additional
		keywords can be specified. Refer to the documentation of `MeasurementWrapper` to
		see how their values will filter through to the various methods of QuBricks.

		:param data: Data from a QuantumSystem.integrate call, or None.
		:type data: numpy.ndarray or None
		:param params: Parameter context to use during this measurement. Parameter types can be anything supported by Parameters.
		:type params: dict
		:param int_kwargs: Keyword arguments to be passed on to any integrator instances, 
			which includes the times and initial states provided to `MeasurementWrapper.integrate`.
		:type int_kwargs: dict
		:param kwargs: Any other keyword arguments not collected explicitly.
		:type kwargs: dict
		'''
		raise NotImplementedError("Measurement.measure has not been implemented.")
	
		
	def result_type(self, *args, **kwargs):
		'''
		This method should return an object suitable for use as the dtype
		argument in a numpy array constructor. Otherwise, no restrictions; other than that it must also
		agree with the data-type returned by `Measurement.measure`.
		
		This method will receive all arguments and keyword arguments passed to
		`iterate_results_init`, where it is used to initialise the storage of 
		measurement results.
		'''
		raise NotImplementedError("Measurement.result_type has not been implemented.")
	
	def result_shape(self, *args, **kwargs):
		'''
		This method should return a tuple describing the shape of the numpy array to be returned
		by Measurement.measure.
		
		This method will receive all arguments and keyword arguments passed to
		`iterate_results_init`, where it is used to initialise the storage of 
		measurement results.
		'''
		raise NotImplementedError("Measurement.result_shape has not been implemented.")
	
#	The following is only needed if you do not want to receive pre-computed integration data.
# 	@property
# 	def is_independent(self):
# 		'''
# 		`True` if this Measurement instance does all required integration internally (and so should
# 		not receive pre-computed integration data). `False` otherwise. The default implementation is 
# 		`False`.
# 		'''
# 		return False

StateOperator

from qubricks import StateOperator

class CustomStateOperator(StateOperator):
	'''
	Refer to the API documentation for `StateOperator` for more information.
	'''
	
	def init(self, **kwargs):
		'''
		This method is called when StateOperator subclasses are
		initialised, which allows subclasses to set themselves up as appropriate.

		:param kwargs: The keyword arguments passed to the StateOperator constructor.
		:type kwargs: dict
		'''
		raise NotImplementedError("StateOperator.init is not implemented.")

	def __call__(self, state, t=0, params={}):
		'''
		StateOperator objects are called on states to effect some desired operation.
		States may be 1-dimensional (as state vectors) or 2-dimensional (as quantum ensembles),
		and each subclass should specify in StateOperator.for_state and StateOperator.for_ensemble
		which kinds of states are supported.
		'''
		raise NotImplementedError("StateOperator.__call__ is not implemented.")
	
	def restrict(self, *indices):
		'''
		This method should return a new StateOperator restricted to the basis states
		with indices `indices`. See Operator.restrict for more information.
		'''
		raise NotImplementedError("StateOperator.restrict is not implemented.")
	
	def connected(self, *indices, **params):
		'''
		This method should return the list of basis state indices which would mix
		with the specified basis state indices `indices` under repeated operation of the
		StateOperator. See Operator.connected for more information.
		'''
		raise NotImplementedError("StateOperator.connected is not implemented.")
	
	def for_state(self):
		'''
		Should be True if the StateOperator supports 1D vector operations; and False otherwise.
		'''
		raise NotImplementedError("StateOperator.for_state is not implemented.")

	def for_ensemble(self):
		'''
		Should be True if the StateOperator supports 2D ensemble operations; and False otherwise.
		'''
		raise NotImplementedError("StateOperator.process_args is not implemented.")
	
	def transform(self, transform_op):
		'''
		This method should transform all future operations on arbitrary input states
		according to the transformation `transform_op`. See Operator.transform
		for more information.
		'''
		raise NotImplementedError("StateOperator.transform is not implemented.")
	
#	The following method is optional	
# 	def collapse(self, *wrt, **params):
# 		'''
# 		This method is a stub to allow subclasses to simplify themselves when
# 		requested. If implemented, and Operators are used, the `collapse` method
# 		should be used on them also. See Operator.collapse for more information.
# 		
# 		Note that unless this method is overridden, no simplification occurs.
# 		'''
# 		return self
	

Basis

from qubricks import Basis

class CustomBasis(Basis):
	'''
	Refer to the API documentation of `Basis` for more details.
	'''
	
	def init(self, **kwargs):
		'''
		This method should do whatever is necessary to prepare the
		Basis instance for use. When this method is called by the 
		Python __init__ method, you can use `Basis.dim` to access
		the raw value of `dim`. If `dim` is necessary to construct 
		the operator, and it is not set, this method should raise an
		exception. All keyword arguments except `dim` and `parameters`
		passed to the `Basis` instance constructor will also be passed to 
		this method. 
		'''
		pass
	
	@property
	def operator(self):
		'''
		This method should return a two dimensional `Operator` object, with basis states as 
		columns. The `Operator` object should use the `Parameters` instance provided by the
		Basis instance. The simplest way to ensure this is to use the `Basis.Operator` method.
		'''
		raise NotImplementedError("Basis operator has not been implemented.")
	
#	The following methods are optional; refer to `Basis` documentation for more information.
# 	def state_info(self, state, params={}):
# 		'''
# 		This method (if implemented) should return a dictionary with more information
# 		about the state provided. There are no further constraints upon what might be
# 		returned.
# 		
# 		:param state: The state about which information should be returned.
# 		:type state: str or iterable
# 		:param params: A dictionary of parameter overrides. (see `parampy.Parameters`)
# 		:type params: dict
# 		'''
# 		return NotImplementedError("Basis.state_info has not been implemented.")
# 
# 	def state_toString(self, state, params={}):
# 		'''
# 		This method (if implemented) should return a string representation of the
# 		provided state, which should then be able to be converted back into the same
# 		state using `Basis.state_fromString`.
# 		
# 		:param state: The state which should be represented as a string.
# 		:type state: iterable
# 		:param params: A dictionary of parameter overrides. (see `parampy.Parameters`)
# 		:type params: dict
# 		'''
# 		raise NotImplementedError("Basis.state_toString has not been implemented.")
# 
# 	def state_fromString(self, string, params={}):
# 		'''
# 		This method (if implemented) should return the state as a numerical array that 
# 		is represented as a string in `string`. Calling `basis.state_toString` should then
# 		return the same (or equivalent) string representation.
# 		
# 		:param string: A string representation of a state.
# 		:type state: str
# 		:param params: A dictionary of parameter overrides. (see `parampy.Parameters`)
# 		:type params: dict
# 		'''
# 		raise NotImplementedError("Basis.state_fromString has not been implemented.")
# 
# 	def state_latex(self, state, params={}):
# 		'''
# 		This method (if implemented) should return string that when compiled by
# 		LaTeX would represent the state.
# 		
# 		:param state: The state which should be represented as a string.
# 		:type state: iterable
# 		:param params: A dictionary of parameter overrides. (see `parampy.Parameters`)
# 		:type params: dict
# 		'''
# 		raise NotImplementedError("Basis.state_latex has not been implemented.")

Integrator

In rare circumstances, you may find it necessary to subclass the Integrator class. Refer to the API documentation for more details on how to do this.

from qubricks import Integrator

class CustomIntegrator(Integrator):
    '''
    Refer to the API documentation for `Integrator` for more information.
    '''
    
    def _integrator(self, f, **kwargs):
        '''
        This method should return the object(s) necessary to perform
        the integration step. `f` is the the function which will return
        the derivative at each step.
        
        :param f: A function with signature f(t,y) which returns the derivative 
            at time `t` for the state `y`. Note that the derivative that is returned
            is that of `_derivative`, but `f` also handles progress reporting.
        :type f: function
        :param kwargs: Any additional keyword arguments passed to the `Integrator`
            constructor.
        :type kwargs: dict
        '''
        pass

    def _integrate(self, integrator, initial, times=None, **kwargs):
        '''
        This method should perform the integration using `integrator`, and
        return a list of two-tuples, each containing
        a time and a corresponding state. The times should be those listed in times,
        which will have been processed into floats.
        
        :param integrator: Whichever value was returned from `_integrator`.
        :type integrator: object
        :param initial: The state at which to start integrating. Will be the type
            returned by `_state_internal2ode`.
        :type initial: object
        :param times: A sequence of times for which to return the state.
        :type times: list of float
        :param kwargs: Additional keyword arguments passed to `Integrator.start`
            and/or `Integrator.extend`.
        :type kwargs: dict
        '''
        pass
    
    def _derivative(self, t, y, dim):
        '''
        This method should return the instantaneous derivative at time `t` 
        with current state `y` with dimensions `dim` (as returned by 
        `_state_internal2ode`. The derivative should be expressed
        in a form understood by the integrator returned by `_integrator`
        as used in `_integrate`.
        
        :param t: The current time.
        :type t: float
        :param y: The current state (in whatever form is returned by the integrator).
        :type y: object
        :param dim: The original dimensions of the state (as returned by `_state_internal2ode`).
        :type dim: object
        '''
        pass
    
    def _state_internal2ode(self, state):
        '''
        This method should return a tuple of a state and its original dimensions in some form.
        The state should be in a form understandable by the integrator returned by `_integrator`,
        and the derivative returned by `_derivative`.
        
        :param state: The state represented as a numpy array. Maybe 1D or 2D.
        :type state: numpy.ndarray
        '''
        pass

    def _state_ode2internal(self, state, dimensions):
        '''
        This method should restore and return the state (currently represented in the form used by the integrator
        returned by `_integrator`) to its representation as a numpy array using the 
        `dimensions` returned by `_state_internal2ode`.
        
        :param state: The state to re-represented as a numpy array.
        :type state: object
        :param dimensions: The dimensions returned by `_state_internal2ode`.
        :type dimensions: object
        '''
        pass

Operator Basics

One class that is worth discussing in more detail is Operator, which is among the most important “bricks” in the QuBricks library. It represents all of the two-dimensional linear operators used in QuBricks. The Operator object is neither directly a symbolic or numeric representation of an operator; but can be used to generate both.

Consider a simple example:

>>> op = Operator([[1,2],[3,4]])
>>> op
<Operator with shape (2,2)>

To generate a matrix representation of this object for inspection, we have two options depending upon whether we want a symbolic or numeric representation.

>>> op() # Numeric representation as a NumPy array
array([[ 1.+0.j, 2.+0.j],
       [ 3.+0.j, 4.+0.j]])
>>> op.symbolic() # Symbolic representation as a SymPy matrix
Matrix([
[1, 2],
[3, 4]])

In this case, there is not much difference, of course, since there are no symbolic parameters used.

Creating an Operator object with named parameters can be done in two ways. Either you must create a dictionary relating parameter names to matrix forms, or you can create a SymPy symbolic matrix. In both cases, one then passes this to the Operator constructor. For example:

>>> op = Operator('B':[[1,0],[0,-1]], 'J':[[0,1],[1,0]])
>>> op.symbolic()
Matrix([
[B, J],
[J, -B]])
>>> op.symbolic(J=10)
Matrix([
[ B, 10],
[10, -B]])
>>> op()
ValueError: Operator requires use of Parameters object; but none specified.

When representing Operator objects symbolically, we can override some parameters and perform parameter substitution. We see that attempting to generate a numeric representation of the Operator object failed, because it did not know how to assign a value to \(B\) and \(J\). Normally, Operator objects will have a reference to a Parameters instance (from python-parameters) passed to it in the constructor phase, for which these parameters can be extracted. This will in most cases be handled for you by QuantumSystem (see QuantumSystem in the API chapters), but for completeness there are two keyword arguments you can pass to Operator instances: parameters, which shold be a reference to an existing Parameters instance, and basis, which should be a reference to an existing Basis object or None (see Basis in the API chapters). For now, let us manually add it for demonstration purposes.

>>> from parameters import Parameters
>>> p = Parameters()
>>> p(B=2,J=1)
< Parameters with 2 definitions >
>>> op = Operator('B':[[1,0],[0,-1]], 'J':[[0,1],[1,0]],parameters=p)
>>> op()
array([[ 2.+0.j,  1.+0.j],
       [ 1.+0.j, -2.+0.j]])
>>> op(J=10,B=1)
array([[  1.+0.j, 10.+0.j],
       [ 10.+0.j, -1.+0.j]])

We see in the above that we can take advantage of temporary parameter overrides for numeric representations too [note that a parameters instance is still necessary for this].

The Parameters instance allows one to have parameters which are functions of one another, which allows for time and/or context dependent operators.

Operator objects support basic arithmetic: addition, subtraction, and multiplication using the standard python syntax. The inverse operation can be performed using the inverse method:

>>> op.inverse()

The Kronecker tensor product can be applied using the tensor method:

>>> op.tensor(other_op)

To apply an Operator object to a vector, you can either use the standard inbuilt multiplication operations, or use the slightly more optimised apply method.

If you are only interested in how certain parameters affect the operator, then to improve performance you can “collapse” the Operator down to only include variables which depend upon those variables.

>>> op.collapse('t',J=1)

The result of the above command would substitute all variables (with a parameter override of \(J=1\)) that do not depend upon \(t\) with their numerical value, and then perform various optimisations to make further substitutions more efficient. This is used, for example, by the integrator.

The last set of key methods of the Operator object are the connected and restrict methods. Operator.connected will return the set of all indicies (of the basis vectors in which the Operator is represented) that are connected by non-zero matrix elements, subject to the provided parameter substitution. Note that this comparison is done with the numerical values of the parameters.

>>> op = Operator('B':[[1,0],[0,-1]], 'J':[[0,1],[1,0]],parameters=p)
>>> op.connected(0)
0,1
>>> op.connected(0, J=0)
0

The restrict method returns a new Operator object which keeps only the entries in the old Operator object which correspond to the basis elements indicated by the indicies.

>>> op = Operator('B':[[1,0],[0,-1]], 'J':[[0,1],[1,0]],parameters=p)
>>> op.restrict(0)
<Operator with shape (1, 1)>
>>> op.symbolic()
Matrix([[B]])

For more detail, please refer to the API documentation for Operator.