The Quantum Exact Simulation Toolkit v4.3.0
Loading...
Searching...
No Matches

Functions for approximate dynamical simulation. More...

Functions

void applyTrotterizedImaginaryTimeEvolution (Qureg qureg, PauliStrSum hamil, qreal tau, int order, int reps, bool permuteTerms)
 
void applyTrotterizedNoisyTimeEvolution (Qureg qureg, PauliStrSum hamil, qreal *damps, PauliStrSum *jumps, int numJumps, qreal time, int order, int reps, bool permuteTerms)
 
void applyTrotterizedUnitaryTimeEvolution (Qureg qureg, PauliStrSum hamil, qreal time, int order, int reps, bool permuteTerms)
 

Detailed Description

Functions for approximate dynamical simulation.

Function Documentation

◆ applyTrotterizedImaginaryTimeEvolution()

void applyTrotterizedImaginaryTimeEvolution ( Qureg qureg,
PauliStrSum hamil,
qreal tau,
int order,
int reps,
bool permuteTerms )

Simulates imaginary-time evolution of qureg for the duration tau under the time-independent Hamiltonian hamil, as approximated by symmetrized Trotterisation of the specified order and number of cycles reps.

Important
This is a non-physical operation and breaks the normalisation of state which can be restored via setQuregToRenormalized().
Formulae

Let \( \hat{H} = \) hamil and \( \tau = \) tau \( \in \mathbb{R} \). This function approximates the action of the non-unitary imaginary-time propagator

\[ \hat{V}(\tau) = \exp \left(- \tau \, \hat{H} \right), \]

as prescribed by Wick rotating (substituting time \( t \) for \( t \rightarrow -\iu \tau \)) the time-independent Schrödinger equation. When qureg is a statevector \( \svpsi \), the resulting state approximates

\[ \approx V(\tau) \svpsi \]

while when qureg is a density matrix \( \dmrho \), the result approximates

\[ \approx V(\tau) \, \dmrho \, V(\tau)^\dagger. \]

‍See applyTrotterizedPauliStrSumGadget() for information about the Trotter method.


Utility

Imaginary-time evolution drives the system toward the (unnormalised) groundstate of the Hamiltonian. Let \( \{ \ket{\phi_i} \} \) and \( \{ \ket{\lambda_i} \} \) be the eigenstates and respective eigenvalues of \( \hat{H} \), which are real due to Hermiticity.

\[ \hat{H} = \sum \limits_i \lambda_i \ket{\phi_i}\bra{\phi_i}, \;\;\;\;\; \lambda_i \in \mathbb{R}. \]

  • When qureg is a statevector \( \svpsi \) and can ergo be expressed in the basis of \( \{ \ket{\phi_i} \} \) as \( \svpsi = \sum_i \alpha_i \ket{\phi_i} \), this function approximates

    \[ \svpsi \, \rightarrow \, \hat{V}(\tau) \svpsi = \sum\limits_i \alpha_i \exp(- \tau \, \lambda_i) \ket{\phi_i}. \]

  • When qureg is a density matrix and is ergo expressible as \( \dmrho = \sum\limits_{ij} \alpha_{ij} \ket{\phi_i}\bra{\phi_j} \), this function effects

    \[ \dmrho \, \rightarrow \, \hat{V}(\tau) \dmrho \hat{V}(\tau)^\dagger = \sum\limits_{ij} \alpha_{ij} \exp(-\tau (\lambda_i + \lambda_j)) \ket{\phi_i}\bra{\phi_j}. \]

As \( \tau \rightarrow \infty \), the resulting unnormalised state approaches statevector \( \svpsi \rightarrow \alpha_0 \exp(-\tau \lambda_0) \ket{\phi_0} \) or density matrix \( \dmrho \rightarrow \alpha_{0,0} \exp(-2 \tau \lambda_0) \ket{\phi_0}\bra{\phi_0} \), where \( \lambda_0 \) is the minimum eigenvalue and \( \ket{\phi_0} \) is the groundstate. Assuming the initial overlap \( \alpha_0 \) is not zero (or exponentially tiny), subsequent renormalisation via setQuregToRenormalized() produces the pure ground-state \( \ket{\phi_0} \) or \( \ket{\phi_0}\bra{\phi_0} \).

