Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Optical activity and dichroism

VeloxChem enables simulations of molecular optical activity by providing access to rotatory strengths, specific rotations, circular dichroism (ECD) spectra, and optical rotatory dispersion (ORD) through density‑functional theory response methods. These chiroptical properties can be obtained either from an eigenvalue‑based TDDFT treatment of excited states or through the complex polarization propagator (CPP) approach, which delivers frequency‑dependent spectra directly comparable to experiment. For systems composed of weakly interacting chromophores, VeloxChem also implements an exciton coupling model to compute ECD spectra from the interactions between local excitations. For further theoretical details, see Norman et al. (2018).

In this section, (S)-methyloxirane is used as the test molecule for ECD rotatory strengths via TDDFT, ECD extinction coefficients via CPP, and ORD specific rotations via CPP. The exciton coupling model is illustrated separately using a helical stack of ethylene molecules.

Loading...

Rotatory strengths

The strength of an ECD band is given by the anisotropy of the decadic molar extinction coefficient

Δϵ(ω)=16πNAln⁡(10)(4πε0)c2π3ℏ∑n>0f(ω;ωn0,γ) ωn0Rn0\Delta\epsilon(\omega) = \frac{ 16\pi N_\mathrm{A} }{ \ln\left(10\right) \left(4\pi\varepsilon_0\right) c^2 } \frac{\pi}{3 \hbar} \sum_{n>0} f(\omega; \omega_{n0},\gamma)\, \omega_{n0} R_{n0}
(1)

where NAN_\mathrm{A} is Avogadro’s constant, ff is the Cauchy distribution, and Rn0R_{n0} is the rotatory strength defined as

Rn0=∑α=x,y,zℑ⟨0∣μ^α∣n⟩⟨n∣m^α∣0⟩=∑α=x,y,z−emeωn0⟨0∣p^α∣n⟩⟨n∣m^α∣0⟩R_{n0} = \sum_{\alpha = x,y,z} \Im \langle 0 | \hat{\mu}_\alpha | n \rangle \langle n | \hat{m}_\alpha | 0\rangle = \sum_{\alpha = x,y,z} \frac{-e}{m_\mathrm{e} \omega_{n0}} \langle 0 | \hat{p}_\alpha | n \rangle \langle n | \hat{m}_\alpha | 0\rangle
(2)

In VeloxChem, the rotatory strength is evaluated in the velocity gauge as given in the second expression, and it is thereby gauge-origin independent.

Python script

import veloxchem as vlx

molecule = vlx.Molecule.read_xyz_string("""10

O             -1.009866880244         1.299407071912         1.951947754409
C              0.031038818173         2.001294498224         1.244371321259
C              0.339506903623         0.795807201271         2.029917688849
C              0.600316751832        -0.544462572436         1.394186918859
H             -0.026285836651         1.930258644799         0.154459855148
H              0.271804379347         2.990968281608         1.641505800108
H              0.794935917667         0.948100885444         3.014899839267
H              0.113110324011        -0.610670447580         0.412641743927
H              1.681576973773        -0.701096452623         1.264343536564
H              0.213803265264        -1.354101659046         2.029564064183""")


basis = vlx.MolecularBasis.read(molecule, "aug-cc-pvdz")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

rsp_drv = vlx.LinearResponseEigenSolver()
rsp_drv.nstates = 8
rsp_results = rsp_drv.compute(molecule, basis, scf_results)
<Figure size 800x500 with 2 Axes>

Text file

@jobs
task: response
@end

@method settings
basis: def2-svpd
xcfun: b3lyp
@end

@response
property: ecd
nstates: 10
nto: yes
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Extinction coefficient

The anisotropy of the decadic molar extinction coefficient can be determined directly from the complex polarization propagator evaluated for mixed electric- and magnetic-dipole operators Jiemchooroj & Norman (2007)

