diff --git a/source/user.rst b/source/user.rst index 62840d0c2..eabb61f0a 100644 --- a/source/user.rst +++ b/source/user.rst @@ -33,4 +33,6 @@ functionalities. user/lo.rst user/sgx.rst user/geomopt.rst + user/gpu.rst user/reference.rst + diff --git a/source/user/gpu.rst b/source/user/gpu.rst new file mode 100644 index 000000000..88695f4b9 --- /dev/null +++ b/source/user/gpu.rst @@ -0,0 +1,154 @@ +.. _user_gpu: + +GPU Acceleration (GPU4PySCF) +**************************** + +*Modules*: :py:mod:`gpu4pyscf` + +.. module:: GPU4PySCF + :synopsis: GPU4PySCF +.. sectionauthor:: Xiaojie Wu . + +Introduction +============ + +Modern GPUs accelerate quantum chemistry calculation significantly, but also have an advantage in cost saving `[1]`_. +Some of basic PySCF modules, such as SCF and DFT, are accelerated with GPU via a plugin package +GPU4PySCF (See the end of this page for the supported functionalities). For the density fitting scheme, +GPU4PySCF on A100-80G can be 1000x faster than PySCF on single-core CPU. The speedup of direct SCF scheme is relatively low. + +.. _[1]: https://arxiv.org/abs/2404.09452 + +Installation +============ +The binary package of GPU4PySCF is released based on the CUDA version. + +.. list-table:: + :widths: 25 25 25 25 + :header-rows: 1 + + * - CUDA version + - GPU4PySCF + - cuTensor + * - CUDA 11.x + - ``pip3 install gpu4pyscf-cuda11x`` + - ``pip3 install cutensor-cu11`` + * - CUDA 12.x + - ``pip3 install gpu4pyscf-cuda12x`` + - ``pip3 install cutensor-cu12`` + +Usage of GPU4PySCF +================== +GPU4PySCF APIs are designed to be compatible with PySCF. When supported, high-level functions and objects are named the same as PySCF. But, GPU4PySCF classes do not directly inherit from PySCF class. +PySCF objects and GPU4PySCF objects can be converted into each other by :func:`to_gpu` and :func:`to_cpu`. In the conversion, the numpy arrays will be converted into cupy array. And the functions will be omitted if they are not supported with GPU acceleration. + +One can use the two modes to accelerate the calculations: directly use GPU4PySCF object:: + + import pyscf + from gpu4pyscf.dft import rks + + atom =''' + O 0.0000000000 -0.0000000000 0.1174000000 + H -0.7570000000 -0.0000000000 -0.4696000000 + H 0.7570000000 0.0000000000 -0.4696000000 + ''' + + mol = pyscf.M(atom=atom, basis='def2-tzvpp') + mf = rks.RKS(mol, xc='LDA').density_fit() + + e_dft = mf.kernel() # compute total energy + print(f"total energy = {e_dft}") + + g = mf.nuc_grad_method() + g_dft = g.kernel() # compute analytical gradient + + h = mf.Hessian() + h_dft = h.kernel() # compute analytical Hessian + +Alternatively, one can convert PySCF object to the corresponding GPU4PySCF object with :func:`to_gpu` since PySCF 2.5.0 :: + + import pyscf + from pyscf.dft import rks + + atom =''' + O 0.0000000000 -0.0000000000 0.1174000000 + H -0.7570000000 -0.0000000000 -0.4696000000 + H 0.7570000000 0.0000000000 -0.4696000000 + ''' + + mol = pyscf.M(atom=atom, basis='def2-tzvpp') + mf = rks.RKS(mol, xc='LDA').density_fit().to_gpu() # move PySCF object to GPU4PySCF object + e_dft = mf.kernel() # compute total energy + + +When the GPU task is done, the GPU4PySCF object can be converted into the corresponding PySCF object via :func:`mf.to_cpu()`. +Then, more sophisticated methods in PySCF can apply. One can also convert the individual CuPy array to numpy array with `Cupy APIs`_. + +.. Cupy APIs: https://docs.cupy.dev/en/stable/user_guide/index.html + +Functionalities supported by GPU4PySCF +====================================== +.. list-table:: + :widths: 25 25 25 25 + :header-rows: 1 + + * - Method + - SCF + - Gradient + - Hessian + * - direct SCF + - O + - GPU + - CPU + * - density fitting + - O + - O + - O + * - LDA + - O + - O + - O + * - GGA + - O + - O + - O + * - mGGA + - O + - O + - O + * - hybrid + - O + - O + - O + * - unrestricted + - O + - O + - O + * - PCM solvent + - GPU + - GPU + - FD + * - SMD solvent + - GPU + - GPU + - FD + * - dispersion correction + - CPU* + - CPU* + - FD + * - nonlocal correlation + - O + - O + - NA + * - ECP + - CPU + - CPU + - CPU + * - MP2 + - GPU + - CPU + - CPU + * - CCSD + - GPU + - CPU + - NA diff --git a/source/user/solvent.rst b/source/user/solvent.rst index 7fd31b57f..d613a0167 100644 --- a/source/user/solvent.rst +++ b/source/user/solvent.rst @@ -8,7 +8,7 @@ Solvation models .. module:: solvent :synopsis: Solvation models and solvent effects -.. sectionauthor:: Qiming Sun . +.. sectionauthor:: Qiming Sun and Xiaojie Wu . Introduction ============ @@ -16,18 +16,18 @@ Introduction Solvation model allows the quantum chemistry calculations to include the interactions between solvents and the quantum solute. Solvents can be treated implicitly, known as continuum solvents, and explicitly. For continuum solvents, -we implemented the ddCOSMO (domain-decomposition COSMO solvation model). For -the explicit solvent environment, we provided the interface to call the -polarizable embedding library CPPE. +we implemented the PCM (polarizable continuum model), ddCOSMO (domain-decomposition COSMO solvation model), +ddPCM (domain-decomposition polarizable continuum model), +and SMD (Solvation Model Density). For the explicit solvent environment, +we provided the interface to call the polarizable embedding library CPPE. +PCM model +========= +PySCF support four types of PCM solvent models, i.e. C-PCM, IEF-PCM, SS(V)PE, and COSMO +(See https://manual.q-chem.com/5.2/Ch12.S2.SS2.html for more detailed descriptions of these methods). +The analytical gradient and semi-analytical Hessian are also supported. -ddCOSMO -======= - -Self-consistent solvents for ground state ------------------------------------------ - -Solvent model can be applied on to an SCF object:: +PCM solvent models can be applied on to an SCF object:: import pyscf mol = pyscf.M(atom=''' @@ -36,7 +36,10 @@ Solvent model can be applied on to an SCF object:: H 0. 0.935 -1.082 H 0. -0.935 -1.082''', basis='6-31g*', verbose=4) - mf = mol.RKS(xc='b3lyp').DDCOSMO().run() + mf = mol.RKS(xc='b3lyp').PCM() + mf.with_solvent.method = 'IEF-PCM' # C-PCM, SS(V)PE, COSMO + mf.with_solvent.eps = 78.3553 # for water + mf.run() In regular MCSCF (CASCI or CASSCF), and post-SCF calculations, the setup for self-consistent solvent is similar:: @@ -48,10 +51,9 @@ self-consistent solvent is similar:: H 0. 0.935 -1.082 H 0. -0.935 -1.082''', basis='6-31g*', verbose=4) - mc = mol.CASCI(4, 4).DDCOSMO().run() - - mp2_model = mol.MP2().DDCOSMO().run() + mc = mol.CASCI(4, 4).PCM().run() + mp2_model = mol.MP2().PCM().run() Solvent for excited states -------------------------- @@ -72,7 +74,7 @@ TDDFT method:: H 0. 0.935 -1.082 H 0. -0.935 -1.082''', basis='6-31g*', verbose=4) - mf = mol.RHF().ddCOSMO().run() + mf = mol.RHF().PCM().run() td = mf.TDA().run() In the diabatic excitations, we would like to let the solvent rapidly responds @@ -81,8 +83,8 @@ should converge to an equilibrium state between solvent and the excited state of the solute. In this scenario, solvent model should be applied on to the excited state methods:: - mf = mol.RHF().ddCOSMO().run() - td = mf.TDA().ddCOSMO() + mf = mol.RHF().PCM().run() + td = mf.TDA().PCM() td.with_solvent.equilibrium_solvation = True td.kernel() @@ -91,7 +93,7 @@ in this case. PySCF by default assumes the slow solvent model for TDDFT. In the complicated procedure which involves for example electronic states from different states (typically in the MCSCF calculations with state-average or -state-specific approximations), PySCF ddCOSMO implementation allows to input a +state-specific approximations), PySCF PCM implementation allows to input a density matrix and freeze the solvent equilibrated against the input density matrix:: @@ -102,7 +104,7 @@ matrix:: H 0. 0.935 -1.082 H 0. -0.935 -1.082''', basis='6-31g*', verbose=4) - mc = mol.CASCI(4, 4).DDCOSMO() + mc = mol.CASCI(4, 4).PCM() mc.fcisolver.nstates = 5 mc.with_solvent.state_id = 1 # Slow solvent wrt the first excited state mc.run() @@ -126,12 +128,12 @@ solute:: H 0. 1.035 -1.082 H 0. -1.035 -1.082''', basis='6-31g*', verbose=4) - mc = mol.CASCI(4, 4).DDCOSMO(dm=scf_dm).run() + mc = mol.CASCI(4, 4).PCM(dm=scf_dm).run() Solvent parameters ------------------ -The default solvent in the ddCOSMO module is water. When studying other types of +The default solvent in the PCM module is water. When studying other types of solvents, you can consider to modify the dielectric parameter ``eps`` using the constants listed below:: @@ -142,203 +144,67 @@ constants listed below:: H 0. 0.935 -1.082 H 0. -0.935 -1.082''', basis='6-31g*', verbose=4) - mf = mol.RHF().ddCOSMO() - mf.with_solvent.eps = 32.613 # methanol + mf = mol.RHF().PCM() + mf.with_solvent.eps = 32.613 # methanol mf.run() These dielectric constants are obtained from https://gaussian.com/scrf/. More dataset can be found in Minnesota Solvent Descriptor Database (https://comp.chem.umn.edu/solvation) -================================== ==================== -Solvent dielectric constant -================================== ==================== -Water 78.3553 -Acetonitrile 35.688 -Methanol 32.613 -Ethanol 24.852 -IsoQuinoline 11.00 -Quinoline 9.16 -Chloroform 4.7113 -DiethylEther 4.2400 -Dichloromethane 8.93 -DiChloroEthane 10.125 -CarbonTetraChloride 2.2280 -Benzene 2.2706 -Toluene 2.3741 -ChloroBenzene 5.6968 -NitroMethane 36.562 -Heptane 1.9113 -CycloHexane 2.0165 -Aniline 6.8882 -Acetone 20.493 -TetraHydroFuran 7.4257 -DiMethylSulfoxide 46.826 -Argon 1.430 -Krypton 1.519 -Xenon 1.706 -n-Octanol 9.8629 -1,1,1-TriChloroEthane 7.0826 -1,1,2-TriChloroEthane 7.1937 -1,2,4-TriMethylBenzene 2.3653 -1,2-DiBromoEthane 4.9313 -1,2-EthaneDiol 40.245 -1,4-Dioxane 2.2099 -1-Bromo-2-MethylPropane 7.7792 -1-BromoOctane 5.0244 -1-BromoPentane 6.269 -1-BromoPropane 8.0496 -1-Butanol 17.332 -1-ChloroHexane 5.9491 -1-ChloroPentane 6.5022 -1-ChloroPropane 8.3548 -1-Decanol 7.5305 -1-FluoroOctane 3.89 -1-Heptanol 11.321 -1-Hexanol 12.51 -1-Hexene 2.0717 -1-Hexyne 2.615 -1-IodoButane 6.173 -1-IodoHexaDecane 3.5338 -1-IodoPentane 5.6973 -1-IodoPropane 6.9626 -1-NitroPropane 23.73 -1-Nonanol 8.5991 -1-Pentanol 15.13 -1-Pentene 1.9905 -1-Propanol 20.524 -2,2,2-TriFluoroEthanol 26.726 -2,2,4-TriMethylPentane 1.9358 -2,4-DiMethylPentane 1.8939 -2,4-DiMethylPyridine 9.4176 -2,6-DiMethylPyridine 7.1735 -2-BromoPropane 9.3610 -2-Butanol 15.944 -2-ChloroButane 8.3930 -2-Heptanone 11.658 -2-Hexanone 14.136 -2-MethoxyEthanol 17.2 -2-Methyl-1-Propanol 16.777 -2-Methyl-2-Propanol 12.47 -2-MethylPentane 1.89 -2-MethylPyridine 9.9533 -2-NitroPropane 25.654 -2-Octanone 9.4678 -2-Pentanone 15.200 -2-Propanol 19.264 -2-Propen-1-ol 19.011 -3-MethylPyridine 11.645 -3-Pentanone 16.78 -4-Heptanone 12.257 -4-Methyl-2-Pentanone 12.887 -4-MethylPyridine 11.957 -5-Nonanone 10.6 -AceticAcid 6.2528 -AcetoPhenone 17.44 -a-ChloroToluene 6.7175 -Anisole 4.2247 -Benzaldehyde 18.220 -BenzoNitrile 25.592 -BenzylAlcohol 12.457 -BromoBenzene 5.3954 -BromoEthane 9.01 -Bromoform 4.2488 -Butanal 13.45 -ButanoicAcid 2.9931 -Butanone 18.246 -ButanoNitrile 24.291 -ButylAmine 4.6178 -ButylEthanoate 4.9941 -CarbonDiSulfide 2.6105 -Cis-1,2-DiMethylCycloHexane 2.06 -Cis-Decalin 2.2139 -CycloHexanone 15.619 -CycloPentane 1.9608 -CycloPentanol 16.989 -CycloPentanone 13.58 -Decalin-mixture 2.196 -DiBromomEthane 7.2273 -DiButylEther 3.0473 -DiEthylAmine 3.5766 -DiEthylSulfide 5.723 -DiIodoMethane 5.32 -DiIsoPropylEther 3.38 -DiMethylDiSulfide 9.6 -DiPhenylEther 3.73 -DiPropylAmine 2.9112 -e-1,2-DiChloroEthene 2.14 -e-2-Pentene 2.051 -EthaneThiol 6.667 -EthylBenzene 2.4339 -EthylEthanoate 5.9867 -EthylMethanoate 8.3310 -EthylPhenylEther 4.1797 -FluoroBenzene 5.42 -Formamide 108.94 -FormicAcid 51.1 -HexanoicAcid 2.6 -IodoBenzene 4.5470 -IodoEthane 7.6177 -IodoMethane 6.8650 -IsoPropylBenzene 2.3712 -m-Cresol 12.44 -Mesitylene 2.2650 -MethylBenzoate 6.7367 -MethylButanoate 5.5607 -MethylCycloHexane 2.024 -MethylEthanoate 6.8615 -MethylMethanoate 8.8377 -MethylPropanoate 6.0777 -m-Xylene 2.3478 -n-ButylBenzene 2.36 -n-Decane 1.9846 -n-Dodecane 2.0060 -n-Hexadecane 2.0402 -n-Hexane 1.8819 -NitroBenzene 34.809 -NitroEthane 28.29 -n-MethylAniline 5.9600 -n-MethylFormamide-mixture 181.56 -n,n-DiMethylAcetamide 37.781 -n,n-DiMethylFormamide 37.219 -n-Nonane 1.9605 -n-Octane 1.9406 -n-Pentadecane 2.0333 -n-Pentane 1.8371 -n-Undecane 1.9910 -o-ChloroToluene 4.6331 -o-Cresol 6.76 -o-DiChloroBenzene 9.9949 -o-NitroToluene 25.669 -o-Xylene 2.5454 -Pentanal 10.0 -PentanoicAcid 2.6924 -PentylAmine 4.2010 -PentylEthanoate 4.7297 -PerFluoroBenzene 2.029 -p-IsoPropylToluene 2.2322 -Propanal 18.5 -PropanoicAcid 3.44 -PropanoNitrile 29.324 -PropylAmine 4.9912 -PropylEthanoate 5.5205 -p-Xylene 2.2705 -Pyridine 12.978 -sec-ButylBenzene 2.3446 -tert-ButylBenzene 2.3447 -TetraChloroEthene 2.268 -TetraHydroThiophene-s,s-dioxide 43.962 -Tetralin 2.771 -Thiophene 2.7270 -Thiophenol 4.2728 -trans-Decalin 2.1781 -TriButylPhosphate 8.1781 -TriChloroEthene 3.422 -TriEthylAmine 2.3832 -Xylene-mixture 2.3879 -z-1,2-DiChloroEthene 9.2 -================================== ==================== +ddCOSMO +======= + +Self-consistent solvents for ground state +----------------------------------------- + +Solvent model can be applied on to an SCF object:: + + import pyscf + mol = pyscf.M(atom=''' + C 0. 0. -0.542 + O 0. 0. 0.677 + H 0. 0.935 -1.082 + H 0. -0.935 -1.082''', + basis='6-31g*', verbose=4) + mf = mol.RKS(xc='b3lyp').DDCOSMO().run() + + +SMD model +========= +SMD model is recommended for computing solvation free energy. The implementation of SMD model +in PySCF is based on IEF-PCM. Other SMx models are not supported yet +(See https://manual.q-chem.com/5.2/Ch12.S2.SS8.html). The source code for CDS contribution +is taken from NWChem (https://github.com/nwchemgit/nwchem/blob/master/src/solvation/mnsol.F) +SMD solvent solvent models can be applied on to an SCF object:: + + import pyscf + mol = pyscf.M(atom=''' + C 0. 0. -0.542 + O 0. 0. 0.677 + H 0. 0.935 -1.082 + H 0. -0.935 -1.082''', + basis='6-31g*', verbose=4) + mf = mol.RKS(xc='b3lyp').SMD() + mf.with_solvent.solvent = 'water' + mf.run() + +The solvent descriptor can be assigned to the ``with_solvent`` directly:: + + import pyscf + mol = pyscf.M(atom=''' + C 0. 0. -0.542 + O 0. 0. 0.677 + H 0. 0.935 -1.082 + H 0. -0.935 -1.082''', + basis='6-31g*', verbose=4) + mf = mol.RKS(xc='b3lyp').SMD() + mf.with_solvent.solvent_descriptors = [1.3843, 1.3766, 0.0, 0.45, 35.06, 13.45, 0.0, 0.0] + mf.run() +The format of solvant names are the same as Minnesota Solvent Descriptor Database +(https://comp.chem.umn.edu/solvation/mnsddb.pdf). One can assign solvent descriptors in the format +``mf.with_solvent.solvent_descriptors = [n, n25, alpha, beta, gamma, epsilon, phi, psi]`` Polarizable embedding =====================