Note degenerate minimum eigenvalues will yield a pure superposition of the corresponding eigenstates, with coefficients informed by the initial, relative populations.

Equivalences
Constraints
  • While the process of imaginary-time evolution is non-unitary (and non-physical), Hermiticity of hamil is still assumed, requiring it contains only real coefficients. Validation will check that hamil is approximately Hermitian, permitting coefficients with imaginary components smaller (in magnitude) than epsilon.

    \[ \max\limits_{i} |c_i| \le \valeps \]

    where the validation epsilon \( \valeps \) can be adjusted with setQuESTValidationEpsilon(). Beware however that imaginary-time evolution under a non-Hermitian Hamiltonian will not necessarily approach the lowest lying eigenstate (the eigenvalues may be non-real) so is likely of limited utility.
  • The tau parameter is necessarily real such that evolution approaches the groundstate (modulo renormalisation). It can generalised to an arbitrary complex number through direct use of applyTrotterizedNonUnitaryPauliStrSumGadget().
  • Simulation is exact such that the effected operation is precisely \( \exp(-\tau \hat{H}) \) only when reps \( \rightarrow \infty \) or all terms in hamil commute with one another.
Example
// pray for a non-zero initial overlap
initRandomPureState(qureg); // works even for density matrices
// minimize then renormalise
qreal tau = 10; // impatient infinity
int order = 4;
int reps = 100;
applyTrotterizedImaginaryTimeEvolution(qureg, hamil, tau, order, reps);
// ground-state (phi_0)
reportQureg(qureg);
// lowest lying eigenvalue (lambda_0)
qreal expec = calcExpecPauliStrSum(qureg, hamil);
reportScalar("expec", expec);
qreal calcExpecPauliStrSum(Qureg qureg, PauliStrSum sum)
qreal setQuregToRenormalized(Qureg qureg)
void initRandomPureState(Qureg qureg)
void reportQureg(Qureg qureg)
Definition qureg.cpp:384
void applyTrotterizedImaginaryTimeEvolution(Qureg qureg, PauliStrSum hamil, qreal tau, int order, int reps, bool permuteTerms)
void reportScalar(const char *label, qcomp num)
Definition types.cpp:55
See also
Parameters
[in,out]quregthe state to modify.
[in]hamilthe Hamiltonian as a a weighted sum of Pauli strings.
[in]tauthe duration over which to simulate imaginary-time evolution.
[in]orderthe order of the Trotter-Suzuki decomposition (e.g. 1, 2, 4, ...).
[in]repsthe number of Trotter repetitions.
[in]permuteTermswhether to randomly reorder Pauli terms at each repetition.
Exceptions
error
  • if qureg or hamil are uninitialised.
  • if hamil contains non-identities on qubits beyond the size of qureg.
  • if hamil is not approximately Hermitian.
  • if order is not 1 nor a positive, even integer.
  • if reps is not a positive integer.
  • if internal allocation needed for term permutation fails.
Author
Tyson Jones
Vasco Ferreira (randomisation)

Definition at line 272 of file trotterisation.cpp.

272 {
273 validate_quregFields(qureg, __func__);
274 validate_pauliStrSumFields(hamil, __func__);
275 validate_pauliStrSumTargets(hamil, qureg, __func__);
276 validate_pauliStrSumIsHermitian(hamil, __func__);
277 validate_trotterParams(order, reps, __func__);
278
279 // exp(-tau H) = exp(x i H) | x=tau*i
280 qcomp angle = qcomp(0, tau);
281 bool onlyLeftApply = false;
282 internal_applyAllTrotterRepetitions(qureg, nullptr, nullptr, 0, hamil, angle, order, reps, onlyLeftApply, permuteTerms, __func__);
283}

◆ applyTrotterizedNoisyTimeEvolution()

