Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

217

218

219

220

221

222

223

224

225

226

227

228

229

230

231

232

233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255

256

257

258

259

260

261

262

263

264

265

266

267

268

269

270

271

272

273

274

275

276

277

278

279

280

281

282

283

284

285

286

287

288

289

290

291

292

293

294

295

296

297

298

299

300

301

302

303

304

305

306

307

308

309

310

311

312

313

314

315

316

317

318

319

320

321

322

323

324

325

326

327

328

329

330

331

332

333

334

335

336

337

338

339

340

341

342

343

344

345

346

347

348

349

350

351

352

353

354

355

356

357

358

359

360

361

362

363

364

365

366

367

368

369

370

371

372

373

374

375

376

377

378

379

380

381

382

383

384

385

386

387

388

389

390

391

392

393

394

395

396

397

398

399

400

401

402

403

404

405

406

407

408

409

410

411

412

413

414

415

416

417

418

419

420

421

422

423

424

425

426

427

428

429

430

431

432

433

434

435

436

437

438

439

440

441

442

443

444

445

446

447

448

449

450

451

452

453

454

# Copyright 2011-2019 Kwant authors. 

# 

# This file is part of Kwant. It is subject to the license terms in the file 

# LICENSE.rst found in the top-level directory of this distribution and at 

# http://kwant-project.org/license. A list of Kwant authors can be found in 

# the file AUTHORS.rst at the top-level directory of this distribution and at 

# http://kwant-project.org/authors. 

 

from keyword import iskeyword 

import functools 

import operator 

import collections 

 

from math import sqrt 

import tinyarray as ta 

import numpy as np 

import sympy 

 

import kwant.lattice 

import kwant.builder 

 

import kwant.continuum 

import kwant.continuum._common 

import kwant.continuum.discretizer 

from kwant.continuum import momentum_operators, position_operators 

 

__all__ = ["to_landau_basis", "discretize_landau"] 

 

coordinate_vectors = dict(zip("xyz", np.eye(3))) 

 

ladder_lower, ladder_raise = sympy.symbols(r"a a^\dagger", commutative=False) 

 

 

def to_landau_basis(hamiltonian, momenta=None): 

r"""Replace two momenta by Landau level ladder operators. 

 

Replaces: 

 

k_0 -> sqrt(B/2) * (a + a^\dagger) 

k_1 -> 1j * sqrt(B/2) * (a - a^\dagger) 

 

Parameters 

---------- 

hamiltonian : str or sympy expression 

The Hamiltonian to which to apply the Landau level transformation. 

momenta : sequence of str (optional) 

The momenta to replace with Landau level ladder operators. If not 

provided, 'k_x' and 'k_y' are used 

 

Returns 

------- 

hamiltonian : sympy expression 

momenta : sequence of sympy atoms 

The momentum operators that have been replaced by ladder operators. 

normal_coordinate : sympy atom 

The remaining position coordinate. May or may not be present 

in 'hamiltonian'. 

""" 

hamiltonian = kwant.continuum.sympify(hamiltonian) 

momenta = _normalize_momenta(momenta) 

normal_coordinate = _find_normal_coordinate(hamiltonian, momenta) 

 

# Substitute ladder operators for Landau level momenta 

B = sympy.symbols("B") 

hamiltonian = hamiltonian.subs( 

{ 

momenta[0]: sympy.sqrt(abs(B) / 2) * (ladder_raise + ladder_lower), 

momenta[1]: sympy.I 

* sympy.sqrt(abs(B) / 2) 

* (ladder_lower - ladder_raise), 

} 

) 

 

return hamiltonian, momenta, normal_coordinate 

 

 

def discretize_landau(hamiltonian, N, momenta=None, grid_spacing=1): 

"""Discretize a Hamiltonian in a basis of Landau levels. 

 

Parameters 

---------- 

hamiltonian : str or sympy expression 

N : positive integer 

The number of Landau levels in the basis. 

momenta : sequence of str (optional) 

The momenta defining the plane perpendicular to the magnetic field. 

If not provided, "k_x" and "k_y" are used. 

grid_spacing : float, default: 1 

The grid spacing to use when discretizing the normal direction 

(parallel to the magnetic field). 

 

Returns 

------- 

builder : `~kwant.builder.Builder` 

'hamiltonian' discretized in a basis of Landau levels in the plane 

defined by 'momenta'. If a third coordinate is present in 'hamiltonian', 

then this also has a translational symmetry in that coordinate direction. 

""" 

 

