diff --git a/kwant/physics/symmetry.py b/kwant/physics/symmetry.py index 900c81cf401244f06c6d045ca44f1bb02f5865f7..9bcfee7e509a49b51f0a16b0e9a13766b2b1b5f8 100644 --- a/kwant/physics/symmetry.py +++ b/kwant/physics/symmetry.py @@ -141,10 +141,11 @@ class DiscreteSymmetry: Returns ------- - broken_symmetry : string or ``None`` - One of "Conservation law", "Time reversal", "Particle-hole", - "Chiral": the symmetry broken by the matrix. If the matrix breaks - more than one symmetry, returns only the first failed check. + broken_symmetries : list + List of strings, the names of symmetries broken by the + matrix: any combination of "Conservation law", "Time reversal", + "Particle-hole", "Chiral". If no symmetries are broken, returns + an empty list. """ # Extra transposes are to enforse sparse dot product in case matrix is # dense. @@ -157,19 +158,22 @@ class DiscreteSymmetry: matrix = new_matrix else: matrix = hstack([matrix, csr_matrix((n, n-m))]) + broken_symmetries = [] if self.projectors is not None: for proj in self.projectors: full = proj.dot(proj.T.conj()) commutator = full.dot(matrix) - (full.T.dot(matrix.T)).T if np.linalg.norm(commutator.data) > 1e-8: - return 'Conservation law' + broken_symmetries.append('Conservation law') + break for symm, conj, sign, name in zip(self[1:], _conj, _signs, _names): if symm is None: continue commutator = symm.T.conj().dot((symm.T.dot(matrix.T)).T) commutator = commutator - sign * cond_conj(matrix, conj) if np.linalg.norm(commutator.data) > 1e-8: - return name + broken_symmetries.append(name) + return broken_symmetries def __getitem__(self, item): return (self.projectors, self.time_reversal, diff --git a/kwant/physics/tests/test_symmetry.py b/kwant/physics/tests/test_symmetry.py index 79f96f09966aea4bf9eb5532b4c8adb4bd7ec82b..1fe225a905bec5fb6970b237e1bb8fb631014d12 100644 --- a/kwant/physics/tests/test_symmetry.py +++ b/kwant/physics/tests/test_symmetry.py @@ -151,22 +151,22 @@ def test_validate(): csr = sparse.csr_matrix sym = DiscreteSymmetry(projectors=[csr(np.array([[1], [0]])), csr(np.array([[0], [1]]))]) - assert sym.validate(csr(np.array([[0], [1]]))) == 'Conservation law' - assert sym.validate(np.array([[1], [0]])) is None - assert sym.validate(np.eye(2)) is None - assert sym.validate(1 - np.eye(2)) == 'Conservation law' + assert sym.validate(csr(np.array([[0], [1]]))) == ['Conservation law'] + assert sym.validate(np.array([[1], [0]])) == [] + assert sym.validate(np.eye(2)) == [] + assert sym.validate(1 - np.eye(2)) == ['Conservation law'] sym = DiscreteSymmetry(particle_hole=sparse.identity(2)) - assert sym.validate(1j * sparse.identity(2)) is None - assert sym.validate(sparse.identity(2)) == 'Particle-hole' + assert sym.validate(1j * sparse.identity(2)) == [] + assert sym.validate(sparse.identity(2)) == ['Particle-hole'] sym = DiscreteSymmetry(time_reversal=sparse.identity(2)) - assert sym.validate(sparse.identity(2)) is None - assert sym.validate(1j * sparse.identity(2)) == 'Time reversal' + assert sym.validate(sparse.identity(2)) == [] + assert sym.validate(1j * sparse.identity(2)) == ['Time reversal'] sym = DiscreteSymmetry(chiral=csr(np.diag((1, -1)))) - assert sym.validate(np.eye(2)) == 'Chiral' - assert sym.validate(1 - np.eye(2)) is None + assert sym.validate(np.eye(2)) == ['Chiral'] + assert sym.validate(1 - np.eye(2)) == [] def random_onsite_hop(n, rng=0): @@ -178,7 +178,13 @@ def random_onsite_hop(n, rng=0): def test_validate_commutator(): - symm_class = ['AI', 'AII', 'D', 'C', 'AIII'] + symm_class = ['AI', 'AII', 'D', 'C', 'AIII', 'BDI'] + sym_dict = {'AI': ['Time reversal'], + 'AII': ['Time reversal'], + 'D': ['Particle-hole'], + 'C': ['Particle-hole'], + 'AIII': ['Chiral'], + 'BDI': ['Time reversal', 'Particle-hole', 'Chiral']} n = 10 rng = 10 for sym in symm_class: @@ -201,6 +207,7 @@ def test_validate_commutator(): disc_symm = DiscreteSymmetry(particle_hole=p_mat, time_reversal=t_mat, chiral=c_mat) - assert disc_symm.validate(h) == None + assert disc_symm.validate(h) == [] a = random_onsite_hop(n, rng=rng)[1] - assert disc_symm.validate(a) in ['Time reversal', 'Particle-hole', 'Chiral'] + for symmetry in disc_symm.validate(a): + assert symmetry in sym_dict[sym] diff --git a/kwant/system.py b/kwant/system.py index dfa28558d0efdb08a7ef8af788299e474e6b2759..4b1edc67a2e848f2f060375c2822cfcd9760d4d1 100644 --- a/kwant/system.py +++ b/kwant/system.py @@ -11,6 +11,7 @@ __all__ = ['System', 'FiniteSystem', 'InfiniteSystem'] import abc +import warnings from copy import copy from . import _system @@ -163,6 +164,17 @@ class FiniteSystem(System, metaclass=abc.ABCMeta): result.leads = new_leads return result + def validate_symmetries(self, args=(), *, params=None): + """Check that the Hamiltonian satisfies discrete symmetries. + + Applies `~kwant.physics.DiscreteSymmetry.validate` to the + Hamiltonian, see its documentation for details on the return + format. + """ + symmetries = self.discrete_symmetry(args=args, params=params) + ham = self.hamiltonian_submatrix(args, sparse=True, params=params) + return symmetries.validate(ham) + class InfiniteSystem(System, metaclass=abc.ABCMeta): """Abstract infinite low-level system. @@ -235,12 +247,19 @@ class InfiniteSystem(System, metaclass=abc.ABCMeta): ham = self.cell_hamiltonian(args, params=params) hop = self.inter_cell_hopping(args, params=params) symmetries = self.discrete_symmetry(args, params=params) - broken = symmetries.validate(ham) - if broken is not None: - raise ValueError("Cell Hamiltonian breaks " + broken.lower()) - broken = symmetries.validate(hop) - if broken is not None: - raise ValueError("Inter-cell hopping breaks " + broken.lower()) + # Check whether each symmetry is broken. + # If a symmetry is broken, it is ignored in the computation. + broken = set(symmetries.validate(ham) + symmetries.validate(hop)) + attribute_names = {'Conservation law': 'projectors', + 'Time reversal': 'time_reversal', + 'Particle-hole': 'particle-hole', + 'Chiral': 'chiral'} + for name in broken: + warnings.warn('Hamiltonian breaks ' + name + + ', ignoring the symmetry in the computation.') + assert name in attribute_names, 'Inconsistent naming of symmetries' + setattr(symmetries, attribute_names[name], None) + shape = ham.shape assert len(shape) == 2 assert shape[0] == shape[1] @@ -269,6 +288,19 @@ class InfiniteSystem(System, metaclass=abc.ABCMeta): return physics.selfenergy(ham, self.inter_cell_hopping(args, params=params)) + def validate_symmetries(self, args=(), *, params=None): + """Check that the Hamiltonian satisfies discrete symmetries. + + Returns `~kwant.physics.DiscreteSymmetry.validate` applied + to the onsite matrix and the hopping. See its documentation for + details on the return format. + """ + symmetries = self.discrete_symmetry(args=args, params=params) + ham = self.cell_hamiltonian(args=args, sparse=True, params=params) + hop = self.inter_cell_hopping(args=args, sparse=True, params=params) + broken = set(symmetries.validate(ham) + symmetries.validate(hop)) + return list(broken) + class PrecalculatedLead: def __init__(self, modes=None, selfenergy=None):