1""" 

2Stream input commands to lammps to perform desired simulations 

3""" 

4from ase.parallel import paropen 

5from ase.calculators.lammps.unitconvert import convert 

6 

7# "End mark" used to indicate that the calculation is done 

8CALCULATION_END_MARK = "__end_of_ase_invoked_calculation__" 

9 

10 

11def lammps_create_atoms(fileobj, parameters, atoms, prismobj): 

12 """Create atoms in lammps with 'create_box' and 'create_atoms' 

13 

14 :param fileobj: open stream for lammps input 

15 :param parameters: dict of all lammps parameters 

16 :type parameters: dict 

17 :param atoms: Atoms object 

18 :type atoms: Atoms 

19 :param prismobj: coordinate transformation between ase and lammps 

20 :type prismobj: Prism 

21 

22 """ 

23 if parameters["verbose"]: 

24 fileobj.write("## Original ase cell\n".encode("utf-8")) 

25 fileobj.write( 

26 "".join( 

27 [ 

28 "# {0:.16} {1:.16} {2:.16}\n".format(*x) 

29 for x in atoms.get_cell() 

30 ] 

31 ).encode("utf-8") 

32 ) 

33 

34 fileobj.write("lattice sc 1.0\n".encode("utf-8")) 

35 

36 # Get cell parameters and convert from ASE units to LAMMPS units 

37 xhi, yhi, zhi, xy, xz, yz = convert(prismobj.get_lammps_prism(), 

38 "distance", "ASE", parameters.units) 

39 

40 if parameters["always_triclinic"] or prismobj.is_skewed(): 

41 fileobj.write( 

42 "region asecell prism 0.0 {0} 0.0 {1} 0.0 {2} ".format( 

43 xhi, yhi, zhi 

44 ).encode("utf-8") 

45 ) 

46 fileobj.write( 

47 "{0} {1} {2} side in units box\n".format(xy, xz, yz).encode( 

48 "utf-8" 

49 ) 

50 ) 

51 else: 

52 fileobj.write( 

53 "region asecell block 0.0 {0} 0.0 {1} 0.0 {2} " 

54 "side in units box\n".format(xhi, yhi, zhi).encode("utf-8") 

55 ) 

56 

57 symbols = atoms.get_chemical_symbols() 

58 try: 

59 # By request, specific atom type ordering 

60 species = parameters["specorder"] 

61 except AttributeError: 

62 # By default, atom types in alphabetic order 

63 species = sorted(set(symbols)) 

64 

65 species_i = {s: i + 1 for i, s in enumerate(species)} 

66 

67 fileobj.write( 

68 "create_box {0} asecell\n" "".format(len(species)).encode("utf-8") 

69 ) 

70 for sym, pos in zip(symbols, atoms.get_positions()): 

71 # Convert position from ASE units to LAMMPS units 

72 pos = convert(pos, "distance", "ASE", parameters.units) 

73 if parameters["verbose"]: 

74 fileobj.write( 

75 "# atom pos in ase cell: {0:.16} {1:.16} {2:.16}\n" 

76 "".format(*tuple(pos)).encode("utf-8") 

77 ) 

78 fileobj.write( 

79 "create_atoms {0} single {1} {2} {3} remap yes units box\n".format( 

80 *((species_i[sym],) + tuple(prismobj.vector_to_lammps(pos))) 

81 ).encode("utf-8") 

82 ) 

83 

84 

85def write_lammps_in(lammps_in, parameters, atoms, prismobj, 

86 lammps_trj=None, lammps_data=None): 

87 """Write a LAMMPS in_ file with run parameters and settings.""" 

88 

89 def write_model_post_and_masses(fileobj, parameters): 

90 # write additional lines needed for some LAMMPS potentials 

91 if 'model_post' in parameters: 

92 mlines = parameters['model_post'] 

93 for ii in range(0, len(mlines)): 

94 fileobj.write(mlines[ii].encode('utf-8')) 

95 

96 if "masses" in parameters: 

97 for mass in parameters["masses"]: 

98 # Note that the variable mass is a string containing 

99 # the type number and value of mass separated by a space 

100 fileobj.write("mass {0} \n".format(mass).encode("utf-8")) 

101 

102 if isinstance(lammps_in, str): 

103 fileobj = paropen(lammps_in, "wb") 

104 close_in_file = True 

105 else: 

106 # Expect lammps_in to be a file-like object 

107 fileobj = lammps_in 

108 close_in_file = False 

109 

110 if parameters["verbose"]: 

111 fileobj.write("# (written by ASE)\n".encode("utf-8")) 

112 

113 # Write variables 

114 fileobj.write( 

115 ( 

116 "clear\n" 

117 'variable dump_file string "{0}"\n' 

118 'variable data_file string "{1}"\n' 

119 ) 

120 .format(lammps_trj, lammps_data) 

121 .encode("utf-8") 

122 ) 

123 

124 if "package" in parameters: 

125 fileobj.write( 

126 ( 

127 "\n".join( 

128 ["package {0}".format(p) for p in parameters["package"]] 

129 ) 

130 + "\n" 

131 ).encode("utf-8") 

132 ) 

