Franck-Condon Vibronic DOS with ADF
***********************************

ADF Franck-Condon Vibrational spectroscopy DOS Spectroscopy

This example uses ADF to model the density of states involved in
charge transfer between two systems with the AH-FC method, here for an
NO₂^🞄 radical and an NO₂⁻ anion.

**Requires:** AMS2026 or later

-[ Related examples ]-

* Quantum ESPRESSO phonon band structure and DOS for BeO

-[ Related tutorials ]-

* Vibrationally resolved spectra with ADF

-[ Related documentation ]-

* AMS FCF module

* ADF vibrationally resolved electronic spectra

[image]


Initial Imports
===============

Import PLAMS components and set up to run jobs in tandem.

   from scm.plams import Settings, AMSJob, Molecule, Atom, FCFJob, config, JobRunner, init, view
   from scm.plams.recipes.fcf_dos import FCFDOS
   from scm.base import ChemicalSystem

   config.default_jobrunner.parallel = JobRunner(parallel=True, maxjobs=2)


Setup Molecules
===============

Create the NO2 molecules using pre-optimized geometries (usually the
geometry optimization step would come first).

   no2_radical = ChemicalSystem()
   no2_radical.add_atom("N", coords=(0.0, 0.0, -0.01857566))
   no2_radical.add_atom("O", coords=(0.0, 1.09915770, -0.49171967))
   no2_radical.add_atom("O", coords=(0.0, -1.09915770, -0.49171967))
   view(no2_radical, width=200, height=200, direction="along_x", guess_bonds=True)

   [30.03|16:39:56] Starting Xvfb...
   [30.03|16:39:57] Xvfb started

   [image: image generated from notebook][image]

   no2_anion = ChemicalSystem()
   no2_anion.add_atom("N", coords=(0.0, 0.0, 0.12041))
   no2_anion.add_atom("O", coords=(0.0, 1.070642, -0.555172))
   no2_anion.add_atom("O", coords=(0.0, -1.070642, -0.555172))
   no2_anion.charge = -1
   view(no2_anion, width=200, height=200, direction="along_x", guess_bonds=True)

   [image: image generated from notebook][image]


Calculate Vibrational Frequencies
=================================

Create the settings objects for the ADF donor/acceptor vibrational
frequencies calculations. Run the calculations in parallel.

   settings_freq = Settings()
   settings_freq.input.adf.symmetry = "NoSym"
   settings_freq.input.adf.basis.type = "DZP"
   settings_freq.input.adf.basis.core = "None"
   settings_freq.input.adf.xc.lda = "SCF VWN"
   settings_freq.input.ams.Task = "SinglePoint"
   settings_freq.input.ams.Properties.NormalModes = "Yes"
   settings_freq.input.adf.title = "Vibrational frequencies"

   settings_freq_radical = settings_freq.copy()
   settings_freq_radical.input.adf.spinpolarization = 1
   settings_freq_radical.input.adf.unrestricted = "Yes"

   settings_freq_anion = settings_freq.copy()

   freq_job_radical = AMSJob(molecule=no2_radical, settings=settings_freq_radical, name="fsradical")
   freq_job_anion = AMSJob(molecule=no2_anion, settings=settings_freq_anion, name="fsanion")

   freq_results = (freq_job_radical.run(), freq_job_anion.run())

   [30.03|16:40:16] JOB fsradical STARTED
   [30.03|16:40:16] JOB fsanion STARTED
   [30.03|16:40:16] JOB fsradical RUNNING
   [30.03|16:40:16] JOB fsanion RUNNING
   [30.03|16:40:24] JOB fsanion FINISHED
   [30.03|16:40:24] JOB fsanion SUCCESSFUL
   [30.03|16:40:25] JOB fsradical FINISHED
   [30.03|16:40:25] JOB fsradical SUCCESSFUL

   rkf = freq_results[0].rkfpath(file="adf")
   !amsspectra {rkf}


Calculate Vibronic Spectra
==========================

Use Frank-Condon jobs to calculate the vibronic spectra.

   def fcf_job(state1, state2, spctype, name):
       settings_fcf = Settings()
       settings_fcf.input.spectrum.type = spctype
       settings_fcf.input.state1 = state1
       settings_fcf.input.state2 = state2
       return FCFJob(inputjob1=state1, inputjob2=state2, settings=settings_fcf, name=name)


   freq_radical = freq_results[0].rkfpath(file="adf")
   freq_anion = freq_results[1].rkfpath(file="adf")

   fc_abs = fcf_job(freq_radical, freq_anion, "absorption", "fcfabs")
   fc_emi = fcf_job(freq_anion, freq_radical, "emission", "fcfemi")

   fc_results = (fc_abs.run(), fc_emi.run())

   [30.03|16:42:14] JOB fcfabs STARTED
   [30.03|16:42:14] JOB fcfemi STARTED
   [30.03|16:42:14] JOB fcfabs RUNNING
   [30.03|16:42:14] JOB fcfemi RUNNING
   [30.03|16:42:16] JOB fcfemi FINISHED
   [30.03|16:42:16] JOB fcfemi SUCCESSFUL
   [30.03|16:42:16] JOB fcfabs FINISHED
   [30.03|16:42:16] JOB fcfabs SUCCESSFUL


