QuantumWalkP1

class QuantumWalkP1(N)[source]

Bases: object

Methods providing a very basic framework of a 1 particle CTQW.

Caution

This object constructs some of the framework of a 1 particle quantum walk, but is not feature complete. For a feature complete construction, see pyCTQW.MPI.Graph() or pyCTQW.MPI.Line() for which this is a base.

Method Summary

QuantumWalkP1.destroy() Destroys the 1 particle quantum walk object, and all associated PETSc matrices/vectors.
QuantumWalkP1.exportState(filename, filetype) Exports the final state of the quantum walk to a file.
QuantumWalkP1.importInitState(filename, filetype) Imports the initial state of the quantum walk from a file.
QuantumWalkP1.marginal(vec)
QuantumWalkP1.plotNodes(filename[, t]) Creates a plot of the node probablities over time.
QuantumWalkP1.propagate(t[, method]) Propagates the quantum walk for time t.
QuantumWalkP1.psiToInit() Copies the state self.psi to self.psi0.
QuantumWalkP1.watch(nodes) Creates a handle that watches node probability during propagation.

Method Details

QuantumWalkP1.destroy()[source]

Destroys the 1 particle quantum walk object, and all associated PETSc matrices/vectors.

QuantumWalkP1.exportState(filename, filetype)[source]

Exports the final state of the quantum walk to a file.

Parameters:
  • filename (str) – path to desired output file.
  • filetype (str) –

    the filetype of the exported adjacency matrix.

    • 'txt' - an \(N\) element column vector in text format.
    • 'bin' - an \(N\) element PETSc binary vector.
QuantumWalkP1.importInitState(filename, filetype)[source]

Imports the initial state of the quantum walk from a file.

Parameters:
  • filename (str) – path to the file containing the input state.
  • filetype (str) –

    the filetype of the imported adjacency matrix.

    • 'txt' - an \(N\) element column vector in text format.
    • 'bin' - an \(N\) element PETSc binary vector.
Returns:

: this creates a PETSc vector containing the initial state, accessed via Graph.psi0.

Return type:

petsc4py.PETSc.Vec()

Warning

The number of elements in the imported input state vector must match the number of nodes the 1P CTQW object is initialized with.

QuantumWalkP1.marginal(vec)[source]
QuantumWalkP1.plotNodes(filename, t=None)[source]

Creates a plot of the node probablities over time.

Parameters:filename (str) – the absolute/relative path to the desired output file.

Important

watch() must be called prior to propagation in order for the probabilities at the specified nodes to be stored.

Note

  • if multiple nodes are watched, they will all be plotted.
  • ensure a file extension is present so that filetype is correctly set (choose one of png, pdf, ps, eps or svg).
  • Timesteps are recorded when propagations occur - to increase the number of timesteps, increase the number of propagations.
QuantumWalkP1.propagate(t, method='chebyshev', **kwargs)[source]

Propagates the quantum walk for time t.

Parameters:
  • t (float) – the timestep over which to propagate the state walk.psi0 via the relationship \(\left|\psi\right\rangle = e^{-iHt}\left|\psi0\right\rangle\).
  • method (str) – the propagation algorithm to use ('chebyshev' (default) or 'krylov').
Keyword Arguments:
 
  • esolver (str) – the eigensolver algorithm to use.

    • 'krylovschur' (default) - Krylov-Schur
    • 'arnoldi' - Arnoldi Method
    • 'lanczos' - Lanczos Method
    • 'power' - Power Iteration/Rayleigh Quotient Iteration
    • 'gd' - Generalized Davidson
    • 'jd' - Jacobi-Davidson,
    • 'lapack' - Uses LAPACK eigenvalue solver routines
    • 'arpack' - only available if SLEPc is compiled with ARPACK linking
  • workType (str) – can be used to set the eigensolver worktype (either 'ncv' or 'mpd'). The default is to let SLEPc decide.

  • workSize (int) – sets the work size if workType is set.

  • tolIn (float) – tolerance of the eigenvalue solver (default 0. (SLEPc decides)).

  • maxIt (int) – maximum number of iterations of the eigenvalue solver (default 0 (SLEPc decides)).

  • verbose (bool) – if True, writes eigensolver information to the console

  • emax_estimate (float) – used to override the calculation of the graphs maximum eigenvalue.

  • emin_estimate (float) – used to override the calculation of the graphs minimum eigenvalue. For a finite graph, emin_estimate=0 is set by default.

    Caution

    • If supplied, the value of emax_estimate (\(\hat{\lambda}_{\max}\)) must satisfy \(\hat{\lambda}_{\max}\geq\lambda_{\max}\), where \(\lambda_{\max}\) is the actual maximum eigenvalue of the graph.
    • Similarly, the value of emin_estimate (\(\hat{\lambda}_{\min}\)) must satisfy \(\hat{\lambda}_{\min}\leq\lambda_{\min}\), where \(\lambda_{\min}\) is the actual minimum eigenvalue of the graph.
    • The greater the difference value of \(|\hat{\lambda}_{\max} -\lambda_{\max}|\) and \(|\hat{\lambda}_{\min} -\lambda_{\min}|\), the longer the convergence time of the chebyshev propagator.

    Note

    • The keyword arguments properties only apply if propagator='chebyshev'
    • For more information regarding these properties,refer to the SLEPc documentation
Returns:

: this creates a PETSc vector containing the propagated state, accessed via the attribute psi, as well as a PETSc vector containing the marginal probabilities, prob.

Return type:

petsc4py.PETSc.Vec()

Note

The EigSolver properties only take effect if method=’chebyshev’.

Fortran interface

This function calls the Fortran function expm(). for the Krylov method, and the function qw_cheby() for the Chebyshev method.

QuantumWalkP1.psiToInit()[source]

Copies the state self.psi to self.psi0.

Example

This can be used to iterate the propagation over discrete timesteps:

for i in range(10):
        walk.propagate(0.01,method='chebyshev')
        walk.psiToInit()
QuantumWalkP1.watch(nodes)[source]

Creates a handle that watches node probability during propagation.

Parameters:nodes (array of ints) – the nodes to watch (e.g. [0,1,4]).
Returns:: creates a handle that can be accessed to retrieve node probabilities for various \(t\)
Return type:pyCTQW.MPI.ctqw.nodeHandle()

Example

>>> walk.watch([0,1,2,3,4])
>>> walk.propagate(5.,method='chebyshev')
>>> timeArray, probArray = walk.handle.getLocalNodes()

Warning

Note that walk.handle attributes are not collective; if running on multiple nodes, only local values will be returned.