if N <= 0: 

raise ValueError("N must be positive") 

 

hamiltonian, momenta, normal_coordinate = to_landau_basis(hamiltonian, momenta) 

 

# Discretize in normal direction and split terms for onsites/hoppings into 

# monomials in ladder operators. 

tb_hamiltonian, _ = kwant.continuum.discretize_symbolic( 

hamiltonian, coords=[normal_coordinate.name] 

) 

tb_hamiltonian = { 

key: kwant.continuum._common.monomials(value, gens=(ladder_lower, ladder_raise)) 

for key, value in tb_hamiltonian.items() 

} 

 

# Replace ladder operator monomials by tuple of integers: 

# e.g. a^\dagger a^2 -> (+1, -2). 

tb_hamiltonian = { 

outer_key: { 

_ladder_term(inner_key): inner_value 

for inner_key, inner_value in outer_value.items() 

} 

for outer_key, outer_value in tb_hamiltonian.items() 

} 

 

# Construct map from LandauLattice HoppingKinds to a sequence of pairs 

# that encode the ladder operators and multiplying expression. 

tb_hoppings = collections.defaultdict(list) 

for outer_key, outer_value in tb_hamiltonian.items(): 

for inner_key, inner_value in outer_value.items(): 

tb_hoppings[(*outer_key, sum(inner_key))] += [(inner_key, inner_value)] 

# Extract the number of orbitals on each site/hopping 

random_element = next(iter(tb_hoppings.values()))[0][1] 

norbs = 1 if isinstance(random_element, sympy.Expr) else random_element.shape[0] 

tb_onsite = tb_hoppings.pop((0, 0), None) 

 

# Construct Builder 

if _has_coordinate(normal_coordinate, hamiltonian): 

sym = kwant.lattice.TranslationalSymmetry([grid_spacing, 0]) 

else: 

sym = kwant.builder.NoSymmetry() 

lat = LandauLattice(grid_spacing, norbs=norbs) 

syst = kwant.Builder(sym) 

 

# Add onsites 

landau_sites = (lat(0, j) for j in range(N)) 

if tb_onsite is None: 

syst[landau_sites] = ta.zeros((norbs, norbs)) 

else: 

syst[landau_sites] = _builder_value( 

tb_onsite, normal_coordinate.name, grid_spacing, is_onsite=True 

) 

 

# Add zero hoppings between adjacent Landau levels. 

# Necessary to be able to use the Landau level builder 

# to populate another builder using builder.fill(). 

syst[kwant.builder.HoppingKind((0, 1), lat)] = ta.zeros((norbs, norbs)) 

 

# Add the hoppings from the Hamiltonian 

for hopping, parts in tb_hoppings.items(): 

syst[kwant.builder.HoppingKind(hopping, lat)] = _builder_value( 

parts, normal_coordinate.name, grid_spacing, is_onsite=False 

) 

 

return syst 

 

 

# This has to subclass lattice so that it will work with TranslationalSymmetry. 

class LandauLattice(kwant.lattice.Monatomic): 

""" 

A `~kwant.lattice.Monatomic` lattice with a Landau level index per site. 

 

Site tags (see `~kwant.builder.SiteFamily`) are pairs of integers, where 

the first integer describes the real space position and the second the 

Landau level index. 

 

The real space Bravais lattice is one dimensional, oriented parallel 

to the magnetic field. The Landau level index represents the harmonic 

oscillator basis states used for the Landau quantization in the plane 

normal to the field. 

 

Parameters 

---------- 

grid_spacing : float 

Real space lattice spacing (parallel to the magnetic field). 

offset : float (optional) 

Displacement of the lattice origin from the real space 

coordinates origin. 

""" 

 

def __init__(self, grid_spacing, offset=None, name="", norbs=None): 

191 ↛ 192line 191 didn't jump to line 192, because the condition on line 191 was never true if offset is not None: 

offset = [offset, 0] 

# The second vector and second coordinate do not matter (they are 

# not used in pos()) 

super().__init__([[grid_spacing, 0], [0, 1]], offset, name, norbs) 

 

def pos(self, tag): 

return ta.array((self.prim_vecs[0, 0] * tag[0] + self.offset[0],)) 

 

def landau_index(self, tag): 

return tag[-1] 

 

 