Δϵ(ω)=16πNAω2ln⁡(10)(4πε0)c2 β(ω)\Delta\epsilon(\omega) = \frac{ 16 \pi N_\mathrm{A} \omega^2 }{ \ln(10) \left(4\pi\varepsilon_0\right) c^2 } \, \beta(\omega)
(3)

where the molecular response property, β(ω)\beta(\omega), is defined as

β(ω)=−13ω(Gxx+Gyy+Gzz)\beta(\omega) = -\frac{1}{3 \omega} (G_{xx} + G_{yy} + G_{zz})
(4)

and

Gαβ=−ℜ⟨⟨μ^α;m^β⟩⟩ωγ=−eωmeℑ⟨⟨p^α;m^β⟩⟩ωγG_{\alpha\beta} = - \Re\langle\langle\hat{\mu}_\alpha;\hat{m}_\beta \rangle\rangle_\omega^\gamma = - \frac{e}{\omega m_e} \Im \langle\langle\hat{p}_\alpha; \hat{m}_\beta \rangle\rangle_\omega^\gamma
(5)

The mixed electric–magnetic dipole tensor, GG, is evaluated in the velocity gauge as given in the second expression. Furthermore, it is complex and calculated with a damping term, ℏγ\hbar \gamma, associated with the inverse finite lifetime of the excited states. The default program setting for this parameter is 0.124 eV (or 0.004556 a.u.).

The resulting values for Δϵ(ω)\Delta \epsilon(\omega) are converted from atomic units to units of L mol−1^{-1} cm−1^{-1} by multiplying with a factor of 10 a0210\, a_0^2.

Python script

import numpy as np
import veloxchem as vlx

molecule = vlx.Molecule.read_xyz_string("""10

O             -1.009866880244         1.299407071912         1.951947754409
C              0.031038818173         2.001294498224         1.244371321259
C              0.339506903623         0.795807201271         2.029917688849
C              0.600316751832        -0.544462572436         1.394186918859
H             -0.026285836651         1.930258644799         0.154459855148
H              0.271804379347         2.990968281608         1.641505800108
H              0.794935917667         0.948100885444         3.014899839267
H              0.113110324011        -0.610670447580         0.412641743927
H              1.681576973773        -0.701096452623         1.264343536564
H              0.213803265264        -1.354101659046         2.029564064183""")

basis = vlx.MolecularBasis.read(molecule, "aug-cc-pvdz")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

cpp_drv = vlx.ComplexResponse()
cpp_drv.frequencies = np.arange(0.207, 0.325, 0.0025)
cpp_drv.property = "ecd"

cpp_results = cpp_drv.compute(molecule, basis, scf_results)
<Figure size 800x500 with 1 Axes>

Text file

@jobs
task: response
@end

@method settings
basis: def2-svpd
xcfun: b3lyp
@end

@response
property: ecd (cpp)
# frequency region (and resolution)
frequencies:  0.207-0.325 (0.0025)
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Optical Rotatory Dispersion

The specific optical rotation can be determined directly from the complex polarization propagator evaluated for mixed electric- and magnetic-dipole operators Norman et al. (2018)

[α]ω=−18000 NAω2π(4πε0)c2M β′(ω)[\alpha]_\omega = \frac{ -18000\, N_\mathrm{A} \omega^2 }{ \pi \left(4\pi\varepsilon_0\right) c^2 M } \, \beta'(\omega)
(6)

where MM is the molar mass and the molecular response property, β′(ω)\beta'(\omega), is defined as

β′(ω)=13ω(Gxx′+Gyy′+Gzz′)\beta'(\omega) = \frac{1}{3 \omega} (G'_{xx} + G'_{yy} + G'_{zz})
(7)

and

Gαβ′=−ℑ⟨⟨μ^α;m^β⟩⟩ωγ=eωmeℜ⟨⟨p^α;m^β⟩⟩ωγG'_{\alpha\beta} = - \Im\langle\langle\hat{\mu}_\alpha;\hat{m}_\beta \rangle\rangle_\omega^\gamma = \frac{e}{\omega m_e} \Re \langle\langle\hat{p}_\alpha; \hat{m}_\beta \rangle\rangle_\omega^\gamma
(8)

