Coverage for /builds/ase/ase/ase/io/gpw.py: 11.86%

1# fmt: off 

2 

3"""Read gpw-file from GPAW.""" 

4import ase.io.ulm as ulm 

5from ase import Atoms 

6from ase.calculators.singlepoint import ( 

7 SinglePointDFTCalculator, 

8 SinglePointKPoint, 

9 all_properties, 

10) 

11from ase.io.trajectory import read_atoms 

12from ase.units import Bohr, Hartree 

13 

14 

15def read_gpw(filename): 

16 try: 

17 reader = ulm.open(filename) 

18 except ulm.InvalidULMFileError: 

19 return read_old_gpw(filename) 

20 

21 atoms = read_atoms(reader.atoms, _try_except=False) 

22 

23 wfs = reader.wave_functions 

24 kpts = wfs.get('kpts') 

25 if kpts is None: 

26 ibzkpts = None 

27 bzkpts = None 

28 bz2ibz = None 

29 else: 

30 ibzkpts = kpts.ibzkpts 

31 bzkpts = kpts.get('bzkpts') 

32 bz2ibz = kpts.get('bz2ibz') 

33 

34 if reader.version >= 3: 

35 efermi = reader.wave_functions.fermi_levels.mean() 

36 else: 

37 efermi = reader.occupations.fermilevel 

38 

39 atoms.calc = SinglePointDFTCalculator( 

40 atoms, 

41 efermi=efermi, 

42 ibzkpts=ibzkpts, 

43 bzkpts=bzkpts, 

44 bz2ibz=bz2ibz, 

45 # New gpw-files may have "non_collinear_magmom(s)" which ASE 

46 # doesn't know: 

47 **{property: value 

48 for property, value in reader.results.asdict().items() 

49 if property in all_properties}) 

50 

51 if kpts is not None: 

52 atoms.calc.kpts = [] 

53 for spin, (eps_kn, f_kn) in enumerate(zip(wfs.eigenvalues, 

54 wfs.occupations)): 

55 for kpt, (weight, eps_n, f_n) in enumerate(zip(kpts.weights, 

56 eps_kn, f_kn)): 

57 atoms.calc.kpts.append( 

58 SinglePointKPoint(weight, spin, kpt, eps_n, f_n)) 

59 reader.close() 

60 

61 return atoms 

62 

63 

64def read_old_gpw(filename): 

65 from gpaw.io.tar import Reader 

66 r = Reader(filename) 

67 positions = r.get('CartesianPositions') * Bohr 

68 numbers = r.get('AtomicNumbers') 

69 cell = r.get('UnitCell') * Bohr 

70 pbc = r.get('BoundaryConditions') 

71 tags = r.get('Tags') 

72 magmoms = r.get('MagneticMoments') 

73 energy = r.get('PotentialEnergy') * Hartree 

74 

75 if r.has_array('CartesianForces'): 

76 forces = r.get('CartesianForces') * Hartree / Bohr 

77 else: 

78 forces = None 

79 

80 atoms = Atoms(positions=positions, 

81 numbers=numbers, 

82 cell=cell, 

83 pbc=pbc) 

84 if tags.any(): 

85 atoms.set_tags(tags) 

86 

87 if magmoms.any(): 

88 atoms.set_initial_magnetic_moments(magmoms) 

89 magmom = magmoms.sum() 

90 else: 

91 magmoms = None 

92 magmom = None 

93 

94 atoms.calc = SinglePointDFTCalculator(atoms, energy=energy, 

95 forces=forces, 

96 magmoms=magmoms, 

97 magmom=magmom) 

98 kpts = [] 

99 if r.has_array('IBZKPoints'): 

100 for w, kpt, eps_n, f_n in zip(r.get('IBZKPointWeights'), 

101 r.get('IBZKPoints'), 

102 r.get('Eigenvalues'), 

103 r.get('OccupationNumbers')): 

104 kpts.append(SinglePointKPoint(w, kpt[0], kpt[1], 

105 eps_n[0], f_n[0])) 

106 atoms.calc.kpts = kpts 

107 

108 return atoms