ADF Fragment Analysis Recipe
****************************

ADF Bonding Analysis EDA

Use the PLAMS ADF fragment recipe to run two fragment calculations
plus a full-system calculation and extract energy decomposition terms.

-[ Related documentation ]-

* ADF fragment recipe details

* ADF molecular fragments

[image]


Initialization
==============

   from scm.plams import Settings, Molecule, AMSJob, Units, view, plot_image_grid
   from scm.plams.recipes.adffragment import ADFFragmentJob


Define the molecules
====================

For convenience we define here two molecules, normally you would read
them from xyz files

   def get_molecule(input_string):
       job = AMSJob.from_input(input_string)
       return job.molecule[""]


   mol1 = get_molecule(
       """
   System
       Atoms
           C      -0.75086900       1.37782400      -2.43303700
           C      -0.05392100       2.51281000      -2.41769100
           H      -1.78964800       1.33942600      -2.09651100
           H      -0.30849400       0.43896500      -2.76734700
           H      -0.49177100       3.45043100      -2.06789100
           H       0.98633900       2.54913500      -2.74329400
       End
   End
   """
   )


   mol2 = get_molecule(
       """
   System
       Atoms
           C       0.14667300      -0.21503500       0.40053800
           C       1.45297400      -0.07836900       0.12424400
           C       2.23119700       1.15868100       0.12912100
           C       1.78331500       2.39701500       0.38779700
           H      -0.48348000       0.63110600       0.67664100
           H      -0.33261900      -1.19332100       0.35411600
           H       2.01546300      -0.97840100      -0.14506700
           H       3.29046200       1.03872500      -0.12139700
           H       2.45728900       3.25301000       0.35150400
           H       0.74193400       2.60120700       0.64028800
       End
   End
   """
   )

   # Molecules could also be read from xyz files
   # mol1 = Molecule("ethene.xyz")
   # mol2 = Molecule("butadiene.xyz")

   plot_image_grid(
       {"1": view(mol1, guess_bonds=True), "2": view(mol2, guess_bonds=True)}, show_labels=False
   );

   [image: image generated from notebook][image]


Setup and run the job
=====================

   common = Settings()  # common settings for all 3 jobs
   common.input.ams.Task = "SinglePoint"
   common.input.adf.basis.type = "DZP"
   common.input.adf.xc.gga = "PBE"
   common.input.adf.symmetry = "NOSYM"

   full = Settings()  # additional settings for full system calculation
   full.input.adf.etsnocv  # empty block
   full.input.adf.print = "etslowdin"

   job = ADFFragmentJob(fragment1=mol1, fragment2=mol2, settings=common, full_settings=full)
   results = job.run()

   [08.04|17:13:19] JOB plamsjob STARTED
   [08.04|17:13:19] JOB plamsjob RUNNING
   [08.04|17:13:19] JOB plamsjob/frag1 STARTED
   [08.04|17:13:19] JOB plamsjob/frag1 RUNNING
   [08.04|17:13:22] JOB plamsjob/frag1 FINISHED
   [08.04|17:13:22] JOB plamsjob/frag1 SUCCESSFUL
   [08.04|17:13:22] JOB plamsjob/frag2 STARTED
   [08.04|17:13:22] JOB plamsjob/frag2 RUNNING
   [08.04|17:13:26] JOB plamsjob/frag2 FINISHED
   [08.04|17:13:26] JOB plamsjob/frag2 SUCCESSFUL
   [08.04|17:13:26] JOB plamsjob/full STARTED
   [08.04|17:13:26] JOB plamsjob/full RUNNING
   [08.04|17:13:38] JOB plamsjob/full FINISHED
   [08.04|17:13:38] JOB plamsjob/full SUCCESSFUL
   [08.04|17:13:38] JOB plamsjob FINISHED
   [08.04|17:13:38] JOB plamsjob SUCCESSFUL


Energy Decomposition Results
============================

Energy terms contributing to the bond energy (with respect to the
fragments):

   import pandas as pd

   bond_energy = results.get_energy(unit="eV")
   decom = results.get_energy_decomposition(unit="eV")

   df = pd.DataFrame(
       [
           {"Term": t, "Energy [eV]": e}
           for t, e in list(decom.items()) + [("Total Bond Energy", bond_energy)]
       ]
   )
   df

Term

Energy [eV]

0

Electrostatic Energy

-0.159357

1

Kinetic Energy

-0.295645

2

Elstat Interaction

0.748751

3

XC Energy

-0.356182

4

E_int

-1.697991

5

E_Pauli

0.187116

6

E_elstat

-99.865788

7

E_orb

-2.454467

8

E_1

-31.722875

9

E_2

-56.338062

10

Total Bond Energy

-0.062434


See also
========


Related documentation
---------------------

* ADF fragment recipe details

* ADF molecular fragments


Python Script
=============

   #!/usr/bin/env python
   # coding: utf-8

   # ## Initialization

   from scm.plams import Settings, Molecule, AMSJob, Units, view, plot_image_grid
   from scm.plams.recipes.adffragment import ADFFragmentJob


   # ## Define the molecules
   # For convenience we define here two molecules, normally you would read them from xyz files


   def get_molecule(input_string):
       job = AMSJob.from_input(input_string)
       return job.molecule[""]


   mol1 = get_molecule(
       """
   System
       Atoms
           C      -0.75086900       1.37782400      -2.43303700
           C      -0.05392100       2.51281000      -2.41769100
           H      -1.78964800       1.33942600      -2.09651100
           H      -0.30849400       0.43896500      -2.76734700
           H      -0.49177100       3.45043100      -2.06789100
           H       0.98633900       2.54913500      -2.74329400
       End
   End
   """
   )


   mol2 = get_molecule(
       """
   System
       Atoms
           C       0.14667300      -0.21503500       0.40053800
           C       1.45297400      -0.07836900       0.12424400
           C       2.23119700       1.15868100       0.12912100
           C       1.78331500       2.39701500       0.38779700
           H      -0.48348000       0.63110600       0.67664100
           H      -0.33261900      -1.19332100       0.35411600
           H       2.01546300      -0.97840100      -0.14506700
           H       3.29046200       1.03872500      -0.12139700
           H       2.45728900       3.25301000       0.35150400
           H       0.74193400       2.60120700       0.64028800
       End
   End
   """
   )

   # Molecules could also be read from xyz files
   # mol1 = Molecule("ethene.xyz")
   # mol2 = Molecule("butadiene.xyz")

   plot_image_grid(
       {"1": view(mol1, guess_bonds=True), "2": view(mol2, guess_bonds=True)}, show_labels=False, save_path="picture1.png"
   )
   # ## Setup and run the job

   common = Settings()  # common settings for all 3 jobs
   common.input.ams.Task = "SinglePoint"
   common.input.adf.basis.type = "DZP"
   common.input.adf.xc.gga = "PBE"
   common.input.adf.symmetry = "NOSYM"

   full = Settings()  # additional settings for full system calculation
   full.input.adf.etsnocv  # empty block
   full.input.adf.print = "etslowdin"

   job = ADFFragmentJob(fragment1=mol1, fragment2=mol2, settings=common, full_settings=full)
   results = job.run()


   # ## Energy Decomposition Results

   # Energy terms contributing to the bond energy (with respect to the fragments):

   import pandas as pd

   bond_energy = results.get_energy(unit="eV")
   decom = results.get_energy_decomposition(unit="eV")

   df = pd.DataFrame(
       [{"Term": t, "Energy [eV]": e} for t, e in list(decom.items()) + [("Total Bond Energy", bond_energy)]]
   )
   print(df)
