RDMFunctional - Maple Help

QuantumChemistry

 RDMFunctional
 compute the ground state energy of a molecule as a functional of the one-electron reduced density matrix

 Calling Sequence RDMFunctional(molecule, options)

Parameters

 molecule - list of lists; each list has 4 elements, the string of an atom's symbol and atom's x, y, and z coordinates options - (optional) equation(s) of the form option = value where option is one of symmetry, unit,  max_memory, xc, nuclear_gradient, populations, conv_tol, diis, diis_space, diis_start_cycle, direct_scf, direct_scf_tol, level_shift_factor, max_cycle, max_rho_cutoff, small_rho_cutoff, grids_atomic_radii, grids_becke_scheme, grids_level, grids_prune, grids_radi_method, grids_radii_adjust, grids_symmetry, excited_states, nstates

Description

 • Reduced Density Matrix Functional Theory (RDMFT) computes the energy of a many-electron atom or molecule as a functional of the one-electron reduced density matrix rather than the many-electron wavefunction ψ(, ..., ${\mathrm{r}}_{\mathrm{N}}$).
 • The command adds a universal correction, based on the 1-electron reduced density matrix (1-RDM) rather than the density alone, to any density functional theory (DFT) functional.
 • A well-known limitation of DFT is its difficulty in predicting the energies and properties of molecules with static correlation, which is important to the computation of charges, van der Waals forces, barrier heights, and bi- and multiradicals.
 • RDMFT's universal correction to DFT allows the orbital occupations to become fractional and thereby, capture strong correlation.
 • Practically, RDMFT has a computational scaling of O(N3) similar to that of DFT.  In RDMFT the energy is minimized by semidefinite programming.
 • RDMFunctional employs "B3LYP" as its default exchange-correlation (xc) functional.  Although "B3LYP" is the default xc functional, the RDMFunctional is currently optimized for semi-local (non-hybrid) functionals such as "PBE" and "SCAN".  The full list of available functionals is provided in the help page Functionals.  Furthermore, the available functionals can be searched with the SearchFunctionals command.
 • On the Windows operating system the RDMFunctional command requires the installation of Microsoft's Windows Subsystem for Linux (WSL).  For Windows 10 (version 2004 and higher) and Windows 11 you can install the WSL by opening the Command Prompt in administrator mode and entering the command: wsl --install -d Ubuntu  For additional details, please refer to: https://learn.microsoft.com/en-us/windows/wsl/install

Outputs

The table of following contents:

 ${t}\left[{\mathrm{e_tot}}\right]$ - float -- total electronic energy of the system ${t}\left[{\mathrm{mo_coeff}}\right]$ - Matrix -- coefficients expressing molecular orbitals (columns) in terms of atomic orbitals (rows) ${t}\left[{\mathrm{mo_occ}}\right]$ - Vector -- molecular orbital occupations ${t}\left[{\mathrm{mo_energy}}\right]$ - Vector -- energies of the molecular orbitals ${t}\left[{\mathrm{mo_symmetry}}\right]$ - Vector -- string labels of the irreducible representations of the molecular orbitals ${t}\left[{\mathrm{group}}\right]$ - string -- name of the molecule's point group symmetry ${t}\left[{\mathrm{aolabels}}\right]$ - Vector -- string label for each atomic orbital consisting of the atomic symbol and the orbital name ${t}\left[{\mathrm{converged}}\right]$ - integer -- 1 or 0, indicating whether the calculation is converged or not ${t}\left[{\mathrm{rdm1}}\right]$ - Matrix -- one-particle reduced density matrix (1-RDM) in the atomic-orbital basis set ${t}\left[{\mathrm{nuclear_gradient}}\right]$ - Matrix -- analytical nuclear gradient ${t}\left[{\mathrm{dipole}}\right]$ - Vector -- dipole moment according to its x, y and z components ${t}\left[{\mathrm{populations}}\right]$ - Matrix -- atomic-orbital populations ${t}\left[{\mathrm{charges}}\right]$ - Vector -- atomic charges from the populations

