Quantum Chemistry Package Tutorial Worksheet
Copyright (c) RDMCHEM LLC 2024
Overview
Load the Package
Definition of a Molecule
Computations
Visualization
Interactive Maplet
Courses
Help
The QuantumChemistry package is an environment for computing and visualizing the quantum energies and properties of many-electron atoms and molecules.
In this tutorial we will teach the user to use the following basic features of the package:
Definition of a molecule's geometry from an interface, an XYZ file, or a web database containing more than 96 million molecules (see section below)
Computation of a molecule's energies and properties from advanced wavefunction, density functional, and reduced density matrix methods (see section below)
Visualization of molecular geometries, densities, and vibrations through 3D interactive plots and animations (see section below)
Launch of an interactive Maplet (see section below) for a graphical user interface to the package commands for defining, computing, and visualizing molecular data
Integration of lessons and curricula in chemistry and physics (see section below) which can be used for formal courses as well as self study
Launch of the Maple help pages for the QuantumChemistry package for detailed information on each command (see section below)
See the What'sNew2024 help page for examples of what's new in the QuantumChemistry Package 2024.
We begin by loading the commands of the Quantum Chemistry package using the Maple with command
withQuantumChemistry;
AOLabels,ActiveSpaceCI,ActiveSpaceSCF,AtomicData,BondAngles,BondDistances,Charges,ChargesPlot,ContractedSchrodinger,CorrelationEnergy,CoupledCluster,DensityFunctional,DensityPlot3D,Dipole,DipolePlot,Energy,ExcitationEnergies,ExcitationSpectra,ExcitationSpectraPlot,ExcitedStateEnergies,ExcitedStateSpins,ExcitonDensityPlot,ExcitonPopulations,ExcitonPopulationsPlot,FullCI,GeometryOptimization,HartreeFock,Interactive,Isotopes,LiteratureSearch,MOCoefficients,MODiagram,MOEnergies,MOIntegrals,MOOccupations,MOOccupationsPlot,MOSymmetries,MP2,MolecularData,MolecularDictionary,MolecularGeometry,NuclearEnergy,NuclearGradient,OscillatorStrengths,Parametric2RDM,PlotMolecule,Populations,Purify2RDM,QuantumComputing,RDM1,RDM2,RTM1,ReadXYZ,Restore,Save,SaveXYZ,SearchBasisSets,SearchFunctionals,SkeletalStructure,SolventDatabase,Thermodynamics,TransitionDipolePlot,TransitionDipoles,TransitionOrbitalPlot,TransitionOrbitals,Variational2RDM,VibrationalModeAnimation,VibrationalModes,Video
We recommend setting the global variable Digits to 15 (double precision) rather than 10 (default):
Digits≔15;
Digits≔15
A molecule's geometry can be defined in the worksheet interface in three ways:
entering the geometry as a Maple list of lists
reading the geometry from an XYZ file
retrieving the geometry from a web database
Maple List of Lists
Let us first enter the geometry of the molecule hydrogen fluoride HF by using a Maple list of lists. Each list has four elements:
a string of the atomic symbol of the atom, i.e. "H" or "F"
the x Cartesian coordinate of the atom
the y Cartesian coordinate of the atom
the z Cartesian coordinate of the atom
The default unit of distance in the package is Angstroms (1 Angstrom= 10−10 m), which we employ in this worksheet.
hf≔H,0,0,−0.48,F,0,0,0.48;
hf≔H,0,0,−0.48000000,F,0,0,0.48000000
Web Database
Second, we can retrieve the geometry of a molecule from a web database containing more than 96 million molecules. The database is especially useful for molecules with more than a few atoms. To perform the search, we use the command MolecularGeometry (the computer must be connected to the web).
hf≔MolecularGeometryhydrogen fluoride;
hf≔F,0,0,0,H,0.43580000,−0.14770000,−0.81880000
XYZ File
Second, we can read the geometry of a molecule from an XYZ file, that is a file containing the Cartesian coordinates of each atom. For example, the XYZ file of hydrogen fluoride has the following format:
2
# Comment Line: Note that we give the total number of atoms in the first line, a comment in the second line, and then the Cartesian coordinates of the atoms in the remaining lines H 0 0 0
F 0 0 0.95
To save an XYZ file, we use the command SaveXYZ to make an XYZ file from the Maple list of lists. To read the resulting XYZ file, we use the command ReadXYZ.
Invoking ReadXYZ without a parameter ReadXYZ() opens a file dialogue for selecting an XYZ file. Invoking ReadXYZ with a string parameter opens the file associated with the string name. If the file is not present in Maple's current working directory, the string must contain the full path to the file or Maple's current working directory must modified by using the currentdir command. For example, the XYZ file for porphin "porphin.xyz" can be read from the data folder of the package. First we construct a full path to the file
file ≔ FileTools:-JoinPathkerneloptstoolboxdir=QuantumChemistry,data,porphin.xyz;
Then we call ReadXYZ to read the XYZ file into the QuantumChemistry package's molecular geometry format
porphin ≔ ReadXYZfile;
porphin≔N,−1.62058500,−1.35841100,−0.00004500,N,1.62042400,−1.35852600,−0.00009500,N,−1.57898000,1.33345500,0.00027900,N,1.57912800,1.33329300,−0.00026600,C,−1.27247000,−2.69668000,0.00024100,C,1.27217700,−2.69676300,0.00031100,C,−3.00281900,−1.21613100,−0.00051500,C,3.00267300,−1.21642300,−0.00028500,C,−0.00016600,−3.27300100,0.00055400,C,−2.51334300,−3.43612200,0.00013200,C,2.51294900,−3.43634000,0.00024600,C,−3.55872200,−2.54069000,−0.00034400,C,3.55842900,−2.54101500,0.00016400,C,−3.63768700,0.02511200,−0.00061100,C,3.63773300,0.02474200,−0.00044500,C,−2.95501700,1.24875700,−0.00024700,C,2.95521500,1.24846000,−0.00035900,C,−1.27258700,2.67294600,0.00036400,C,1.27286200,2.67274100,0.00017600,C,0.00017500,3.25342200,0.00036700,C,−3.55030400,2.58382300,−0.00002800,C,3.55059000,2.58346400,−0.00045100,C,−2.50983200,3.46576200,0.00046900,C,2.51016700,3.46547700,0.00019800,H,−0.00017100,−4.35835700,0.00090700,H,−1.06231600,−0.49153600,0.00028100,H,1.06235700,−0.49148500,0.00004000,H,−2.57314000,−4.51439200,0.00038500,H,2.57252600,−4.51461900,0.00038400,H,−4.61490100,−2.76464200,−0.00062200,H,4.61460100,−2.76500800,0.00017400,H,−4.72223800,0.03505700,−0.00090700,H,4.72228200,0.03449200,−0.00050400,H,0.00026800,4.34141600,0.00063100,H,−4.60947300,2.79979000,−0.00015400,H,4.60975700,2.79940900,−0.00066000,H,−2.55314700,4.54611800,0.00079700,H,2.55354400,4.54583900,0.00050700
Computations of molecules can be performed in one of two ways:
Maple commands based on the energy or property to be computed
Maple commands based on the quantum chemistry method to be employed
Energies/Properties
The QuantumChemistry package contains commands to compute directly the following energies and properties:
energy
correlation energy
nuclear energy
molecular orbital coefficients
molecular orbital occupations
molecular orbital energies
orbital symmetries
atomic-orbital populations and charges
dipole moments
nuclear gradients
excited-state energies
excited-state spins
excitation spectra
transition dipole moments
oscillator strengths
natural transition orbitals
thermodynamic properties
1-electron reduced density matrix (1-RDM)
2-electron reduced density matrix (2-RDM)
1-electron reduced transition matrices (1-RTMs)
Unless otherwise noted, all energies are in atomic units (a.u.).
The energy of hydrogen fluoride HF, for example, can be directly computed from the Hartree-Fock method from the command Energy
Energyhf,method=HartreeFock;
−98.57249547
Each energy/property command can take all of the options, specified by optional keywords, available for the given method. For instance, we can specify a larger basis set with the basis keyword
Energyhf,method=HartreeFock,basis=cc-pvdz;
−100.01804664
We can select a different method for computing the energy
Energyhf,method=DensityFunctional,basis=cc-pvdz;
−100.39764610
The nuclear energy from the repulsion of the nuclei in the molecule can be computed with the command NuclearEnergy
NuclearEnergyhf;
5.07069443
The energies of the molecular orbitals can be computed with the command MOEnergies
MOEnergieshf,method=HartreeFock;
−25.90202575−1.46443647−0.57855098−0.46350983−0.463509830.60611202
The occupations of the molecular orbitals can be computed with the command MOOccupations
MOOccupationshf,method=HartreeFock;
2.000000002.000000002.000000002.000000002.000000000.
The occupations being 2 (filled) or 0 (unfilled) reflect that the Hartree-Fock method does not correlate the orbitals describing the electrons. Using the command MOOccupations with a correlated method like the variational 2-RDM method generates fractional occupations, that is occupations between 0 and 2
MOOccupationshf,method=Variational2RDM, active=6,4;
2.000000002.000000001.999252851.999252841.974011660.02748258
The 1-electron reduced density matrix (1-RDM) can be calculated from the variational 2-RDM method with the command RDM1
dm1≔RDM1hf,method=Variational2RDM,active=6,4;
dm1≔1.973724200.0.−0.023653360.1.999252840.0.0.0.1.999252850.−0.023653360.0.0.02777004
The eigenvalues of the 1-RDM are the occupations of its eigenvectors, known as natural orbitals. It is the occupations of the natural orbitals that are returned by the command MOOccupations. To confirm that this is true for the variational 2-RDM method, we can calculate the eigenvalues of the 1-RDM with the Eigenvalues command in the LinearAlgebra package. First, we ensure that dm1 is a symmetric matrix which will have real eigenvalues,
dm1≔Matrixdm1,shape=symmetric;
and then we compute the eigenvalues
eigs≔LinearAlgebra:-Eigenvaluesdm1;
eigs≔0.027482581.974011661.999252841.99925285
The eigenvalues reveal the number of electrons in each natural orbital. The deviation of the eigenvalues from 0 and 2 reveals the presence of electron correlation. With a knowledge of the 1-RDM other important one-electron properties can be computed.
The Energy command is especially useful for generating potential energy curves (and surfaces) with just a few lines of code. For example, the following five lines generate the potential energy of hydrogen fluoride from the Hartree-Fock method with the internuclear distance ranging from 0.8 Angstroms to 2.0 Angstroms in 0.05 increments
rdata≔seqr,r=0.8..2.0,0.05:
molecules≔seqH,0,0,0,F,0,0,i,i=0.8..2.0,0.05:
edata_hf≔mapEnergy,molecules,method=HartreeFock:
pes_hf≔splinerdata,edata_hf,R,1:
plot_hf≔plotpes_hf,R=0.8..2.0,axes=BOXED,font=ROMAN,16,labelfont=ROMAN,16,labels=R,E,color=blue;
A similar point plot can be generated from the full solution of the Schrodinger equation in the given basis set, known as full configuration interaction
edata_ci≔mapEnergy,molecules,method=ActiveSpaceCI,active=3,4,5,6:
pes_hf≔seqrdatai,edata_cii,i=1..nopsrdata:
plot_ci≔plots:-pointplotpes_hf,axes=BOXED,font=ROMAN,16,labelfont=ROMAN,16,symbolsize=16,labels=R,E,color=black;
Finally, a plot can be generated from the variational 2-RDM method in which the correlated 2-RDM is directly computed without computing the many-electron wavefunction
edata_2rdm≔mapEnergy,molecules,method=Variational2RDM,symmetry=true,conditions=DQ,active=3,4,5,6:
pes_2rdm≔splinerdata,edata_2rdm,R,1:
plot_2rdm≔plotpes_2rdm,R=0.8..2.0,axes=BOXED,font=ROMAN,16,labelfont=ROMAN,16,labels=R,E,color=red;
The three plots can be readily displayed together to reveal that the energies from the variational 2-RDM calculations exactly match those from the full configuration interaction. .
plots:-displayplot_ci,plot_hf,plot_2rdm;
In contrast to the 2-RDM and full configuration methods, the Hartree-Fock method yields a potential energy curve whose minimum energy occurs at an internuclear distance R that is too short and whose curvature at the minimum energy is too large. Optimization of a molecule's geometry to reach a first-order critical point such as a minimum or a transition state can be performed with the command GeometryOptimization. Furthermore, the vibrational frequencies of a molecule in the normal-mode approximation can be computed with the command VibrationalModes. For example, we can optimize the internuclear distance of the hydrogen fluoride molecule by the Hartree-Fock method
hf_hf,data_hf≔GeometryOptimizationhf,HartreeFock: hf_hf;
F,−0.00376378,0.00127561,0.00707155,H,0.43956378,−0.14897561,−0.82587155
Because the fluorine is located at the origin, the bond distance can be computed as follows
bonddistance_hf ≔ BondDistanceshf_hf,1,2;
bonddistance_hf≔0.95546270
Once the optimized geometry is known, the vibrational frequency of the bond can be determined to be 4459 cm−1 from
emode,vmode≔VibrationalModeshf_hf,HartreeFock;
emode,vmode≔4458.50363636,0.10414645−0.03529699−0.19567489−0.452153250.153242390.84952527
Similarly, we can find optimize the internuclear distance by the variational 2-RDM method, using the D and Q N-representability conditions and freezing the core orbitals,
hf_2rdm,data_2rdm≔GeometryOptimizationhf,Variational2RDM,conditions=DQ,active=3,4,5,6: hf_2rdm;
F,−0.01249964,0.00423634,0.02348487,H,0.44829964,−0.15193634,−0.84228487
In this case the equilibrium bond distance, the distance at which the minimum energy occurs, is
bonddistance_2rdm≔BondDistanceshf_2rdm,1,2;
bonddistance_2rdm≔0.99311789
The vibrational frequency of the bond can be determined in cm−1 from
emode,vmode≔VibrationalModeshf_2rdm,Variational2RDM,conditions=DQ,active=3,4,5,6;
emode,vmode≔3780.92884067,0.10408256−0.03519865−0.19572664−0.451876320.152815440.84974949
As in the plots above, we observe that the treatment of electron correlation by the variational 2-RDM method causes the equilibrium bond length to increase to 0.993 from 0.955 Angstroms and the vibrational frequency in wavenumbers to decrease significantly from 4459 cm−1.
Methods
The QuantumChemistry package can solve the Schrodinger equation of atoms and molecules by several different methods. The methods can be employed in the Energy and other commands of the previous section, and they can also be called directly. Supported methods include:
Hartree-Fock (HF) methods
Density Functional Theory (DFT) methods
RDM Functional Theory (RDMFT) method
Variational Two-electron Reduced Density Matrix (2-RDM) method
Parametric Two-electron Reduced Density Matrix (2-RDM) method
Anti-Hermitian Contracted Schrödinger Equation (CSE) method
Full Configuration Interaction (FCI) method
Second-order Many-body Perturbation Theory (MP2) method
Coupled Cluster Singles-Doubles (CCSD) method
Complete-Active-Space Configuration Interaction (CAS-CI) method
Complete-Active-Space Self-Consistent-Field (CAS-SCF) method
Time-dependent Hartree-Fock (TDHF) method (via the HartreeFock command)
Configuration-Interaction Singles (CIS) method (via the HartreeFock command)
Time-dependent Density Functional Theory (TDDFT) method (via the DensityFunctional command)
Tamm-Dancoff Approximation (TDA) method (via the DensityFunctional command)
The Hartree-Fock method, for example, can be invoked by the HartreeFock command
data≔HartreeFockhf, basis=dz;
As described in the help page for the HartreeFock command, which can be opened by clicking on the hyperlink, each method has the molecule's geometry as a Maple list of lists as its first argument. Additional keyword arguments can be provided to control the charge, spin, symmetry, and basis set as well as technical details of the solution method. Complete information about these options is provided on the help page. Running a method's command returns a Maple table containing entries that correspond to generated data such as the total electronic energy, the coefficients of the molecular orbitals in terms of the atomic orbitals, the electron occupations of the molecular orbitals, the energies of the orbitals, the symmetries of the orbitals, the 1- and 2-electron reduced density matrices, the name of the point-group symmetry, and a parameter indicating convergence (set to one) or non-convergence (set to zero). The total energy from the data table can be accessed as follows:
datae_tot
−100.02155522
Similarly, the energies of the molecular orbitals can be retrieved as follows
datamo_energy
Rerunning the command takes a negligible amount of CPU time because the result is stored in a cache table.
The variational 2-RDM method, as another example, can be invoked by the Variational2RDM command
data≔Variational2RDMhf,basis=cc-pvdz, active=3,6,conditions=DQ:
The basis set of the molecule can be changed using the basis keyword. The basis keyword supports a wide range of basis sets including mixed basis sets in which each atom type has its own basis set as well as effective core potential (ECP) basis sets for heavy atoms (see Basis for details).
data1≔HartreeFockhf,basis=aug-cc-pvdz:
The charge of the molecule can be changed using the charge and spin keywords
data2≔HartreeFockhf,basis=aug-cc-pvdz,charge=1,spin=1:
From the two energies we can estimate the ionization energy, the energy to remove an electron as 0.531 a.u.
ie≔data2e_tot−data1e_tot
ie≔0.53141295
which (at least at this level of approximation) is slightly greater than the energy to ionize the hydrogen atom.
Density functional theory (DFT) can be invoked by the DensityFunctional command
data1≔DensityFunctionalhf,basis=aug-cc-pvdz:
As with the Hartree-Fock method the charge of the molecule can be changed using the charge and spin keywords
data2≔DensityFunctionalhf,basis=aug-cc-pvdz,charge=1,spin=1:
From DFT with the default B3LYP functional the energy to remove an electron is 0.596 a.u.
ie≔0.59619956
Second-order many-body perturbation theory (MP2) can be run with the MP2 command
data1≔MP2hf,basis=aug-cc-pvdz:
As with the Hartree-Fock and the DFT method the molecule can be treated with the +1 charge
data2≔MP2hf,basis=aug-cc-pvdz,charge=1,spin=1:
From MP2 the energy to remove an electron is 0.552 atomic units of energy
ie≔0.55160015
MP2 in the aug-cc-pVDZ basis set predicts a lower ionization energy for hydrogen fluoride than DFT with the B3LYP functional.
The parametric 2-RDM method can be run with the Parametric2RDM command
data1≔Parametric2RDMhf,basis=aug-cc-pvdz:
data2≔Parametric2RDMhf,basis=aug-cc-pvdz,charge=1,spin=1:
From parametric 2-RDM the energy to remove an electron is 0.586 atomic units of energy
ie≔0.58562861
Parametric 2-RDM in the aug-cc-pVDZ basis set predicts a slightly lower ionization energy for hydrogen fluoride than DFT with the B3LYP functional.
The anti-Hermitian contracted Schrödinger equation (ACSE) method can be run with the ContractedSchrodinger command
data1≔ContractedSchrodingerhf,basis=aug-cc-pvdz, active=6,4, active_cse=10,10:
As with Hartree-Fock, DFT, and the Parametric2RDM method the molecule can be treated with the +1 charge
data2≔ContractedSchrodingerhf,basis=aug-cc-pvdz,charge=1,spin=1, active=5,4, active_cse=9,10:
From the ACSE the energy to remove an electron is 0.539 atomic units of energy
ie≔0.53884898
The ACSE ionization energy is higher than that from the Hartree-Fock method but lower than that from the parametric 2-RDM method.
The coupled cluster method with single and double excitations can be run with the CoupledCluster command
data1≔CoupledClusterhf,basis=aug-cc-pvdz:
data2≔CoupledClusterhf,basis=aug-cc-pvdz,charge=1,spin=1:
From coupled cluster the energy to remove an electron is 0.585 atomic units of energy
ie≔0.58533448
The coupled cluster singles-doubles result is close to that from the parametric 2-RDM method.
The excited spectra of hydrogen fluoride can be computed from time-dependent density functional theory (TDDFT) with the ExcitationSpectra command
ex_table ≔ ExcitationSpectrahf, method=DensityFunctional, basis=aug-cc-pvdz,showtable=true:
State
Energy
Wavelength
Spin
Oscillator
1
8.83272265eV
140.36917303nm
Triplet
0.03753967
8.83272575eV
140.36912368nm
0.03753965
3
9.25846857eV
133.91436875nm
Singlet
0.03516694
4
9.25847182eV
133.91432172nm
0.03516693
5
11.99880830eV
103.33042610nm
0.40290203
6
12.38162341eV
100.13565532nm
0.01137906
7
12.38162658eV
100.13562965nm
0.01137904
8
12.64008197eV
98.08812769nm
0.00602209
9
12.64008532eV
98.08810170nm
10
13.18586379eV
94.02811934nm
0.13597576
11
13.21703694eV
93.80634854nm
0.00272417
12
14.32133712eV
86.57305972nm
0.00177742
The ExcitationSpectraPlot command generates an excitation spectra from the computed data
ExcitationSpectraPlothf,ex_table;
Other available methods include the coupled cluster method with single, double, and perturbative triple excitations [CCSD(T)]. For example, the ground-state CCSD(T) energy of hydrogen fluoride in the default basis set can be computed with the CoupledCluster command
data_cc≔CoupledClusterhf, ccsdt=true;
data_cc≔tablepopulations=1.999140761.951820471.832892451.980277671.411609490.82425915,converged=1,charges=−0.175740850.17574085,dipole=0.52012820−0.17628025−0.97723949,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.60007479,rdm1=1.99999925−0.000016854.9502559810−60.0.0.00006529−0.000016851.996983820.008387990.0.−0.003796984.9502559810−60.008387991.971030310.0.−0.021119980.0.0.1.999403970.0.0.0.0.0.1.999403970.0.00006529−0.00379698−0.021119980.0.0.03317869,mo_occ=1.999999711.999505631.999403971.999403971.968745180.03294154,e_corr=−0.02757932,e_tot_ccsdt=−98.60007479,mo_coeff=−1.002278990.232468970.0.0.006288670.080994570.00884878−1.026227010.0.−0.13908643−0.52734083−0.002837850.05398513−0.85059038−0.247400340.31832859−0.383847280.00096179−0.01829647−0.353938570.92195357−0.107886950.130092340.00533187−0.10142961−0.38887465−0.29798438−0.598089600.721188980.005369430.009184530.0.0.561109541.06349216,e_tot_mp2=−98.59082940
The core 1s orbital of fluorine can be frozen to an occupation of one with the frozen keyword
data_ccf≔CoupledClusterhf,ccsdt=true, ccsdfrozen=1;
data_ccf≔tablepopulations=1.999140761.951820471.832892451.980277671.411609490.82425915,converged=1,charges=−0.175740850.17574085,dipole=0.52012820−0.17628025−0.97723949,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.60007479,rdm1=1.99999925−0.000016854.9502559810−60.0.0.00006529−0.000016851.996983820.008387990.0.−0.003796984.9502559810−60.008387991.971030310.0.−0.021119980.0.0.1.999403970.0.0.0.0.0.1.999403970.0.00006529−0.00379698−0.021119980.0.0.03317869,mo_occ=1.999999711.999505631.999403971.999403971.968745180.03294154,e_corr=−0.02757932,e_tot_ccsdt=−98.60007479,mo_coeff=−1.002278990.232468970.0.0.006288670.080994570.00884878−1.026227010.0.−0.13908643−0.52734083−0.002837850.05398513−0.85059038−0.247400340.31832859−0.383847280.00096179−0.01829647−0.353938570.92195357−0.107886950.130092340.00533187−0.10142961−0.38887465−0.29798438−0.598089600.721188980.005369430.009184530.0.0.561109541.06349216,e_tot_mp2=−98.59082940
Similarly, we can perform a variational 2-RDM calculation with the Variational2RDM command. The 1s orbital of fluorine can be frozen by specifying the active (non-frozen) orbitals with the active keyword
data_2rdmf≔Variational2RDMhf,active=2,3,4,5,6;
data_2rdmf≔tablepopulations=1.999125871.950911991.833695011.980352171.414493160.82142167,converged=1,charges=−0.178578190.17857833,dipole=0.52844734−0.17909975−0.99286985,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.60038641,group=C1,rdm1=1.997032250.008416890.0.−0.002774270.008416891.971259980.0.−0.017364940.0.1.999383980.0.0.0.0.1.999383980.−0.00277427−0.017364940.0.0.03293968,active_orbitals=2..6,mo_occ=2.000000001.999567151.999383981.999383981.968884350.03278041,e_corr=−0.02789094,mo_coeff=0.99472636−0.262834490.0.−0.005517430.081312470.022413951.025650510.0.0.13682957−0.528649110.00129454−0.05487058−0.018410730.88564777−0.31891137−0.38324583−0.000438740.01859655−0.985616360.061897280.108084460.12988850−0.002432240.103093200.167992350.460213820.599184540.72005893−0.00552213−0.008646660.0.−0.559281351.06445845,mo_coeff_canonical=0.99472636−0.25006317−0.080042140.0.0.082375940.022413950.942406830.420890820.0.−0.533733510.001294540.03851262−0.324652140.018174570.88565265−0.38040816−0.00043874−0.013052580.110030110.985599820.062160100.12892677−0.00243224−0.072359190.60997056−0.168115070.460169000.71472739−0.005522130.15381697−0.528650890.0.1.06902327
A full configuration interaction calculation can be performed with the FullCI command
data_ci≔FullCIhf;
data_ci≔tablepopulations=1.999140761.951820451.832892471.980277671.411609560.82425909,charges=−0.175740910.17574091,dipole=0.52012838−0.17628032−0.97723983,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.60007479,rdm1=1.99999925−0.000016854.9503280710−60.0.0.00006529−0.000016851.996983830.008387980.0.−0.003796964.9503280710−60.008387981.971030320.0.−0.021119900.0.0.1.999403960.0.0.0.0.0.1.999403960.0.00006529−0.00379696−0.021119900.0.0.03317868,mo_occ=1.999999711.999505631.999403961.999403961.968745200.03294154,ci_coeff=0.991596700.0.−0.012005030.001591270.000030990.−0.017263210.0.0.0.0.0.−0.017263210.0.0.−0.012005030.0.−0.11594211−0.02996891−0.000099940.001591270.0.−0.02996891−0.024641600.000458120.000030990.0.−0.000099940.00045812−0.00039491,e_corr=−0.02757932,mo_coeff=−1.00227898−0.232469030.0.0.006288690.080994580.008848711.026226990.0.−0.13908651−0.52734086−0.00283784−0.05398511−0.087555750.881501510.31832861−0.383847270.000961790.01829647−0.98743945−0.01530444−0.107886960.130092340.005331860.101429570.131519310.47193310−0.598089640.721188950.00536943−0.009184450.0.0.561109501.06349218
We can freeze orbitals in the configuration interaction by using the command ActiveSpaceCI with the active keyword
data_cif≔ActiveSpaceCIhf,active=2,3,4,5,6;
data_cif≔tablepopulations=1.999131631.951826651.832900911.980276621.411645160.82421904,charges=−0.175780960.17578096,dipole=0.52021054−0.17630816−0.97739420,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.60005428,rdm1=1.999505860.0.0.0.0.1.999401670.0.0.0.0.1.999401670.0.0.0.0.1.968738150.0.0.0.0.0.03295264,active_orbitals=2..6,mo_occ=2.000000001.999505861.999401671.999401671.968738150.03295264,ci_coeff=0.991727620.0.0.0.0.−0.125023700.0.0.0.0.−0.017296320.0.0.0.0.−0.017296320.0.0.0.0.−0.01571841,e_corr=−0.02755881,spin_squared=0,mo_coeff=0.994726360.25006317−0.078721310.014480990.0.082375940.02241395−0.942406830.41394543−0.076146360.−0.533733510.00129454−0.03851262−0.316006750.07660979−0.88565265−0.38040816−0.000438740.013052580.286526320.94942946−0.062160100.12892677−0.002432240.072359190.56949015−0.27569501−0.460169000.71472739−0.00552213−0.15381697−0.519927280.095642000.1.06902327
Note that the energy from configuration interaction (CI) and the variational 2-RDM method are extremely close
data_2rdmfe_tot−data_cife_tot
−0.00033213
Finally, we can also perform active-space calculations, that is calculations with a set of orbitals designated as active (not frozen to an occupation of 1 or 0), where we optimize the active orbitals through orbital rotations. These types of calculations can be performed with complete-active-space self-consistent-field (CASSCF) with the ActiveSpaceSCF command or with the variational 2-RDM method with the Variational2RDM command. For example, orbitals 3 and 6 can be selected as active. For example, with the variational 2-RDM method we have
data_scf≔Variational2RDMhf,basis=6-31g,active=3,6,casscf=true;
Molecular geometries, densities, and vibrations can be visualized through 3D interactive plots and animations. See the What'sNew2020 help page for new plots for excited-state spectra.
Molecular geometries
To illustrate the PlotMolecule command, we examine three molecules: (1) hydrogen fluoride, (2) 1,3-diboromobenzene, and (3) tylenol
hydrogen fluoride
As shown previously, the geometry can be retrieved from the web database
hf≔MolecularGeometryhydrogenfluoride;
The resulting geometry as a Maple can then be passed to PlotMolecule
PlotMoleculehf;
The keyword option model allows us to choose a "ballandstick" (default) or "spacefilling" model
PlotMoleculehf,model=spacefilling;
1,3-dibromobenzene
The geometry can be retrieved from the web database
bromobenzene≔MolecularGeometry1,3-dibromobenzene;
bromobenzene≔Br,−2.84530000,−1.23120000,0,Br,2.84540000,−1.23100000,0,C,0.00010000,−0.98450000,0.00010000,C,−1.20800000,−0.28710000,0,C,1.20800000,−0.28700000,−0.00010000,C,−1.20810000,1.10770000,−0.00010000,C,1.20790000,1.10780000,0.00010000,C,−0.00020000,1.80520000,0,H,0.00020000,−2.07230000,0.00010000,H,−2.14060000,1.66640000,−0.00010000,H,2.14030000,1.66660000,0.00010000,H,−0.00020000,2.89130000,0
We can also retrieve the skeletal structure with the command SkeletalStructure
SkeletalStructure1,3-dibromobenzene;
PlotMoleculebromobenzene;
tylenol
tylenol≔MolecularGeometrytylenol;
tylenol≔O,3.85090000,0.45160000,0.00120000,O,−2.59990000,1.40410000,−0.00180000,N,−1.57050000,−0.71710000,0.00010000,C,−0.20660000,−0.42310000,−0.00020000,C,0.22050000,0.90470000,0.00040000,C,0.72980000,−1.45700000,−0.00070000,C,1.58410000,1.19860000,0.00020000,C,2.09330000,−1.16290000,−0.00070000,C,2.52040000,0.16480000,−0.00030000,C,−2.64850000,0.17820000,0.00090000,C,−3.97350000,−0.54200000,0.00100000,H,−0.44360000,1.75770000,0.00120000,H,0.41130000,−2.49630000,−0.00100000,H,−1.80100000,−1.70860000,0.00010000,H,1.90530000,2.23700000,0.00090000,H,2.81800000,−1.97260000,−0.00080000,H,−4.06550000,−1.14630000,−0.90580000,H,−4.79040000,0.18440000,0.02880000,H,−4.04450000,−1.18860000,0.88020000,H,3.96500000,1.41760000,0.00170000
SkeletalStructuretylenol;
Or its dictionary entry with the command MolecularDictionary
MolecularDictionarytylenol;
Acetaminophen: Paracetamol is a member of the class of phenols that is 4-aminophenol in which one of the hydrogens attached to the amino group has been replaced by an acetyl group. It has a role as a cyclooxygenase 2 inhibitor, a cyclooxygenase 1 inhibitor, a non-narcotic analgesic, an antipyretic, a non-steroidal anti-inflammatory drug, a cyclooxygenase 3 inhibitor, a xenobiotic, an environmental contaminant, a human blood serum metabolite and a hepatotoxic agent. It is a member of phenols and a member of acetamides. It derives from a 4-aminophenol. -ChEBI Acetaminophen: Acetaminophen is a p-aminophenol derivative with analgesic and antipyretic activities. Although the exact mechanism through which acetaminophen exert its effects has yet to be fully determined, acetaminophen may inhibit the nitric oxide (NO) pathway mediated by a variety of neurotransmitter receptors including N-methyl-D-aspartate (NMDA) and substance P, resulting in elevation of the pain threshold. The antipyretic activity may result from inhibition of prostaglandin synthesis and release in the central nervous system (CNS) and prostaglandin-mediated effects on the heat-regulating center in the anterior hypothalamus. -NCI Thesaurus (NCIt)
The resulting geometry as a Maple list of lists can then be passed to PlotMolecule
PlotMoleculetylenol;
Molecular densities
To illustrate the DensityPlot3D command, we examine two molecules: (1) hydrogen fluoride and (2) tylenol
We can compute the orbitals of hydrogen fluoride from the HartreeFock command
data≔HartreeFockhf;
data≔tablepopulations=1.999090641.946540421.839455071.981559061.433267890.80008693,converged=1,charges=−0.199913070.19991307,dipole=0.58797708−0.19927539−1.10471692,aolabels=0 F 1s0 F 2s0 F 2px0 F 2py0 F 2pz1 H 1s,e_tot=−98.57249547,group=C1,rdm1=2.11683774−0.494108980.03528596−0.01195901−0.06629679−0.00328525−0.494108982.13156421−0.200639060.067999980.37696940−0.155339830.03528596−0.200639061.783189680.073480690.407352670.35508878−0.011959010.067999980.073480691.97509615−0.13805872−0.12034560−0.066296790.376969400.40735267−0.138058721.23464808−0.66715625−0.00328525−0.155339830.35508878−0.12034560−0.667156250.60632385,mo_energy=−25.90202575−1.46443647−0.57855098−0.46350983−0.463509830.60611202,mo_occ=2.000000002.000000002.000000002.000000002.000000000.,mo_coeff=0.99472636−0.25006317−0.080042140.0.0.082375940.022413950.942406830.420890820.0.−0.533733510.001294540.03851262−0.324652140.499919400.73129442−0.38040816−0.00043874−0.013052580.110030110.85891086−0.487383850.12892677−0.00243224−0.072359190.609970560.111142820.477143020.71472739−0.005522130.15381697−0.528650890.0.1.06902327,mo_symmetry=AAAAAA
The electron density of the molecule can be displayed with the DensityPlot3D command
DensityPlot3Dhf,data;
The density of a given orbital can be selected with the keyword orbitalindex
DensityPlot3Dhf,data,orbitalindex=6;
We can compute the orbitals of tylenol from the HartreeFock command
data≔HartreeFocktylenol;
DensityPlot3Dtylenol,data;
DensityPlot3Dtylenol,data,orbitalindex=36;
Vibrational modes
To illustrate the VibrationalModeAnimation command, we examine the water molecule
water≔MolecularGeometrywater;
water≔O,0,0,0,H,0.27740000,0.89290000,0.25440000,H,0.60680000,−0.23830000,−0.71690000
After obtaining an initial geometry from the web database, we perform a geometry optimization with the Hartree-Fock method
water,data≔GeometryOptimizationwater, HartreeFock: water;
O,−0.01931358,−0.01429454,0.01010464,H,0.28821797,0.89605410,0.24592050,H,0.61529561,−0.22715956,−0.71852514
The vibrational modes are then computed
emodes,vmodes≔VibrationalModeswater,HartreeFock, basis=cc-pvdz;
We can finally animate the lowest vibrational mode of water (at this level of theory and in this basis set) with the VibrationalModeAnimation command
VibrationalModeAnimationwater,emodes,vmodes,1;
An interactive Maplet can be launched with the Interactive command for a graphical user interface to the package
Interactive;
A video tutorial for the Interactive command is available from the Video command
VideoInteractive.mp4;
Lessons and curricula in chemistry and physics are directly integrated into the QuantumChemistry package. Curricula are available for chemistry courses in general chemistry, physical chemistry (quantum mechanics), physical chemistry (statistical thermodynamics), computational chemistry, and graduate-level quantum chemistry as well as physics courses in quantum mechanics and statistical thermodynamics. The lessons and curricula, which are accessible through the help pages of the package, employ advanced features of the package to facilitate the learning and understanding of both introductory and advanced concepts in the study of atoms and molecules. They can be used for formal courses as well as self study.
Maple help pages for the QuantumChemistry package provide detailed information on each command . To access the main QuantumChemistry help page, you can issue the command
help QuantumChemistry/Overview;
Other commands are also accessible in this fashion. For example, to obtain the HartreeFock help page, you can issue the command
help HartreeFock;
Download Help Document
