1C> \ingroup nwpwxc 2C> @{ 3C> 4C> \file nwpwxc_c_b97.F 5C> The B97 family of correlation functionals 6C> 7C> @} 8#include "nwpwxcP_xc_b97.h" 9C> 10C> \ingroup nwpwxc_priv 11C> @{ 12C> 13C> \brief Evaluate the B97 family of correlation functionals 14C> 15C> This code evaluates correlation functionals from the 16C> B97 family of functionals [1,2]. 17C> 18C> ### References ### 19C> 20C> [1] A.D. Becke, "Density-functional thermochemistry. V. Systematic 21C> optimization of exchange-correlation functionals", J. Chem. Phys. 22C> 107 (1997) 8554-8560, DOI: 23C> <a href="https://doi.org/10.1063/1.475007"> 24C> 10.1063/1.475007</a>. 25C> 26C> [2] S. Grimme, "Semiempirical GGA-type density functional constructed 27C> with a long-range dispersion correction", J. Comput. Chem. 27 28C> (2006) 1787-1799, DOI: 29C> <a href="https://doi.org/10.1002/jcc.20495"> 30C> 10.1002/jcc.20495</a>. 31C> 32 Subroutine nwpwxc_c_b97(param,tol_rho,ipol,nq,wght,rho,rgamma, 33 & func,Amat,Cmat) 34c 35c $Id$ 36c 37 implicit none 38c 39#include "nwpwxc_param.fh" 40c 41c Input and other parameters 42c 43 double precision param(*)!< [Input] Parameters of functional as 44 !< defined in [1]: 45 !< - param(1): \f$m\f$ of Eqs.(20). 46 !< - param(2): \f$C_{C\sigma\sigma,0}\f$ 47 !< - param(3): \f$C_{C\alpha\beta,0}\f$ 48 !< - param(4): \f$C_{C\sigma\sigma,1}\f$ 49 !< - param(5): \f$C_{C\alpha\beta,1}\f$ 50 !< - param(6): \f$C_{C\sigma\sigma,2}\f$ 51 !< - param(7): \f$C_{C\alpha\beta,2}\f$ 52 !< - param(8): \f$C_{C\sigma\sigma,3}\f$ 53 !< - param(9): \f$C_{C\alpha\beta,3}\f$ 54 !< - param(10): \f$C_{C\sigma\sigma,4}\f$ 55 !< - param(11): \f$C_{C\alpha\beta,4}\f$ 56 double precision tol_rho !< [Input] The lower limit on the density 57 integer ipol !< [Input] The number of spin channels 58 integer nq !< [Input] The number of points 59 double precision wght !< [Input] The weight of the functional 60c 61c Charge Density 62c 63 double precision rho(nq,*) !< [Input] The density 64c 65c Charge Density Gradient 66c 67 double precision rgamma(nq,*) !< [Input] The norm of the density gradients 68c 69c Sampling Matrices for the XC Potential 70c 71 double precision func(nq) !< [Output] The value of the functional 72 double precision Amat(nq,*) !< [Output] The derivative wrt rho 73 double precision Cmat(nq,*) !< [Output] The derivative wrt rgamma 74c 75c Local variables 76c 77 integer i 78 double precision rho_a(0:1) 79 double precision rho_b(0:1) 80 double precision FC(0:_FXC_NUMDERI) 81c 82c Code 83c 84 if (ipol.eq.1) then 85 do i = 1, nq 86 rho_a(0) = rho(i,R_T)*0.5d0 87 rho_b(0) = rho_a(0) 88 rho_a(1) = rgamma(i,G_TT)*0.25d0 89 rho_b(1) = rho_a(1) 90 if (rho_a(0).gt.tol_rho) then 91 call nwpwxcp_c_b97(rho_a,rho_b,wght,tol_rho,FC,param) 92 func(i) = func(i) + FC(_FXC_E)*wght 93 Amat(i,D1_RA) = Amat(i,D1_RA) + FC(_FXC_RA)*wght 94 Cmat(i,D1_GAA) = Cmat(i,D1_GAA) + FC(_FXC_GAA)*wght 95 endif 96 enddo 97 else 98 do i = 1, nq 99 rho_a(0) = rho(i,R_A) 100 rho_b(0) = rho(i,R_B) 101 rho_a(1) = rgamma(i,G_AA) 102 rho_b(1) = rgamma(i,G_BB) 103 if (rho_a(0).gt.tol_rho.or.rho_b(0).gt.tol_rho) then 104 call nwpwxcp_c_b97(rho_a,rho_b,wght,tol_rho,FC,param) 105 func(i) = func(i) + FC(_FXC_E)*wght 106 Amat(i,D1_RA) = Amat(i,D1_RA) + FC(_FXC_RA)*wght 107 Amat(i,D1_RB) = Amat(i,D1_RB) + FC(_FXC_RB)*wght 108 Cmat(i,D1_GAA) = Cmat(i,D1_GAA) + FC(_FXC_GAA)*wght 109 Cmat(i,D1_GBB) = Cmat(i,D1_GBB) + FC(_FXC_GBB)*wght 110 endif 111 enddo 112 endif 113c 114 end 115c 116 Subroutine nwpwxc_c_b97_d2(param,tol_rho,ipol,nq,wght,rho,rgamma, 117 & func,Amat,Amat2,Cmat,Cmat2) 118c 119c $Id$ 120c 121 implicit none 122c 123#include "nwpwxc_param.fh" 124c 125c Input and other parameters 126c 127 double precision param(*)!< [Input] Parameters of functional as 128 !< defined in [1]: 129 !< - param(1): \f$m\f$ of Eqs.(20). 130 !< - param(2): \f$C_{C\sigma\sigma,0}\f$ 131 !< - param(3): \f$C_{C\alpha\beta,0}\f$ 132 !< - param(4): \f$C_{C\sigma\sigma,1}\f$ 133 !< - param(5): \f$C_{C\alpha\beta,1}\f$ 134 !< - param(6): \f$C_{C\sigma\sigma,2}\f$ 135 !< - param(7): \f$C_{C\alpha\beta,2}\f$ 136 !< - param(8): \f$C_{C\sigma\sigma,3}\f$ 137 !< - param(9): \f$C_{C\alpha\beta,3}\f$ 138 !< - param(10): \f$C_{C\sigma\sigma,4}\f$ 139 !< - param(11): \f$C_{C\alpha\beta,4}\f$ 140 double precision tol_rho !< [Input] The lower limit on the density 141 integer ipol !< [Input] The number of spin channels 142 integer nq !< [Input] The number of points 143 double precision wght !< [Input] The weight of the functional 144c 145c Charge Density 146c 147 double precision rho(nq,*) !< [Input] The density 148c 149c Charge Density Gradient 150c 151 double precision rgamma(nq,*) !< [Input] The norm of the density gradients 152c 153c Sampling Matrices for the XC Potential 154c 155 double precision func(nq) !< [Output] The value of the functional 156 double precision Amat(nq,*) !< [Output] The derivative wrt rho 157 double precision Cmat(nq,*) !< [Output] The derivative wrt rgamma 158 double precision Amat2(nq,*) !< [Output] The 2nd derivative wrt rho 159 double precision Cmat2(nq,*) !< [Output] The 2nd derivative wrt 160 !< rho and rgamma 161c 162c Local variables 163c 164 integer i 165 double precision rho_a(0:1) 166 double precision rho_b(0:1) 167 double precision FC(0:_FXC_NUMDERI) 168c 169c Code 170c 171 if (ipol.eq.1) then 172 do i = 1, nq 173 rho_a(0) = rho(i,R_T)*0.5d0 174 rho_b(0) = rho_a(0) 175 rho_a(1) = rgamma(i,G_TT)*0.25d0 176 rho_b(1) = rho_a(1) 177 if (rho_a(0).gt.tol_rho) then 178 call nwpwxcp_c_b97(rho_a,rho_b,wght,tol_rho,FC,param) 179 func(i) = func(i) + FC(_FXC_E)*wght 180 Amat(i,D1_RA) = Amat(i,D1_RA) + FC(_FXC_RA)*wght 181 Cmat(i,D1_GAA) = Cmat(i,D1_GAA) + FC(_FXC_GAA)*wght 182 Amat2(i,D2_RA_RA) = Amat2(i,D2_RA_RA) + FC(_FXC_RARA)*wght 183 Amat2(i,D2_RA_RB) = Amat2(i,D2_RA_RB) + FC(_FXC_RARB)*wght 184 Cmat2(i,D2_RA_GAA) = Cmat2(i,D2_RA_GAA) 185 & + FC(_FXC_RAGAA)*wght 186 Cmat2(i,D2_RA_GBB) = Cmat2(i,D2_RA_GBB) 187 & + FC(_FXC_RAGBB)*wght 188 Cmat2(i,D2_GAA_GAA) = Cmat2(i,D2_GAA_GAA) 189 & + FC(_FXC_GAAGAA)*wght 190 Cmat2(i,D2_GAA_GBB) = Cmat2(i,D2_GAA_GBB) 191 & + FC(_FXC_GAAGBB)*wght 192 endif 193 enddo 194 else 195 do i = 1, nq 196 rho_a(0) = rho(i,R_A) 197 rho_b(0) = rho(i,R_B) 198 rho_a(1) = rgamma(i,G_AA) 199 rho_b(1) = rgamma(i,G_BB) 200 if (rho_a(0).gt.tol_rho.or.rho_b(0).gt.tol_rho) then 201 call nwpwxcp_c_b97(rho_a,rho_b,wght,tol_rho,FC,param) 202 func(i) = func(i) + FC(_FXC_E)*wght 203 Amat(i,D1_RA) = Amat(i,D1_RA) + FC(_FXC_RA)*wght 204 Amat(i,D1_RB) = Amat(i,D1_RB) + FC(_FXC_RB)*wght 205 Cmat(i,D1_GAA) = Cmat(i,D1_GAA) + FC(_FXC_GAA)*wght 206 Cmat(i,D1_GBB) = Cmat(i,D1_GBB) + FC(_FXC_GBB)*wght 207 Amat2(i,D2_RA_RA) = Amat2(i,D2_RA_RA) + FC(_FXC_RARA)*wght 208 Amat2(i,D2_RA_RB) = Amat2(i,D2_RA_RB) + FC(_FXC_RARB)*wght 209 Amat2(i,D2_RB_RB) = Amat2(i,D2_RB_RB) + FC(_FXC_RBRB)*wght 210 Cmat2(i,D2_RA_GAA) = Cmat2(i,D2_RA_GAA) 211 & + FC(_FXC_RAGAA)*wght 212 Cmat2(i,D2_RA_GBB) = Cmat2(i,D2_RA_GBB) 213 & + FC(_FXC_RAGBB)*wght 214 Cmat2(i,D2_RB_GAA) = Cmat2(i,D2_RB_GAA) 215 & + FC(_FXC_RBGAA)*wght 216 Cmat2(i,D2_RB_GBB) = Cmat2(i,D2_RB_GBB) 217 & + FC(_FXC_RBGBB)*wght 218 Cmat2(i,D2_GAA_GAA) = Cmat2(i,D2_GAA_GAA) 219 & + FC(_FXC_GAAGAA)*wght 220 Cmat2(i,D2_GAA_GBB) = Cmat2(i,D2_GAA_GBB) 221 & + FC(_FXC_GAAGBB)*wght 222 Cmat2(i,D2_GBB_GBB) = Cmat2(i,D2_GBB_GBB) 223 & + FC(_FXC_GBBGBB)*wght 224 endif 225 enddo 226 endif 227c 228 end 229c 230 subroutine nwpwxcp_c_pwlda(ra,rb,FCLDA) 231 implicit none 232c 233 double precision ra 234 double precision rb 235 double precision FCLDA(0:_FCLDA_ELEMENTS) 236c 237 double precision ec 238 double precision rho(2) 239 double precision Amat(2) 240 double precision Amat2(3) 241c 242 rho(R_A) = ra 243 rho(R_B) = rb 244 ec = 0.0d0 245 Amat(D1_RA) = 0.0d0 246 Amat(D1_RB) = 0.0d0 247 Amat2(D2_RA_RA) = 0.0d0 248 Amat2(D2_RA_RB) = 0.0d0 249 Amat2(D2_RB_RB) = 0.0d0 250c 251 call nwpwxc_c_pw91lda_d2(1.0d-20,2,1,1.0d0,rho,ec,Amat,Amat2) 252c 253 FCLDA(_FXC_E) = ec 254 FCLDA(_FXC_RA) = Amat(D1_RA) 255 FCLDA(_FXC_RB) = Amat(D1_RB) 256 FCLDA(_FXC_RARA) = Amat2(D2_RA_RA) 257 FCLDA(_FXC_RARB) = Amat2(D2_RA_RB) 258 FCLDA(_FXC_RBRB) = Amat2(D2_RB_RB) 259c 260 end 261C> 262C> @} 263