Refactor ABINIT testcases

CedricTravelletti · Feb 7, 2020 · 36a1e94 · 36a1e94
1 parent 5c00d55
commit 36a1e94
Show file tree

Hide file tree

Showing 16 changed files with 198 additions and 278 deletions.
diff --git a/test/bzmesh.jl b/test/bzmesh.jl
@@ -51,8 +51,8 @@ end
     test_reduction(silicon, [ 2,  3,  4])
     test_reduction(silicon, [ 9, 11, 13])
 
-    test_reduction(manganese, [ 2,  3,  2])
-    test_reduction(manganese, [ 3,  3,  3])
-    test_reduction(manganese, [ 2,  3,  4])
-    test_reduction(manganese, [ 9, 11, 13])
+    test_reduction(magnesium, [ 2,  3,  2])
+    test_reduction(magnesium, [ 3,  3,  3])
+    test_reduction(magnesium, [ 2,  3,  4])
+    test_reduction(magnesium, [ 9, 11, 13])
 end
diff --git a/test/compute_density.jl b/test/compute_density.jl
@@ -84,7 +84,7 @@ include("testcases.jl")
 
     test_full_vs_irreducible(silicon, [3, 3, 3], Ecut=5, tol=1e-6)
     test_full_vs_irreducible(silicon, [2, 3, 4], Ecut=5, tol=1e-6)
-    test_full_vs_irreducible(manganese, [2, 3, 4], Ecut=5, tol=1e-6)
+    test_full_vs_irreducible(magnesium, [2, 3, 4], Ecut=5, tol=1e-6)
     #
     # That's pretty expensive:
     # test_full_vs_irreducible([4, 4, 4], Ecut=5, tol=1e-6)

diff --git a/test/occupation.jl b/test/occupation.jl
@@ -66,7 +66,7 @@ end
 end
 
 @testset "Smearing for metals" begin
-    testcase = manganese
+    testcase = magnesium
     Ecut = 5
     fft_size = [15, 15, 15]
     kgrid_size = [2, 3, 4]

diff --git a/test/testcases.jl b/test/testcases.jl
@@ -44,7 +44,7 @@ silicon = (
        ]),
 )
 
-manganese = (
+magnesium = (
     lattice = [-3.0179389205999998 -3.0179389205999998 0.0000000000000000;
                -5.2272235447000002 5.2272235447000002 0.0000000000000000;
                0.0000000000000000 0.0000000000000000 -9.7736219469000005],

diff --git a/test/testcases_ABINIT/graphite_E20_k8_LDA_Methfessel_Paxton/generate.jl b/test/testcases_ABINIT/graphite_E20_k8_LDA_Methfessel_Paxton/generate.jl
@@ -1,18 +1,22 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using PyCall
+using DFTK
+import DFTK.units
 
-structure = build_graphite_structure()
-pspfile = joinpath(abidata.__path__[end], "hgh_pseudos/6c.4.hgh")
+abinitpseudos = [joinpath(pyimport("abipy.data").__path__[end], "hgh_pseudos/6c.4.hgh")]
+pspmap = Dict(6 => "hgh/lda/c-q4", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[8, 8, 8], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=20,        # Hartree
-    nband=12,       # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc=1,          # Teter LDA
-    occopt=6,       # Methfessel and Paxton, Hermite polynomial degree 2
-    tsmear=0.01,    # Hartree
-)
-infile.extra = Dict("pspmap" => Dict(6 => "hgh/lda/c-q4", ), )
+# Note: This is not exactly the minimum-energy structure
+a = 1.228units.Ǎ
+b = 2.12695839units.Ǎ
+c = 7units.Ǎ
+lattice = [[a a 0]; [-b b 0]; [0 0 c]]
+C = Element(6, load_psp(pspmap[6]))
+atoms = [C => [[0, 0, 1/4], [0, 0, 3/4], [1/3, 2/3, 1/4], [2/3, 1/3, 3/4]]]
 
-run_ABINIT_scf(infile, @__DIR__)
+model = model_dft(lattice, :lda_xc_teter93, atoms, temperature=0.01,
+                  smearing=Smearing.MethfesselPaxton2())
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=20, kgrid=[8, 8, 8], n_bands=12, tol=1e-10,
+                    iscf=3) # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/magnesium_E15_k3_GGA_Methfessel_Paxton/generate.jl b/test/testcases_ABINIT/magnesium_E15_k3_GGA_Methfessel_Paxton/generate.jl
@@ -1,18 +1,15 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using PyCall
+using DFTK
 
-structure = build_magnesium_structure()
-pspfile = joinpath(abidata.__path__[end], "hgh_pseudos/12mg.2.hgh")
+abinitpseudos = [joinpath(pyimport("abipy.data").__path__[end], "hgh_pseudos/12mg.2.hgh")]
+pspmap = Dict(12 => "hgh/lda/mg-q2", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[3, 3, 3], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=15,        # Hartree
-    nband=10,       # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc="-101130",  # PBE X and PBE C
-    occopt=6,       # Methfessel and Paxton, Hermite polynomial degree 2
-    tsmear=0.01,    # Hartree
-)
-infile.extra = Dict("pspmap" => Dict(12 => "hgh/lda/mg-q2", ), )
+atoms = [Element(12, load_psp(pspmap[12])) => magnesium.positions]
+model = model_dft(magnesium.lattice, [:gga_x_pbe, :gga_c_pbe], atoms,
+                  temperature=0.01, smearing=Smearing.MethfesselPaxton2())
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=15, kgrid=[3, 3, 3], n_bands=10, tol=1e-10,
+                    iscf=3)  # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/magnesium_E15_k3_LDA_Fermi_Dirac/generate.jl b/test/testcases_ABINIT/magnesium_E15_k3_LDA_Fermi_Dirac/generate.jl
@@ -1,19 +1,15 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using PyCall
+using DFTK
 
