Coverage for ase / md / logger.py: 76.36%

1# fmt: off 

2 

3"""Logging for molecular dynamics.""" 

4import weakref 

5from typing import IO, Any 

6 

7from ase import Atoms, units 

8from ase.parallel import world 

9from ase.utils import IOContext 

10 

11 

12class MDLogger(IOContext): 

13 """Class for logging molecular dynamics simulations. 

14 

15 Parameters 

16 ---------- 

17 dyn : :class:`~ase.md.md.MolecularDynamics` 

18 ASE :class:`~ase.md.md.MolecularDynamics` object. 

19 Only a weak reference is kept. 

20 atoms : :class:`~ase.Atoms` 

21 ASE :class:`~ase.Atoms` object. 

22 logfile : IO | str 

23 File name or open file, "-" meaning standard output. 

24 stress : bool, default: :obj:`False` 

25 If :obj:`True`, the logger writes stress. 

26 peratom : bool, default: :obj:`False` 

27 If :obj:`True`, the logger writes energies per atom. 

28 mode : str, default: ``'a'`` 

29 How the file is opened if logfile is a filename. 

30 

31 """ 

32 

33 def __init__( 

34 self, 

35 dyn: Any, # not fully annotated so far to avoid a circular import 

36 atoms: Atoms, 

37 logfile: IO | str, 

38 header: bool = True, 

39 stress: bool = False, 

40 peratom: bool = False, 

41 mode: str = "a", 

42 comm=world, 

43 ): 

44 self.dyn = weakref.proxy(dyn) if hasattr(dyn, "get_time") else None 

45 self.atoms = atoms 

46 global_natoms = atoms.get_global_number_of_atoms() 

47 self.logfile = self.openfile(file=logfile, mode=mode, comm=comm) 

48 self.stress = stress 

49 self.peratom = peratom 

50 if self.dyn is not None: 

51 self.hdr = "%-9s " % ("Time[ps]",) 

52 self.fmt = "%-10.4f " 

53 else: 

54 self.hdr = "" 

55 self.fmt = "" 

56 if self.peratom: 

57 self.hdr += "%12s %12s %12s %6s" % ("Etot/N[eV]", "Epot/N[eV]", 

58 "Ekin/N[eV]", "T[K]") 

59 self.fmt += "%12.4f %12.4f %12.4f %6.1f" 

60 else: 

61 self.hdr += "%12s %12s %12s %6s" % ("Etot[eV]", "Epot[eV]", 

62 "Ekin[eV]", "T[K]") 

63 # Choose a sensible number of decimals 

64 if global_natoms <= 100: 

65 digits = 4 

66 elif global_natoms <= 1000: 

67 digits = 3 

68 elif global_natoms <= 10000: 

69 digits = 2 

70 else: 

71 digits = 1 

72 self.fmt += 3 * ("%%12.%df " % (digits,)) + " %6.1f" 

73 if self.stress: 

74 self.hdr += (' ---------------------- stress [GPa] ' 

75 '-----------------------') 

76 self.fmt += 6 * " %10.3f" 

77 self.fmt += "\n" 

78 if header: 

79 self.logfile.write(self.hdr + "\n") 

80 

81 def __del__(self): 

82 self.close() 

83 

84 def __call__(self): 

85 epot = self.atoms.get_potential_energy() 

86 ekin = self.atoms.get_kinetic_energy() 

87 temp = self.atoms.get_temperature() 

88 global_natoms = self.atoms.get_global_number_of_atoms() 

89 if self.peratom: 

90 epot /= global_natoms 

91 ekin /= global_natoms 

92 if self.dyn is not None: 

93 t = self.dyn.get_time() / (1000 * units.fs) 

94 dat = (t,) 

95 else: 

96 dat = () 

97 dat += (epot + ekin, epot, ekin, temp) 

98 if self.stress: 

99 dat += tuple(self.atoms.get_stress( 

100 include_ideal_gas=True) / units.GPa) 

101 self.logfile.write(self.fmt % dat) 

102 self.logfile.flush()