void applyTrotterizedNoisyTimeEvolution ( Qureg qureg,
PauliStrSum hamil,
qreal * damps,
PauliStrSum * jumps,
int numJumps,
qreal time,
int order,
int reps,
bool permuteTerms )
Warning
This function has not yet been unit tested and may contain bugs. Please use with caution!

Simulates open dynamics of qureg as per the Lindblad master equation, under the time-independent Hamiltonian hamil and jump operators jumps with corresponding damping rates damps, with evolution approximated by symmetrized Trotterisation of the specified order and number of cycles reps.

Note the ordering of all passed PauliStrSum (through functions like sortPauliStrSumMagnitude()) will affect that of the internally created super-propagator and ergo the Trotter accuracy. This is overridden by passing permuteTerms=true, whereby the super-propagator order is randomised every Trotter repetition.

Formulae

Let \( \rho = \) qureg, \( \hat{H} = \) hamil, \( t = \) time, and denote the \( i \)-th element of damps and jumps as \( \gamma_i \) and \( \hat{J}_i \) respectively. The Lindblad master equation prescribes that \( \rho \) time-evolves according to

\[ \frac{\mathrm{d}}{\mathrm{d}t} \rho = -\iu [\hat{H}, \rho] + \sum\limits_i \gamma_i \left( \hat{J}_i \rho \hat{J}_i^\dagger - \frac{1}{2} \left\{ \hat{J}_i^\dagger \hat{J}_i, \rho \right\} \right). \]

This function works by building a superoperator of the right-hand-side which acts upon the space of linearised \(\rho\),

\[ \boldsymbol{L} = -\iu \left( \hat{\id} \otimes \hat{H} - \hat{H}^* \otimes \hat{\id} \right) + \sum\limits_i \gamma_i \left( \hat{J}_i^* \otimes \hat{J}_i - \frac{1}{2} \hat{\id} \otimes (\hat{J}^\dagger J_i) - \frac{1}{2} (\hat{J}^\dagger J_i)^* \otimes \hat{\id} \right), \]

as a non-Hermitian weighted sum of Pauli strings (a PauliStrSum). The superoperator \( \boldsymbol{L} \) informs a superpropagator which exactly solves evolution as:

\[ \ket{\rho(t)} = \exp\left( t \boldsymbol{L} \right) \ket{\rho(0)}. \]

This function approximates the superpropagator \( \exp\left( t \boldsymbol{L} \right) \) using a higher-order symmetrized Suzuki-Trotter decomposition, as informed by parameters order and reps.

‍See applyTrotterizedPauliStrSumGadget() for information about the Trotter method.


Utility

This function simulates time evolution of an open system, where the jump operators model interactions with the environment. This can capture sophisticated decoherence processes of the quantum state which are untenable to model as discrete operations with functions like mixKrausMap(). This function also proves useful for preparing realistic, physical input states to quantum metrological circuits, or the general high-performance simulation of digital time evolution of condensed matter systems.

Equivalences
  • When numJumps = 0, evolution is unitary and the Lindblad master equation simplifes to the Liouville–von Neumann equation, which is equivalently (and more efficiently) simulated via applyTrotterizedUnitaryTimeEvolution().
Constraints
  • Each damping rate in damps is expected to be a zero or positive number, in order for evolution to be trace preserving. Validation will assert that each damping rate \( \gamma_i \) satisfies

    \[ \min\limits_{i} \gamma_i \ge - \valeps \]

    where the validation epsilon \( \valeps \) can be adjusted with setQuESTValidationEpsilon(). Non-trace-preserving, negative damping rates can be simulated by disabling numerical validation via setQuESTValidationEpsilon(0).
  • The time parameter is necessarily real, and cannot be generalised to imaginary or complex like in other functions. Generalisation is trivially numerically possible, but has no established physical meaning and so is not exposed in the API. Please open an issue on Github for advice on complex-time simulation.
  • Simulation is exact only when reps \( \rightarrow \infty \) or all terms in the superoperator \( \boldsymbol{L} \) incidentally commute with one another, and otherwise incorporates Trotter error. Unlike for unitary evolution, Trotter error does break normalisation of the state and so this function is generally non-trace-preserving. In theory, normalisation can be restored with setQuregToRenormalized() though noticable norm-breaking indicates evolution was inaccurate, and should instead be repeated with increased order or reps parameters.
  • The function instantiates superoperator \( \boldsymbol{L} \) above as a temporary PauliStrSum, incurring a memory and time overhead which grows quadratically with the number of terms in hamil, plus quadratically with the number in each jump operator. These overheads may prove prohibitively costly for PauliStrSum containing very many terms.
