Coverage for /builds/kinetik161/ase/ase/calculators/lammps/inputwriter.py: 76.00%
100 statements
« prev ^ index » next coverage.py v7.2.7, created at 2023-12-10 11:04 +0000
1"""
2Stream input commands to lammps to perform desired simulations
3"""
4from ase.calculators.lammps.unitconvert import convert
5from ase.parallel import paropen
7# "End mark" used to indicate that the calculation is done
8CALCULATION_END_MARK = "__end_of_ase_invoked_calculation__"
11def lammps_create_atoms(fileobj, parameters, atoms, prismobj):
12 """Create atoms in lammps with 'create_box' and 'create_atoms'
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
22 """
23 if parameters["verbose"]:
24 fileobj.write("## Original ase cell\n")
25 fileobj.write(
26 "".join(
27 [
28 "# {:.16} {:.16} {:.16}\n".format(*x)
29 for x in atoms.get_cell()
30 ]
31 )
32 )
34 fileobj.write("lattice sc 1.0\n")
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)
40 if parameters["always_triclinic"] or prismobj.is_skewed():
41 fileobj.write(
42 "region asecell prism 0.0 {} 0.0 {} 0.0 {} ".format(
43 xhi, yhi, zhi
44 )
45 )
46 fileobj.write(
47 f"{xy} {xz} {yz} side in units box\n"
48 )
49 else:
50 fileobj.write(
51 "region asecell block 0.0 {} 0.0 {} 0.0 {} "
52 "side in units box\n".format(xhi, yhi, zhi)
53 )
55 symbols = atoms.get_chemical_symbols()
56 try:
57 # By request, specific atom type ordering
58 species = parameters["specorder"]
59 except AttributeError:
60 # By default, atom types in alphabetic order
61 species = sorted(set(symbols))
63 species_i = {s: i + 1 for i, s in enumerate(species)}
65 fileobj.write(
66 "create_box {} asecell\n" "".format(len(species))
67 )
68 for sym, pos in zip(symbols, atoms.get_positions()):
69 # Convert position from ASE units to LAMMPS units
70 pos = convert(pos, "distance", "ASE", parameters.units)
71 if parameters["verbose"]:
72 fileobj.write(
73 "# atom pos in ase cell: {:.16} {:.16} {:.16}\n"
74 "".format(*tuple(pos))
75 )
76 fileobj.write(
77 "create_atoms {} single {} {} {} remap yes units box\n".format(
78 *((species_i[sym],) + tuple(prismobj.vector_to_lammps(pos)))
79 )
80 )
83def write_lammps_in(lammps_in, parameters, atoms, prismobj,
84 lammps_trj=None, lammps_data=None):
85 """Write a LAMMPS in_ file with run parameters and settings."""
87 def write_model_post_and_masses(fileobj, parameters):
88 # write additional lines needed for some LAMMPS potentials
89 if 'model_post' in parameters:
90 mlines = parameters['model_post']
91 for ii in range(0, len(mlines)):
92 fileobj.write(mlines[ii])
94 if "masses" in parameters:
95 for mass in parameters["masses"]:
96 # Note that the variable mass is a string containing
97 # the type number and value of mass separated by a space
98 fileobj.write(f"mass {mass} \n")
100 if isinstance(lammps_in, str):
101 fileobj = paropen(lammps_in, "w")
102 close_in_file = True
103 else:
104 # Expect lammps_in to be a file-like object
105 fileobj = lammps_in
106 close_in_file = False
108 if parameters["verbose"]:
109 fileobj.write("# (written by ASE)\n")
111 # Write variables
112 fileobj.write(
113 (
114 "clear\n"
115 'variable dump_file string "{}"\n'
116 'variable data_file string "{}"\n'
117 ).format(lammps_trj, lammps_data)
118 )
120 if "package" in parameters:
121 fileobj.write(
122 "\n".join(
123 [f"package {p}" for p in parameters["package"]]
124 )
125 + "\n"
126 )
128 # setup styles except 'pair_style'
129 for style_type in ("atom", "bond", "angle",
130 "dihedral", "improper", "kspace"):
131 style = style_type + "_style"
132 if style in parameters:
133 fileobj.write(
134 '{} {} \n'.format(
135 style,
136 parameters[style]))
138 # write initialization lines needed for some LAMMPS potentials
139 if 'model_init' in parameters:
140 mlines = parameters['model_init']
141 for ii in range(0, len(mlines)):
142 fileobj.write(mlines[ii])
144 # write units
145 if 'units' in parameters:
146 units_line = 'units ' + parameters['units'] + '\n'
147 fileobj.write(units_line)
148 else:
149 fileobj.write('units metal\n')
151 pbc = atoms.get_pbc()
152 if "boundary" in parameters:
153 fileobj.write(
154 "boundary {} \n".format(parameters["boundary"])
155 )
156 else:
157 fileobj.write(
158 "boundary {} {} {} \n".format(
159 *tuple("sp"[int(x)] for x in pbc)
160 )
161 )
162 # Prior to version 22Dec2022, `box tilt large` is necessary to run systems
163 # with large tilts. Since version 22Dec2022, this command is ignored, and
164 # systems with large tilts can be run by default.
165 # https://docs.lammps.org/Commands_removed.html#box-command
166 # This command does not affect the efficiency for systems with small tilts
167 # and therefore worth written always.
168 fileobj.write("box tilt large \n")
169 fileobj.write("atom_modify sort 0 0.0 \n")
170 for key in ("neighbor", "newton"):
171 if key in parameters:
172 fileobj.write(
173 f"{key} {parameters[key]} \n"
174 )
175 fileobj.write("\n")
177 # write the simulation box and the atoms
178 if not lammps_data:
179 lammps_create_atoms(fileobj, parameters, atoms, prismobj)
180 # or simply refer to the data-file
181 else:
182 fileobj.write(f"read_data {lammps_data}\n")
184 # Write interaction stuff
185 fileobj.write("\n### interactions\n")
186 if "kim_interactions" in parameters:
187 fileobj.write(
188 "{}\n".format(
189 parameters["kim_interactions"]))
190 write_model_post_and_masses(fileobj, parameters)
192 elif ("pair_style" in parameters) and ("pair_coeff" in parameters):
193 pair_style = parameters["pair_style"]
194 fileobj.write(f"pair_style {pair_style} \n")
195 for pair_coeff in parameters["pair_coeff"]:
196 fileobj.write(
197 "pair_coeff {} \n" "".format(pair_coeff)
198 )
199 write_model_post_and_masses(fileobj, parameters)
201 else:
202 # simple default parameters
203 # that should always make the LAMMPS calculation run
204 fileobj.write(
205 "pair_style lj/cut 2.5 \n"
206 "pair_coeff * * 1 1 \n"
207 "mass * 1.0 \n"
208 )
210 if "group" in parameters:
211 fileobj.write(
212 "\n".join([f"group {p}" for p in parameters["group"]])
213 + "\n"
214 )
216 fileobj.write("\n### run\n" "fix fix_nve all nve\n")
218 if "fix" in parameters:
219 fileobj.write(
220 "\n".join([f"fix {p}" for p in parameters["fix"]])
221 + "\n"
222 )
224 fileobj.write(
225 "dump dump_all all custom {1} {0} id type x y z vx vy vz "
226 "fx fy fz\n"
227 "".format(lammps_trj, parameters["dump_period"])
228 )
229 fileobj.write(
230 "thermo_style custom {}\n"
231 "thermo_modify flush yes format float %23.16g\n"
232 "thermo 1\n".format(" ".join(parameters["thermo_args"]))
233 )
235 if "timestep" in parameters:
236 fileobj.write(
237 "timestep {}\n".format(parameters["timestep"])
238 )
240 if "minimize" in parameters:
241 fileobj.write(
242 "minimize {}\n".format(parameters["minimize"])
243 )
244 if "run" in parameters:
245 fileobj.write("run {}\n".format(parameters["run"]))
246 if not (("minimize" in parameters) or ("run" in parameters)):
247 fileobj.write("run 0\n")
249 fileobj.write(
250 f'print "{CALCULATION_END_MARK}" \n'
251 )
252 # Force LAMMPS to flush log
253 fileobj.write("log /dev/stdout\n")
255 fileobj.flush()
256 if close_in_file:
257 fileobj.close()