First steps with PLAMS¶
We will start with the simple example of a water molecule geometry optimization using ADF as computational engine. In this script we will show how to:
- create a
- specify the input for ADF through the
- create and run the
- read from the calculation’s
# H2O_opti.py # Geometry optimization of a water molecule using ADF # Create the molecule: mol = Molecule() mol.add_atom(Atom(symbol='O', coords=(0,0,0))) mol.add_atom(Atom(symbol='H', coords=(1,0,0))) mol.add_atom(Atom(symbol='H', coords=(0,1,0))) # Initialize the settings for the ADF calculation: sett = Settings() sett.input.basis.type = 'DZ' sett.input.basis.core = 'None' sett.input.NumericalQuality = 'Basic' sett.input.XC.GGA = 'PBE' sett.input.geometry # Create and run the job: job = ADFJob(molecule=mol, settings=sett, name='water_GO') job.run() # Fetch and print some results: energy = job.results.readkf('Energy', 'Bond Energy') opt_mol = job.results.get_molecule('Geometry', 'xyz') bond_angle = opt_mol.atoms.angle(opt_mol.atoms, opt_mol.atoms) print ('== Water optimization Results ==') print ('Bonding energy: %f kcal/mol' % Units.convert(energy, 'au', 'kcal/mol')) print ('Bond angle: %f degree' % Units.convert(bond_angle, 'rad', 'degree'))
Running the script¶
Create a working directory
plams_tutorial, move into it, and save the script
H2O_opti.py. To execute the script type in a terminal:
The calculation should take just a few seconds and you should obtain the following result:
[15:23:50] PLAMS working folder: /Users/username/plams_tutorial/plams.98957 [15:23:50] Job water_GO started [15:23:53] Job water_GO finished with status 'successful' == Water optimization Results == Bonding energy: -299.288756164 kcal/mol Bond angle: 107.746310525 degree
Every time you run a PLAMS script, a uniquely named working directory is created (
This folder will contain one subdirectory per job (in our case, one folder named
water_GO). Each job folder contains the job’s input and results files.
An instance of a Molecule class can be created by either
initializing an empty molecule and manually adding the atoms
# Create the molecule: mol = Molecule() mol.add_atom(Atom(symbol='O', coords=(0,0,0))) mol.add_atom(Atom(symbol='H', coords=(1,0,0))) mol.add_atom(Atom(symbol='H', coords=(0,1,0)))
importing the atomic coordinates from an external file (supported formats:
mol = Molecule('molecules/H2O.xyz')
The Molecule class provides several methods for basic molecular manipulations, such as rotate, translate. Furthermore, the list of Atoms contained in the corresponding Molecule object can be directly accessed and manipulated. See the Molecule section of the PLAMS documentation for the complete list of methods.
The Settings class extends the regular Python dictionary by many useful features.
In the following we show how the
Settings object is used to define the input for an
ADFJob. The reader is strongly encouraged to consult the Settings section of the PLAMS documentation for a detailed explanation of the Settings class.
One important aspect of the PLAMS library is that PLAMS knows (almost) nothing about the input options of the various computational engines.
In other words, you need to know the structure of the ASCII input file of the computational engine in order to set up the corresponding
Although input files for the various computational engines (ADF, BAND, and DFTB) use different sets of keywords, they share a similar structure – they consist of blocks and subblocks which in turn contain keys and values.
To give an example, the
Settings object in
# Initialize the settings for the ADF calculation: sett = Settings() sett.input.basis.type = 'DZ' sett.input.basis.core = 'None' sett.input.NumericalQuality = 'Basic' sett.input.XC.GGA = 'PBE' sett.input.geometry
will translate into the following ADF input file when supplied to an
Basis Core None Type DZ End Geometry End Numericalquality Basic Xc Gga PBE End
With the exception of the atomic coordinates (these will be picked up from the Molecule object; see above) all desired input options must be specified in the Settings object.
See the ADF Suite section of the PLAMS documentation for a comprehensive description of how the settings object is converted into an ASCII input file.
You can get the computational engine’s input file via the job’s get_input() method. For example:
job = ADFJob(molecule=mol, settings=sett, name='water_GO') print(job.get_input())
Creating and running the Job¶
Job class is the most important object in PLAMS as it takes care of the preparing, executing and collecting the results of a Job.
For a detailed description of the Job class, see the Jobs section of the PLAMS documentation.
PLAMS ships three ready-to-use Job subclasses, ADFJob, BANDJob, and DFTBJob. The base class Job can readily be extended to be interfaced with other computational engines.
Creating and running the
ADFJob of our example is straightforward:
# Create and run the job: job = ADFJob(molecule=mol, settings=sett, name='water_GO') job.run()
Each job has its own subdirectory (with the same name as the job) in the working folder.
After the successful completion of a job, the results can be accessed via the
Results object associated with the
There are two ways of retrieving data from the job results:
The recommended approach is to read the binary results file (in KF format) by using readkf. The readkf method will access the appropriate binary file (TAPE21 for ADF, RUNKF for BAND, and RKF for DFTB) and returns the requested data in a python built-in type (int, float, bool, ...). If the requested data is an array, it returns a python list.
# Fetch and print some results: energy = job.results.readkf('Energy', 'Bond Energy')
The two arguments are the respective “Section” and “Variable” names in the binary KF result file.
You can check the content of a binary KF file with the kfbrowser tool
Most results in the binary KF file are provided in atomic units. The Units utility can be used to convert the results into the desired units:
energy_in_kcal_mol = Units.convert(energy_in_au, 'au', 'kcal/mol')
The order of the atomic coordinates in the KF result file might differ from the initially defined order due to the internal machinery of ADF, BAND, or DFTB. It is therefore advised to extract the molecular geometry from the binary KF result file via the get_molecule method:
opt_mol = job.results.get_molecule('Geometry', 'xyz')
The two arguments are again the “Section” and “Variable” in the KF file.