1#!/usr/bin/env python 2 3''' 4User defined XC functional 5 6See also 7* The parser parse_xc function implemented in pyscf.dft.libxc 8* Example 24-define_xc_functional.py to input a functional which is not 9 provided in Libxc or XcFun library. 10''' 11 12from pyscf import gto 13from pyscf import dft 14 15mol = gto.M( 16 atom = ''' 17 O 0. 0. 0. 18 H 0. -0.757 0.587 19 H 0. 0.757 0.587 ''', 20 basis = 'ccpvdz') 21 22# 23# DFT can parse the custom XC functional, following the rules: 24# * The given functional description must be a one-line string. 25# * The functional description is case-insensitive. 26# * The functional description string has two parts, separated by ",". The 27# first part describes the exchange functional, the second is the correlation 28# functional. 29# - If "," was not presented in string, the entire string is treated as a 30# compound XC functional (including both X and C functionals, such as b3lyp). 31# - To input only X functional (without C functional), leave the second part 32# blank. E.g. description='slater,' for pure LDA functional. 33# - To neglect the X functional (only input C functional), leave the first 34# part blank. E.g. description=',vwn' means pure VWN functional 35# - If compound XC functional is specified, no matter whehter it is in the X 36# part (the string in front of comma) or the C part (the string behind 37# comma), both X and C functionals of the compound XC functional will be 38# used. 39# * The functional name can be placed in arbitrary order. Two names needs to 40# be separated by operators + or -. Blank spaces are ignored. 41# NOTE the parser only reads operators + - *. / is not supported. 42# * A functional name can have at most one factor. If the factor is not 43# given, it is set to 1. Compound functional can be scaled as a unit. For 44# example '0.5*b3lyp' is equivalent to 45# 'HF*0.1 + .04*LDA + .36*B88, .405*LYP + .095*VWN' 46# * String "HF" stands for exact exchange (HF K matrix). It is allowed to 47# put "HF" in C (correlation) functional part. 48# * String "RSH" means range-separated operator. Its format is 49# RSH(omega, alpha, beta). Another way to input RSH is to use keywords 50# SR_HF and LR_HF: "SR_HF(0.1) * alpha_plus_beta" and "LR_HF(0.1) * 51# alpha" where the number in parenthesis is the value of omega. 52# * Be careful with the libxc convention on GGA functional, in which the LDA 53# contribution has been included. 54 55mf = dft.RKS(mol) 56# B3LYP can be constructed 57mf.xc = 'HF*0.2 + .08*LDA + .72*B88, .81*LYP + .19*VWN' 58e1 = mf.kernel() 59print('E = %.15g ref = -76.3832244350081' % e1) 60 61# 62# No correlation functional 63# 64mf.xc = '.2*HF + .08*LDA + .72*B88' 65e1 = mf.kernel() 66print('E = %.15g ref = -75.9807850596666' % e1) 67 68# 69# If not given, the factor for each functional equals 1 by default. 70# 71mf = dft.RKS(mol) 72mf.xc = 'b88,lyp' 73e1 = mf.kernel() 74 75mf = dft.RKS(mol) 76mf.xc = 'b88*1,lyp*1' 77eref = mf.kernel() 78print('%.15g == %.15g' % (e1, eref)) 79 80# 81# Compound functionals can be used as part of the definition 82# 83mf = dft.RKS(mol) 84mf.xc = '0.5*b3lyp, 0.5*lyp' 85e1 = mf.kernel() 86print('E = %.15g ref = -71.9508340443282' % e1) 87 88# Compound XC functional can be presented in the C part (the string behind 89# comma). Both X and C functionals of the compound XC functional will be used. 90# Compound XC functional can be scaled as a unit. 91# 92mf = dft.RKS(mol) 93mf.xc = '0.5*b3lyp, 0.5*b3p86' 94e1 = mf.kernel() 95 96mf = dft.RKS(mol) 97mf.xc = '0.5*b3lyp + 0.5*b3p86' 98e2 = mf.kernel() 99print('E1 = %.15g E2 = %.15g ref = -76.3923625924023' % (e1, e2)) 100 101# 102# More examples of customized functionals. NOTE These customized functionals 103# are presented for the purpose of demonstrating the feature of the XC input. 104# They are not reported in any literature. DO NOT use them in the actual 105# calculations. 106# 107# Half HF exchange plus half B3LYP plus half VWN functional 108mf.xc = '.5*HF+.5*B3LYP,VWN*.5' 109 110# "-" to subtract one functional from another 111mf.xc = 'B88 - SLATER*.5' 112 113# The functional below gives omega = 0.33, alpha = 0.6 * 0.65 = 0.39 114# beta = -0.46 * 0.6 + 0.4 * 0.2(from HF of B3P86) = -0.196 115mf.xc = '0.6*CAM_B3LYP+0.4*B3P86' 116 117# The input XC description does not depend on the order of functionals. The 118# functional below is the same to the functional above 119mf.xc = '0.4*B3P86+0.6*CAM_B3LYP' 120 121# Use SR_HF/LR_HF keywords to input range-separated functionals. 122# When exact HF exchange is presented, it is split into SR and LR parts 123# alpha = 0.8(from HF) + 0.22 124# beta = 0.5 + 0.8(from HF) 125mf.xc = '0.5*SR-HF(0.3) + .8*HF + .22*LR_HF' 126 127# RSH is another keyword to input range-separated functionals 128mf.xc = '0.5*RSH(0.3,2.04,0.56) + 0.5*BP86' 129 130# A shorthand to input 'PBE,PBE', which is a compound functional. Note the 131# shorthand input is different to the two examples 'PBE,' and ',PBE' below. 132mf.xc = 'PBE' 133 134# Exchange part only 135mf.xc = 'PBE,' 136 137# Correlation part only 138mf.xc = ',PBE' 139 140# When defining custom functionals, compound functional will affect both the 141# exchange and correlation parts 142mf.xc = 'PBE + SLATER*.5' 143 144# The above input is equivalent to 145mf.xc = 'PBE + SLATER*.5, PBE' 146 147