Example
// |+><+|
1 IIX
2 IYI
3 ZZZ
)");
// |0><0|
0.5 I
0.5 Z
)");
// |1><0|
0.5 X
-0.5i Y
)");
// "noisiness"
qreal damps[] = {.3, .4};
PauliStrSum jumps[] = {jump1, jump2};
int numJumps = 2;
reportScalar("initial energy", calcExpecPauliStrSum(qureg, hamil));
// time and accuracy
qreal time = 0.5;
int order = 4;
int reps = 100;
applyTrotterizedNoisyTimeEvolution(qureg, hamil, damps, jumps, numJumps, time, order, reps);
reportScalar("final energy", calcExpecPauliStrSum(qureg, hamil));
void initPlusState(Qureg qureg)
PauliStrSum createInlinePauliStrSum(const char *str)
Definition paulis.cpp:198
Qureg createDensityQureg(int numQubits)
Definition qureg.cpp:297
void applyTrotterizedNoisyTimeEvolution(Qureg qureg, PauliStrSum hamil, qreal *damps, PauliStrSum *jumps, int numJumps, qreal time, int order, int reps, bool permuteTerms)
Definition qureg.h:49
See also
Parameters
[in,out]quregthe density-matrix state to evolve and modify.
[in]hamilthe Hamiltonian of the qubit system (excludes any environment).
[in]dampsthe damping rates of each jump operator in jumps.
[in]jumpsthe jump operators specified as PauliStrSum.
[in]numJumpsthe length of list jumps (and damps).
[in]timethe duration through which to evolve the state.
[in]orderthe order of the Trotter-Suzuki decomposition (e.g. 1, 2, 4, ...).
[in]repsthe number of Trotter repetitions.
[in]permuteTermswhether to randomly reorder Pauli terms at each repetition.
Exceptions
error
  • if qureg, hamil or any element of jumps are uninitialised.
  • if qureg is not a density matrix.
  • if hamil or any element of jumps contains non-identities on qubits beyond the size of qureg.
  • if hamil is not approximately Hermitian.
  • if numJumps is negative.
  • if any element of damps is not approximately positive.
  • if the total number of Lindbladian superoperator terms overflows the qindex type.
  • if all Lindbladian superoperator terms cannot simultaneously fit into CPU memory.
  • if memory allocation of the Lindbladian superoperator terms unexpectedly fails.
  • if order is not 1 nor a positive, even integer.
  • if reps is not a positive integer.
  • if internal allocation needed for term permutation fails.
Author
Tyson Jones
Vasco Ferreira (randomisation)

Definition at line 295 of file trotterisation.cpp.

