ASE calculators¶
calorine provides two ASE calculators for NEP calculations, one that uses the GPU implementation and one that uses the CPU implementation of NEP. For smaller calculations the CPU calculators is usually more performant. For very large simulations and for comparison the GPU calculator can be useful as well.
GPU calculator¶
- class calorine.calculators.GPUNEP(model_filename, directory=None, label='GPUNEP', atoms=None, command='gpumd')[source]¶
This class provides an ASE calculator for NEP calculations with GPUMD.
This calculator writes files that are input to the
gpumd
executable. It is thus likely to be slow if many calculations are to be performed.- Parameters:
model_filename (str) – Path to file in
nep.txt
format with model parameters.directory (str) – Directory to run GPUMD in. If None, a temporary directory will be created and removed once the calculations are finished. If specified, the directory will not be deleted. In the latter case, it is advisable to do no more than one calculation with this calculator (unless you know exactly what you are doing).
label (str) – Label for this calculator.
atoms (Atoms) – Atoms to attach to this calculator.
command (str) – Command to run GPUMD with. Default:
gpumd
prefix (str) – Filename (excluding file ending) for run script. Default:
to_run
Example
>>> calc = GPUNEP('nep.txt') >>> atoms.calc = calc >>> atoms.get_potential_energy()
-
command:
Optional
[str
] = 'gpumd'¶ Command used to start calculation
- run_custom_md(parameters, return_last_atoms=False, only_prepare=False)[source]¶
Run a custom MD simulation.
- Parameters:
parameters (
List
[Tuple
[str
,Any
]]) –Parameters to be specified in the run.in file. The potential keyword is set automatically, all other keywords need to be set via this argument. Example:
[('dump_thermo', 100), ('dump_position', 1000), ('velocity', 300), ('time_step', 1), ('ensemble', ['nvt_ber', 300, 300, 100]), ('run', 10000)]
return_last_atoms (
bool
) – IfTrue
the last saved snapshot will be returned.only_prepare (
bool
) –- If
True
the necessary input files will be written but theMD run will not be executed.
- If
- Returns:
The last snapshot if
return_last_atoms
isTrue
.
CPU calculator¶
- class calorine.calculators.CPUNEP(model_filename, atoms=None, label=None, debug=False)[source]¶
This class provides an ASE calculator for nep_cpu, the in-memory CPU implementation of GPUMD.
- Parameters:
model_filename (str) – Path to file in
nep.txt
format with model parametersatoms (Atoms) – Atoms to attach the calculator to
label (str) – Label for this calclator
debug (bool, optional) – Flag to toggle debug mode. Prints GPUMD output. Defaults to False.
- Raises:
FileNotFoundError – Raises
FileNotFoundError
ifmodel_filename
does not point to a valid file.ValueError – Raises
ValueError
atoms are not defined when trying to get energies and forces.
Example
>>> calc = CPUNEP('nep.txt') >>> atoms.calc = calc >>> atoms.get_potential_energy()