-structure = build_magnesium_structure()
-pspfile = joinpath(abidata.__path__[end], "hgh_pseudos/12mg.2.hgh")
+abinitpseudos = [joinpath(pyimport("abipy.data").__path__[end], "hgh_pseudos/12mg.2.hgh")]
+pspmap = Dict(12 => "hgh/lda/mg-q2", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[3, 3, 3], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=15,        # Hartree
-    nband=10,       # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc="-020",     # Teter1993 LDA reparametrisation
-    # occopt=6,     # Methfessel and Paxton, Hermite polynomial degree 2
-    occopt=3,       # Fermi-Dirac
-    tsmear=0.01,    # Hartree
-)
-infile.extra = Dict("pspmap" => Dict(12 => "hgh/lda/mg-q2", ), )
+atoms = [Element(12, load_psp(pspmap[12])) => magnesium.positions]
+model = model_dft(magnesium.lattice, :lda_xc_teter93, atoms,
+                  temperature=0.01, smearing=Smearing.FermiDirac())
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=15, kgrid=[3, 3, 3], n_bands=10, tol=1e-10,
+                    iscf=3)  # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/magnesium_E25_k5_GGA_Methfessel_Paxton/generate.jl b/test/testcases_ABINIT/magnesium_E25_k5_GGA_Methfessel_Paxton/generate.jl
@@ -1,18 +1,15 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using PyCall
+using DFTK
 
