Typo and bugfix

typo in GGA derivatives

Higher-level julia interface and test

Ignore vim temporaries

Make testdata local

patch version 0.1.2

.gitignore

@@ -4,3 +4,8 @@
 deps/usr
 deps/deps.jl
 deps/build.log
+
+# vim temporaries
+*~
+.*.swp
+.*.swo

Project.toml

@@ -1,7 +1,7 @@
 name = "Libxc"
 uuid = "66e17ffc-8502-11e9-23b5-c9248d0eb96d"
 authors = ["Jason Eu <[email protected]>"]
-version = "0.1.0"
+version = "0.1.2"
 
 [compat]
 julia = "1"

README.md

@@ -8,7 +8,7 @@ Binary built by using [BinaryBuilder](https://github.com/JuliaPackaging/BinaryBu
 Now it is registered in [JuliaRegisties](https://github.com/JuliaRegistries/General), thus can be installed by running:
 
 ```sh
-(v1.1) pkg> add Libxcs
+(v1.1) pkg> add Libxc
 ```
 
 When wrapping new libxc version, first run `update_functionals_info.jl` to update functionals' info

src/Functional.jl

@@ -0,0 +1,90 @@
+"""Return the list of available libxc functionals as strings"""
+function available_functionals()
+    n_xc = ccall((:xc_number_of_functionals, libxc), Cint, ())
+    max_string_length = ccall((:xc_maximum_name_length, libxc), Cint, ())
+
+    funcnames = Vector{String}(undef, n_xc)
+    for i in 1:n_xc
+        funcnames[i] = "\0"^(max_string_length + 2)
+    end
+    ccall((:xc_available_functional_names, libxc), Cvoid, (Ptr{Ptr{UInt8}}, ), funcnames)
+    [string(split(funcnames[i], "\0")[1]) for i in 1:n_xc]
+end
+
+@enum FunctionalKind begin
+    functional_exchange             = 0
+    functional_correlation          = 1
+    functional_exchange_correlation = 2
+    functional_kinetic              = 3
+end
+
+@enum FunctionalFamily begin
+    family_unknown  = -1
+    family_lda      = 1
+    family_gga      = 2
+    family_mggai    = 4
+    family_lca      = 8
+    family_oep      = 16
+    family_hyb_gga  = 32
+    family_hyb_mgga = 64
+end
+
+mutable struct Functional
+    number::Int
+    name::String
+    kind::FunctionalKind
+    family::FunctionalFamily
+    n_spin::Int
+
+    # Pointer holding the LibXC representation of this functional
+    pointer::Ptr{XCFuncType}
+end
+
+
+"""
+    Functional(identifier; n_spin::Integer = 1)
+
+Construct a Functional from a libxc `identifier` and the number
+of spins `n_spin` to consider. `
+"""
+function Functional(identifier; n_spin::Integer = 1)
+    if n_spin != 1 && n_spin != 2
+        error("n_spin needs to be 1 or 2")
+    end
+
+    number = ccall((:xc_functional_get_number, libxc), Cint, (Cstring, ),
+                   string(identifier))
+    if number == -1
+        error("Functional $identifier is not known.")
+    end
+
+    function pointer_cleanup(ptr::Ptr{XCFuncType})
+        if ptr != C_NULL
+            xc_func_end(ptr)
+            xc_func_free(ptr)
+        end
+    end
+
+    pointer = xc_func_alloc()
+    ret = xc_func_init(pointer, number, n_spin)
+    if ret != 0
+        error("Something went wrong initialising the functional")
+    end
+
+    try
+        funcinfo = xc_func_get_info(pointer)
+        kind = xc_func_info_get_kind(funcinfo)
+        family = xc_func_info_get_family(funcinfo)
+        flags = xc_func_info_get_flags(funcinfo)
+        # TODO Extract references ....
+
+        # Make functional and attach finalizer for cleaning up the pointer
+        func = Functional(number, string(identifier), FunctionalKind(kind),
+                          FunctionalFamily(family), n_spin, pointer)
+        finalizer(cls -> pointer_cleanup(cls.pointer), func)
+        return func
+    catch
+        pointer_cleanup(pointer)
+        rethrow()
+    end
+end

src/Libxc.jl

@@ -14,11 +14,15 @@ function __init__()
     check_deps()
 end
 
-export XCFuncType
-
 include("version.jl")
 include("xc.jl")
 include("functionals_info.jl")
 
+export Functional
+include("Functional.jl")
+
+export evaluate_lda!
+export evaluate_gga!
+include("evaluate.jl")
 
 end # module

src/evaluate.jl

@@ -0,0 +1,70 @@
+function evaluate_lda!(func::Functional, ρ::Array{Float64}; kwargs...)
+    @assert func.family == family_lda
+    n_p = Int(length(ρ) / func.n_spin)
+
+    sizes = Dict(:E => n_p, :Vρ => 2n_p, :V2ρ2 => 3n_p, :V3ρ3 => 4n_p)
+    if func.n_spin == 1
+        for key in keys(sizes)
+            sizes[key] = n_p
+        end
+    end
+
+    dargs = Dict(kwargs)
+    for key in keys(dargs)
+        if !haskey(sizes, key)
+            throw(ArgumentError("Unknown keyword argument: $(string(key))"))
+        end
+        if length(dargs[key]) != sizes[key]
+            throw(DimensionMismatch("Length of keyword argument $(string(key)) == " *
+                                    "$(length(dargs[key])) does not agree with expected " *
+                                    "value $(sizes[key])"))
+        end
+    end
+
+    ccall((:xc_lda, libxc), Cvoid, (Ptr{XCFuncType}, Cint, Ptr{Float64}, Ptr{Float64},
+                                    Ptr{Float64}, Ptr{Float64}, Ptr{Float64}),
+          func.pointer, n_p, ρ, get(dargs, :E, C_NULL),
+          get(dargs, :Vρ, C_NULL), get(dargs, :V2ρ2, C_NULL), get(dargs, :V3ρ3, C_NULL)
+    )
+end
+
+
+function evaluate_gga!(func::Functional, ρ::Array{Float64}, σ::Array{Float64}; kwargs...)
+    @assert func.family == family_gga
+    n_p = Int(length(ρ) / func.n_spin)
+
+    n_σ = 3n_p
+    sizes = Dict(:E => n_p, :Vρ => 2n_p, :Vσ => 3n_p,
+                 :V2ρ2 => 3n_p, :V2ρσ => 6n_p, :V2σ2 => 6n_p,
+                 :V3ρ3 => 4n_p, :V3ρ2σ => 9n_p, :V3ρσ2 => 12n_p, :V3σ3 => 10n_p)
+    if func.n_spin == 1
+        n_σ = n_p
+        for key in keys(sizes)
+            sizes[key] = n_p
+        end
+    end
+
+    dargs = Dict(kwargs)
+    for key in keys(dargs)
+        if !haskey(sizes, key)
+            throw(ArgumentError("Unknown keyword argument: $(string(key))"))
+        end
+        if length(dargs[key]) != sizes[key]
+            throw(DimensionMismatch("Length of keyword argument $(string(key)) == " *
+                                    "$(length(dargs[key])) does not agree with expected " *
+                                    "value $(sizes[key])"))
+        end
+    end
+
+    ccall((:xc_gga, libxc), Cvoid, (Ptr{XCFuncType}, Cint, Ptr{Float64}, Ptr{Float64},
+                                    Ptr{Float64}, Ptr{Float64}, Ptr{Float64},
+                                    Ptr{Float64}, Ptr{Float64}, Ptr{Float64},
+                                    Ptr{Float64}, Ptr{Float64}, Ptr{Float64},
+                                    Ptr{Float64}),
+          func.pointer, n_p, ρ, σ,
+          get(dargs, :E, C_NULL),    get(dargs, :Vρ, C_NULL),    get(dargs, :Vσ, C_NULL),
+          get(dargs, :V2ρ2, C_NULL), get(dargs, :V2ρσ, C_NULL),  get(dargs, :V2σ2, C_NULL),
+          get(dargs, :V3ρ3, C_NULL), get(dargs, :V3ρ2σ, C_NULL), get(dargs, :V3ρσ2, C_NULL),
+          get(dargs, :V3σ3, C_NULL)
+    )
+end

src/version.jl

@@ -1,11 +1,8 @@
 function xc_version()
-    a = Ref{Cint}(0)
-    b = Ref{Cint}(0)
-    c = Ref{Cint}(0)
-    ccall( (:xc_version, libxc), Cvoid, (Ref{Cint}, Ref{Cint}, Ref{Cint}), a, b, c)
-    return Base.convert(Int64, a[]),
-           Base.convert(Int64, b[]),
-           Base.convert(Int64, c[])
+    varray = zeros(Cint, 3)
+    ccall((:xc_version, libxc), Cvoid, (Ptr{Cint}, Ptr{Cint}, Ptr{Cint}),
+          pointer(varray, 1), pointer(varray, 2), pointer(varray, 3))
+    VersionNumber(varray[1], varray[2], varray[3])
 end
 
 function xc_version_string()

src/xc.jl

@@ -109,13 +109,13 @@ mutable struct XCFuncType
 end
 
 xc_func_alloc() = ccall( (:xc_func_alloc, libxc), Ptr{XCFuncType}, () )
-xc_func_init(p::Ptr{XCFuncType}, functional::Int, nspin::Int) =
+xc_func_init(p::Ptr{XCFuncType}, functional::Integer, nspin::Integer) =
     ccall( (:xc_func_init, libxc), Cint, (Ptr{XCFuncType}, Cint, Cint), p, functional, nspin)
 xc_func_end(p::Ptr{XCFuncType}) = ccall( (:xc_func_end, libxc), Cvoid, (Ptr{XCFuncType},), p)
 xc_func_free(p::Ptr{XCFuncType}) = ccall( (:xc_func_free, libxc), Cvoid, (Ptr{XCFuncType},), p)
 xc_func_get_info(p::Ptr{XCFuncType}) = ccall( (:xc_func_get_info, libxc), Ptr{XCFuncInfoType}, (Ref{XCFuncType},), p )
 xc_func_set_dens_threshold(p::Ptr{XCFuncType}, dens_threshold::Float64) = ccall( (:xc_func_set_dens_threshold, libxc), Cvoid, (Ptr{XCFuncType}, Float64), p, dens_threshold)
-xc_func_set_ext_params(p::Ptr{XCFuncType}, ext_params::Ptr{Float64}) = ccall( (:xc_func_set_dens_threshold, libxc), Cvoid, (Ptr{XCFuncType}, Ptr{Float64}), p, ext_params)
+xc_func_set_ext_params(p::Ptr{XCFuncType}, ext_params::Ptr{Float64}) = ccall( (:xc_func_set_ext_params, libxc), Cvoid, (Ptr{XCFuncType}, Ptr{Float64}), p, ext_params)
 
 ######################################################################
 ##  LDA-xc
@@ -142,7 +142,7 @@ xc_gga!(p::Ptr{XCFuncType}, np::Int,
         rho::Array{Float64, 1}, sigma::Array{Float64, 1},
         zk::Array{Float64, 1}, vrho::Array{Float64, 1}, vsigma::Array{Float64, 1},
         v2rho2::Array{Float64, 1}, v2rhosigma::Array{Float64, 1}, v2sigma2::Array{Float64, 1},
-        v3rho3::Array{Float64, 1}, v3rhosigma::Array{Float64, 1}, v3rhosigma2::Array{Float64, 1}, v3sigma3::Array{Float64, 1}) =
+        v3rho3::Array{Float64, 1}, v3rho2sigma::Array{Float64, 1}, v3rhosigma2::Array{Float64, 1}, v3sigma3::Array{Float64, 1}) =
     ccall( (:xc_gga, libxc), Cvoid, (Ptr{XCFuncType}, Cint,
                                      Ptr{Float64}, Ptr{Float64},
                                      Ptr{Float64}, Ptr{Float64}, Ptr{Float64},
@@ -152,7 +152,7 @@ xc_gga!(p::Ptr{XCFuncType}, np::Int,
                                      rho, sigma,
                                      zk, vrho, vsigma,
                                      v2rho2, v2rhosigma, v2sigma2,
-                                     v3rho3, v3rhosigma, v3rhosigma2, v3sigma3)
+                                     v3rho3, v3rho2sigma, v3rhosigma2, v3sigma3)
 
 xc_gga_exc!(p::Ptr{XCFuncType}, np::Int, rho::Array{Float64, 1}, sigma::Array{Float64, 1},
             zk::Array{Float64, 1}) =

test/runtests.jl

@@ -2,7 +2,7 @@ using Test
 using Libxc
 
 @testset "Version" begin
-    @test isa(Libxc.xc_version(), Tuple{Int64, Int64, Int64})
+    @test isa(Libxc.xc_version(), VersionNumber)
     @test isa(Libxc.xc_version_string(), String)
 end
 
@@ -27,10 +27,11 @@ end
     Libxc.xc_func_free(reff)
 end
 
-rho = [0.1, 0.2, 0.3, 0.4, 0.5]
-sigma = [0.2, 0.3, 0.4, 0.5, 0.6]
-result = zeros(Float64, 5)
 @testset "XCFuncType xc calc" begin
+    rho = [0.1, 0.2, 0.3, 0.4, 0.5]
+    sigma = [0.2, 0.3, 0.4, 0.5, 0.6]
+    result = zeros(Float64, 5)
+
     ptr = Libxc.xc_func_alloc()
     Libxc.xc_func_init(ptr, Libxc.LDA_X, 1)
     Libxc.xc_lda_exc!(ptr, 5, rho, result)
@@ -45,3 +46,17 @@ result = zeros(Float64, 5)
     Libxc.xc_func_end(ptr)
     Libxc.xc_func_free(ptr)
 end
+
+@testset "High-level interface" begin
+    rho = [0.1, 0.2, 0.3, 0.4, 0.5]
+    sigma = [0.2, 0.3, 0.4, 0.5, 0.6]
+    result = zeros(Float64, 5)
+
+    func = Libxc.Functional(:lda_x)
+    Libxc.evaluate_lda!(func, rho, E=result)
+    @test result ≈ [-0.342809, -0.431912, -0.494416, -0.544175, -0.586194] atol=1e-5
+
+    func = Libxc.Functional(:gga_x_pbe)
+    Libxc.evaluate_gga!(func, rho, sigma, E=result)
+    @test result ≈ [-0.452598, -0.478878, -0.520674, -0.561428, -0.598661] atol=1e-5
+end