-
Notifications
You must be signed in to change notification settings - Fork 83
/
xc_fallback.jl
74 lines (63 loc) · 2.31 KB
/
xc_fallback.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
using Test
using Libxc
using DFTK: xc_fallback!, energy_per_particle
@testset "Fallback LDA" begin
for func_name in (:lda_x, :lda_c_vwn, :lda_c_pw)
func = Libxc.Functional(func_name)
# Create reference
ρ = abs.(randn(100))
Eref = similar(ρ)
Vref = similar(ρ)
V2ref = similar(ρ)
Libxc.evaluate!(func, rho=ρ, zk=Eref, vrho=Vref, v2rho2=V2ref)
@testset "$func_name" begin
# Compute in DFTK in elevated precision
ρbig = Array{BigFloat}(ρ)
E = similar(ρbig)
V = similar(ρbig)
V2 = similar(ρbig)
xc_fallback!(func, Val(:lda), ρbig; zk=E, vrho=V, v2rho2=V2)
@test E ≈ Eref atol=5e-15
@test V ≈ Vref atol=5e-15
@test V2 ≈ V2ref atol=5e-15
# TODO Check floating-point type consistency:
# ε = energy_per_particle(Val(func_name), rand(Float32))
# @test Float32 == eltype(ε)
end
end
end
@testset "Fallback GGA" begin
gga_fallback = (
:gga_x_pbe, :gga_c_pbe, :gga_x_pbe_r,
:gga_x_pbe_sol, :gga_c_pbe_sol,
:gga_x_pbefe, :gga_c_pbefe,
:gga_x_xpbe, :gga_c_xpbe,
:gga_x_pbe_mol, :gga_c_pbe_mol,
:gga_x_apbe, :gga_c_apbe,
)
for func_name in gga_fallback
func = Libxc.Functional(func_name)
# Create reference
ρ = abs.(randn(100))
σ = abs.(randn(100))
Eref = similar(ρ)
Vρref = similar(ρ)
Vσref = similar(ρ)
Libxc.evaluate!(func, rho=ρ, sigma=σ, zk=Eref, vrho=Vρref, vsigma=Vσref)
@testset "$func_name" begin
# Compute in DFTK in elevated precision
ρbig = Array{BigFloat}(ρ)
σbig = Array{BigFloat}(σ)
E = similar(ρbig)
Vρ = similar(ρbig)
Vσ = similar(ρbig)
xc_fallback!(func, Val(:gga), ρbig; sigma=σbig, zk=E, vrho=Vρ, vsigma=Vσ)
@test E ≈ Eref atol=5e-15
@test Vρ ≈ Vρref atol=5e-15
@test Vσ ≈ Vσref atol=5e-15
# TODO Check floating-point type consistency:
# ε = energy_per_particle(Val(func_name), rand(Float32))
# @test Float32 == eltype(ε)
end
end
end