pymatgen.io.nwchem module

class NwInput(mol, tasks, directives=None, geometry_options=(u’units’, u’angstroms’), symmetry_options=None, memory_options=None)[source]

Bases: monty.json.MSONable

An object representing a Nwchem input file, which is essentially a list of tasks on a particular molecule.

Parameters:
  • mol – Input molecule. If molecule is a single string, it is used as a direct input to the geometry section of the Gaussian input file.
  • tasks – List of NwTasks.
  • directives – List of root level directives as tuple. E.g., [(“start”, “water”), (“print”, “high”)]
  • geometry_options – Additional list of options to be supplied to the geometry. E.g., [“units”, “angstroms”, “noautoz”]. Defaults to (“units”, “angstroms”).
  • symmetry_options – Addition list of option to be supplied to the symmetry. E.g. [“c1”] to turn off the symmetry
  • memory_options – Memory controlling options. str. E.g “total 1000 mb stack 400 mb”
as_dict()[source]
classmethod from_dict(d)[source]
classmethod from_file(filename)[source]

Read an NwInput from a file. Currently tested to work with files generated from this class itself.

Parameters:filename – Filename to parse.
Returns:NwInput object
classmethod from_string(string_input)[source]

Read an NwInput from a string. Currently tested to work with files generated from this class itself.

Parameters:string_input – string_input to parse.
Returns:NwInput object
molecule

Returns molecule associated with this GaussianInput.

write_file(filename)[source]
exception NwInputError[source]

Bases: exceptions.Exception

Error class for NwInput.

class NwOutput(filename)[source]

Bases: object

A Nwchem output file parser. Very basic for now - supports only dft and only parses energies and geometries. Please note that Nwchem typically outputs energies in either au or kJ/mol. All energies are converted to eV in the parser.

Parameters:filename – Filename to read.
get_excitation_spectrum(width=0.1, npoints=2000)[source]

Generate an excitation spectra from the singlet roots of TDDFT calculations.

Parameters:
  • width (float) – Width for Gaussian smearing.
  • npoints (int) – Number of energy points. More points => smoother curve.
Returns:

(ExcitationSpectrum) which can be plotted using

pymatgen.vis.plotters.SpectrumPlotter.

parse_tddft()[source]

Parses TDDFT roots. Adapted from nw_spectrum.py script.

Returns:
{
“singlet”: [
{
“energy”: float, “osc_strength: float

}

], “triplet”: [

{
“energy”: float

}

]

}

class NwTask(charge, spin_multiplicity, basis_set, basis_set_option=u’cartesian’, title=None, theory=u’dft’, operation=u’optimize’, theory_directives=None, alternate_directives=None)[source]

Bases: monty.json.MSONable

Base task for Nwchem.

Very flexible arguments to support many types of potential setups. Users should use more friendly static methods unless they need the flexibility.

Parameters:
  • charge – Charge of the molecule. If None, charge on molecule is used. Defaults to None. This allows the input file to be set a charge independently from the molecule itself.
  • spin_multiplicity – Spin multiplicity of molecule. Defaults to None, which means that the spin multiplicity is set to 1 if the molecule has no unpaired electrons and to 2 if there are unpaired electrons.
  • basis_set – The basis set used for the task as a dict. E.g., {“C”: “6-311++G**”, “H”: “6-31++G**”}.
  • basis_set_option – cartesian (default) | spherical,
  • title – Title for the task. Defaults to None, which means a title based on the theory and operation of the task is autogenerated.
  • theory – The theory used for the task. Defaults to “dft”.
  • operation – The operation for the task. Defaults to “optimize”.
  • theory_directives – A dict of theory directives. For example, if you are running dft calculations, you may specify the exchange correlation functional using {“xc”: “b3lyp”}.
  • alternate_directives – A dict of alternate directives. For example, to perform cosmo calculations and dielectric constant of 78, you’d supply {‘cosmo’: {“dielectric”: 78}}.
as_dict()[source]
classmethod dft_task(mol, xc=u’b3lyp’, **kwargs)[source]

