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# Copyright 2011-2019 Kwant authors. 

2# 

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

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

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

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

7# https://kwant-project.org/authors. 

8 

9from keyword import iskeyword 

10from collections import defaultdict 

11import itertools 

12import warnings 

13 

14import numpy as np 

15import tinyarray as ta 

16 

17import sympy 

18from sympy.matrices.matrices import MatrixBase 

19from sympy.utilities.lambdify import lambdastr 

20from sympy.printing.lambdarepr import LambdaPrinter 

21from sympy.printing.precedence import precedence 

22from sympy.core.function import AppliedUndef 

23 

24from .. import builder, lattice, system 

25from .. import KwantDeprecationWarning 

26from .._common import reraise_warnings 

27from ._common import (sympify, gcd, position_operators, momentum_operators, 

28 monomials) 

29 

30 

31__all__ = ['discretize'] 

32 

33 

34_wf = sympy.Function('_internal_unique_name', commutative=False) 

35_momentum_operators = {s.name: s for s in momentum_operators} 

36_position_operators = {s.name: s for s in position_operators} 

37_displacements = {s: sympy.Symbol('a_{}'.format(s)) for s in 'xyz'} 

38 

39 

40class _DiscretizedBuilder(builder.Builder): 

41 """A builder that is made from a discretized model and knows how to 

42 pretty-print itself.""" 

43 

44 def __init__(self, coords, lattice, symmetry=None, **kwargs): 

45 super().__init__(symmetry, **kwargs) 

46 self.lattice = lattice 

47 self._coords = coords 

48 

49 def __str__(self): 

50 result = [] 

51 

52 sv = list(s for s in self.site_value_pairs()) 

53 if len(sv) != 1: 53 ↛ 54line 53 didn't jump to line 54, because the condition on line 53 was never true

54 raise ValueError("Cannot pretty-print _DiscretizedBuilder: " 

55 "must contain a single site.") 

56 site, site_value = sv[0] 

57 if any(e != 0 for e in site.tag): 57 ↛ 58line 57 didn't jump to line 58, because the condition on line 57 was never true

58 raise ValueError("Cannot pretty-print _DiscretizedBuilder: " 

59 "site must be located at origin.") 

60 

61 result.extend(["# Discrete coordinates: ", 

62 " ".join(self._coords), 

63 "\n\n"]) 

64 

65 for key, val in itertools.chain(self.site_value_pairs(), 

66 self.hopping_value_pairs()): 

67 if isinstance(key, system.Site): 

68 result.append("# Onsite element:\n") 

69 else: 

70 a, b = key 

71 assert a is site 

72 result.extend(["# Hopping from ", 

73 str(tuple(b.tag)), 

74 ":\n"]) 

75 result.append(val._source if callable(val) else repr(val)) 

76 result.append('\n\n') 

77 

78 result.pop() 

79 

80 return "".join(result) 

81 

82 # Enable human-readable rendering in Jupyter notebooks: 

83 def _repr_html_(self): 

84 return self.__str__() 

85 

86 

87################ Interface functions 

88 

89 

90def discretize(hamiltonian, coords=None, *, grid=None, locals=None, 

91 grid_spacing=None): 