133 

134 # setup styles except 'pair_style' 

135 for style_type in ("atom", "bond", "angle", 

136 "dihedral", "improper", "kspace"): 

137 style = style_type + "_style" 

138 if style in parameters: 

139 fileobj.write('{} {} \n'.format(style, parameters[style]).encode("utf-8")) 

140 

141 # write initialization lines needed for some LAMMPS potentials 

142 if 'model_init' in parameters: 

143 mlines = parameters['model_init'] 

144 for ii in range(0, len(mlines)): 

145 fileobj.write(mlines[ii].encode('utf-8')) 

146 

147 # write units 

148 if 'units' in parameters: 

149 units_line = 'units ' + parameters['units'] + '\n' 

150 fileobj.write(units_line.encode('utf-8')) 

151 else: 

152 fileobj.write('units metal\n'.encode('utf-8')) 

153 

154 pbc = atoms.get_pbc() 

155 if "boundary" in parameters: 

156 fileobj.write( 

157 "boundary {0} \n".format(parameters["boundary"]).encode("utf-8") 

158 ) 

159 else: 

160 fileobj.write( 

161 "boundary {0} {1} {2} \n".format( 

162 *tuple("sp"[int(x)] for x in pbc) 

163 ).encode("utf-8") 

164 ) 

165 fileobj.write("atom_modify sort 0 0.0 \n".encode("utf-8")) 

166 for key in ("neighbor", "newton"): 

167 if key in parameters: 

168 fileobj.write( 

169 "{0} {1} \n".format(key, parameters[key]).encode("utf-8") 

170 ) 

171 fileobj.write("\n".encode("utf-8")) 

172 

173 # write the simulation box and the atoms 

174 if not lammps_data: 

175 lammps_create_atoms(fileobj, parameters, atoms, prismobj) 

176 # or simply refer to the data-file 

177 else: 

178 fileobj.write("read_data {0}\n".format(lammps_data).encode("utf-8")) 

179 

180 # Write interaction stuff 

181 fileobj.write("\n### interactions\n".encode("utf-8")) 

182 if "kim_interactions" in parameters: 

183 fileobj.write("{}\n".format(parameters["kim_interactions"]).encode("utf-8")) 

184 write_model_post_and_masses(fileobj, parameters) 

185 

186 elif ("pair_style" in parameters) and ("pair_coeff" in parameters): 

187 pair_style = parameters["pair_style"] 

188 fileobj.write("pair_style {0} \n".format(pair_style).encode("utf-8")) 

189 for pair_coeff in parameters["pair_coeff"]: 

190 fileobj.write( 

191 "pair_coeff {0} \n" "".format(pair_coeff).encode("utf-8") 

192 ) 

193 write_model_post_and_masses(fileobj, parameters) 

194 

195 else: 

196 # simple default parameters 

197 # that should always make the LAMMPS calculation run 

198 fileobj.write( 

199 "pair_style lj/cut 2.5 \n" 

200 "pair_coeff * * 1 1 \n" 

201 "mass * 1.0 \n".encode("utf-8") 

202 ) 

203 

204 if "group" in parameters: 

205 fileobj.write( 

206 ( 

207 "\n".join(["group {0}".format(p) for p in parameters["group"]]) 

208 + "\n" 

209 ).encode("utf-8") 

210 ) 

211 

212 fileobj.write("\n### run\n" "fix fix_nve all nve\n".encode("utf-8")) 

213 

214 if "fix" in parameters: 

215 fileobj.write( 

216 ( 

217 "\n".join(["fix {0}".format(p) for p in parameters["fix"]]) 

218 + "\n" 

219 ).encode("utf-8") 

220 ) 

221 

222 fileobj.write( 

223 "dump dump_all all custom {1} {0} id type x y z vx vy vz " 

224 "fx fy fz\n" 

225 "".format(lammps_trj, parameters["dump_period"]).encode("utf-8") 

226 ) 

227 fileobj.write( 

228 "thermo_style custom {0}\n" 

229 "thermo_modify flush yes format float %23.16g\n" 

230 "thermo 1\n".format(" ".join(parameters["thermo_args"])).encode( 

231 "utf-8" 

232 ) 

233 ) 

234 

235 if "timestep" in parameters: 

236 fileobj.write( 

237 "timestep {0}\n".format(parameters["timestep"]).encode("utf-8") 

238 ) 

239 

240 if "minimize" in parameters: 

241 fileobj.write( 

242 "minimize {0}\n".format(parameters["minimize"]).encode("utf-8") 

243 ) 

244 if "run" in parameters: 

245 fileobj.write("run {0}\n".format(parameters["run"]).encode("utf-8")) 

246 if not (("minimize" in parameters) or ("run" in parameters)): 

247 fileobj.write("run 0\n".encode("utf-8")) 

248 

249 fileobj.write( 

250 'print "{0}" \n'.format(CALCULATION_END_MARK).encode("utf-8") 

251 ) 

252 # Force LAMMPS to flush log 

253 fileobj.write("log /dev/stdout\n".encode("utf-8")) 

254 

255 fileobj.flush() 

256 if close_in_file: 

257 fileobj.close()