Calculate Density of States
===========================

Calculate the DOS by computing the overlap of the absorption and
emission FCF spectra.

   job = FCFDOS(fc_results[0].kfpath(), fc_results[1].kfpath(), 10000.0, 10000.0)
   dos = job.dos()

   print(f"The density of states is {dos:.8e}")

   The density of states is 1.30090295e-08


See also
========


Related examples
----------------

* Quantum ESPRESSO phonon band structure and DOS for BeO


Related tutorials
-----------------

* Vibrationally resolved spectra with ADF


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

* AMS FCF module

* ADF vibrationally resolved electronic spectra


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

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

   # ## Initial Imports
   # Import PLAMS components and set up to run jobs in tandem.

   from scm.plams import Settings, AMSJob, Molecule, Atom, FCFJob, config, JobRunner, init, view
   from scm.plams.recipes.fcf_dos import FCFDOS
   from scm.base import ChemicalSystem


   config.default_jobrunner.parallel = JobRunner(parallel=True, maxjobs=2)


   # ## Setup Molecules
   # Create the NO<sub>2</sub> molecules using pre-optimized geometries (usually the geometry optimization step would come first).

   no2_radical = ChemicalSystem()
   no2_radical.add_atom("N", coords=(0.0, 0.0, -0.01857566))
   no2_radical.add_atom("O", coords=(0.0, 1.09915770, -0.49171967))
   no2_radical.add_atom("O", coords=(0.0, -1.09915770, -0.49171967))
   view(no2_radical, width=200, height=200, direction="along_x", guess_bonds=True, picture_path="picture1.png")


   no2_anion = ChemicalSystem()
   no2_anion.add_atom("N", coords=(0.0, 0.0, 0.12041))
   no2_anion.add_atom("O", coords=(0.0, 1.070642, -0.555172))
   no2_anion.add_atom("O", coords=(0.0, -1.070642, -0.555172))
   no2_anion.charge = -1
   view(no2_anion, width=200, height=200, direction="along_x", guess_bonds=True, picture_path="picture2.png")


   # ## Calculate Vibrational Frequencies
   # Create the settings objects for the ADF donor/acceptor vibrational frequencies calculations. Run the calculations in parallel.

   settings_freq = Settings()
   settings_freq.input.adf.symmetry = "NoSym"
   settings_freq.input.adf.basis.type = "DZP"
   settings_freq.input.adf.basis.core = "None"
   settings_freq.input.adf.xc.lda = "SCF VWN"
   settings_freq.input.ams.Task = "SinglePoint"
   settings_freq.input.ams.Properties.NormalModes = "Yes"
   settings_freq.input.adf.title = "Vibrational frequencies"

   settings_freq_radical = settings_freq.copy()
   settings_freq_radical.input.adf.spinpolarization = 1
   settings_freq_radical.input.adf.unrestricted = "Yes"

   settings_freq_anion = settings_freq.copy()


   freq_job_radical = AMSJob(molecule=no2_radical, settings=settings_freq_radical, name="fsradical")
   freq_job_anion = AMSJob(molecule=no2_anion, settings=settings_freq_anion, name="fsanion")

   freq_results = (freq_job_radical.run(), freq_job_anion.run())


   rkf = freq_results[0].rkfpath(file="adf")
   # get_ipython().system('amsspectra {rkf}')


   # ## Calculate Vibronic Spectra
   # Use Frank-Condon jobs to calculate the vibronic spectra.


   def fcf_job(state1, state2, spctype, name):
       settings_fcf = Settings()
       settings_fcf.input.spectrum.type = spctype
       settings_fcf.input.state1 = state1
       settings_fcf.input.state2 = state2
       return FCFJob(inputjob1=state1, inputjob2=state2, settings=settings_fcf, name=name)


   freq_radical = freq_results[0].rkfpath(file="adf")
   freq_anion = freq_results[1].rkfpath(file="adf")

   fc_abs = fcf_job(freq_radical, freq_anion, "absorption", "fcfabs")
   fc_emi = fcf_job(freq_anion, freq_radical, "emission", "fcfemi")

   fc_results = (fc_abs.run(), fc_emi.run())


   # ## Calculate Density of States
   #
   # Calculate the DOS by computing the overlap of the absorption and emission FCF spectra.

   job = FCFDOS(fc_results[0].kfpath(), fc_results[1].kfpath(), 10000.0, 10000.0)
   dos = job.dos()


   print(f"The density of states is {dos:.8e}")