92 """Construct a tight-binding model from a continuum Hamiltonian. 

93 

94 If necessary, the given Hamiltonian is sympified using 

95 `kwant.continuum.sympify`. It is then discretized symbolically and turned 

96 into a `~kwant.builder.Builder` instance that may be used with 

97 `~kwant.builder.Builder.fill`. 

98 

99 This is a convenience function that is equivalent to first calling 

100 `~kwant.continuum.discretize_symbolic` and feeding its result into 

101 `~kwant.continuum.build_discretized`. 

102 

103 .. warning:: 

104 This function uses ``eval`` (because it calls ``sympy.sympify``), and 

105 thus should not be used on unsanitized input. 

106 

107 Parameters 

108 ---------- 

109 hamiltonian : str or SymPy expression 

110 Symbolic representation of a continuous Hamiltonian. It is 

111 converted to a SymPy expression using `kwant.continuum.sympify`. 

112 coords : sequence of strings, optional 

113 The coordinates for which momentum operators will be treated as 

114 differential operators. May contain only "x", "y" and "z" and must be 

115 sorted. If not provided, `coords` will be obtained from the input 

116 Hamiltonian by reading the present coordinates and momentum operators. 

117 grid : scalar or kwant.lattice.Monatomic instance, optional 

118 Lattice that will be used as a discretization grid. It must have 

119 orthogonal primitive vectors. If a scalar value is given, a lattice 

120 with the appropriate grid spacing will be generated. If not provided, 

121 a lattice with grid spacing 1 in all directions will be generated. 

122 locals : dict, optional 

123 Additional namespace entries for `~kwant.continuum.sympify`. May be 

124 used to simplify input of matrices or modify input before proceeding 

125 further. For example: 

126 ``locals={'k': 'k_x + I * k_y'}`` or 

127 ``locals={'sigma_plus': [[0, 2], [0, 0]]}``. 

128 grid_spacing : int or float, optional 

129 (deprecated) Spacing of the discretization grid. If unset the default 

130 value will be 1. Cannot be used together with ``grid``. 

131 

132 Returns 

133 ------- 

134 model : `~kwant.builder.Builder` 

135 The translationally symmetric builder that corresponds to the provided 

136 Hamiltonian. This builder instance belongs to a subclass of the 

137 standard builder that may be printed to obtain the source code of the 

138 value functions. It also holds the discretization lattice (a 

139 `~kwant.lattice.Monatomic` instance with lattice constant 

140 `grid_spacing`) in the ``lattice`` attribute. 

141 """ 

142 tb, coords = discretize_symbolic(hamiltonian, coords, locals=locals) 

143 return build_discretized(tb, coords, grid=grid, grid_spacing=grid_spacing) 

144 

145 

146def discretize_symbolic(hamiltonian, coords=None, *, locals=None): 

147 """Discretize a continuous Hamiltonian into a tight-binding representation. 

148 

149 If necessary, the given Hamiltonian is sympified using 

150 `kwant.continuum.sympify`. It is then discretized symbolically. 

151 

152 The two return values may be used directly as the first two arguments for 

153 `~kwant.continuum.build_discretized`. 

154 

155 .. warning:: 

156 This function uses ``eval`` (because it calls ``sympy.sympify``), and 

157 thus should not be used on unsanitized input. 

158 

159 Parameters 

160 ---------- 

161 hamiltonian : str or SymPy expression 

162 Symbolic representation of a continuous Hamiltonian. It is 

163 converted to a SymPy expression using `kwant.continuum.sympify`. 

164 coords : sequence of strings, optional 

165 The coordinates for which momentum operators will be treated as 

166 differential operators. May contain only "x", "y" and "z" and must be 

167 sorted. If not provided, `coords` will be obtained from the input 

168 Hamiltonian by reading the present coordinates and momentum operators. 

169 locals : dict, optional 

170 Additional namespace entries for `~kwant.continuum.sympify`. May be 

171 used to simplify input of matrices or modify input before proceeding 

172 further. For example: 

173 ``locals={'k': 'k_x + I * k_y'}`` or 

174 ``locals={'sigma_plus': [[0, 2], [0, 0]]}``. 

175 

176 Returns 

177 ------- 

178 tb_hamiltonian : dict 

179 Keys are tuples of integers; the offsets of the hoppings ((0, 0, 0) for 

180 the onsite). Values are symbolic expressions for the hoppings/onsite. 

181 coords : list of strings 

182 The coordinates that have been discretized. 

183 """ 

184 with reraise_warnings(): 

185 hamiltonian = sympify(hamiltonian, locals) 

186 

187 atoms_names = [s.name for s in hamiltonian.atoms(sympy.Symbol)] 

188 if any(s in ('a_x', 'a_y', 'a_z') for s in atoms_names): 188 ↛ 189line 188 didn't jump to line 189, because the condition on line 188 was never true

189 raise TypeError("'a_x', 'a_y' and 'a_z' are symbols used internally " 

190 "to represent grid spacings; please use a different " 

191 "symbol.") 

192 

193 hamiltonian = sympy.expand(hamiltonian) 

194 if coords is None: 

195 used_momenta = set(_momentum_operators) & set(atoms_names) 

196 coords = {k[-1] for k in used_momenta} 