A class method for quickly creating DFT tasks with optional cosmo parameter .

Parameters:
  • mol – Input molecule
  • xc – Exchange correlation to use.
  • **kwargs – Any of the other kwargs supported by NwTask. Note the theory is always “dft” for a dft task.
classmethod esp_task(mol, **kwargs)[source]

A class method for quickly creating ESP tasks with RESP charge fitting.

Parameters:
  • mol – Input molecule
  • **kwargs – Any of the other kwargs supported by NwTask. Note the theory is always “dft” for a dft task.
classmethod from_dict(d)[source]
classmethod from_molecule(mol, theory, charge=None, spin_multiplicity=None, basis_set=u‘6-31g’, basis_set_option=u’cartesian’, title=None, operation=u’optimize’, theory_directives=None, alternate_directives=None)[source]

Very flexible arguments to support many types of potential setups. Users should use more friendly static methods unless they need the flexibility.

Parameters:
  • mol – Input molecule
  • charge – Charge of the molecule. If None, charge on molecule is used. Defaults to None. This allows the input file to be set a charge independently from the molecule itself.
  • spin_multiplicity – Spin multiplicity of molecule. Defaults to None, which means that the spin multiplicity is set to 1 if the molecule has no unpaired electrons and to 2 if there are unpaired electrons.
  • basis_set – The basis set to be used as string or a dict. E.g., {“C”: “6-311++G**”, “H”: “6-31++G**”} or “6-31G”. If string, same basis set is used for all elements.
  • basis_set_option – cartesian (default) | spherical,
  • title – Title for the task. Defaults to None, which means a title based on the theory and operation of the task is autogenerated.
  • theory – The theory used for the task. Defaults to “dft”.
  • operation – The operation for the task. Defaults to “optimize”.
  • theory_directives – A dict of theory directives. For example, if you are running dft calculations, you may specify the exchange correlation functional using {“xc”: “b3lyp”}.
  • alternate_directives – A dict of alternate directives. For example, to perform cosmo calculations with DFT, you’d supply {‘cosmo’: “cosmo”}.
operations = {u’‘: u’dummy’, u’frequencies’: u’Compute second derivatives and print out an analysis of molecular vibrations.’, u’energy’: u’Evaluate the single point energy.’, u’vscf’: u’Compute anharmonic contributions to the vibrational modes.’, u’freq’: u’Same as frequencies.’, u’dynamics’: u’Perform classical molecular dynamics.’, u’hessian’: u’Compute second derivatives.’, u’saddle’: u’Conduct a search for a transition state (or saddle point).’, u’gradient’: u’Evaluate the derivative of the energy with respect to nuclear coordinates.’, u’optimize’: u’Minimize the energy by varying the molecular structure.’, u’thermodynamics’: u’Perform multi-configuration thermodynamic integration using classical MD.’, u’property’: u’Calculate the properties for the wave function.’}
theories = {u’md’: u’Classical molecular dynamics simulation’, u’ccsd(t)’: u’Coupled-cluster linearized triples approximation’, u’ccsd’: u’Coupled-cluster single and double excitations’, u’esp’: u’ESP’, u’dft’: u’DFT’, u’ccsd+t(ccsd)’: u’Fourth order triples contribution’, u’scf’: u’Hartree-Fock’, u’sodft’: u’Spin-Orbit DFT’, u’selci’: u’Selected CI with perturbation correction’, u’band’: u’Pseudopotential plane-wave DFT for solids using NWPW’, u’direct_mp2’: u’MP2 using a full-direct algorithm’, u’mp2’: u’MP2 using a semi-direct algorithm’, u’g3gn’: u’some description’, u’tddft’: u’Time Dependent DFT’, u’rimp2’: u’MP2 using the RI approximation’, u’pspw’: u’Pseudopotential plane-wave DFT for molecules and insulating solids using NWPW’, u’mcscf’: u’Multiconfiguration SCF’, u’tce’: u’Tensor Contraction Engine’}