The calculation of molecular polarizabilities in VeloxChem is based on linear‑response theory, where the induced dipole moment is evaluated as a function of an external electric field to obtain the polarizability tensor. Starting from the ground state DFT reference, the frequency‑dependent response equations are solved to provide both static and dynamic polarizabilities, which serve as the foundation for describing optical properties such as refractive indices and light‑scattering behavior.
The linear electric-dipole polarizabilty is determined from the linear response function
where μ^α is the electric-dipole operator along Cartesian coordinate α and ω is the optical angular frequency.
At optical frequencies well separated from transition frequencies in the system, the polarizability is real, whereas in near-resonance or resonance regions of the spectrum, it becomes complex. For further theoretical details, see Norman et al. (2018).
In the nonresonant region, where the driving field lies far from any electronic transition, the molecular polarizability varies smoothly with frequency, and the response is dominated by elastic Rayleigh scattering arising from the real part of the polarizability tensor, while the imaginary component remains vanishingly small due to the absence of absorption.
The frequencies of the perturbing electric field is specified as a list in the Python script input file and as a frequency region with a frequency point separation in parenthesis in the text file input.
In the resonant region, the molecular polarizability exhibits strong frequency dependence as the probing field approaches an electronic transition, leading to large dispersive features in the real part and the emergence of an imaginary component directly associated with absorption and radiative damping processes.
Python script
import numpy as np
import veloxchem as vlx
molecule = vlx.Molecule.read_molecule_string("""
C 0.67759997 0.00000000 0.00000000
C -0.67759997 0.00000000 0.00000000
H 1.21655197 0.92414474 0.00000000
H 1.21655197 -0.92414474 0.00000000
H -1.21655197 -0.92414474 0.00000000
H -1.21655197 0.92414474 0.00000000
""")
basis = vlx.MolecularBasis.read(molecule, "def2-svpd")
scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)
crs = vlx.ComplexResponse()
# to be implemented (together with a plot function)
#crs.property = "polarizability"
crs.a_operator = "electric dipole"
crs.b_operator = "electric dipole"
crs.a_components = ["x"]
crs.b_components = ["x"]
crs.frequencies = np.arange(0.0, 0.4, 0.0025)
crs_results = crs.compute(molecule, basis, scf_results)
Self Consistent Field Driver Setup
====================================
Wave Function Model : Spin-Restricted Kohn-Sham
Initial Guess Model : Superposition of Atomic Densities
Convergence Accelerator : Two Level Direct Inversion of Iterative Subspace
Max. Number of Iterations : 50
Max. Number of Error Vectors : 10
Convergence Threshold : 1.0e-06
ERI Screening Threshold : 1.0e-12
Linear Dependence Threshold : 1.0e-06
Exchange-Correlation Functional : B3LYP
Molecular Grid Level : 4
* Info * Using the B3LYP functional.
P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch., J. Phys. Chem. 98, 11623 (1994)
* Info * Using the Libxc library (v7.0.0).
S. Lehtola, C. Steigemann, M. J.T. Oliveira, and M. A.L. Marques., SoftwareX 7, 1–5 (2018)
* Info * Using the following algorithm for XC numerical integration.
J. Kussmann, H. Laqua and C. Ochsenfeld, J. Chem. Theory Comput. 2021, 17, 1512-1521
* Info * Starting Reduced Basis SCF calculation...
* Info * ...done. SCF energy in reduced basis set: -77.937969432013 a.u. Time: 0.32 sec.
Iter. | Kohn-Sham Energy | Energy Change | Gradient Norm | Max. Gradient | Density Change
--------------------------------------------------------------------------------------------
1 -78.532278436460 0.0000000000 0.14021310 0.01417127 0.00000000
2 -78.532278549051 -0.0000001126 0.14259745 0.01533326 0.12560857
3 -78.535117485840 -0.0028389368 0.00098432 0.00010919 0.05517977
4 -78.535117641211 -0.0000001554 0.00018577 0.00002292 0.00165399
5 -78.535117646117 -0.0000000049 0.00003883 0.00000353 0.00021718
6 -78.535117646301 -0.0000000002 0.00000213 0.00000022 0.00004693
7 -78.535117646302 -0.0000000000 0.00000022 0.00000003 0.00000175
*** SCF converged in 7 iterations. Time: 8.23 sec.
Spin-Restricted Kohn-Sham:
--------------------------
Total Energy : -78.5351176463 a.u.
Electronic Energy : -111.8437512523 a.u.
Nuclear Repulsion Energy : 33.3086336060 a.u.
------------------------------------
Gradient Norm : 0.0000002150 a.u.
Ground State Information
------------------------
Charge of Molecule : 0.0
Multiplicity (2S+1) : 1
Magnetic Quantum Number (M_S) : 0.0
Complex Response Solver Setup
===============================
Number of Frequencies : 160
Max. Number of Iterations : 150
Convergence Threshold : 1.0e-04
ERI Screening Threshold : 1.0e-12
Exchange-Correlation Functional : B3LYP
Molecular Grid Level : 4
* Info * Using the B3LYP functional.
P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch., J. Phys. Chem. 98, 11623 (1994)
* Info * Using the Libxc library (v7.0.0).
S. Lehtola, C. Steigemann, M. J.T. Oliveira, and M. A.L. Marques., SoftwareX 7, 1–5 (2018)
* Info * Using the following algorithm for XC numerical integration.
J. Kussmann, H. Laqua and C. Ochsenfeld, J. Chem. Theory Comput. 2021, 17, 1512-1521
* Info * 5 gerade trial vectors in reduced space
* Info * 6 ungerade trial vectors in reduced space
*** Iteration: 1 * Residuals (Max,Min): 5.52e-01 and 3.58e-01
* Info * 10 gerade trial vectors in reduced space
* Info * 13 ungerade trial vectors in reduced space
*** Iteration: 2 * Residuals (Max,Min): 2.50e-02 and 2.20e-02
* Info * 17 gerade trial vectors in reduced space
* Info * 19 ungerade trial vectors in reduced space
*** Iteration: 3 * Residuals (Max,Min): 1.17e-03 and 1.04e-03
* Info * 22 gerade trial vectors in reduced space
* Info * 25 ungerade trial vectors in reduced space
*** Iteration: 4 * Residuals (Max,Min): 5.64e-05 and 3.70e-05
*** Complex response converged in 4 iterations. Time: 30.10 sec