Coverage for /builds/ase/ase/ase/calculators/nwchem.py: 71.74%
46 statements
« prev ^ index » next coverage.py v7.5.3, created at 2025-08-02 00:12 +0000
1# fmt: off
3"""This module defines an ASE interface to NWchem
5https://nwchemgit.github.io
6"""
7import os
8from pathlib import Path
10import numpy as np
12from ase import io
13from ase.calculators.calculator import FileIOCalculator
14from ase.spectrum.band_structure import BandStructure
15from ase.units import Hartree
18class NWChem(FileIOCalculator):
19 implemented_properties = ['energy', 'free_energy',
20 'forces', 'stress', 'dipole']
21 _legacy_default_command = 'nwchem PREFIX.nwi > PREFIX.nwo'
22 accepts_bandpath_keyword = True
23 discard_results_on_any_change = True
25 fileio_rules = FileIOCalculator.ruleset(
26 extend_argv=['{prefix}.nwi'],
27 stdout_name='{prefix}.nwo')
29 def __init__(self, restart=None,
30 ignore_bad_restart_file=FileIOCalculator._deprecated,
31 label='nwchem', atoms=None, command=None, **kwargs):
32 """
33 NWChem keywords are specified using (potentially nested)
34 dictionaries. Consider the following input file block::
36 dft
37 odft
38 mult 2
39 convergence energy 1e-9 density 1e-7 gradient 5e-6
40 end
42 This can be generated by the NWChem calculator by using the
43 following settings:
44 >>> from ase.calculators.nwchem import NWChem
45 >>> calc = NWChem(dft={'odft': None,
46 ... 'mult': 2,
47 ... 'convergence': {'energy': 1e-9,
48 ... 'density': 1e-7,
49 ... 'gradient': 5e-6,
50 ... },
51 ... },
52 ... )
54 In addition, the calculator supports several special keywords:
56 theory: str
57 Which NWChem module should be used to calculate the
58 energies and forces. Supported values are ``'dft'``,
59 ``'scf'``, ``'mp2'``, ``'ccsd'``, ``'tce'``, ``'tddft'``,
60 ``'pspw'``, ``'band'``, and ``'paw'``. If not provided, the
61 calculator will attempt to guess which theory to use based
62 on the keywords provided by the user.
63 xc: str
64 The exchange-correlation functional to use. Only relevant
65 for DFT calculations.
66 task: str
67 What type of calculation is to be performed, e.g.
68 ``'energy'``, ``'gradient'``, ``'optimize'``, etc. When
69 using ``'SocketIOCalculator'``, ``task`` should be set
70 to ``'optimize'``. In most other circumstances, ``task``
71 should not be set manually.
72 basis: str or dict
73 Which basis set to use for gaussian-type orbital
74 calculations. Set to a string to use the same basis for all
75 elements. To use a different basis for different elements,
76 provide a dict of the form:
78 >>> calc = NWChem(...,
79 ... basis={'O': '3-21G',
80 ... 'Si': '6-31g'})
82 basispar: str
83 Additional keywords to put in the NWChem ``basis`` block,
84 e.g. ``'rel'`` for relativistic bases.
85 symmetry: int or str
86 The point group (for gaussian-type orbital calculations) or
87 space group (for plane-wave calculations) of the system.
88 Supports both group names (e.g. ``'c2v'``, ``'Fm3m'``) and
89 numbers (e.g. ``225``).
90 autosym: bool
91 Whether NWChem should automatically determine the symmetry
92 of the structure (defaults to ``False``).
93 center: bool
94 Whether NWChem should automatically center the structure
95 (defaults to ``False``). Enable at your own risk.
96 autoz: bool
97 Whether NWChem should automatically construct a Z-matrix
98 for your molecular system (defaults to ``False``).
99 geompar: str
100 Additional keywords to put in the NWChem `geometry` block,
101 e.g. ``'nucleus finite'`` for gaussian-shaped nuclear
102 charges. Do not set ``'autosym'``, ``'center'``, or
103 ``'autoz'`` in this way; instead, use the appropriate
104 keyword described above for these settings.
105 set: dict
106 Used to manually create or modify entries in the NWChem
107 rtdb. For example, the following settings enable
108 pseudopotential filtering for plane-wave calculations::
110 set nwpw:kbpp_ray .true.
111 set nwpw:kbpp_filter .true.
113 These settings are generated by the NWChem calculator by
114 passing the arguments:
116 >>> calc = NWChem(...,
117 >>> set={'nwpw:kbpp_ray': True,
118 >>> 'nwpw:kbpp_filter': True})
120 kpts: (int, int, int), or dict
121 Indicates which k-point mesh to use. Supported syntax is
122 similar to that of GPAW. Implies ``theory='band'``.
123 bandpath: BandPath object
124 The band path to use for a band structure calculation.
125 Implies ``theory='band'``.
126 pretasks: list of dict
127 Tasks used to produce a better initial guess
128 for the wavefunction.
129 These task typically use a cheaper level of theory
130 or smaller basis set (but not both).
131 The output energy and forces should remain unchanged
132 regardless of the number of tasks or their parameters,
133 but the runtime may be significantly improved.
135 For example, a MP2 calculation preceded by guesses at the
136 DFT and HF levels would be
138 >>> calc = NWChem(theory='mp2', basis='aug-cc-pvdz',
139 >>> pretasks=[
140 >>> {'dft': {'xc': 'hfexch'},
141 >>> 'set': {'lindep:n_dep': 0}},
142 >>> {'theory': 'scf', 'set': {'lindep:n_dep': 0}},
143 >>> ])
145 Each dictionary could contain any of the other parameters,
146 except those which pertain to global configurations
147 (e.g., geometry details, scratch dir).
149 The default basis set is that of the final step in the calculation,
150 or that of the previous step that which defines a basis set.
151 For example, all steps in the example will use aug-cc-pvdz
152 because the last step is the only one which defines a basis.
154 Steps which change basis set must use the same theory.
155 The following specification would perform SCF using the 3-21G
156 basis set first, then B3LYP//3-21g, and then B3LYP//6-31G(2df,p)
158 >>> calc = NWChem(theory='dft', xc='b3lyp', basis='6-31g(2df,p)',
159 >>> pretasks=[
160 >>> {'theory': 'scf', 'basis': '3-21g',
161 >>> 'set': {'lindep:n_dep': 0}},
162 >>> {'dft': {'xc': 'b3lyp'}},
163 >>> ])
165 The :code:`'set': {'lindep:n_dep': 0}` option is highly suggested
166 as a way to avoid errors relating to symmetry changes between tasks.
168 The calculator will configure appropriate options for saving
169 and loading intermediate wavefunctions, and
170 place an "ignore" task directive between each step so that
171 convergence errors in intermediate steps do not halt execution.
172 """
173 FileIOCalculator.__init__(self, restart, ignore_bad_restart_file,
174 label, atoms, command, **kwargs)
175 if self.prefix is None:
176 self.prefix = 'nwchem'
177 self.calc = None
179 def input_filename(self):
180 return f'{self.prefix}.nwi'
182 def output_filename(self):
183 return f'{self.prefix}.nwo'
185 def write_input(self, atoms, properties=None, system_changes=None):
186 FileIOCalculator.write_input(self, atoms, properties, system_changes)
188 # Prepare perm and scratch directories
189 perm = os.path.abspath(self.parameters.get('perm', self.label))
190 scratch = os.path.abspath(self.parameters.get('scratch', self.label))
191 os.makedirs(perm, exist_ok=True)
192 os.makedirs(scratch, exist_ok=True)
194 io.write(Path(self.directory) / self.input_filename(),
195 atoms, properties=properties,
196 label=self.label, **self.parameters)
198 def read_results(self):
199 output = io.read(Path(self.directory) / self.output_filename())
200 self.calc = output.calc
201 self.results = output.calc.results
203 def band_structure(self):
204 self.calculate()
205 perm = self.parameters.get('perm', self.label)
206 if self.calc.get_spin_polarized():
207 alpha = np.loadtxt(os.path.join(perm, self.label + '.alpha_band'))
208 beta = np.loadtxt(os.path.join(perm, self.label + '.beta_band'))
209 energies = np.array([alpha[:, 1:], beta[:, 1:]]) * Hartree
210 else:
211 data = np.loadtxt(os.path.join(perm,
212 self.label + '.restricted_band'))
213 energies = data[np.newaxis, :, 1:] * Hartree
214 eref = self.calc.get_fermi_level()
215 if eref is None:
216 eref = 0.
217 return BandStructure(self.parameters.bandpath, energies, eref)