1% 2% $Id$ 3% 4\label{sec:ecp} 5\def\ell{l} 6Effective core potentials (ECPs) are a useful means of replacing the core 7electrons in a calculation with an effective potential, thereby eliminating 8the need for the core basis functions, which usually require a large set of 9Gaussians to describe them. In addition to replacing the core, they may be 10used to represent relativistic effects, which are largely confined to the 11core. In this context, both the scalar (spin-free) relativistic effects and 12spin-orbit (spin-dependent) relativistic effects may be included in 13effective potentials. NWChem has the facility to use both, and these are 14described in the next two sections. 15 16A brief recapitulation of the development of RECPs is given here, following 17Pacios and Christiansen\footnote{l.~F.~Pacios and P.~A.~Christiansen, 18J.~Chem.~Phys.~{\bf 82}, 2664 (1985)}. The process can be viewed as starting 19from an atomic Dirac-Hartree-Fock calculation, done in {\it jj} coupling, 20and producing relativistic effective potentials (REPs) for each $\ell$ and 21$j$ value, $U^{\rm REP}_{\ell j}$. From these, a local potential is 22extracted, which for example contains the Coulomb potential of the core 23electrons balanced by the part of the nuclear attraction which cancels the 24core electron charge. The residue is expressed in a semi-local form, 25\begin{equation} 26U^{\rm REP} = U^{\rm REP}_{LJ}(r) + \sum_{\ell=0}^{L-1} 27\sum_{j=|\ell-1/2}^{\ell+1/2} \left[ U^{\rm REP}_{\ell j}(r) - 28U^{\rm REP}_{LJ}(r) \right] \sum_m | \ell j m \rangle \langle \ell j m | 29\end{equation} 30where $L$ is one larger than the maximum angular momentum in the atom. 31The scalar potential is obtained by averaging the REPs for each $j$ for a 32given $\ell$ to give an averaged relativistic effective potential, or AREP, 33\begin{equation} 34U^{\rm AREP}_\ell(r) = \frac{1}{2\ell+1} \left[ \ell U^{\rm REP}_{\ell-1/2}(r) 35+ (\ell+1) U^{\rm REP}_{\ell+1/2}(r) \right]. 36\end{equation} 37These are summed into the full potential. 38%\begin{equation} 39%U^{\rm AREP} = U^{\rm AREP}_{L}(r) + \sum_{\ell=0}^{L-1} 40%\left[ U^{\rm AREP}_{\ell}(r) - U^{\rm AREP}_{L}(r) 41%\sum_m | \ell m \rangle \langle \ell m |. 42%\end{equation} 43 44The spin-orbit potential is obtained from the difference between the REPs 45for the two $j$ values for a given \ell, and may be 46represented in terms of an effective spin-orbit operator, 47\begin{equation} 48H^{\rm SO} = {\bf s} \cdot \sum_{\ell=1}^{L-1} \frac{2}{2\ell+1} 49\Delta U^{\rm REP}_{\ell} \sum_{mm'} 50| \ell m \rangle \langle \ell m | \hat\ell | \ell m' \rangle \langle \ell m' |. 51\end{equation} 52where 53\begin{equation} 54\Delta U^{\rm REP}_{\ell} = U^{\rm REP}_{\ell+1/2}(r) 55 - U^{\rm REP}_{\ell-1/2}(r). 56\end{equation} 57The spin-orbit integrals generated by NWChem are the integrals over the sum, 58including the factor of $2/(2\ell+1)$, so that they may be treated as an 59effective spin-orbit operator without further factors introduced. 60 61The effective potentials, both scalar and spin-orbit, are fitted to 62Gaussians with the form 63\[ 64 r^2U_l(r) = \sum_{k} A_{lk} r^{n_{lk}} e^{-B_{lk}r^{2}} 65\] 66where $A_{lk}$ is the contraction coefficient, $n_{lk}$ is the 67exponent of the ``r'' term (r-exponent), and $B_{lk}$ is the Gaussian 68exponent. The $n_{lk}$ is shifted by 2, in accordance with most of the ECP 69literature and implementations, i.e., an $n_{lk} = 0$ implies 70$r^{-2}$. The current implementation allows $n_{lk}$ values 71of only 0, 1, or 2. 72 73\section{Scalar ECPs} 74\label{sec:scalar_ecp} 75 76The optional directive \verb+ECP+ allows the user to describe an effective core 77potential (ECP) in terms of contracted Gaussian functions as given above. 78Potentials using these functions must be specified explicitly by user input 79in the \verb+ECP+ directive. This directive has essentially the same form 80and properties as the standard \verb+BASIS+ directive, except for essential 81differences required for ECPs. Because of this, the ECP is treated 82internally as a basis set. The form of the input for the 83\verb+ECP+ directive is as follows: 84 85% [spherical || cartesian default cartesian] 86% [segment || nosegment default segment] 87 88\begin{verbatim} 89 ECP [<string name default "ecp basis">] \ 90 [print || noprint default print] 91 92 <string tag> library [<string tag_in_lib>] \ 93 <string standard_set> [file <filename>] \ 94 [except <string tag list>] 95 96 <string tag> [nelec] <integer number_of_electrons_replaced> 97 98 ... 99 100 <string tag> <string shell_type> 101 <real r-exponent> <real Gaussian-exponent> <real list_of_coefficients> 102 ... 103 104 END 105\end{verbatim} 106 107ECPs are automatically segmented, even if general contractions are input. 108The projection operators defined in an ECP are spherical by default, so 109there is no need to include the \verb+CARTESIAN+ or \verb+SPHERICAL+ keyword 110as there is for a standard basis set. ECPs are associated with centers in 111geometries through tags or names of centers. These tags must match in the 112same manner as for basis sets the tags in a \verb+GEOMETRY+ and 113\verb+ECP+ directives, and are limited to sixteen (16) characters. 114Each center with the same tag will have the same ECP. By default, the 115input module prints each ECP that it encounters. The \verb+NOPRINT+ 116option can be used to disable printing. There can be only one active 117ECP, even though several may exist in the input deck. The ECP modules 118load ``ecp basis'' inputs along with any ``ao basis'' inputs present. 119ECPs may be used in both energy and gradient calculations. 120 121ECPs are named in the same fashion as geometries or regular basis 122sets, with the default name being \verb+"ecp basis"+. It should be 123clear from the above discussion on geometries and database entries how 124indirection is supported. All directives that are in common with the 125standard Gaussian basis set input have the same function and syntax. 126 127As for regular basis sets, ECPs may be obtained from the standard library. 128The names of the sets of ECPs available in the standard 129library (their coverage is described in Appendix \ref{sec:knownbasis}) are 130\begin{itemize} 131\item \verb,"Hay-Wadt MB (n+1) ECP", 132\item \verb,"Hay-Wadt VDZ (n+1) ECP", 133\item \verb+"LANL2DZ ECP"+ 134\item \verb+"SBKJC VDZ ECP"+ 135\item \verb+"Stuttgart RLC ECP"+ 136\item \verb+"Stuttgart RSC ECP"+ 137\item \verb+"CRENBL ECP"+ 138\item \verb+"CRENBS ECP"+ 139\end{itemize} 140 141The keyword \verb+nelec+ allows the user to specify the number of core 142electrons replaced by the ECP. Additional input lines define the 143specific coefficients and exponents. The variable \verb+<shell_type>+ 144is used to specify the components of the ECP. The keyword \verb+ul+ 145entered for \verb+<shell_type>+ denotes the local part of the ECP. 146This is equivalent to the highest angular momentum functions specified 147in the literature for most ECPs. The standard entries (\verb+s, p, d+, 148etc.) for \verb+shell_type+ specify the angular momentum projector 149onto the local function. The shell type label of \verb+s+ indicates 150the \verb+ul-s+ projector input, \verb+p+ indicates the \verb+ul-p+, 151etc. 152 153For example, the Christiansen, Ross and Ermler ARECPs are available in 154the standard basis set libary named \verb+{crenbl_ecp}+. To perform a 155calculation on uranyl (UO$_2^{2+}$) with all-electron oxygen 156(aug-cc-pvdz basis), and uranium with an ARECP and using the 157corresponding basis the following input can be used 158\begin{verbatim} 159 geometry 160 U 0 0 0 161 O 0 0 1.65 162 O 0 0 -1.65 163 end 164 basis 165 U library crenbl_ecp 166 O library aug-cc-pvdz 167 end 168 ecp 169 U library crenbl_ecp 170 end 171\end{verbatim} 172 173The following is an example of explicit input of an ECP for H$_2$CO. 174It defines an ECP for the carbon and oxygen atoms in the molecule. 175 176% \centerline{{\bf H$_2$CO }} 177 178\begin{verbatim} 179 ecp 180 C nelec 2 # ecp replaces 2 electrons on C 181 C ul # d 182 1 80.0000000 -1.60000000 183 1 30.0000000 -0.40000000 184 2 0.5498205 -0.03990210 185 C s # s - d 186 0 0.7374760 0.63810832 187 0 135.2354832 11.00916230 188 2 8.5605569 20.13797020 189 C p # p - d 190 2 10.6863587 -3.24684280 191 2 23.4979897 0.78505765 192 O nelec 2 # ecp replaces 2 electrons on O 193 O ul # d 194 1 80.0000000 -1.60000000 195 1 30.0000000 -0.40000000 196 2 1.0953760 -0.06623814 197 O s # s - d 198 0 0.9212952 0.39552179 199 0 28.6481971 2.51654843 200 2 9.3033500 17.04478500 201 O p # p - s 202 2 52.3427019 27.97790770 203 2 30.7220233 -16.49630500 204 end 205\end{verbatim} 206 207Various ECPs without a local function are available, including those of 208the Stuttgart group. For those, no "ul" part needs to be defined. To 209define the absence of the local potential, simply specify one contraction 210with a zero coefficient: 211 212\begin{verbatim} 213 <string tag> ul 214 2 1.00000 0.00000 215\end{verbatim} 216 217\section{Spin-orbit ECPs} 218\label{sec:spinorb_ecp} 219 220The Spin-orbit ECPs can be used with the Density Functional Approach, but 221one has to run the calculations without symmetry. Note: when a Hartree-Fock 222method is specified the spin-orbit input will be ignored. 223 224Spin-orbit ECPs are fitted in precisely the same functional form as the 225scalar RECPs and have the same properties, with the exception that there is 226no local potential ul, no $s$ potential and no effective charge has to be 227defined. Spin-orbit potentials are 228specified in the same way as ECPs except that the directive \verb+SO+ is 229used instead of \verb+ECP+. Note that there currently are no spin-orbit 230ECPs defined in the standard NWChem library. The \verb+SO+ 231directive is as follows: 232 233\begin{verbatim} 234 SO [<string name default "so basis">] \ 235 [print || noprint default print] 236 237 <string tag> library [<string tag_in_lib>] \ 238 <string standard_set> [file <filename>] 239 [except <string tag list>] 240 ... 241 242 <string tag> <string shell_type> 243 <real r-exponent> <real Gaussian-exponent> <real list_of_coefficients> 244 ... 245 246 END 247\end{verbatim} 248 249Note: in the literature the coefficients of the spin-orbit potentials are NOT 250always defined in the same manner. The NWChem code assumes that the spin-orbit 251potential defined in the input is of the form: 252\begin{equation} 253\Delta U^{\rm NWChem}_{\ell} = \frac{2}{2\ell+1} \Delta U_{\ell} 254\end{equation} 255For example, in the literature the Stuttgart potentials are defined as 256$\Delta U_{\ell}$ and, hence, have to be multiplied by $2/(2{\ell}+1)$. On the 257other hand, the CRENBL potentials in the published papers are defined as 258$\frac{\ell}{2\ell+1} \Delta U_{\ell}$ and, hence, have to be multiplied by 259$2/{\ell}$ (Warning: on the CRENBL website the spin-orbit potentials already 260have been corrected with the $2/{\ell}$ factor). 261 262%RJH: Align six comment symbols (#) with the other two.--fmr. 263%Generally contracted ECP basis sets are not in wide use, but the 264%functionality has been included in NWChem for applications where they 265%might be useful. 266 267