298 {
299 validate_quregFields(qureg, __func__);
300 validate_quregIsDensityMatrix(qureg, __func__);
301 validate_pauliStrSumFields(hamil, __func__);
302 validate_pauliStrSumTargets(hamil, qureg, __func__);
303 validate_pauliStrSumIsHermitian(hamil, __func__);
304 validate_trotterParams(order, reps, __func__);
305 validate_lindbladJumpOps(jumps, numJumps, qureg, __func__);
306 validate_lindbladDampingRates(damps, numJumps, __func__);
307
308 qindex numSuperTerms = internal_getNumTotalSuperPropagatorTerms(hamil, jumps, numJumps); // 0 indicates overflow
309 validate_numLindbladSuperPropagatorTerms(numSuperTerms, __func__);
310
311 // validate memory allocations for all super-propagator terms
312 vector<PauliStr> superStrings;
313 vector<qcomp> superCoeffs;
314 auto callbackString = [&]() { validate_tempListAllocSucceeded(false, numSuperTerms, sizeof(PauliStr), __func__); };
315 auto callbackCoeff = [&]() { validate_tempListAllocSucceeded(false, numSuperTerms, sizeof(qcomp), __func__); };
316 util_tryAllocVector(superStrings, numSuperTerms, callbackString);
317 util_tryAllocVector(superCoeffs, numSuperTerms, callbackCoeff);
318
319 qindex superTermInd = 0;
320
321 // collect -i[H,rho] terms
322 for (qindex n=0; n<hamil.numTerms; n++) {
323 PauliStr oldStr = hamil.strings[n];
324 qcomp oldCoeff = hamil.coeffs[n];
325
326 // term of -i Id (x) H
327 superStrings[superTermInd] = oldStr;
328 superCoeffs [superTermInd] = -1_i * oldCoeff;
329 superTermInd++;
330
331 // term of i conj(H) (x) I
332 superStrings[superTermInd] = paulis_getShiftedPauliStr(oldStr, qureg.numQubits);
333 superCoeffs [superTermInd] = 1_i * paulis_getSignOfPauliStrConj(oldStr) * std::conj(oldCoeff);
334 superTermInd++;
335 }
336
337 // below we bind superStrings/Coeffs to a spoofed PauliStrSum to pass to paulis functions
338 PauliStrSum temp;
339 int flagForDebugSafety = -1;
340 temp.isApproxHermitian = &flagForDebugSafety;
341
342 // collect jump terms
343 for (int n=0; n<numJumps; n++) {
344
345 // damp conj(J) (x) J
346 temp.strings = &superStrings[superTermInd];
347 temp.coeffs = &superCoeffs[superTermInd];
348 temp.numTerms = jumps[n].numTerms * jumps[n].numTerms;
349 superTermInd += temp.numTerms;
350 paulis_setPauliStrSumToScaledTensorProdOfConjWithSelf(temp, damps[n], jumps[n], qureg.numQubits);
351
352 // -damp/2 I (x) (adj(J) . J)
353 temp.strings = &superStrings[superTermInd];
354 temp.coeffs = &superCoeffs[superTermInd];
355 temp.numTerms = paulis_getNumTermsInPauliStrSumProdOfAdjointWithSelf(jumps[n]);
356 superTermInd += temp.numTerms;
357 paulis_setPauliStrSumToScaledProdOfAdjointWithSelf(temp, -damps[n]/2, jumps[n]);
358
359 // -damp/2 conj(adj(J) . J) (x) I = conj(above) when damp is real
360 PauliStrSum temp2;
361 temp2.strings = &superStrings[superTermInd];
362 temp2.coeffs = &superCoeffs[superTermInd];
363 temp2.numTerms = temp.numTerms;
364 superTermInd += temp2.numTerms;
365 paulis_setPauliStrSumToShiftedConj(temp2, temp, qureg.numQubits);
366 }
367
368 // defensively check we didn't write too few (or too many, though that'd segfault
369 // above) Lindblad terms, in case the above code changes when jump ops are generalised
370 if (superTermInd != numSuperTerms)
371 error_unexpectedNumLindbladSuperpropTerms();
372
373 // pass superpropagator terms as temporary PauliStrSum
374 PauliStrSum superSum;
375 superSum.numTerms = numSuperTerms;
376 superSum.strings = superStrings.data();
377 superSum.coeffs = superCoeffs.data();
378 superSum.isApproxHermitian = nullptr; // will not be queried
379
380 // effect exp(t S) = exp(x i S) | x=-i*time, left-multiplying only
381 qcomp angle = qcomp(0, -time);
382 bool onlyLeftApply = true;
383 internal_applyAllTrotterRepetitions(qureg, nullptr, nullptr, 0, superSum, angle, order, reps, onlyLeftApply, permuteTerms, __func__);
384}

◆ applyTrotterizedUnitaryTimeEvolution()

void applyTrotterizedUnitaryTimeEvolution ( Qureg qureg,
PauliStrSum hamil,
qreal time,
int order,
int reps,
bool permuteTerms )