197 else: 

198 coords = list(coords) 

199 if coords != sorted(coords): 199 ↛ 200line 199 didn't jump to line 200, because the condition on line 199 was never true

200 raise ValueError("The argument 'coords' must be sorted.") 

201 if any(c not in 'xyz' for c in coords): 201 ↛ 202line 201 didn't jump to line 202, because the condition on line 201 was never true

202 raise ValueError("The argument 'coords' may only contain " 

203 "'x', 'y', or 'z'.") 

204 

205 coords = sorted(coords) 

206 

207 if len(coords) == 0: 207 ↛ 208line 207 didn't jump to line 208, because the condition on line 207 was never true

208 raise ValueError("Failed to read any discrete coordinates. This is " 

209 "probably due to a lack of momentum operators in " 

210 "your input. You can use the 'coords' " 

211 "parameter to provide them.") 

212 

213 onsite_zeros = (0,) * len(coords) 

214 

215 if not isinstance(hamiltonian, MatrixBase): 

216 hamiltonian = sympy.Matrix([hamiltonian]) 

217 _input_format = 'expression' 

218 else: 

219 _input_format = 'matrix' 

220 

221 shape = hamiltonian.shape 

222 tb = defaultdict(lambda: sympy.zeros(*shape)) 

223 tb[onsite_zeros] = sympy.zeros(*shape) 

224 

225 for (i, j), expression in np.ndenumerate(hamiltonian): 

226 hoppings = _discretize_expression(expression, coords) 

227 

228 for offset, hop in hoppings.items(): 

229 tb[offset][i, j] += hop 

230 

231 # do not include Hermitian conjugates of hoppings 