The mixed electric–magnetic dipole tensor, G′G', is evaluated in the velocity gauge as given in the second expression, and is thereby gauge-origin independent. It is complex and calculated with a damping term, ℏγ\hbar \gamma, associated with the inverse finite lifetime of the excited states. The default program setting for this parameter is 0.124 eV (or 0.004556 a.u.).

The resulting values for [α]ω[\alpha]_\omega are reported in units of deg dm−1^{-1} (g cm−3^{-3})−1^{-1}.

Python script

import numpy as np
import veloxchem as vlx

molecule = vlx.Molecule.read_xyz_string("""10

O             -1.009866880244         1.299407071912         1.951947754409
C              0.031038818173         2.001294498224         1.244371321259
C              0.339506903623         0.795807201271         2.029917688849
C              0.600316751832        -0.544462572436         1.394186918859
H             -0.026285836651         1.930258644799         0.154459855148
H              0.271804379347         2.990968281608         1.641505800108
H              0.794935917667         0.948100885444         3.014899839267
H              0.113110324011        -0.610670447580         0.412641743927
H              1.681576973773        -0.701096452623         1.264343536564
H              0.213803265264        -1.354101659046         2.029564064183""")

basis = vlx.MolecularBasis.read(molecule, "aug-cc-pvdz")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

cpp_drv = vlx.ComplexResponse()
cpp_drv.frequencies = np.arange(0.207, 0.325, 0.0025)
cpp_drv.property = "ord"

cpp_results = cpp_drv.compute(molecule, basis, scf_results)
<Figure size 800x500 with 1 Axes>

Text file

@jobs
task: response
@end

@method settings
basis: def2-svpd
xcfun: b3lyp
@end

@response
property: ord (cpp)
# frequency region (and resolution)
frequencies: 0.207-0.325 (0.0025)
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Exciton coupling model

VeloxChem implements the exciton coupling model to determine circular dichroism spectra. For further theoretical details about this model, see Li et al. (2017).

Loading...

Figure: A helical stack of 10 ethylene molecules (or fragments).

Python script

import veloxchem as vlx

# the molecule object can be built from a molecular string listing all fragments 
# in the helical stack and all atoms in each fragment
molecule = vlx.Molecule.read_str(mol_str)
basis = vlx.MolecularBasis.read(molecule, "def2-svp")

exciton_drv = vlx.ExcitonModelDriver()

exciton_drv.update_settings({"fragments": "10", 'atoms_per_fragment': "6"})

exciton_drv.xcfun_label = "b3lyp"
exciton_drv.nstates = 1
exciton_drv.ct_nocc = 0
exciton_drv.ct_nvir = 0

exciton_results = exciton_drv.compute(molecule, basis)
<Figure size 700x400 with 1 Axes>

Text file

@jobs
task: exciton
@end

@method settings
xcfun: b3lyp
basis: cc-pvdz
@end

@exciton
fragments: 40
atoms_per_fragment: 55
nstates: 5
ct_nocc: 0
ct_nvir: 0
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
# XYZ coordinates for 40 x 55 atoms
...
@end
References
  1. Norman, P., Ruud, K., & Saue, T. (2018). Principles and practices of molecular properties. John Wiley & Sons, Ltd.
  2. Jiemchooroj, A., & Norman, P. (2007). Electronic circular dichroism spectra from the complex polarization propagator. J. Chem. Phys., 126(13), 134102. 10.1063/1.2716660
  3. Li, X., Parrish, R. M., Liu, F., Kokkila Schumacher, S. I. L., & Martinez, T. J. (2017). An Ab Initio Exciton Model Including Charge-Transfer Excited States. J. Chem. Theory Comput., 13(8), 3493–3504. 10.1021/acs.jctc.7b00171
Properties
Polarizability
Properties
Vibrational spectroscopies