1""" 

2This is the implementation of the exciting I/O functions 

3The functions are called with read write using the format "exciting" 

4 

5""" 

6 

7import numpy as np 

8import xml.etree.ElementTree as ET 

9from ase.atoms import Atoms 

10from ase.units import Bohr 

11from ase.utils import writer 

12from xml.dom import minidom 

13 

14 

15def read_exciting(fileobj, index=-1): 

16 """Reads structure from exiting xml file. 

17 

18 Parameters 

19 ---------- 

20 fileobj: file object 

21 File handle from which data should be read. 

22 

23 Other parameters 

24 ---------------- 

25 index: integer -1 

26 Not used in this implementation. 

27 """ 

28 

29 # Parse file into element tree 

30 doc = ET.parse(fileobj) 

31 root = doc.getroot() 

32 speciesnodes = root.find('structure').iter('species') 

33 symbols = [] 

34 positions = [] 

35 basevects = [] 

36 atoms = None 

37 # Collect data from tree 

38 for speciesnode in speciesnodes: 

39 symbol = speciesnode.get('speciesfile').split('.')[0] 

40 natoms = speciesnode.iter('atom') 

41 for atom in natoms: 

42 x, y, z = atom.get('coord').split() 

43 positions.append([float(x), float(y), float(z)]) 

44 symbols.append(symbol) 

45 # scale unit cell accorting to scaling attributes 

46 if 'scale' in doc.find('structure/crystal').attrib: 

47 scale = float(str(doc.find('structure/crystal').attrib['scale'])) 

48 else: 

49 scale = 1 

50 

51 if 'stretch' in doc.find('structure/crystal').attrib: 

52 a, b, c = doc.find('structure/crystal').attrib['stretch'].text.split() 

53 stretch = np.array([float(a), float(b), float(c)]) 

54 else: 

55 stretch = np.array([1.0, 1.0, 1.0]) 

56 basevectsn = root.findall('structure/crystal/basevect') 

57 for basevect in basevectsn: 

58 x, y, z = basevect.text.split() 

59 basevects.append(np.array([float(x) * Bohr * stretch[0], 

60 float(y) * Bohr * stretch[1], 

61 float(z) * Bohr * stretch[2] 

62 ]) * scale) 

63 atoms = Atoms(symbols=symbols, cell=basevects) 

64 

65 atoms.set_scaled_positions(positions) 

66 if 'molecule' in root.find('structure').attrib.keys(): 

67 if root.find('structure').attrib['molecule']: 

68 atoms.set_pbc(False) 

69 else: 

70 atoms.set_pbc(True) 

71 

72 return atoms 

73 

74 

75@writer 

76def write_exciting(fileobj, images): 

77 """writes exciting input structure in XML 

78 

79 Parameters 

80 ---------- 

81 filename : str 

82 Name of file to which data should be written. 

83 images : Atom Object or List of Atoms objects 

84 This function will write the first Atoms object to file. 

85 

86 Returns 

87 ------- 

88 """ 

89 root = atoms2etree(images) 

90 rough_string = ET.tostring(root, 'utf-8') 

91 reparsed = minidom.parseString(rough_string) 

92 pretty = reparsed.toprettyxml(indent="\t") 

93 fileobj.write(pretty) 

94 

95 

96def atoms2etree(images): 

97 """This function creates the XML DOM corresponding 

98 to the structure for use in write and calculator 

99 

100 Parameters 

101 ---------- 

102 

103 images : Atom Object or List of Atoms objects 

104 

105 Returns 

106 ------- 

107 root : etree object 

108 Element tree of exciting input file containing the structure 

109 """ 

110 if not isinstance(images, (list, tuple)): 

111 images = [images] 

112 

113 root = ET.Element('input') 

114 root.set( 

115 '{http://www.w3.org/2001/XMLSchema-instance}noNamespaceSchemaLocation', 

116 'http://xml.exciting-code.org/excitinginput.xsd') 

117 

118 title = ET.SubElement(root, 'title') 

119 title.text = '' 

120 structure = ET.SubElement(root, 'structure') 

121 crystal = ET.SubElement(structure, 'crystal') 

122 atoms = images[0] 

123 for vec in atoms.cell: 

124 basevect = ET.SubElement(crystal, 'basevect') 

125 basevect.text = '%.14f %.14f %.14f' % tuple(vec / Bohr) 

126 

127 oldsymbol = '' 

128 oldrmt = -1 

129 newrmt = -1 

130 scaled = atoms.get_scaled_positions() 

131 for aindex, symbol in enumerate(atoms.get_chemical_symbols()): 

132 if 'rmt' in atoms.arrays: 

133 newrmt = atoms.get_array('rmt')[aindex] / Bohr 

134 if symbol != oldsymbol or newrmt != oldrmt: 

135 speciesnode = ET.SubElement(structure, 'species', 

136 speciesfile='%s.xml' % symbol, 

137 chemicalSymbol=symbol) 

138 oldsymbol = symbol 

139 if 'rmt' in atoms.arrays: 

140 oldrmt = atoms.get_array('rmt')[aindex] / Bohr 

141 if oldrmt > 0: 

142 speciesnode.attrib['rmt'] = '%.4f' % oldrmt 

143 

144 atom = ET.SubElement(speciesnode, 'atom', 

145 coord='%.14f %.14f %.14f' % tuple(scaled[aindex])) 

146 if 'momenta' in atoms.arrays: 

147 atom.attrib['bfcmt'] = '%.14f %.14f %.14f' % tuple( 

148 atoms.get_array('mommenta')[aindex]) 

149 

150 return root