232 wanted_hoppings = sorted(list(tb))[len(list(tb)) // 2:] 

233 tb = {k: v for k, v in tb.items() if k in wanted_hoppings} 

234 

235 if _input_format == 'expression': 

236 tb = {k: v[0, 0] for k, v in tb.items()} 

237 

238 return tb, coords 

239 

240 

241def build_discretized(tb_hamiltonian, coords, *, grid=None, locals=None, 

242 grid_spacing=None): 

243 """Create a template builder from a symbolic tight-binding Hamiltonian. 

244 

245 The provided symbolic tight-binding Hamiltonian is put on a (hyper) square 

246 lattice and turned into Python functions. These functions are used to 

247 create a `~kwant.builder.Builder` instance that may be used with 

248 `~kwant.builder.Builder.fill` to construct a system of a desired shape. 

249 

250 The return values of `~kwant.continuum.discretize_symbolic` may be used 

251 directly for the first two arguments of this function. 

252 

253 .. warning:: 

254 This function uses ``eval`` (because it calls ``sympy.sympify``), and 

255 thus should not be used on unsanitized input. 

256 

257 Parameters 

258 ---------- 

259 tb_hamiltonian : dict 

260 Keys must be the offsets of the hoppings, represented by tuples of 

261 integers ((0, 0, 0) for onsite). Values must be symbolic expressions 

262 for the hoppings/onsite or expressions that can by sympified with 

263 `kwant.continuum.sympify`. 

264 coords : sequence of strings 

265 The coordinates for which momentum operators will be treated as 

266 differential operators. May contain only "x", "y" and "z" and must be 

267 sorted. 

268 grid : scalar or kwant.lattice.Monatomic instance, optional 

269 Lattice that will be used as a discretization grid. It must have 

270 orthogonal primitive vectors. If a scalar value is given, a lattice 

271 with the appropriate grid spacing will be generated. If not provided, 

272 a lattice with grid spacing 1 in all directions will be generated. 

273 locals : dict, optional 

274 Additional namespace entries for `~kwant.continuum.sympify`. May be 

275 used to simplify input of matrices or modify input before proceeding 

276 further. For example: 

277 ``locals={'k': 'k_x + I * k_y'}`` or 

278 ``locals={'sigma_plus': [[0, 2], [0, 0]]}``. 

279 grid_spacing : int or float, optional 

280 (deprecated) Spacing of the discretization grid. If not provided, 

281 the default value will be 1. Cannot be used together with ``grid``. 

282 

283 Returns 

284 ------- 

285 model : `~kwant.builder.Builder` 

286 The translationally symmetric builder that corresponds to the provided 

287 Hamiltonian. This builder instance belongs to a subclass of the 

288 standard builder that may be printed to obtain the source code of the 

289 value functions. It also holds the discretization lattice (a 

290 `~kwant.lattice.Monatomic` instance with lattice constant 

291 `grid_spacing`) in the ``lattice`` attribute. 

292 """ 

293 # check already available constraints (grid will be check later) 

294 if len(coords) == 0: 294 ↛ 295line 294 didn't jump to line 295, because the condition on line 294 was never true

295 raise ValueError('Discrete coordinates cannot be empty.') 

296 

297 if grid_spacing is not None: # TODO: remove when we remove 'grid_spacing' 

298 warnings.warn('The "grid_spacing" parameter is deprecated. Use ' 

299 '"grid" instead.', KwantDeprecationWarning, stacklevel=3) 

300 if grid is None and grid_spacing is None: 

301 grid = 1 # default case 

302 elif grid is None: # TODO: remove when we remove 'grid_spacing' 

303 grid = grid_spacing 

304 elif grid_spacing is not None: 304 ↛ 305line 304 didn't jump to line 305, because the condition on line 304 was never true

305 raise ValueError('"grid_spacing" and "grid" are mutually exclusive.') 

306 

307 coords = list(coords) 

308 grid_dim = len(coords) 

309 

310 if coords != sorted(coords): 310 ↛ 311line 310 didn't jump to line 311, because the condition on line 310 was never true

311 raise ValueError("The argument 'coords' must be sorted.") 

312 

313 # run sympifcation on hamiltonian values 

314 with reraise_warnings(): 

315 for k, v in tb_hamiltonian.items(): 

316 tb_hamiltonian[k] = sympify(v, locals) 

317 

318 # generate grid if required, check constraints if provided 

319 random_element = next(iter(tb_hamiltonian.values())) 

320 norbs = (1 if isinstance(random_element, sympy.Expr) 

321 else random_element.shape[0]) 

322 

323 if np.isscalar(grid): 

324 lat = lattice.Monatomic(grid * np.eye(grid_dim), norbs=norbs) 

325 else: 

326 lat = grid 

327 

328 # check grid constraints 

329 is_diagonal = lambda m: np.allclose(m, np.diag(np.diagonal(m))) 

330 if not (lat.prim_vecs.shape[0] == grid_dim and 

331 is_diagonal(lat.prim_vecs)): 

332 raise ValueError('"grid" has to be an orthogonal lattice ' 

333 'of dimension matching number of "coords".') 

334 

335 if (lat.norbs is not None) and (lat.norbs != norbs): 

336 raise ValueError( 

337 'Number of lattice orbitals does not match the number ' 

338 'of orbitals in the Hamiltonian.' 

339 ) 

340 

341 # continue with building the template 

342 tb = {} 

343 for n, (offset, hopping) in enumerate(tb_hamiltonian.items()): 

344 onsite = all(i == 0 for i in offset) 

345 

346 if onsite: 

347 name = 'onsite' 

348 else: 

349 name = 'hopping_{}'.format(n) 

350 

351 tb[offset] = _builder_value(hopping, coords, np.diag(lat.prim_vecs), 

352 onsite, name) 

353 

354 onsite_zeros = (0,) * grid_dim 

355 onsite = tb.pop(onsite_zeros) 

356 # 'delta' parameter to HoppingKind is the negative of the 'hopping offset' 

357 hoppings = {builder.HoppingKind(tuple(-i for i in d), lat): val 

358 for d, val in tb.items()} 

359 

360 syst = _DiscretizedBuilder( 

361 coords, lat, lattice.TranslationalSymmetry(*lat.prim_vecs) 

362 ) 

363 syst[lat(*onsite_zeros)] = onsite 

364 for hop, val in hoppings.items(): 

365 syst[hop] = val 

366 

367 return syst 

368 

369 

370def _differentiate(expression, coordinate_name): 

371 """Calculate derivative of an expression for given coordinate. 

372 

373 Parameters: 

374 ----------- 

375 expression : sympy.Expr 

376 Sympy expression containing function to be derivated. 

377 coordinate_name : string 

378 Coordinate over which derivative is calculated. 

379 

380 Returns 

381 ------- 

382 sympy.Expr 

383 """ 

384 assert coordinate_name in 'xyz' 

385 coordinate = _position_operators[coordinate_name] 

386 h = _displacements[coordinate_name] 

387 

388 expr1 = expression.subs(coordinate, coordinate + h) 

389 expr2 = expression.subs(coordinate, coordinate - h) 

390 

391 return (expr1 - expr2) / (2 * h) 

392 

393 

394def _discretize_summand(summand, coords): 

395 """Discretize a product of factors. 

396 

397 Parameters 

398 ---------- 

399 summand : sympy.Expr 

400 coords : sequence of strings 

401 Must be a subset of ``{'x', 'y', 'z'}``. 

402 

403 Returns 

404 ------- 

405 sympy.Expr 

406 """ 

407 assert not isinstance(summand, sympy.Add), "Input should be one summand." 

408 momenta = ['k_{}'.format(s) for s in coords] 

409 

410 factors = reversed(summand.as_ordered_factors()) 

411 result = 1 

412 for factor in factors: 

413 symbol, exponent = factor.as_base_exp() 

414 if isinstance(symbol, sympy.Symbol) and symbol.name in momenta: 

415 for i in range(exponent): 

416 coordinate = symbol.name[-1] 

417 # apply momentum as differential operator '-i d/dx' 

418 result = -sympy.I * _differentiate(result, coordinate) 

419 else: 

420 result = factor * result 

421 

422 return result 

423 

424 

425def _discretize_expression(expression, coords): 

426 """Discretize an expression into a discrete (tight-binding) representation. 

427 

428 Parameters 

429 ---------- 

430 expression : sympy.Expr 

431 coords : sequence of strings 

432 Must be a subset of ``{'x', 'y', 'z'}``. 

433 

434 Returns 

435 ------- 

436 dict 

437 Keys are tuples of integers; the offsets of the hoppings 

438 ((0, 0, 0) for the onsite). Values are symbolic expressions 

439 for the hoppings/onsite. 

440 """ 

441 def _read_offset(wf): 

442 # e.g. wf(x + h, y, z + h) -> (1, 0, 1) 

443 assert wf.func == _wf 

444 

445 offset = [] 

446 for c, arg in zip(coords, wf.args): 

447 coefficients = arg.as_coefficients_dict() 

448 assert coefficients[_position_operators[c]] == 1 

449 

450 ai = _displacements[c] 

451 offset.append(coefficients.pop(ai, 0)) 

452 return tuple(offset) 

453 

454 def _extract_hoppings(expr): 

455 """Read hoppings and perform shortening operation.""" 

456 expr = sympy.expand(expr) 

457 summands = [e.as_ordered_factors() for e in expr.as_ordered_terms()] 

458 

459 offset = [_read_offset(s[-1]) for s in summands] 

460 coeffs = [sympy.Mul(*s[:-1]) for s in summands] 

461 offset = np.array(offset, dtype=int) 

462 # rescale the offsets for each coordinate by their greatest 

463 # common divisor across the summands. e.g: 

464 # wf(x+2h) + wf(x+4h) --> wf(x+h) + wf(x+2h) and a_x //= 2 

465 subs = {} 

466 for i, xi in enumerate(coords): 

467 factor = int(gcd(*offset[:, i])) 

468 if factor < 1: 

469 continue 

470 offset[:, i] //= factor 

471 subs[_displacements[xi]] = _displacements[xi] / factor 

472 # apply the rescaling to the hoppings 

473 output = defaultdict(lambda: sympy.Integer(0)) 

474 for n, c in enumerate(coeffs): 

475 output[tuple(offset[n].tolist())] += c.subs(subs) 

476 return dict(output) 

477 

478 # if there are no momenta in the expression, then it is an onsite 

479 atoms_names = [s.name for s in expression.atoms(sympy.Symbol)] 

480 if not set(_momentum_operators) & set(atoms_names): 

481 n = len(coords) 

482 return {(0,) * n: expression} 

483 

484 # make sure we have list of summands 

485 summands = expression.as_ordered_terms() 

486 

487 # discretize every summand 

488 coordinates = tuple(_position_operators[s] for s in coords) 

489 wf = _wf(*coordinates) 

490 

491 discrete_expression = defaultdict(int) 

492 for summand in summands: 

493 summand = _discretize_summand(summand * wf, coords) 

494 hops = _extract_hoppings(summand) 

495 for k, v in hops.items(): 

496 discrete_expression[k] += v 

497 

498 return dict(discrete_expression) 

499 

500 

501################ string processing 

502 

503class _NumericPrinter(LambdaPrinter): 

504 

505 def __init__(self): 

506 if 'allow_unknown_functions' in LambdaPrinter._default_settings: 506 ↛ 510line 506 didn't jump to line 510, because the condition on line 506 was never false

507 settings = {'allow_unknown_functions': True} 

508 else: 

509 # We're on Sympy without "allow_unknown_functions" setting 

510 settings = {} 

511 

512 LambdaPrinter.__init__(self, settings=settings) 

513 

514 self.known_functions = {} 

515 self.known_constants = {'pi': 'pi', 'Pi': 'pi', 'I': 'I'} 

516 

517 def _print_ImaginaryUnit(self, expr): 

518 # prevent sympy from printing 'I' for imaginary unit 

519 return "1j" 

520 

521 def _print_Pow(self, expr): 

522 # copied from sympy's StrPrinter with the code paths 

523 # to print 'sqrt' removed. 

524 PREC = precedence(expr) 

525 

526 if expr.is_commutative and expr.exp is -sympy.S.One: 526 ↛ 527line 526 didn't jump to line 527, because the condition on line 526 was never true

527 return '1/%s' % self.parenthesize(expr.base, PREC) 

528 

529 e = self.parenthesize(expr.exp, PREC) 

530 if (self.printmethod == '_sympyrepr' and 530 ↛ 532line 530 didn't jump to line 532, because the condition on line 530 was never true

531 expr.exp.is_Rational and expr.exp.q != 1): 

532 if e.startswith('(Rational'): 

533 return '%s**%s' % (self.parenthesize(expr.base, PREC), e[1:-1]) 

534 return '%s**%s' % (self.parenthesize(expr.base, PREC), e) 

535 

536 

537def _print_sympy(expr): 

538 return lambdastr((), expr, printer=_NumericPrinter)[len('lambda : '):] 

539 

540 

541def _return_string(expr, coords): 

542 """Process a sympy expression into an evaluatable Python return statement. 

543 

544 Parameters 

545 ---------- 

546 expr : sympy.Expr 

547 

548 Returns 

549 ------- 

550 output : string 

551 A return string that can be used to assemble a Kwant value function. 

552 map_func_calls : dict 

553 mapping of function calls to assigned constants. 

554 const_symbols : sequance of sympy.Symbol 

555 All constants that appear in the expression. 

556 _cache: dict 

557 mapping of cache symbols to cached matrices. 

558 """ 

559 _cache = {} 

560 def cache(x): 

561 s = sympy.symbols('_cache_{}'.format(len(_cache))) 

562 _cache[str(s)] = ta.array(x.tolist(), complex) 

563 return s 

564 

565 blacklisted = set(coords) | {'site', 'site1', 'site2'} 

566 const_symbols = {s for s in expr.atoms(sympy.Symbol) 

567 if s.name not in blacklisted} 

568 

569 # functions will be evaluated within the function body and the 

570 # result assigned to a symbol '_const_<n>', so we replace all 

571 # function calls by these symbols in the return statement. 

572 map_func_calls = expr.atoms(AppliedUndef, sympy.Function) 

573 map_func_calls = {s: sympy.symbols('_const_{}'.format(n)) 

574 for n, s in enumerate(map_func_calls)} 

575 

576 expr = expr.subs(map_func_calls) 

577 

578 if isinstance(expr, MatrixBase): 

579 # express matrix return values in terms of sums of known matrices, 

580 # which will be assigned to '_cache_n' in the function body. 

581 mons = monomials(expr, expr.atoms(sympy.Symbol)) 

582 mons = {k: cache(v) for k, v in mons.items()} 

583 mons = ["{} * {}".format(_print_sympy(k), _print_sympy(v)) 

584 for k, v in mons.items()] 

585 output = " + ".join(mons) 

586 else: 

587 output = _print_sympy(expr) 

588 

589 return 'return {}'.format(output), map_func_calls, const_symbols, _cache 

590 

591 

592def _assign_symbols(map_func_calls, coords, onsite): 

593 """Generate a series of assignments. 

594 

595 Parameters 

596 ---------- 

597 map_func_calls : dict 

598 mapping of function calls to assigned constants. 

599 coords : sequence of strings 

600 If left as None coordinates will not be read from a site. 

601 onsite : bool 

602 True if function is called for onsite, false for hoppings 

603 

604 Returns 

605 ------- 

606 assignments : list of strings 

607 List of lines used for including in a function. 

608 """ 

609 lines = [] 

610 

611 if coords: 

612 site = 'site' if onsite else 'site1' 

613 args = ', '.join(coords), site 

614 lines.append('({}, ) = {}.pos'.format(*args)) 

615 

616 for k, v in map_func_calls.items(): 

617 lines.append("{} = {}".format(v, _print_sympy(k))) 

618 

619 return lines 

620 

621 

622def _builder_value(expr, coords, grid_spacing, onsite, 

623 name='_anonymous_func'): 

624 """Generate a builder value from a sympy expression. 

625 

626 Parameters 

627 ---------- 

628 expr : sympy.Expr or sympy.matrix 

629 Expr that from which value function will be generated. 

630 coords : sequence of strings 

631 List of coodinates present in the system. 

632 grid_spacing : sequence of scalars 

633 Lattice spacing of the system in each coordinate. 

634 

635 Returns 

636 ------- 

637 `expr` transformed into an object that can be used as a 

638 `kwant.builder.Builder` value. Either a numerical value 

639 (``tinyarray.array`` instance or complex number) or a value function. In 

640 the case of a function, the source code is available in its `_source` 

641 attribute. 

642 """ 

643 

644 expr = expr.subs({_displacements[c]: grid_spacing[n] 

645 for n, c in enumerate(coords)}) 

646 return_string, map_func_calls, const_symbols, _cache = _return_string( 

647 expr, coords=coords) 

648 

649 # first check if value function needs to read coordinates 

650 atoms_names = {s.name for s in expr.atoms(sympy.Symbol)} 

651 if not set(coords) & atoms_names: 

652 coords = None 

653 

654 # constants and functions in the sympy input will be passed 

655 # as arguments to the value function 

656 arg_names = set.union({s.name for s in const_symbols}, 

657 {str(k.func) for k in map_func_calls}) 

658 

659 # check if all argument names are valid python identifiers 

660 for arg_name in arg_names: 

661 if not (arg_name.isidentifier() and not iskeyword(arg_name)): 

662 raise ValueError("Invalid name in used symbols: {}\n" 

663 "Names of symbols used in Hamiltonian " 

664 "must be valid Python identifiers and " 

665 "may not be keywords".format(arg_name)) 

666 

667 arg_names = ', '.join(sorted(arg_names)) 

668 

669 if (not arg_names) and (coords is None): 

670 # we can just use a constant value instead of a value function 

671 if isinstance(expr, sympy.MatrixBase): 

672 return ta.array(expr.tolist(), complex) 

673 else: 

674 return complex(expr) 

675 

676 lines = _assign_symbols(map_func_calls, onsite=onsite, coords=coords) 

677 lines.append(return_string) 

678 

679 separator = '\n ' 

680 # 'site_string' is tightly coupled to the symbols used in '_assign_symbol' 

681 site_string = 'site' if onsite else 'site1, site2' 

682 if arg_names: 

683 header_str = 'def {}({}, {}):' 

684 header = header_str.format(name, site_string, arg_names) 

685 else: 

686 header = 'def {}({}):'.format(name, site_string) 

687 func_code = separator.join([header] + list(lines)) 

688 

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

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

691 # of some function. 

692 namespace = {'pi': np.pi, 'I': 1j} 

693 namespace.update(_cache) 

694 

695 source = [] 

696 for k, v in _cache.items(): 

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

698 source.append(func_code) 

699 

700 exec(func_code, namespace) 

701 f = namespace[name] 

702 f._source = "".join(source) 

703 

704 return f