Options

 • basis = string -- name of the basis set.  See Basis for a list of available basis sets.  Default is "sto-3g".
 • spin = nonnegint -- twice the total spin S (= 2S). Default is 0.
 • charge = nonnegint -- net charge of the molecule. Default is 0.
 • unit = "Angstrom" or "Bohr". Default is "Angstrom".
 • max_memory = posint -- allowed memory in MB.
 • xc = string -- The full list of available functionals is provided in the help page Functionals.  Default functional is "B3LYP".
 • solvent = string -- name of the solvent.  See Solvent.  Default is "None".
 • ghost = list of lists -- each list has the string of an atom's symbol and the atom's x, y, and z coordinates.  See Ghost Atoms.
 • nuclear_gradient = boolean -- option to return the analytical nuclear gradient if available. Default is false.
 • active = list -- optional list [N,r] where r is the maximum number of orbitals allowed to be fractional and N is the number of electrons populating those r orbitals.  By default, all orbitals are allowed to become fractional.
 • populations = string -- atomic-orbital population analysis: "Mulliken" and "Mulliken/meta-Lowdin". Default is "Mulliken".
 • conv_tol = float -- converge threshold. Default is ${10}^{-10}.$
 • diis = boolean -- whether to employ diis. Default is true.
 • diis_space = posint -- diis's space size. By default, 8 Fock matrices and errors vector are stored.
 • diis_start_cycle = posint -- the step to start diis. Default is 1.
 • direct_scf = boolean -- direct SCF is used by default.
 • direct_scf_tol = float -- direct SCF cutoff threshold. Default is ${10}^{-13}.$
 • level_shift_factor = float/int -- level shift (in au) for virtual space. Default is $0.$
 • max_cycle = posint -- max number of iterations. Default is 50.
 • max_rho_cutoff = float -- ${10}^{-7}.$
 • small_rho_cutoff = float -- $0.$
 • grids_becke_scheme = string -- "gen_grid.original_becke" (default) or "gen_grid.stratmann".
 • grids_level = (0-9) -- big number for large mesh grids. Default is 3.
 • grids_prune = "gen_grid.nwchem_prune" (default), "gen_grid.sg1_prune", "gen_grid.treutler_prune", or "None".
 • grids_symmetry = boolean -- true or false to symmetrize mesh grids.

References

 1 D. Gibney, J.-N. Boyn, and D. A. Mazziotti, Phys. Rev. Lett. 131, 243003 (2023). "Universal Generalization of Density Functional Theory for Static Correlation"
 2 D. Gibney, J.-N. Boyn, and D. A. Mazziotti, J. Phys. Chem. Lett. 13, 1382–1388 (2022). "Density Functional Theory Transformed into a One-Electron Reduced-Density-Matrix Functional Theory for the Capture of Static Correlation"
 3 D. Gibney, J.-N. Boyn, and D. A. Mazziotti, J. Chem. Theory Comput. 18, 6600–6607 (2022). "Comparison of Density-Matrix Corrections to Density Functional Theory"
 4 D. Gibney, J.-N. Boyn, and D. A. Mazziotti, J. Phys. Chem. Lett. 12, 385–391 (2021). "Toward a Resolution of the Static Correlation Problem in Density Functional Theory from Semidefinite Programming"
 5 R. G. Parr and W. Yang, Density Functional Theory of Atoms and Molecules (Oxford University Press, New York, 1989).

Examples

While the generalized DFT is applicable to a wide range of molecules, we demonstrate the command RDMFunctional with the dissociation of the diatomic molecule N2.  First, we load the QuantumChemistry package

 > $\mathrm{with}\left(\mathrm{QuantumChemistry}\right);$
 $\left[{\mathrm{AOLabels}}{,}{\mathrm{ActiveSpaceCI}}{,}{\mathrm{ActiveSpaceSCF}}{,}{\mathrm{AtomicData}}{,}{\mathrm{BondAngles}}{,}{\mathrm{BondDistances}}{,}{\mathrm{Charges}}{,}{\mathrm{ChargesPlot}}{,}{\mathrm{Chat}}{,}{\mathrm{ContractedSchrodinger}}{,}{\mathrm{CorrelationEnergy}}{,}{\mathrm{CoupledCluster}}{,}{\mathrm{DensityFunctional}}{,}{\mathrm{DensityPlot3D}}{,}{\mathrm{Dipole}}{,}{\mathrm{DipolePlot}}{,}{\mathrm{Energy}}{,}{\mathrm{ExcitationEnergies}}{,}{\mathrm{ExcitationSpectra}}{,}{\mathrm{ExcitationSpectraPlot}}{,}{\mathrm{ExcitedStateEnergies}}{,}{\mathrm{ExcitedStateSpins}}{,}{\mathrm{ExcitonDensityPlot}}{,}{\mathrm{ExcitonPopulations}}{,}{\mathrm{ExcitonPopulationsPlot}}{,}{\mathrm{FullCI}}{,}{\mathrm{GeometryOptimization}}{,}{\mathrm{HartreeFock}}{,}{\mathrm{Interactive}}{,}{\mathrm{Isotopes}}{,}{\mathrm{LiteratureSearch}}{,}{\mathrm{MOCoefficients}}{,}{\mathrm{MODiagram}}{,}{\mathrm{MOEnergies}}{,}{\mathrm{MOIntegrals}}{,}{\mathrm{MOOccupations}}{,}{\mathrm{MOOccupationsPlot}}{,}{\mathrm{MOSymmetries}}{,}{\mathrm{MP2}}{,}{\mathrm{MolecularData}}{,}{\mathrm{MolecularDictionary}}{,}{\mathrm{MolecularGeometry}}{,}{\mathrm{NuclearEnergy}}{,}{\mathrm{NuclearGradient}}{,}{\mathrm{OscillatorStrengths}}{,}{\mathrm{Parametric2RDM}}{,}{\mathrm{PlotMolecule}}{,}{\mathrm{Populations}}{,}{\mathrm{Purify2RDM}}{,}{\mathrm{QuantumComputing}}{,}{\mathrm{RDM1}}{,}{\mathrm{RDM2}}{,}{\mathrm{RDMFunctional}}{,}{\mathrm{RTM1}}{,}{\mathrm{ReadXYZ}}{,}{\mathrm{Restore}}{,}{\mathrm{Save}}{,}{\mathrm{SaveXYZ}}{,}{\mathrm{SearchBasisSets}}{,}{\mathrm{SearchFunctionals}}{,}{\mathrm{SkeletalStructure}}{,}{\mathrm{SolventDatabase}}{,}{\mathrm{Thermodynamics}}{,}{\mathrm{TransitionDipolePlot}}{,}{\mathrm{TransitionDipoles}}{,}{\mathrm{TransitionOrbitalPlot}}{,}{\mathrm{TransitionOrbitals}}{,}{\mathrm{Variational2RDM}}{,}{\mathrm{VibrationalModeAnimation}}{,}{\mathrm{VibrationalModes}}{,}{\mathrm{Video}}\right]$ (1)

