Bases: object
Methods providing a very basic framework of a 2 particle CTQW.
Caution
This object constructs some of the framework of a 2 particle quantum walk, but is not feature complete. For a feature complete construction, see pyCTQW.MPI.Graph2P() or pyCTQW.MPI.Line2P() for which this is a base.
QuantumWalkP2.destroy() | Destroys the 2 particle quantum walk object, and all associated PETSc matrices/vectors. |
QuantumWalkP2.exportPartialTrace(filename, ...) | Exports the partial trace of the density matrix to a file. |
QuantumWalkP2.exportState(filename, filetype) | Exports the final state psi of the quantum walk to a file. |
QuantumWalkP2.importInitState(filename, filetype) | Imports the initial state of the quantum walk from a file. |
QuantumWalkP2.marginal(vec) | |
QuantumWalkP2.plotEntanglement(filename) | Creates a plot of the entanglement over time. |
QuantumWalkP2.plotNode(filename, node[, t]) | Creates a plot of the marginal probablities on a specified node over time. |
QuantumWalkP2.plotNodes(filename[, p, t]) | Creates a plot of the node probablities over time. |
QuantumWalkP2.propagate(t[, method]) | Propagates the quantum walk for time t. |
QuantumWalkP2.psiToInit() | Copies the state self.psi to self.psi0. |
QuantumWalkP2.watch(nodes[, watchtype]) | Creates a handle that watches either marginal probability or walker entanglement during propagation. |
Destroys the 2 particle quantum walk object, and all associated PETSc matrices/vectors.
Exports the partial trace of the density matrix to a file.
Parameters: |
|
---|
For example, when p=1, \(\rho_2 = \text{Tr}_1(\rho_{12})\), where \(\rho_{12} = |\psi(t)\rangle\langle\psi(t)|\), is exported.
Fortran interface
This function calls the Fortran function partial_trace_mat() for bin filetypes, and the function partial_trace_array() otherwise.
Exports the final state psi of the quantum walk to a file.
Parameters: |
|
---|
Imports the initial state of the quantum walk from a file.
Parameters: |
|
---|---|
Returns: | : this creates a PETSc vector containing the initial statespace, accessed via self.psi0. |
Return type: | petsc4py.PETSc.Vec() |
Warning
The number of elements in the imported input statespace must equal \(N^2\) (either as a \(N\times N\) text file or a \(N^2\) PETSc binary vector), where \(N\) is the number of nodes the 2P CTQW object is initialized with.
Creates a plot of the entanglement 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 entanglement to be stored during propagation.
Note
Creates a plot of the marginal probablities on a specified node over time.
Parameters: |
|
---|
Important
watch() must be called prior to propagation in order for the probabilities at the specified nodes to be stored.
Note
Creates a plot of the node probablities over time.
Parameters: |
|
---|
Important
watch() must be called prior to propagation in order for the probabilities at the specified nodes to be stored.
Note
Propagates the quantum walk for time t.
Parameters: |
|
---|---|
Keyword Arguments: | |
|
|
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.
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()
Creates a handle that watches either marginal probability or walker entanglement during propagation.
Parameters: |
|
---|---|
Keyword Arguments: | |
: If ``watchtype=’entanglement’``, EigSolver keywords can also be passed; – for more details of the available EigSolver properties, see propagate(). |
|
Returns: | : * if watchtype='prob', creates a handle that can be accessed to retrieve marginal node probabilities for various \(t\) * if watchtype='entanglment', creates a handle that can be accessed to retrieve entanglement values for various \(t\) |
Return type: |
|
Examples
>>> walk.watch(None,watchtype='entanglement')
>>> walk.propagate(5.,method='chebyshev')
>>> timeArray, entArray = walk.entanglementHandle.getEntanglement()
Note
>>> walk.watch([0,1,4])
>>> walk.propagate(2.,method='chebyshev')
>>> timeArray, probXArray, probYArray = walk.handle.getLocalNode(4,p=2)
Warning
Note that walk.handle attributes are not collective; if running on multiple nodes, only local values will be returned.