def _builder_value(terms, normal_coordinate, grid_spacing, is_onsite): 

"""Construct an onsite/hopping value function from a list of terms 

 

Parameters 

---------- 

terms : list 

Each element is a pair (ladder_term, sympy expression/matrix). 

ladder_term is a tuple of integers that encodes a string of 

Landau raising/lowering operators and the sympy expression 

is the rest 

normal_coordinate : str 

grid_spacing : float 

The grid spacing in the normal direction 

is_onsite : bool 

True if we are constructing an onsite value function 

""" 

ladder_term_symbols = [sympy.Symbol(_ladder_term_name(lt)) for lt, _ in terms] 

# Construct a single expression from the terms, multiplying each part 

# by the placeholder that represents the prefactor from the ladder operator term. 

(ladder, (_, part)), *rest = zip(ladder_term_symbols, terms) 

expr = ladder * part 

for ladder, (_, part) in rest: 

expr += ladder * part 

expr = expr.subs( 

{kwant.continuum.discretizer._displacements[normal_coordinate]: grid_spacing} 

) 

# Construct the return string and temporary variable names 

# for function calls. 

return_string, map_func_calls, const_symbols, _cache = kwant.continuum.discretizer._return_string( 

expr, coords=[normal_coordinate] 

) 

 

# Remove the ladder term placeholders, as these are not parameters 

const_symbols = set(const_symbols) 

for ladder_term in ladder_term_symbols: 

const_symbols.discard(ladder_term) 

 

# Construct the argument part of signature. Arguments 

# consist of any constants and functions in the return string. 

arg_names = set.union( 

{s.name for s in const_symbols}, {str(k.func) for k in map_func_calls} 

) 

arg_names.discard("Abs") # Abs function is not a parameter 

for arg_name in arg_names: 

248 ↛ 249line 248 didn't jump to line 249, because the condition on line 248 was never true if not (arg_name.isidentifier() and not iskeyword(arg_name)): 

raise ValueError( 

"Invalid name in used symbols: {}\n" 

"Names of symbols used in Hamiltonian " 

"must be valid Python identifiers and " 

"may not be keywords".format(arg_name) 

) 

arg_names = ", ".join(sorted(arg_names)) 

# Construct site part of the function signature 

site_string = "from_site" if is_onsite else "to_site, from_site" 

# Construct function signature 

if arg_names: 

function_header = "def _({}, {}):".format(site_string, arg_names) 

else: 

function_header = "def _({}):".format(site_string) 

# Construct function body 

function_body = [] 

if "B" not in arg_names: 

# B is not a parameter for terms with no ladder operators but we still 

# need something to pass to _evaluate_ladder_term 

function_body.append("B = +1") 

function_body.extend( 

[ 

"{}, = from_site.pos".format(normal_coordinate), 

"_ll_index = from_site.family.landau_index(from_site.tag)", 

] 

) 

# To get the correct hopping if B < 0, we need to set the Hermitian 

# conjugate of the ladder operator matrix element, which swaps the 

# from_site and to_site Landau level indices. 

if not is_onsite: 

function_body.extend( 

["if B < 0:", 

" _ll_index = to_site.family.landau_index(to_site.tag)" 

] 

) 

function_body.extend( 

"{} = _evaluate_ladder_term({}, _ll_index, B)".format(symb.name, lt) 

for symb, (lt, _) in zip(ladder_term_symbols, terms) 

) 

function_body.extend( 

"{} = {}".format(v, kwant.continuum.discretizer._print_sympy(k)) 

for k, v in map_func_calls.items() 

) 

function_body.append(return_string) 

func_code = "\n ".join([function_header] + function_body) 

# Add "I" to namespace just in case sympy again would miss to replace it 

# with Python's 1j as it was the case with SymPy 1.2 when "I" was argument 

# of some function. 

namespace = { 

"pi": np.pi, 

"I": 1j, 

"_evaluate_ladder_term": _evaluate_ladder_term, 

"Abs": abs, 

} 

namespace.update(_cache) 

# Construct full source, including cached arrays 

source = [] 

for k, v in _cache.items(): 

source.append("{} = (\n{})\n".format(k, repr(np.array(v)))) 

source.append(func_code) 

exec(func_code, namespace) 

f = namespace["_"] 

f._source = "".join(source) 

 

return f 

 

 

def _ladder_term(operator_string): 