We select a set of bond distances from the roots of the sixth-order Chebyshev polynomial

 >
 ${\mathrm{bond_distances}}{≔}\left[{1.01921472}{,}{1.16853039}{,}{1.44442977}{,}{1.80490968}{,}{2.19509032}{,}{2.55557023}{,}{2.83146961}{,}{2.98078528}\right]$ (2)

and define a list of molecular geometries with each geometry corresponding to one of the bond distances

 >
 ${\mathrm{molecules}}{≔}\left[\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{1.01921472}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{1.16853039}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{1.44442977}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{1.80490968}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{2.19509032}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{2.55557023}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{2.83146961}\right]\right]{,}\left[\left[{"N"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"N"}{,}{0}{,}{0}{,}{2.98078528}\right]\right]\right]$ (3)

The energies for each geometry may be then readily computed from DensityFunctional with the Energy  command in the Quantum Chemistry package,

 >
 ${\mathrm{energies_dft}}{≔}\left[{-109.50786339}{,}{-109.52331040}{,}{-109.37588586}{,}{-109.17741963}{,}{-109.03358249}{,}{-108.95146495}{,}{-108.90997138}{,}{-108.89291096}\right]$ (4)

and using polynomial interpolation, we generate a polynomial in the bond distance R

 >
 ${\mathrm{pes_dft}}{≔}{-}{0.14339035}{}{{R}}^{{7}}{+}{2.20942302}{}{{R}}^{{6}}{-}{14.44738600}{}{{R}}^{{5}}{+}{51.96535307}{}{{R}}^{{4}}{-}{110.97891955}{}{{R}}^{{3}}{+}{140.38786698}{}{{R}}^{{2}}{-}{96.50876625}{}{R}{-}{81.97856683}$ (5)

Similarly, we can use the RDMFunctional method in the Energy command

 >
 ${\mathrm{energies_rdm}}{≔}\left[{-109.50786339}{,}{-109.52331040}{,}{-109.37588586}{,}{-109.18887170}{,}{-109.10062254}{,}{-109.06672079}{,}{-109.05806429}{,}{-109.05565802}\right]$ (6)

and generate a polynomial in R

 >
 ${\mathrm{pes_rdm}}{≔}{0.00403016}{}{{R}}^{{7}}{+}{0.13029220}{}{{R}}^{{6}}{-}{2.19462464}{}{{R}}^{{5}}{+}{12.96198119}{}{{R}}^{{4}}{-}{38.72108443}{}{{R}}^{{3}}{+}{62.54064536}{}{{R}}^{{2}}{-}{51.32917602}{}{R}{-}{92.88749381}$ (7)

The potential energy curves from DensityFunctional (red) and RDMFunctional (blue) can be plotted together

 >

While the DFT dissociation curve (red) rises too quickly after 1.8$⟦Å⟧$due to its incorrect treatment of static correlation, the RDMFunctional dissociation curve (blue) exhibits the correct behavior, leveling off as the bond breaks.

As a second example, we compute the energies and properties of hydrogen fluoride at a stretched geometry

 >
 ${\mathrm{molecule}}{≔}\left[\left[{"H"}{,}{0}{,}{0}{,}{0}\right]{,}\left[{"F"}{,}{0}{,}{0}{,}{2.00000000}\right]\right]$ (8)
 >
 ${\mathrm{data}}{≔}{table}{}\left(\left[{\mathrm{charges}}{=}\left[\right]{,}{\mathrm{mo_coeff}}{=}\left[\right]{,}{\mathrm{mo_symmetry}}{=}\left[\right]{,}{\mathrm{populations}}{=}\left[\right]{,}{\mathrm{mo_occ}}{=}\left[\right]{,}{\mathrm{converged}}{=}{1}{,}{\mathrm{group}}{=}{"C1"}{,}{\mathrm{aolabels}}{=}\left[\right]{,}{\mathrm{mo_energy}}{=}\left[\right]{,}{\mathrm{dipole}}{=}\left[\right]{,}{\mathrm{e_tot}}{=}{-98.77669656}{,}{\mathrm{rdm1}}{=}\left[\right]\right]\right)$ (9)
 >