Coverage for gamdpy/configuration/make_lattice.py: 96%

26 statements  

« prev     ^ index     » next       coverage.py v7.4.4, created at 2025-06-14 15:25 +0200

1def make_lattice(unit_cell: dict, cells: list = None, rho=None) -> tuple: 

2 """ Returns a configuration of a crystal lattice. 

3 The lattice is constructed by replicating the unit cell in all directions. 

4 The `unit_cell` is a dictonary with `fractional_coordinates` for particles, and 

5 the `lattice_constants` as a list of lengths of the unit cell in all directions. 

6 The `cells` are the number of unit cells in each direction. 

7 

8 Returns a list of positions of the atoms in the lattice, and the box vector of the lattice. 

9 The returned positions are in a box of size -L/2 to L/2 for each direction. 

10 

11 Example 

12 ------- 

13 

14 >>> import gamdpy as gp 

15 >>> positions, box_vector = gp.configuration.make_lattice(gp.unit_cells.FCC, cells=[8, 8, 8], rho=1.0) 

16 >>> configuration = gp.Configuration(D=3) 

17 >>> configuration['r'] = positions 

18 >>> configuration.simbox = gp.Orthorhombic(configuration.D, box_vector) 

19 

20 See also 

21 -------- 

22 

23 :meth:`gamdpy.Configuration.make_lattice` 

24 

25 """ 

26 import numpy as np 

27 pos = unit_cell["fractional_coordinates"] 

28 lat = unit_cell["lattice_constants"] 

29 particles_in_unit_cell = len(pos) 

30 spatial_dimension = len(pos[0]) 

31 if cells is None: 

32 cells = [1] * spatial_dimension 

33 number_of_cells = np.prod(cells) 

34 positions = np.zeros( 

35 shape=(particles_in_unit_cell * number_of_cells, spatial_dimension), 

36 dtype=np.float64, 

37 ) 

38 for cell_index in range(number_of_cells): 

39 cell_coordinates = np.array( 

40 np.unravel_index(cell_index, cells), dtype=np.float64 

41 ) 

42 for particle_index in range(particles_in_unit_cell): 

43 positions[cell_index * particles_in_unit_cell + particle_index] = ( 

44 pos[particle_index] + cell_coordinates 

45 ) 

46 positions *= lat 

47 box_vector = np.array(lat) * np.array(cells) 

48 if rho is not None: 

49 box_volume = np.prod(box_vector) 

50 number_of_particles = len(positions) 

51 volume_per_particle = box_volume / number_of_particles 

52 target_volume_per_particle = 1.0 / rho 

53 scale_factor = (target_volume_per_particle / volume_per_particle) ** (1.0 / 3.0) 

54 positions *= scale_factor 

55 box_vector *= scale_factor 

56 

57 # Center the box (-L/2 to L/2) 

58 positions -= box_vector/2.0 

59 

60 return positions, box_vector