-
Notifications
You must be signed in to change notification settings - Fork 624
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
HEOM Parity support #2261
base: master
Are you sure you want to change the base?
HEOM Parity support #2261
Conversation
… created an utils class
@mcditoos Thanks more making the PR. Would you mind reverting all the unrelated changes? We don't mind small clean ups, but big clean ups should rather go into separate PRs and it should be clear why they are being made. Once you're ready for this to be reviewed, remove the draft status and ping me on the PR. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thank you for removing all the unnecessary changes. The parity support changes are much clearer now.
I've request some small changes. In addition we will need a test for the new feature and ideally for the feature to be mentioned in the HEOM guide in the documentation (perhaps also with a small example).
Co-authored-by: Simon Cross <hodgestar+github@gmail.com>
Co-authored-by: Simon Cross <hodgestar+github@gmail.com>
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for the changes. Next round of reviewing done.
qutip/solver/heom/bofin_solvers.py
Outdated
sign1 = (-1) ** (n_excite + 1 - self.odd_parity) | ||
sign2 = (-1)**(n_excite_before_m + self.odd_parity) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could you put these lines back where they were so that there is a small block of code for the calculation of each sign instead of interleaving the two calculations?
n_excite_before_m = sum(he_fermionic_n[:k]) | ||
sign2 = (-1) ** (n_excite_before_m) | ||
sign1 = (-1) ** (n_excite + 1 - self.odd_parity) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Same here.
qutip/solver/heom/bofin_solvers.py
Outdated
@@ -1281,6 +1293,7 @@ class _GatherHEOMRHS: | |||
nhe : int | |||
The number of ADOs in the hierarchy. | |||
""" | |||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could you again undo these unrelated changes?
@@ -54,8 +54,8 @@ def assert_raises_steady_state_time_dependent(hsolver): | |||
with pytest.raises(ValueError) as err: | |||
hsolver.steady_state() | |||
assert str(err.value) == ( | |||
"A steady state cannot be determined for a time-dependent" | |||
" system" | |||
"A steady state cannot be determined for a time-dependent" |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This can stay -- it's a valid indentation fix.
@@ -232,6 +232,7 @@ class DrudeLorentzPureDephasingModel: | |||
""" Analytic Drude-Lorentz pure-dephasing model for testing the HEOM | |||
solver. | |||
""" | |||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Remove these too.
@@ -9,7 +9,7 @@ | |||
|
|||
from qutip import ( | |||
basis, destroy, expect, liouvillian, qeye, sigmax, sigmaz, | |||
tensor, Qobj, QobjEvo | |||
tensor, Qobj, QobjEvo, fcreate, fdestroy |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
fcreate
appears to be unused. Could you remove it?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Sure, I shouldn't have included it
resultHEOMPade = HEOMSolver(L, [bath1, bath2], depth, odd_parity=True) | ||
rhoss, _ = resultHEOMPade.steady_state(use_mkl=False) | ||
rhoss = rhoss.full() | ||
expected_odd = np.diag([-0.18472, 0.68472, 0.68472, -0.18472]) | ||
expected = np.diag([0.10623, 0.39376, 0.39376, 0.10623]) | ||
assert np.isclose(rhoss, expected_odd,atol=1e-5).all() | ||
resultHEOMPade = HEOMSolver(L, [bath1, bath2], depth, odd_parity=False) | ||
rhoss, _ = resultHEOMPade.steady_state(use_mkl=False) | ||
rhoss = rhoss.full() | ||
assert np.isclose(rhoss, expected,atol=1e-5).all() |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
resultHEOMPade = HEOMSolver(L, [bath1, bath2], depth, odd_parity=True) | |
rhoss, _ = resultHEOMPade.steady_state(use_mkl=False) | |
rhoss = rhoss.full() | |
expected_odd = np.diag([-0.18472, 0.68472, 0.68472, -0.18472]) | |
expected = np.diag([0.10623, 0.39376, 0.39376, 0.10623]) | |
assert np.isclose(rhoss, expected_odd,atol=1e-5).all() | |
resultHEOMPade = HEOMSolver(L, [bath1, bath2], depth, odd_parity=False) | |
rhoss, _ = resultHEOMPade.steady_state(use_mkl=False) | |
rhoss = rhoss.full() | |
assert np.isclose(rhoss, expected,atol=1e-5).all() | |
solver = HEOMSolver(L, [bath1, bath2], depth, odd_parity=True) | |
rhoss, _ = solver.steady_state(use_mkl=False) | |
rhoss = rhoss.full() | |
expected_odd = np.diag([-0.18472, 0.68472, 0.68472, -0.18472]) | |
np.testing.assert_allclose(rhoss, expected_odd, atol=1e-5) | |
solver = HEOMSolver(L, [bath1, bath2], depth, odd_parity=False) | |
rhoss, _ = solver.steady_state(use_mkl=False) | |
rhoss = rhoss.full() | |
expected = np.diag([0.10623, 0.39376, 0.39376, 0.10623]) | |
np.testing.assert_allclose(rhoss, expected, atol=1e-5) |
|
||
|
||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
There appear to be too many blank lines here.
Q=d_2, gamma=2 * Gamma, w=W, mu=mu, T=1 / beta, Nk=Nk, | ||
tag="Lead 2") | ||
resultHEOMPade = HEOMSolver(L, [bath1, bath2], depth, odd_parity=True) | ||
rhoss, _ = resultHEOMPade.steady_state(use_mkl=False) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could you explain why you needed to pass use_mkl=False
here? I think we would expect it to automatically do something sensible for users in all cases and we don't appear to have had to use it in any of the other HEOM solver tests yet.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I used used_mkl=False, because if I don't for the test fails for python 3.10 with no cython (only in that setup it passes in all others), the test fails with Exception: Zero pivot, numerical factorization or iterative refinement problem. I found this old thread #791 suggesting that sometimes mkl solvers give some trouble so I decided to set it to false
This is normally prohibited in 'normal' electronic systems, where | ||
charge is conserved, but is useful for the construction of | ||
electron spectrum (under linear response assumptions). | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Simon asked us to help distill and clarify this a little bit; how about:
" For fermionic baths only. Default is "False". "Parity" refers to the parity of the initial system state used with the HEOM. An example of an odd parity state is one made from applying an odd number of fermionic creation operators to a physical density operator. Physical systems have even parity, but allowing the generalization to odd-parity states allows one to calculate useful physical quantities like the system power spectrum or density of states. The form of the HEOM differs depending on the parity of the initial system state, so if this option is set to "True", a different RHS is constructed, which can then be used with a system state of odd parity. "
Description
This PR adds odd parity support to HEOM's fermionic baths, the user can choose whether to use odd or even parity according to his needs when construncting the right hand side of HEOM's equations
Related issues or PRs
Heom tutorial for this feature qutip/qutip-tutorials#74