-structure = build_magnesium_structure()
-pspfile = joinpath(abidata.__path__[end], "hgh_pseudos/12mg.2.hgh")
+abinitpseudos = [joinpath(pyimport("abipy.data").__path__[end], "hgh_pseudos/12mg.2.hgh")]
+pspmap = Dict(12 => "hgh/lda/mg-q2", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[5, 5, 5], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=25,        # Hartree
-    nband=10,       # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc="-101130",  # PBE X and PBE C
-    occopt=6,       # Methfessel and Paxton, Hermite polynomial degree 2
-    tsmear=0.01,    # Hartree
-)
-infile.extra = Dict("pspmap" => Dict(12 => "hgh/lda/mg-q2", ), )
+atoms = [Element(12, load_psp(pspmap[12])) => magnesium.positions]
+model = model_dft(magnesium.lattice, [:gga_x_pbe, :gga_c_pbe], atoms,
+                  temperature=0.01, smearing=Smearing.MethfesselPaxton2())
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=25, kgrid=[5, 5, 5], n_bands=10, tol=1e-10,
+                    iscf=3)  # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/runtests.jl b/test/testcases_ABINIT/runtests.jl
@@ -1,24 +1,76 @@
-include("testcases.jl")
-
-test_folder(Float64, "silicon_E10_k3_LDA", test_tol=5e-8)
-test_folder(Float64, "silicon_E15_k4_LDA", n_ignored=1, test_tol=5e-8)
-test_folder(Float64, "silicon_E25_k3_LDA", n_ignored=3, test_tol=1e-8)
-test_folder(Float64, "silicon_E15_k4_GGA")
-test_folder(Float64, "silicon_E25_k3_GGA")
-test_folder(Float64, "silicon_E25_k4_GGA", n_ignored=1, test_tol=1e-7)
-test_folder(Float64, "magnesium_E15_k3_LDA_Fermi_Dirac", n_ignored=2)
-test_folder(Float64, "magnesium_E15_k3_GGA_Methfessel_Paxton", n_ignored=1)
-test_folder(Float64, "magnesium_E25_k5_GGA_Methfessel_Paxton", n_ignored=2)
-test_folder(Float64, "graphite_E20_k8_LDA_Methfessel_Paxton", scf_tol=7e-5, test_tol=7e-4,
-            n_ignored=2)
-
-# for file in readdir()
-#     if isdir(file) && isfile(joinpath(file, "generate.jl"))
-#         @testset "$file" begin
-#             if isfile(joinpath(file, "extra.jld")) && isfile(joinpath(file, "out_GSR.nc"))
-#                 test_folder(Float64, file)
-#             end
-#         end
-#         run_folder(file)
-#     end
-# end
+using Test
+using DFTK
+
+function load_reference(folder::EtsfFolder)
+    n_kpoints = size(folder.gsr["reduced_coordinates_of_kpoints"], 2)
+    bands = Vector{Vector{Float64}}(undef, n_kpoints)
+    for ik in 1:n_kpoints
+        bands[ik] = Vector(folder.eig["Eigenvalues"][:, ik, 1])
+    end
+
+    energies = Dict{Symbol, Float64}(
+                     :Ewald => folder.gsr["e_ewald"][:],
+                     :PspCorrection => folder.gsr["e_corepsp"][:],
+                     :PotXC => folder.gsr["e_xc"][:],
+                     :Kinetic => folder.gsr["e_kinetic"][:],
+                     :PotHartree => folder.gsr["e_hartree"][:],
+                     :PotLocal => folder.gsr["e_localpsp"][:],
+                     :PotNonLocal => folder.gsr["e_nonlocalpsp"][:],
+               )
+
+    (energies=energies, bands=bands)
+end
+
+
+function test_folder(T, folder; scf_tol=1e-8, n_ignored=0, test_tol=1e-6)
+    @testset "$folder" begin
+        etsf = EtsfFolder(folder)
+
+        basis = load_basis(T, etsf)
+        atoms = load_atoms(T, etsf)
+        ref = load_reference(etsf)
+        n_bands = length(ref.bands[1])
+
+        ham = Hamiltonian(basis, guess_density(basis))
+        scfres = self_consistent_field(ham, n_bands, tol=scf_tol)
+
+        println("etot    ", sum(values(scfres.energies)) - sum(values(ref.energies)))
+        for ik in 1:length(basis.kpoints)
+            @test eltype(scfres.orben[ik]) == T
+            @test eltype(scfres.Psi[ik]) == Complex{T}
+            println(ik, "  ", abs.(scfres.orben[ik] - ref.bands[ik]))
+        end
+        for ik in 1:length(basis.kpoints)
+            # Ignore last few bands, because these eigenvalues are hardest to converge
+            # and typically a bit random and unstable in the LOBPCG
+            diff = abs.(scfres.orben[ik] - ref.bands[ik])
+            @test maximum(diff[1:n_bands - n_ignored]) < test_tol
+        end
+        for (key, value) in pairs(scfres.energies)
+            if haskey(ref.energies, key)
+                @test value ≈ ref.energies[key] atol=5test_tol
+            end
+        end
+        @test sum(values(scfres.energies)) ≈ sum(values(ref.energies)) atol=test_tol
+    end
+end
+
+
+function main()
+    test_folder(Float64, "silicon_E10_k3_LDA", test_tol=5e-8)
+    test_folder(Float64, "silicon_E15_k4_LDA", n_ignored=1, test_tol=5e-8)
+    test_folder(Float64, "silicon_E25_k3_LDA", n_ignored=3, test_tol=1e-8)
+    test_folder(Float64, "silicon_E15_k4_GGA")
+    test_folder(Float64, "silicon_E25_k3_GGA")
+    test_folder(Float64, "silicon_E25_k4_GGA", n_ignored=1, test_tol=1e-7)
+    test_folder(Float64, "magnesium_E15_k3_LDA_Fermi_Dirac", n_ignored=2)
+    test_folder(Float64, "magnesium_E15_k3_GGA_Methfessel_Paxton", n_ignored=1)
+    test_folder(Float64, "magnesium_E25_k5_GGA_Methfessel_Paxton", n_ignored=2)
+    test_folder(Float64, "graphite_E20_k8_LDA_Methfessel_Paxton", scf_tol=7e-5, test_tol=7e-4,
+                n_ignored=2)
+end
+
+if abspath(PROGRAM_FILE) == @__FILE__
+    main()
+end
+
diff --git a/test/testcases_ABINIT/silicon_E10_k3_LDA/generate.jl b/test/testcases_ABINIT/silicon_E10_k3_LDA/generate.jl
@@ -1,16 +1,13 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using DFTK
 
