Coverage for /builds/ase/ase/ase/calculators/calculator.py : 81.82%
Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1import os
2import copy
3import subprocess
4from math import pi, sqrt
5from pathlib import Path
6from typing import Union, Optional, List, Set, Dict, Any
7import warnings
8from abc import abstractmethod
10import numpy as np
12from ase.cell import Cell
13from ase.outputs import Properties, all_outputs
14from ase.utils import jsonable
15from ase.calculators.abc import GetPropertiesMixin
17from .names import names
20class CalculatorError(RuntimeError):
21 """Base class of error types related to ASE calculators."""
24class CalculatorSetupError(CalculatorError):
25 """Calculation cannot be performed with the given parameters.
27 Reasons to raise this errors are:
28 * The calculator is not properly configured
29 (missing executable, environment variables, ...)
30 * The given atoms object is not supported
31 * Calculator parameters are unsupported
33 Typically raised before a calculation."""
36class EnvironmentError(CalculatorSetupError):
37 """Raised if calculator is not properly set up with ASE.
39 May be missing an executable or environment variables."""
42class InputError(CalculatorSetupError):
43 """Raised if inputs given to the calculator were incorrect.
45 Bad input keywords or values, or missing pseudopotentials.
47 This may be raised before or during calculation, depending on
48 when the problem is detected."""
51class CalculationFailed(CalculatorError):
52 """Calculation failed unexpectedly.
54 Reasons to raise this error are:
55 * Calculation did not converge
56 * Calculation ran out of memory
57 * Segmentation fault or other abnormal termination
58 * Arithmetic trouble (singular matrices, NaN, ...)
60 Typically raised during calculation."""
63class SCFError(CalculationFailed):
64 """SCF loop did not converge."""
67class ReadError(CalculatorError):
68 """Unexpected irrecoverable error while reading calculation results."""
71class PropertyNotImplementedError(NotImplementedError):
72 """Raised if a calculator does not implement the requested property."""
75class PropertyNotPresent(CalculatorError):
76 """Requested property is missing.
78 Maybe it was never calculated, or for some reason was not extracted
79 with the rest of the results, without being a fatal ReadError."""
82def compare_atoms(atoms1, atoms2, tol=1e-15, excluded_properties=None):
83 """Check for system changes since last calculation. Properties in
84 ``excluded_properties`` are not checked."""
85 if atoms1 is None:
86 system_changes = all_changes[:]
87 else:
88 system_changes = []
90 properties_to_check = set(all_changes)
91 if excluded_properties:
92 properties_to_check -= set(excluded_properties)
94 # Check properties that aren't in Atoms.arrays but are attributes of
95 # Atoms objects
96 for prop in ['cell', 'pbc']:
97 if prop in properties_to_check:
98 properties_to_check.remove(prop)
99 if not equal(getattr(atoms1, prop), getattr(atoms2, prop),
100 atol=tol):
101 system_changes.append(prop)
103 arrays1 = set(atoms1.arrays)
104 arrays2 = set(atoms2.arrays)
106 # Add any properties that are only in atoms1.arrays or only in
107 # atoms2.arrays (and aren't excluded). Note that if, e.g. arrays1 has
108 # `initial_charges` which is merely zeros and arrays2 does not have
109 # this array, we'll still assume that the system has changed. However,
110 # this should only occur rarely.
111 system_changes += properties_to_check & (arrays1 ^ arrays2)
113 # Finally, check all of the non-excluded properties shared by the atoms
114 # arrays.
115 for prop in properties_to_check & arrays1 & arrays2:
116 if not equal(atoms1.arrays[prop], atoms2.arrays[prop], atol=tol):
117 system_changes.append(prop)
119 return system_changes
122all_properties = ['energy', 'forces', 'stress', 'stresses', 'dipole',
123 'charges', 'magmom', 'magmoms', 'free_energy', 'energies',
124 'dielectric_tensor', 'born_effective_charges', 'polarization']
127all_changes = ['positions', 'numbers', 'cell', 'pbc',
128 'initial_charges', 'initial_magmoms']
131special = {'cp2k': 'CP2K',
132 'demonnano': 'DemonNano',
133 'dftd3': 'DFTD3',
134 'dmol': 'DMol3',
135 'eam': 'EAM',
136 'elk': 'ELK',
137 'emt': 'EMT',
138 'exciting': 'ExcitingGroundStateCalculator',
139 'crystal': 'CRYSTAL',
140 'ff': 'ForceField',
141 'gamess_us': 'GAMESSUS',
142 'gulp': 'GULP',
143 'kim': 'KIM',
144 'lammpsrun': 'LAMMPS',
145 'lammpslib': 'LAMMPSlib',
146 'lj': 'LennardJones',
147 'mopac': 'MOPAC',
148 'morse': 'MorsePotential',
149 'nwchem': 'NWChem',
150 'openmx': 'OpenMX',
151 'orca': 'ORCA',
152 'qchem': 'QChem',
153 'tip3p': 'TIP3P',
154 'tip4p': 'TIP4P'}
157external_calculators = {}
160def register_calculator_class(name, cls):
161 """ Add the class into the database. """
162 assert name not in external_calculators
163 external_calculators[name] = cls
164 names.append(name)
165 names.sort()
168def get_calculator_class(name):
169 """Return calculator class."""
170 if name == 'asap':
171 from asap3 import EMT as Calculator
172 elif name == 'gpaw':
173 from gpaw import GPAW as Calculator
174 elif name == 'hotbit':
175 from hotbit import Calculator
176 elif name == 'vasp2':
177 from ase.calculators.vasp import Vasp2 as Calculator
178 elif name == 'ace':
179 from ase.calculators.acemolecule import ACE as Calculator
180 elif name == 'Psi4':
181 from ase.calculators.psi4 import Psi4 as Calculator
182 elif name in external_calculators:
183 Calculator = external_calculators[name]
184 else:
185 classname = special.get(name, name.title())
186 module = __import__('ase.calculators.' + name, {}, None, [classname])
187 Calculator = getattr(module, classname)
188 return Calculator
191def equal(a, b, tol=None, rtol=None, atol=None):
192 """ndarray-enabled comparison function."""
193 # XXX Known bugs:
194 # * Comparing cell objects (pbc not part of array representation)
195 # * Infinite recursion for cyclic dicts
196 # * Can of worms is open
197 if tol is not None:
198 msg = 'Use `equal(a, b, rtol=..., atol=...)` instead of `tol=...`'
199 warnings.warn(msg, DeprecationWarning)
200 assert rtol is None and atol is None, \
201 'Do not use deprecated `tol` with `atol` and/or `rtol`'
202 rtol = tol
203 atol = tol
205 a_is_dict = isinstance(a, dict)
206 b_is_dict = isinstance(b, dict)
207 if a_is_dict or b_is_dict:
208 # Check that both a and b are dicts
209 if not (a_is_dict and b_is_dict):
210 return False
211 if a.keys() != b.keys():
212 return False
213 return all(equal(a[key], b[key], rtol=rtol, atol=atol) for key in a)
215 if np.shape(a) != np.shape(b):
216 return False
218 if rtol is None and atol is None:
219 return np.array_equal(a, b)
221 if rtol is None:
222 rtol = 0
223 if atol is None:
224 atol = 0
226 return np.allclose(a, b, rtol=rtol, atol=atol)
229def kptdensity2monkhorstpack(atoms, kptdensity=3.5, even=True):
230 """Convert k-point density to Monkhorst-Pack grid size.
232 atoms: Atoms object
233 Contains unit cell and information about boundary conditions.
234 kptdensity: float
235 Required k-point density. Default value is 3.5 point per Ang^-1.
236 even: bool
237 Round up to even numbers.
238 """
240 recipcell = atoms.cell.reciprocal()
241 kpts = []
242 for i in range(3):
243 if atoms.pbc[i]:
244 k = 2 * pi * sqrt((recipcell[i]**2).sum()) * kptdensity
245 if even:
246 kpts.append(2 * int(np.ceil(k / 2)))
247 else:
248 kpts.append(int(np.ceil(k)))
249 else:
250 kpts.append(1)
251 return np.array(kpts)
254def kpts2mp(atoms, kpts, even=False):
255 if kpts is None:
256 return np.array([1, 1, 1])
257 if isinstance(kpts, (float, int)):
258 return kptdensity2monkhorstpack(atoms, kpts, even)
259 else:
260 return kpts
263def kpts2sizeandoffsets(size=None, density=None, gamma=None, even=None,
264 atoms=None):
265 """Helper function for selecting k-points.
267 Use either size or density.
269 size: 3 ints
270 Number of k-points.
271 density: float
272 K-point density in units of k-points per Ang^-1.
273 gamma: None or bool
274 Should the Gamma-point be included? Yes / no / don't care:
275 True / False / None.
276 even: None or bool
277 Should the number of k-points be even? Yes / no / don't care:
278 True / False / None.
279 atoms: Atoms object
280 Needed for calculating k-point density.
282 """
284 if size is not None and density is not None:
285 raise ValueError('Cannot specify k-point mesh size and '
286 'density simultaneously')
287 elif density is not None and atoms is None:
288 raise ValueError('Cannot set k-points from "density" unless '
289 'Atoms are provided (need BZ dimensions).')
291 if size is None:
292 if density is None:
293 size = [1, 1, 1]
294 else:
295 size = kptdensity2monkhorstpack(atoms, density, None)
297 # Not using the rounding from kptdensity2monkhorstpack as it doesn't do
298 # rounding to odd numbers
299 if even is not None:
300 size = np.array(size)
301 remainder = size % 2
302 if even:
303 size += remainder
304 else: # Round up to odd numbers
305 size += (1 - remainder)
307 offsets = [0, 0, 0]
308 if atoms is None:
309 pbc = [True, True, True]
310 else:
311 pbc = atoms.pbc
313 if gamma is not None:
314 for i, s in enumerate(size):
315 if pbc[i] and s % 2 != bool(gamma):
316 offsets[i] = 0.5 / s
318 return size, offsets
321@jsonable('kpoints')
322class KPoints:
323 def __init__(self, kpts=None):
324 if kpts is None:
325 kpts = np.zeros((1, 3))
326 self.kpts = kpts
328 def todict(self):
329 return vars(self)
332def kpts2kpts(kpts, atoms=None):
333 from ase.dft.kpoints import monkhorst_pack
335 if kpts is None:
336 return KPoints()
338 if hasattr(kpts, 'kpts'):
339 return kpts
341 if isinstance(kpts, dict):
342 if 'kpts' in kpts:
343 return KPoints(kpts['kpts'])
344 if 'path' in kpts:
345 cell = Cell.ascell(atoms.cell)
346 return cell.bandpath(pbc=atoms.pbc, **kpts)
347 size, offsets = kpts2sizeandoffsets(atoms=atoms, **kpts)
348 return KPoints(monkhorst_pack(size) + offsets)
350 if isinstance(kpts[0], int):
351 return KPoints(monkhorst_pack(kpts))
353 return KPoints(np.array(kpts))
356def kpts2ndarray(kpts, atoms=None):
357 """Convert kpts keyword to 2-d ndarray of scaled k-points."""
358 return kpts2kpts(kpts, atoms=atoms).kpts
361class EigenvalOccupationMixin:
362 """Define 'eigenvalues' and 'occupations' properties on class.
364 eigenvalues and occupations will be arrays of shape (spin, kpts, nbands).
366 Classes must implement the old-fashioned get_eigenvalues and
367 get_occupations methods."""
368 # We should maybe deprecate this and rely on the new
369 # Properties object for eigenvalues/occupations.
371 @property
372 def eigenvalues(self):
373 return self._propwrapper().eigenvalues
375 @property
376 def occupations(self):
377 return self._propwrapper().occupations
379 def _propwrapper(self):
380 from ase.calculator.singlepoint import OutputPropertyWrapper
381 return OutputPropertyWrapper(self)
384class Parameters(dict):
385 """Dictionary for parameters.
387 Special feature: If param is a Parameters instance, then param.xc
388 is a shorthand for param['xc'].
389 """
391 def __getattr__(self, key):
392 if key not in self:
393 return dict.__getattribute__(self, key)
394 return self[key]
396 def __setattr__(self, key, value):
397 self[key] = value
399 @classmethod
400 def read(cls, filename):
401 """Read parameters from file."""
402 # We use ast to evaluate literals, avoiding eval()
403 # for security reasons.
404 import ast
405 with open(filename) as fd:
406 txt = fd.read().strip()
407 assert txt.startswith('dict(')
408 assert txt.endswith(')')
409 txt = txt[5:-1]
411 # The tostring() representation "dict(...)" is not actually
412 # a literal, so we manually parse that along with the other
413 # formatting that we did manually:
414 dct = {}
415 for line in txt.splitlines():
416 key, val = line.split('=', 1)
417 key = key.strip()
418 val = val.strip()
419 if val[-1] == ',':
420 val = val[:-1]
421 dct[key] = ast.literal_eval(val)
423 parameters = cls(dct)
424 return parameters
426 def tostring(self):
427 keys = sorted(self)
428 return 'dict(' + ',\n '.join(
429 '{}={!r}'.format(key, self[key]) for key in keys) + ')\n'
431 def write(self, filename):
432 Path(filename).write_text(self.tostring())
435class BaseCalculator(GetPropertiesMixin):
436 implemented_properties: List[str] = []
437 'Properties calculator can handle (energy, forces, ...)'
439 # Placeholder object for deprecated arguments. Let deprecated keywords
440 # default to _deprecated and then issue a warning if the user passed
441 # any other object (such as None).
442 _deprecated = object()
444 def __init__(self, parameters=None, use_cache=True):
445 if parameters is None:
446 parameters = {}
448 self.parameters = dict(parameters)
449 self.atoms = None
450 self.results = {}
451 self.use_cache = use_cache
453 def calculate_properties(self, atoms, properties):
454 """This method is experimental; currently for internal use."""
455 for name in properties:
456 if name not in all_outputs:
457 raise ValueError(f'No such property: {name}')
459 # We ignore system changes for now.
460 self.calculate(atoms, properties, system_changes=all_changes)
462 props = self.export_properties()
464 for name in properties:
465 if name not in props:
466 raise PropertyNotPresent(name)
467 return props
469 @abstractmethod
470 def calculate(self, atoms, properties, system_changes):
471 ...
473 def check_state(self, atoms, tol=1e-15):
474 """Check for any system changes since last calculation."""
475 if self.use_cache:
476 return compare_atoms(self.atoms, atoms, tol=tol)
477 else:
478 return all_changes
480 def get_property(self, name, atoms=None, allow_calculation=True):
481 if name not in self.implemented_properties:
482 raise PropertyNotImplementedError('{} property not implemented'
483 .format(name))
485 if atoms is None:
486 atoms = self.atoms
487 system_changes = []
488 else:
489 system_changes = self.check_state(atoms)
491 if system_changes:
492 self.atoms = None
493 self.results = {}
495 if name not in self.results:
496 if not allow_calculation:
497 return None
499 if self.use_cache:
500 self.atoms = atoms.copy()
502 self.calculate(atoms, [name], system_changes)
504 if name not in self.results:
505 # For some reason the calculator was not able to do what we want,
506 # and that is OK.
507 raise PropertyNotImplementedError('{} not present in this '
508 'calculation'.format(name))
510 result = self.results[name]
511 if isinstance(result, np.ndarray):
512 result = result.copy()
513 return result
515 def calculation_required(self, atoms, properties):
516 assert not isinstance(properties, str)
517 system_changes = self.check_state(atoms)
518 if system_changes:
519 return True
520 for name in properties:
521 if name not in self.results:
522 return True
523 return False
525 def export_properties(self):
526 return Properties(self.results)
528 def _get_name(self) -> str: # child class can override this
529 return self.__class__.__name__.lower()
531 @property
532 def name(self) -> str:
533 return self._get_name()
536class Calculator(BaseCalculator):
537 """Base-class for all ASE calculators.
539 A calculator must raise PropertyNotImplementedError if asked for a
540 property that it can't calculate. So, if calculation of the
541 stress tensor has not been implemented, get_stress(atoms) should
542 raise PropertyNotImplementedError. This can be achieved simply by not
543 including the string 'stress' in the list implemented_properties
544 which is a class member. These are the names of the standard
545 properties: 'energy', 'forces', 'stress', 'dipole', 'charges',
546 'magmom' and 'magmoms'.
547 """
549 default_parameters: Dict[str, Any] = {}
550 'Default parameters'
552 ignored_changes: Set[str] = set()
553 'Properties of Atoms which we ignore for the purposes of cache '
554 'invalidation with check_state().'
556 discard_results_on_any_change = False
557 'Whether we purge the results following any change in the set() method. '
558 'Most (file I/O) calculators will probably want this.'
560 def __init__(self, restart=None,
561 ignore_bad_restart_file=BaseCalculator._deprecated,
562 label=None, atoms=None, directory='.',
563 **kwargs):
564 """Basic calculator implementation.
566 restart: str
567 Prefix for restart file. May contain a directory. Default
568 is None: don't restart.
569 ignore_bad_restart_file: bool
570 Deprecated, please do not use.
571 Passing more than one positional argument to Calculator()
572 is deprecated and will stop working in the future.
573 Ignore broken or missing restart file. By default, it is an
574 error if the restart file is missing or broken.
575 directory: str or PurePath
576 Working directory in which to read and write files and
577 perform calculations.
578 label: str
579 Name used for all files. Not supported by all calculators.
580 May contain a directory, but please use the directory parameter
581 for that instead.
582 atoms: Atoms object
583 Optional Atoms object to which the calculator will be
584 attached. When restarting, atoms will get its positions and
585 unit-cell updated from file.
586 """
587 self.atoms = None # copy of atoms object from last calculation
588 self.results = {} # calculated properties (energy, forces, ...)
589 self.parameters = None # calculational parameters
590 self._directory = None # Initialize
592 if ignore_bad_restart_file is self._deprecated:
593 ignore_bad_restart_file = False
594 else:
595 warnings.warn(FutureWarning(
596 'The keyword "ignore_bad_restart_file" is deprecated and '
597 'will be removed in a future version of ASE. Passing more '
598 'than one positional argument to Calculator is also '
599 'deprecated and will stop functioning in the future. '
600 'Please pass arguments by keyword (key=value) except '
601 'optionally the "restart" keyword.'
602 ))
604 if restart is not None:
605 try:
606 self.read(restart) # read parameters, atoms and results
607 except ReadError:
608 if ignore_bad_restart_file:
609 self.reset()
610 else:
611 raise
613 self.directory = directory
614 self.prefix = None
615 if label is not None:
616 if self.directory == '.' and '/' in label:
617 # We specified directory in label, and nothing in the diretory
618 # key
619 self.label = label
620 elif '/' not in label:
621 # We specified our directory in the directory keyword
622 # or not at all
623 self.label = '/'.join((self.directory, label))
624 else:
625 raise ValueError('Directory redundantly specified though '
626 'directory="{}" and label="{}". '
627 'Please omit "/" in label.'
628 .format(self.directory, label))
630 if self.parameters is None:
631 # Use default parameters if they were not read from file:
632 self.parameters = self.get_default_parameters()
634 if atoms is not None:
635 atoms.calc = self
636 if self.atoms is not None:
637 # Atoms were read from file. Update atoms:
638 if not (equal(atoms.numbers, self.atoms.numbers) and
639 (atoms.pbc == self.atoms.pbc).all()):
640 raise CalculatorError('Atoms not compatible with file')
641 atoms.positions = self.atoms.positions
642 atoms.cell = self.atoms.cell
644 self.set(**kwargs)
646 if not hasattr(self, 'get_spin_polarized'):
647 self.get_spin_polarized = self._deprecated_get_spin_polarized
648 # XXX We are very naughty and do not call super constructor!
650 # For historical reasons we have a particular caching protocol.
651 # We disable the superclass' optional cache.
652 self.use_cache = False
654 @property
655 def directory(self) -> str:
656 return self._directory
658 @directory.setter
659 def directory(self, directory: Union[str, os.PathLike]):
660 self._directory = str(Path(directory)) # Normalize path.
662 @property
663 def label(self):
664 if self.directory == '.':
665 return self.prefix
667 # Generally, label ~ directory/prefix
668 #
669 # We use '/' rather than os.pathsep because
670 # 1) directory/prefix does not represent any actual path
671 # 2) We want the same string to work the same on all platforms
672 if self.prefix is None:
673 return self.directory + '/'
675 return '{}/{}'.format(self.directory, self.prefix)
677 @label.setter
678 def label(self, label):
679 if label is None:
680 self.directory = '.'
681 self.prefix = None
682 return
684 tokens = label.rsplit('/', 1)
685 if len(tokens) == 2:
686 directory, prefix = tokens
687 else:
688 assert len(tokens) == 1
689 directory = '.'
690 prefix = tokens[0]
691 if prefix == '':
692 prefix = None
693 self.directory = directory
694 self.prefix = prefix
696 def set_label(self, label):
697 """Set label and convert label to directory and prefix.
699 Examples:
701 * label='abc': (directory='.', prefix='abc')
702 * label='dir1/abc': (directory='dir1', prefix='abc')
703 * label=None: (directory='.', prefix=None)
704 """
705 self.label = label
707 def get_default_parameters(self):
708 return Parameters(copy.deepcopy(self.default_parameters))
710 def todict(self, skip_default=True):
711 defaults = self.get_default_parameters()
712 dct = {}
713 for key, value in self.parameters.items():
714 if hasattr(value, 'todict'):
715 value = value.todict()
716 if skip_default:
717 default = defaults.get(key, '_no_default_')
718 if default != '_no_default_' and equal(value, default):
719 continue
720 dct[key] = value
721 return dct
723 def reset(self):
724 """Clear all information from old calculation."""
726 self.atoms = None
727 self.results = {}
729 def read(self, label):
730 """Read atoms, parameters and calculated properties from output file.
732 Read result from self.label file. Raise ReadError if the file
733 is not there. If the file is corrupted or contains an error
734 message from the calculation, a ReadError should also be
735 raised. In case of succes, these attributes must set:
737 atoms: Atoms object
738 The state of the atoms from last calculation.
739 parameters: Parameters object
740 The parameter dictionary.
741 results: dict
742 Calculated properties like energy and forces.
744 The FileIOCalculator.read() method will typically read atoms
745 and parameters and get the results dict by calling the
746 read_results() method."""
748 self.set_label(label)
750 def get_atoms(self):
751 if self.atoms is None:
752 raise ValueError('Calculator has no atoms')
753 atoms = self.atoms.copy()
754 atoms.calc = self
755 return atoms
757 @classmethod
758 def read_atoms(cls, restart, **kwargs):
759 return cls(restart=restart, label=restart, **kwargs).get_atoms()
761 def set(self, **kwargs):
762 """Set parameters like set(key1=value1, key2=value2, ...).
764 A dictionary containing the parameters that have been changed
765 is returned.
767 Subclasses must implement a set() method that will look at the
768 chaneged parameters and decide if a call to reset() is needed.
769 If the changed parameters are harmless, like a change in
770 verbosity, then there is no need to call reset().
772 The special keyword 'parameters' can be used to read
773 parameters from a file."""
775 if 'parameters' in kwargs:
776 filename = kwargs.pop('parameters')
777 parameters = Parameters.read(filename)
778 parameters.update(kwargs)
779 kwargs = parameters
781 changed_parameters = {}
783 for key, value in kwargs.items():
784 oldvalue = self.parameters.get(key)
785 if key not in self.parameters or not equal(value, oldvalue):
786 changed_parameters[key] = value
787 self.parameters[key] = value
789 if self.discard_results_on_any_change and changed_parameters:
790 self.reset()
791 return changed_parameters
793 def check_state(self, atoms, tol=1e-15):
794 """Check for any system changes since last calculation."""
795 return compare_atoms(self.atoms, atoms, tol=tol,
796 excluded_properties=set(self.ignored_changes))
798 def calculate(self, atoms=None, properties=['energy'],
799 system_changes=all_changes):
800 """Do the calculation.
802 properties: list of str
803 List of what needs to be calculated. Can be any combination
804 of 'energy', 'forces', 'stress', 'dipole', 'charges', 'magmom'
805 and 'magmoms'.
806 system_changes: list of str
807 List of what has changed since last calculation. Can be
808 any combination of these six: 'positions', 'numbers', 'cell',
809 'pbc', 'initial_charges' and 'initial_magmoms'.
811 Subclasses need to implement this, but can ignore properties
812 and system_changes if they want. Calculated properties should
813 be inserted into results dictionary like shown in this dummy
814 example::
816 self.results = {'energy': 0.0,
817 'forces': np.zeros((len(atoms), 3)),
818 'stress': np.zeros(6),
819 'dipole': np.zeros(3),
820 'charges': np.zeros(len(atoms)),
821 'magmom': 0.0,
822 'magmoms': np.zeros(len(atoms))}
824 The subclass implementation should first call this
825 implementation to set the atoms attribute and create any missing
826 directories.
827 """
828 if atoms is not None:
829 self.atoms = atoms.copy()
830 if not os.path.isdir(self._directory):
831 try:
832 os.makedirs(self._directory)
833 except FileExistsError as e:
834 # We can only end up here in case of a race condition if
835 # multiple Calculators are running concurrently *and* use the
836 # same _directory, which cannot be expected to work anyway.
837 msg = ('Concurrent use of directory ' + self._directory +
838 'by multiple Calculator instances detected. Please '
839 'use one directory per instance.')
840 raise RuntimeError(msg) from e
842 def calculate_numerical_forces(self, atoms, d=0.001):
843 """Calculate numerical forces using finite difference.
845 All atoms will be displaced by +d and -d in all directions."""
846 from ase.calculators.test import numeric_forces
847 return numeric_forces(atoms, d=d)
849 def calculate_numerical_stress(self, atoms, d=1e-6, voigt=True):
850 """Calculate numerical stress using finite difference."""
851 from ase.calculators.test import numeric_stress
852 return numeric_stress(atoms, d=d, voigt=voigt)
854 def _deprecated_get_spin_polarized(self):
855 msg = ('This calculator does not implement get_spin_polarized(). '
856 'In the future, calc.get_spin_polarized() will work only on '
857 'calculator classes that explicitly implement this method or '
858 'inherit the method via specialized subclasses.')
859 warnings.warn(msg, FutureWarning)
860 return False
862 def band_structure(self):
863 """Create band-structure object for plotting."""
864 from ase.spectrum.band_structure import get_band_structure
865 # XXX This calculator is supposed to just have done a band structure
866 # calculation, but the calculator may not have the correct Fermi level
867 # if it updated the Fermi level after changing k-points.
868 # This will be a problem with some calculators (currently GPAW), and
869 # the user would have to override this by providing the Fermi level
870 # from the selfconsistent calculation.
871 return get_band_structure(calc=self)
874class FileIOCalculator(Calculator):
875 """Base class for calculators that write/read input/output files."""
877 command: Optional[str] = None
878 'Command used to start calculation'
880 def __init__(self, restart=None,
881 ignore_bad_restart_file=Calculator._deprecated,
882 label=None, atoms=None, command=None, **kwargs):
883 """File-IO calculator.
885 command: str
886 Command used to start calculation.
887 """
889 Calculator.__init__(self, restart, ignore_bad_restart_file, label,
890 atoms, **kwargs)
892 if command is not None:
893 self.command = command
894 else:
895 name = 'ASE_' + self.name.upper() + '_COMMAND'
896 self.command = os.environ.get(name, self.command)
898 def calculate(self, atoms=None, properties=['energy'],
899 system_changes=all_changes):
900 Calculator.calculate(self, atoms, properties, system_changes)
901 self.write_input(self.atoms, properties, system_changes)
902 self.execute()
903 self.read_results()
905 def execute(self):
906 if self.command is None:
907 raise CalculatorSetupError(
908 'Please set ${} environment variable '
909 .format('ASE_' + self.name.upper() + '_COMMAND') +
910 'or supply the command keyword')
911 command = self.command
912 if 'PREFIX' in command:
913 command = command.replace('PREFIX', self.prefix)
915 try:
916 proc = subprocess.Popen(command, shell=True, cwd=self.directory)
917 except OSError as err:
918 # Actually this may never happen with shell=True, since
919 # probably the shell launches successfully. But we soon want
920 # to allow calling the subprocess directly, and then this
921 # distinction (failed to launch vs failed to run) is useful.
922 msg = 'Failed to execute "{}"'.format(command)
923 raise EnvironmentError(msg) from err
925 errorcode = proc.wait()
927 if errorcode:
928 path = os.path.abspath(self.directory)
929 msg = ('Calculator "{}" failed with command "{}" failed in '
930 '{} with error code {}'.format(self.name, command,
931 path, errorcode))
932 raise CalculationFailed(msg)
934 def write_input(self, atoms, properties=None, system_changes=None):
935 """Write input file(s).
937 Call this method first in subclasses so that directories are
938 created automatically."""
940 absdir = os.path.abspath(self.directory)
941 if absdir != os.curdir and not os.path.isdir(self.directory):
942 os.makedirs(self.directory)
944 def read_results(self):
945 """Read energy, forces, ... from output file(s)."""
946 pass