r"""Return a tuple of integers representing a string of ladder operators 

 

Parameters 

---------- 

operator_string : Sympy expression 

Monomial in ladder operators, e.g. a^\dagger a^2 a^\dagger. 

 

Returns 

------- 

ladder_term : tuple of int 

e.g. a^\dagger a^2 -> (+1, -2) 

""" 

ret = [] 

for factor in operator_string.as_ordered_factors(): 

ladder_op, exponent = factor.as_base_exp() 

if ladder_op == ladder_lower: 

sign = -1 

elif ladder_op == ladder_raise: 

sign = +1 

else: 

sign = 0 

ret.append(sign * int(exponent)) 

return tuple(ret) 

 

 

def _ladder_term_name(ladder_term): 

""" 

Parameters 

---------- 

ladder_term : tuple of int 

 

Returns 

------- 

ladder_term_name : str 

""" 

 

def ns(i): 

if i >= 0: 

return str(i) 

else: 

return "_" + str(-i) 

 

return "_ladder_{}".format("_".join(map(ns, ladder_term))) 

 

 

def _evaluate_ladder_term(ladder_term, n, B): 

r"""Evaluates the prefactor for a ladder operator on a landau level. 

 

Example: a^\dagger a^2 -> (n - 1) * sqrt(n) 

 

Parameters 

---------- 

ladder_term : tuple of int 

Represents a string of ladder operators. Positive 

integers represent powers of the raising operator, 

negative integers powers of the lowering operator. 

n : non-negative int 

Landau level index on which to act with ladder_term. 

B : float 

Magnetic field with sign 

 

Returns 

------- 

ladder_term_prefactor : float 

""" 

assert n >= 0 

# For negative B we swap a and a^\dagger. 

if B < 0: 

ladder_term = tuple(-i for i in ladder_term) 

ret = 1 

for m in reversed(ladder_term): 

if m > 0: 

factors = range(n + 1, n + m + 1) 

elif m < 0: 

factors = range(n + m + 1, n + 1) 

if n == 0: 

return 0 # a|0> = 0 

else: 

factors = (1,) 

ret *= sqrt(functools.reduce(operator.mul, factors)) 

n += m 

return ret 

 

 

def _normalize_momenta(momenta=None): 

"""Return Landau level momenta as Sympy atoms 

 

Parameters 

---------- 

momenta : None or pair of int or pair of str 

The momenta to choose. If None then 'k_x' and 'k_y' 

are chosen. If integers, then these are the indices 

of the momenta: 0 → k_x, 1 → k_y, 2 → k_z. If strings, 

then these name the momenta. 

 

Returns 

------- 

The specified momenta as sympy atoms. 

""" 

 

# Specify which momenta to substitute for the Landau level basis. 

if momenta is None: 

# Use k_x and k_y by default 

momenta = momentum_operators[:2] 

else: 

422 ↛ 423line 422 didn't jump to line 423, because the condition on line 422 was never true if len(momenta) != 2: 

raise ValueError("Two momenta must be specified.") 

 

k_names = [k.name for k in momentum_operators] 

 

427 ↛ exit,   427 ↛ 4302 missed branches: 1) line 428 didn't finish the list comprehension on line 428, 2) line 427 didn't jump to line 430, because the condition on line 427 was never true if all([type(i) is int for i in momenta]) and all( 

[i >= 0 and i < 3 for i in momenta] 

): 

momenta = [momentum_operators[i] for i in momenta] 

431 ↛ 438line 431 didn't jump to line 438, because the condition on line 431 was never false elif all([isinstance(momentum, str) for momentum in momenta]) and all( 

[momentum in k_names for momentum in momenta] 

): 

momenta = [ 

momentum_operators[k_names.index(momentum)] for momentum in momenta 

] 

else: 

raise ValueError("Momenta must all be integers or strings.") 

 

return tuple(momenta) 

 

 

def _find_normal_coordinate(hamiltonian, momenta): 

444 ↛ exitline 445 didn't finish the generator expression on line 445 discrete_momentum = next( 

momentum for momentum in momentum_operators if momentum not in momenta 

) 

normal_coordinate = position_operators[momentum_operators.index(discrete_momentum)] 

return normal_coordinate 

 

 

def _has_coordinate(coord, expr): 

momentum = momentum_operators[position_operators.index(coord)] 

atoms = set(expr.atoms()) 

return coord in atoms or momentum in atoms