-structure = build_silicon_structure()
-pspfile = joinpath(@__DIR__, "Si-q4-pade.abinit.hgh")
+abinitpseudos = [joinpath(@__DIR__, "Si-q4-pade.abinit.hgh")]
+pspmap = Dict(14 => "hgh/lda/si-q4", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[3, 3, 3], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=10,        # Hartree
-    nband=6,        # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc=1           # LDA_XC_TETER93
-)
-infile.extra = Dict("pspmap" => Dict(14 => "hgh/lda/si-q4", ), )
+atoms = [Element(14, load_psp(pspmap[14])) => silicon.positions]
+model = model_dft(silicon.lattice, :lda_xc_teter93, atoms)
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=10, kgrid=[3, 3, 3], n_bands=6, tol=1e-10,
+                    iscf=3)  # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/silicon_E15_k4_GGA/generate.jl b/test/testcases_ABINIT/silicon_E15_k4_GGA/generate.jl
@@ -1,16 +1,13 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using DFTK
 
-structure = build_silicon_structure()
-pspfile = joinpath(@__DIR__, "Si-q4-pbe.abinit.hgh")
+abinitpseudos = [joinpath(@__DIR__, "Si-q4-pbe.abinit.hgh")]
+pspmap = Dict(14 => "hgh/pbe/si-q4", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[4, 4, 4], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=15,        # Hartree
-    nband=6,        # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc=11          # PBE C and X
-)
-infile.extra = Dict("pspmap" => Dict(14 => "hgh/pbe/si-q4", ), )
+atoms = [Element(14, load_psp(pspmap[14])) => silicon.positions]
+model = model_dft(silicon.lattice, [:gga_c_pbe, :gga_x_pbe], atoms)
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=15, kgrid=[4, 4, 4], n_bands=6, tol=1e-10,
+                    iscf=3) # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/silicon_E15_k4_LDA/generate.jl b/test/testcases_ABINIT/silicon_E15_k4_LDA/generate.jl
@@ -1,16 +1,13 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using DFTK
 
-structure = build_silicon_structure()
-pspfile = joinpath(@__DIR__, "Si-q4-pade.abinit.hgh")
+abinitpseudos = [joinpath(@__DIR__, "Si-q4-pade.abinit.hgh")]
+pspmap = Dict(14 => "hgh/lda/si-q4", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[4, 4, 4], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=15,        # Hartree
-    nband=6,        # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc=1           # LDA_XC_TETER93
-)
-infile.extra = Dict("pspmap" => Dict(14 => "hgh/lda/si-q4", ), )
+atoms = [Element(14, load_psp(pspmap[14])) => silicon.positions]
+model = model_dft(silicon.lattice, :lda_xc_teter93, atoms)
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=15, kgrid=[4, 4, 4], n_bands=6, tol=1e-10,
+                    iscf=3) # Use Anderson mixing instead of minimisation
diff --git a/test/testcases_ABINIT/silicon_E25_k3_GGA/generate.jl b/test/testcases_ABINIT/silicon_E25_k3_GGA/generate.jl
@@ -1,16 +1,13 @@
-include("../testcases.jl")
+include("../../testcases.jl")
+using DFTK
 
-structure = build_silicon_structure()
-pspfile = joinpath(@__DIR__, "Si-q4-pbe.abinit.hgh")
+abinitpseudos = [joinpath(@__DIR__, "Si-q4-pbe.abinit.hgh")]
+pspmap = Dict(14 => "hgh/pbe/si-q4", )
 
-infile = abilab.AbinitInput(structure=structure, pseudos=abidata.pseudos(pspfile))
-infile.set_kmesh(ngkpt=[3, 3, 3], shiftk=[0, 0, 0])
-infile.set_vars(
-    ecut=25,        # Hartree
-    nband=10,       # Number of bands
-    tolvrs=1e-10,   # General tolerance settings
-    ixc="-101130",  # PBE C and X
-)
-infile.extra = Dict("pspmap" => Dict(14 => "hgh/pbe/si-q4", ), )
+atoms = [Element(14, load_psp(pspmap[14])) => silicon.positions]
+model = model_dft(silicon.lattice, [:gga_c_pbe, :gga_x_pbe], atoms)
 
-run_ABINIT_scf(infile, @__DIR__)
+DFTK.run_abinit_scf(model, @__DIR__;
+                    abinitpseudos=abinitpseudos, pspmap=pspmap,
+                    Ecut=25, kgrid=[3, 3, 3], n_bands=10, tol=1e-10,
+                    iscf=3) # Use Anderson mixing instead of minimisation