1# XCFun, an arbitrary order exchange-correlation library 2# Copyright (C) 2020 Ulf Ekström and contributors. 3# 4# This file is part of XCFun. 5# 6# This Source Code Form is subject to the terms of the Mozilla Public 7# License, v. 2.0. If a copy of the MPL was not distributed with this 8# file, You can obtain one at http://mozilla.org/MPL/2.0/. 9# 10# For information on the complete list of contributors to the 11# XCFun library, see: <https://xcfun.readthedocs.io/> 12 13import pytest 14 15import numpy 16from numpy.testing import assert_allclose 17 18import xcfun 19 20 21def test_splash(): 22 23 version = xcfun.xcfun_version() 24 splash = xcfun.xcfun_splash() 25 26 print("\n Sample use of python interface to xcfun") 27 28 print("\n\nXCFun version ", version) 29 print(splash) 30 31 32def test_lda_old(): 33 fun_xc_name = 'LDA' 34 fun_xc_weight = 1.0 35 fun_xc = xcfun.Functional({fun_xc_name: fun_xc_weight}) 36 37 n = numpy.array([1.0]) 38 result = fun_xc.eval_potential_n(n) 39 40 energy = result[0, 0] 41 potential = result[0, 1] 42 43 print("Functional used: ", fun_xc_name, " weight ", fun_xc_weight) 44 print(" density: ", n) 45 print(" energy: ", energy) 46 print(" potential: ", potential) 47 48 assert_allclose(energy, -0.8101513) 49 assert_allclose(potential, -1.06468341) 50 51 52@pytest.fixture 53def lda_fun(): 54 return xcfun.Functional({'LDA': 1.0}) 55 56 57@pytest.fixture 58def pbe_fun(): 59 return xcfun.Functional({'PBE': 1.0}) 60 61 62def test_type(lda_fun, pbe_fun): 63 assert lda_fun.type == 0 64 assert pbe_fun.type == 1 65 66 67@pytest.fixture 68def dens(): 69 return numpy.arange(0.1, 5.0, 0.5) 70 71 72@pytest.fixture 73def spindens(dens): 74 return 0.1 * numpy.sin(dens) 75 76 77@pytest.fixture 78def densgrad(dens): 79 dg = numpy.zeros((dens.shape[0], 3)) 80 dg[:, 0] = 0.2 * numpy.cos(dens) 81 dg[:, 1] = 0.3 * numpy.sin(dens) 82 dg[:, 2] = 0.1 * numpy.sin(dens)**2 83 return dg 84 85 86@pytest.fixture 87def denshess(dens): 88 dh = numpy.zeros((dens.shape[0], 6)) 89 dh[:, 0] = 0.5 90 dh[:, 1] = 0.1 * dens 91 dh[:, 2] = numpy.cos(dens)**2 92 dh[:, 3] = numpy.sin(dens) 93 dh[:, 4] = 0.6 * numpy.cos(dens)**2 94 dh[:, 5] = 0.3 * numpy.sin(dens) 95 return dh 96 97 98@pytest.fixture 99def refen_lda(): 100 return numpy.array([ 101 -0.03962059, -0.41418026, -0.91826807, -1.50300217, -2.14968647, -2.84774973, -3.59023064, -4.37213299, 102 -5.18966423, -6.03982917 103 ]) 104 105 106@pytest.fixture 107def refpot_lda(): 108 return numpy.array([ 109 -0.51789018, -0.90610372, -1.09730105, -1.23582171, -1.34751238, -1.44246858, -1.52581304, -1.60054497, 110 -1.66858921, -1.73126329 111 ]) 112 113 114def test_lda_energy_restr(lda_fun, dens, refen_lda): 115 out = lda_fun.eval_energy_n(dens) 116 assert_allclose(out, refen_lda) 117 118 119def test_lda_potential_restr(lda_fun, dens, refen_lda, refpot_lda): 120 out = lda_fun.eval_potential_n(dens) 121 assert_allclose(out[:, 0], refen_lda) 122 assert_allclose(out[:, 1], refpot_lda) 123 124 125def test_lda_energy_unrestr_closedshell(lda_fun, dens, refen_lda): 126 rho = numpy.zeros((dens.shape[0], 2)) 127 rho[:, 0] = 0.5 * dens[:] 128 rho[:, 1] = 0.5 * dens[:] 129 130 out = lda_fun.eval_energy_ab(rho) 131 assert_allclose(out, refen_lda) 132 133 134def test_lda_potential_unrestr_closedshell(lda_fun, dens, refen_lda, refpot_lda): 135 rho = numpy.zeros((dens.shape[0], 2)) 136 rho[:, 0] = 0.5 * dens[:] 137 rho[:, 1] = 0.5 * dens[:] 138 139 out = lda_fun.eval_potential_ab(rho) 140 assert_allclose(out[:, 0], refen_lda) 141 assert_allclose(out[:, 1], refpot_lda) 142 assert_allclose(out[:, 2], refpot_lda) 143 144 145def test_lda_energy_unrestr_openshell(lda_fun, dens, spindens): 146 rho = numpy.zeros((dens.shape[0], 2)) 147 rho[:, 0] = 0.5 * dens[:] + spindens[:] 148 rho[:, 1] = 0.5 * dens[:] - spindens[:] 149 150 out = lda_fun.eval_energy_ab(rho) 151 152 refen_lda_openshell = numpy.array([ 153 -0.03985427, -0.41666061, -0.92248897, -1.50718908, -2.15231233, -2.84856692, -3.59023539, -4.3726223, 154 -5.19120401, -6.04193876 155 ]) 156 157 assert_allclose(out, refen_lda_openshell) 158 159 160def test_lda_potential_unrestr_openshell(lda_fun, dens, spindens): 161 rho = numpy.zeros((dens.shape[0], 2)) 162 rho[:, 0] = 0.5 * dens[:] + spindens[:] 163 rho[:, 1] = 0.5 * dens[:] - spindens[:] 164 165 out = lda_fun.eval_potential_ab(rho) 166 167 refpot_lda_alpha = numpy.array([ 168 -0.53994997, -0.94756428, -1.14234107, -1.27610516, -1.37715554, -1.4581231, -1.52695501, -1.58940021, 169 -1.64952872, -1.70973604 170 ]) 171 refpot_lda_beta = numpy.array([ 172 -0.49296826, -0.85942489, -1.04738939, -1.1922106, -1.31627816, -1.42641378, -1.52466911, -1.61151607, 173 -1.68716914, -1.75220302 174 ]) 175 176 assert_allclose(out[:, 1], refpot_lda_alpha) 177 assert_allclose(out[:, 2], refpot_lda_beta) 178 179 outen = lda_fun.eval_energy_ab(rho) 180 assert_allclose(out[:, 0], outen) 181 182 183@pytest.fixture 184def refen_pbe(): 185 return numpy.array([ 186 -0.04051725, -0.41397904, -0.91782841, -1.50231227, -2.14873903, -2.84654021, -3.58875595, -4.37039097, 187 -5.18765335, -6.03754832 188 ]) 189 190 191@pytest.fixture 192def refpot_pbe(): 193 return numpy.array([ 194 -0.48524858, -0.90474127, -1.09660442, -1.23524903, -1.34697602, -1.44193844, -1.52527933, -1.60000832, 195 -1.66805094, -1.73072351 196 ]) 197 198 199def test_pbe_energy_restr(pbe_fun, dens, densgrad, refen_pbe): 200 out = pbe_fun.eval_energy_n(dens, densgrad) 201 assert_allclose(out, refen_pbe) 202 203 204def test_pbe_potential_restr(pbe_fun, dens, densgrad, denshess, refen_pbe, refpot_pbe): 205 out = pbe_fun.eval_potential_n(dens, densgrad, denshess) 206 assert_allclose(out[:, 0], refen_pbe) 207 assert_allclose(out[:, 1], refpot_pbe) 208 209 210def test_pbe_energy_unrestr_closedshell(pbe_fun, dens, densgrad, refen_pbe): 211 rho = numpy.zeros((dens.shape[0], 2)) 212 rho[:, 0] = 0.5 * dens[:] 213 rho[:, 1] = 0.5 * dens[:] 214 215 rhograd = numpy.zeros((dens.shape[0], 3, 2)) 216 rhograd[:, 0:3, 0] = 0.5 * densgrad[:, 0:3] 217 rhograd[:, 0:3, 1] = 0.5 * densgrad[:, 0:3] 218 219 out = pbe_fun.eval_energy_ab(rho, rhograd) 220 assert_allclose(out, refen_pbe) 221 222 223def test_pbe_potential_unrestr_closedshell(pbe_fun, dens, densgrad, denshess, refen_pbe, refpot_pbe): 224 rho = numpy.zeros((dens.shape[0], 2)) 225 rho[:, 0] = 0.5 * dens[:] 226 rho[:, 1] = 0.5 * dens[:] 227 228 rhograd = numpy.zeros((dens.shape[0], 3, 2)) 229 rhograd[:, 0:3, 0] = 0.5 * densgrad[:, 0:3] 230 rhograd[:, 0:3, 1] = 0.5 * densgrad[:, 0:3] 231 232 rhohess = numpy.zeros((dens.shape[0], 6, 2)) 233 rhohess[:, 0:6, 0] = 0.5 * denshess[:, 0:6] 234 rhohess[:, 0:6, 1] = 0.5 * denshess[:, 0:6] 235 236 out = pbe_fun.eval_potential_ab(rho, rhograd, rhohess) 237 assert_allclose(out[:, 0], refen_pbe) 238 assert_allclose(out[:, 1], refpot_pbe) 239 assert_allclose(out[:, 2], refpot_pbe) 240 241 242def test_pbe_energy_unrestr_openshell(pbe_fun, dens, spindens, densgrad, refen_pbe): 243 rho = numpy.zeros((dens.shape[0], 2)) 244 rho[:, 0] = dens[:] + spindens[:] 245 rho[:, 1] = dens[:] - spindens[:] 246 247 rhograd = numpy.zeros((dens.shape[0], 3, 2)) 248 rhograd[:, 0:3, 0] = 0.4 * densgrad[:, 0:3] 249 rhograd[:, 0:3, 1] = 0.6 * densgrad[:, 0:3] 250 251 out = pbe_fun.eval_energy_ab(rho, rhograd) 252 253 refen_pbe_openshell = numpy.array([ 254 -0.09848243, -1.03067727, -2.28717791, -3.74476386, -5.35679716, -7.09762348, -8.95059771, -10.90338775, 255 -12.94601448, -15.07029572 256 ]) 257 258 assert_allclose(out, refen_pbe_openshell) 259 260 261def test_pbe_potential_unrestr_closedshell(pbe_fun, dens, spindens, densgrad, denshess, refen_pbe, refpot_pbe): 262 rho = numpy.zeros((dens.shape[0], 2)) 263 rho[:, 0] = dens[:] + spindens[:] 264 rho[:, 1] = dens[:] - spindens[:] 265 266 rhograd = numpy.zeros((dens.shape[0], 3, 2)) 267 rhograd[:, 0:3, 0] = 0.4 * densgrad[:, 0:3] 268 rhograd[:, 0:3, 1] = 0.6 * densgrad[:, 0:3] 269 270 rhohess = numpy.zeros((dens.shape[0], 6, 2)) 271 rhohess[:, 0:6, 0] = 0.7 * denshess[:, 0:6] 272 rhohess[:, 0:6, 1] = 0.2 * denshess[:, 0:6] 273 274 out = pbe_fun.eval_potential_ab(rho, rhograd, rhohess) 275 276 refpot_pbe_alpha = numpy.array([ 277 -0.62768401, -1.15201996, -1.39499424, -1.56595725, -1.69939803, -1.80975816, -1.90530719, -1.99151265, 278 -2.07206904, -2.1492557 279 ]) 280 refpot_pbe_beta = numpy.array([ 281 -0.64069084, -1.10200908, -1.33808862, -1.51488446, -1.66214897, -1.79045119, -1.90425337, -2.00565175, 282 -2.09589401, -2.17609771 283 ]) 284 assert_allclose(out[:, 1], refpot_pbe_alpha) 285 assert_allclose(out[:, 2], refpot_pbe_beta) 286 287 outen = pbe_fun.eval_energy_ab(rho, rhograd) 288 assert_allclose(out[:, 0], outen) 289 290 291def test_raises_eval_energy_n(lda_fun, pbe_fun): 292 # correct shape for all arguments (LDA) 293 lda_fun.eval_energy_n(numpy.zeros((10, ))) 294 295 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 296 lda_fun.eval_energy_n(numpy.zeros((10, 2))) 297 298 # correct shape for all arguments (PBE) 299 pbe_fun.eval_energy_n(numpy.zeros((10, )), numpy.zeros((10, 3))) 300 301 with pytest.raises(xcfun.XCFunException, match="Density gradient required"): 302 pbe_fun.eval_energy_n(numpy.zeros((10, ))) 303 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 304 pbe_fun.eval_energy_n(numpy.zeros((10, )), numpy.zeros((10, ))) 305 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 306 pbe_fun.eval_energy_n(numpy.zeros((10, )), numpy.zeros((10, 4))) 307 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 308 pbe_fun.eval_energy_n(numpy.zeros((10, )), numpy.zeros((20, 3))) 309 310 311def test_raises_eval_potential_n(lda_fun, pbe_fun): 312 # correct shape for all arguments (LDA) 313 lda_fun.eval_potential_n(numpy.zeros((10, ))) 314 315 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 316 lda_fun.eval_potential_n(numpy.zeros((10, 2))) 317 318 # correct shape for all arguments (PBE) 319 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 3)), numpy.zeros((10, 6))) 320 321 with pytest.raises(xcfun.XCFunException, match="Density gradient and Hessian required"): 322 pbe_fun.eval_potential_n(numpy.zeros((10, ))) 323 with pytest.raises(xcfun.XCFunException, match="Density gradient and Hessian required"): 324 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 3))) 325 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 326 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, )), numpy.zeros((10, 6))) 327 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 328 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 4)), numpy.zeros((10, 6))) 329 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 330 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((20, 3)), numpy.zeros((10, 6))) 331 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 332 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 3)), numpy.zeros((10, ))) 333 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 334 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 3)), numpy.zeros((10, 3))) 335 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 336 pbe_fun.eval_potential_n(numpy.zeros((10, )), numpy.zeros((10, 3)), numpy.zeros((20, 6))) 337 338 339def test_raises_eval_energy_ab(lda_fun, pbe_fun): 340 # correct shape for all arguments (LDA) 341 lda_fun.eval_energy_ab(numpy.zeros((10, 2))) 342 343 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 344 lda_fun.eval_energy_ab(numpy.zeros((10, ))) 345 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 346 lda_fun.eval_energy_ab(numpy.zeros((10, 4))) 347 348 # correct shape for all arguments (PBE) 349 pbe_fun.eval_energy_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2))) 350 351 with pytest.raises(xcfun.XCFunException, match="Density gradient required"): 352 pbe_fun.eval_energy_ab(numpy.zeros((10, 2))) 353 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 354 pbe_fun.eval_energy_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3))) 355 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 356 pbe_fun.eval_energy_ab(numpy.zeros((10, 2)), numpy.zeros((10, 4, 2))) 357 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 358 pbe_fun.eval_energy_ab(numpy.zeros((10, 2)), numpy.zeros((20, 3, 2))) 359 360 361def test_raises_eval_potential_ab(lda_fun, pbe_fun): 362 # correct shape for all arguments (LDA) 363 lda_fun.eval_potential_ab(numpy.zeros((10, 2))) 364 365 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 366 lda_fun.eval_potential_ab(numpy.zeros((10, ))) 367 with pytest.raises(xcfun.XCFunException, match="Wrong shape of density argument"): 368 lda_fun.eval_potential_ab(numpy.zeros((10, 4))) 369 370 # correct shape for all arguments (PBE) 371 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2)), numpy.zeros((10, 6, 2))) 372 373 with pytest.raises(xcfun.XCFunException, match="Density gradient and Hessian required"): 374 pbe_fun.eval_potential_ab(numpy.zeros((10, 2))) 375 with pytest.raises(xcfun.XCFunException, match="Density gradient and Hessian required"): 376 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2))) 377 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 378 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3)), numpy.zeros((10, 6, 2))) 379 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 380 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 4, 2)), numpy.zeros((10, 6, 2))) 381 with pytest.raises(xcfun.XCFunException, match="Wrong shape of densgrad argument"): 382 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((20, 3, 2)), numpy.zeros((10, 6, 2))) 383 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 384 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2)), numpy.zeros((10, 6))) 385 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 386 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2)), numpy.zeros((10, 3, 2))) 387 with pytest.raises(xcfun.XCFunException, match="Wrong shape of denshess argument"): 388 pbe_fun.eval_potential_ab(numpy.zeros((10, 2)), numpy.zeros((10, 3, 2)), numpy.zeros((20, 6, 2))) 389