Coverage for /builds/debichem-team/python-ase/ase/io/gpw.py: 11.86%

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

2import ase.io.ulm as ulm 

3from ase import Atoms 

4from ase.calculators.singlepoint import ( 

5 SinglePointDFTCalculator, 

6 SinglePointKPoint, 

7 all_properties, 

8) 

9from ase.io.trajectory import read_atoms 

10from ase.units import Bohr, Hartree 

11 

12 

13def read_gpw(filename): 

14 try: 

15 reader = ulm.open(filename) 

16 except ulm.InvalidULMFileError: 

17 return read_old_gpw(filename) 

18 

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

20 

21 wfs = reader.wave_functions 

22 kpts = wfs.get('kpts') 

23 if kpts is None: 

24 ibzkpts = None 

25 bzkpts = None 

26 bz2ibz = None 

27 else: 

28 ibzkpts = kpts.ibzkpts 

29 bzkpts = kpts.get('bzkpts') 

30 bz2ibz = kpts.get('bz2ibz') 

31 

32 if reader.version >= 3: 

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

34 else: 

35 efermi = reader.occupations.fermilevel 

36 

37 atoms.calc = SinglePointDFTCalculator( 

38 atoms, 

39 efermi=efermi, 

40 ibzkpts=ibzkpts, 

41 bzkpts=bzkpts, 

42 bz2ibz=bz2ibz, 

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

44 # doesn't know: 

45 **{property: value 

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

47 if property in all_properties}) 

48 

49 if kpts is not None: 

50 atoms.calc.kpts = [] 

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

52 wfs.occupations)): 

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

54 eps_kn, f_kn)): 

55 atoms.calc.kpts.append( 

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

57 reader.close() 

58 

59 return atoms 

60 

61 

62def read_old_gpw(filename): 

63 from gpaw.io.tar import Reader 

64 r = Reader(filename) 

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

66 numbers = r.get('AtomicNumbers') 

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

68 pbc = r.get('BoundaryConditions') 

69 tags = r.get('Tags') 

70 magmoms = r.get('MagneticMoments') 

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

72 

73 if r.has_array('CartesianForces'): 

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

75 else: 

76 forces = None 

77 

78 atoms = Atoms(positions=positions, 

79 numbers=numbers, 

80 cell=cell, 

81 pbc=pbc) 

82 if tags.any(): 

83 atoms.set_tags(tags) 

84 

85 if magmoms.any(): 

86 atoms.set_initial_magnetic_moments(magmoms) 

87 magmom = magmoms.sum() 

88 else: 

89 magmoms = None 

90 magmom = None 

91 

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

93 forces=forces, 

94 magmoms=magmoms, 

95 magmom=magmom) 

96 kpts = [] 

97 if r.has_array('IBZKPoints'): 

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

99 r.get('IBZKPoints'), 

100 r.get('Eigenvalues'), 

101 r.get('OccupationNumbers')): 

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

103 eps_n[0], f_n[0])) 

104 atoms.calc.kpts = kpts 

105 

106 return atoms