Unitarily time evolves qureg for the duration time under the time-independent Hamiltonian hamil, as approximated by symmetrized Trotterisation of the specified order and number of cycles reps.

Formulae

Let \( \hat{H} = \) hamil and \( t = \) time \( \in \mathbb{R} \). This function approximates the action of the unitary-time evolution operator/propagator

\[ \hat{U}(t) = \exp \left(- \iu \, t \, \hat{H} \right), \]

as solves the time-independent Schrödinger equation. When qureg is a statevector \( \svpsi \), the resulting state approximates

\[ \approx U(t) \svpsi \]

while when qureg is a density matrix \( \dmrho \), the result approximates

\[ \approx U(t) \, \dmrho \, U(t)^\dagger. \]

‍See applyTrotterizedPauliStrSumGadget() for information about the Trotter method.


Equivalences
Constraints
  • Unitarity requires that hamil is Hermitian and ergo contains only real coefficients. Validation will check that hamil is approximately Hermitian, permitting coefficients with imaginary components smaller (in magnitude) than epsilon.

    \[ \max\limits_{i} |c_i| \le \valeps \]

    where the validation epsilon \( \valeps \) can be adjusted with setQuESTValidationEpsilon(). The imaginary components of the Hamiltonian are considered during simulation.
  • The time parameter is necessarily real to retain unitarity. It can be substituted for a strictly imaginary scalar to perform imaginary-time evolution (as per Wick rotation \( t \rightarrow - \iu \tau \)) via applyTrotterizedImaginaryTimeEvolution(), or generalised to an arbitrary complex number through direct use of applyTrotterizedNonUnitaryPauliStrSumGadget().
  • The simulated system is closed with dynamics described fully by the Hamiltonian hamil. Open or otherwise noisy system dynamics can be simulated with applyTrotterizedNoisyTimeEvolution().
  • Simulation is exact such that the effected operation is precisely \( \exp(-\iu t \hat{H}) \) only when reps \( \rightarrow \infty \) or all terms in hamil commute with one another. Conveniently, Trotter error does not break normalisation of the state since the approximating circuit remains unitary.
Example
1 ZZI
2 IZZ
3 ZIZ
1.5 XII
2.5 IXI
3.5 IIX
)");
qreal time = 0.8 * hbar;
int order = 4;
int reps = 20;
applyTrotterizedUnitaryTimeEvolution(qureg, hamil, time, order, reps);
void applyTrotterizedUnitaryTimeEvolution(Qureg qureg, PauliStrSum hamil, qreal time, int order, int reps, bool permuteTerms)
See also
Parameters
[in,out]quregthe state to modify.
[in]hamilthe Hamiltonian as a a weighted sum of Pauli strings.
[in]timethe duration over which to simulate evolution.
[in]orderthe order of the Trotter-Suzuki decomposition (e.g. 1, 2, 4, ...).
[in]repsthe number of Trotter repetitions.
[in]permuteTermswhether to randomly reorder Pauli terms at each repetition.
Exceptions
error
  • if qureg or hamil are uninitialised.
  • if hamil contains non-identities on qubits beyond the size of qureg.
  • if hamil is not approximately Hermitian.
  • if order is not 1 nor a positive, even integer.
  • if reps is not a positive integer.
  • if internal allocation needed for term permutation fails.
Author
Tyson Jones
Vasco Ferreira (randomisation)

Definition at line 259 of file trotterisation.cpp.

259 {
260 validate_quregFields(qureg, __func__);
261 validate_pauliStrSumFields(hamil, __func__);
262 validate_pauliStrSumTargets(hamil, qureg, __func__);
263 validate_pauliStrSumIsHermitian(hamil, __func__);
264 validate_trotterParams(order, reps, __func__);
265
266 // exp(-i t H) = exp(x i H) | x=-t
267 qcomp angle = - time;
268 bool onlyLeftApply = false;
269 internal_applyAllTrotterRepetitions(qureg, nullptr, nullptr, 0, hamil, angle, order, reps, onlyLeftApply, permuteTerms, __func__);
270}