JuliaMolSim / DFTK.jl

Density-functional toolkit
https://docs.dftk.org
MIT License
442 stars 89 forks source link

Custom pseudopotential identifier fails on irrfbz_path #968

Open niklasschmitz opened 8 months ago

niklasschmitz commented 8 months ago

Currently (as of DFTK v0.6.18) there is an error when one tries to give a custom pseudopotential a new identifier (to avoid confusion with the default parameters) and then tries to extract the irrfbz_path from the model.

using DFTK

function create_silicon_model(psp)
    a = 10.26  # Silicon lattice constant in Bohr
    lattice = a / 2 * [[0 1 1.];
                        [1 0 1.];
                        [1 1 0.]]
    Si = ElementPsp(:Si; psp)
    atoms     = [Si, Si]
    positions = [ones(3)/8, -ones(3)/8]
    model_LDA(lattice, atoms, positions; temperature=1e-2)
end

psp_hgh  = load_psp("hgh/lda/si-q4.hgh")

irrfbz_path(create_silicon_model(psp_hgh))  # works as expected

psp_custom = PspHgh{Float64}(4, 0.44, [-7.33610297, 0.0, 0.0, 0.0], 1, [0.42273813, 0.48427842], [[5.90692831 -1.26189397; -1.26189397 3.25819622], [2.72701346;;]], "my_custom_identifier", "my custom name")

irrfbz_path(create_silicon_model(psp_custom))  # throws an error

The full stack trace:

julia> psp_hgh  = load_psp("hgh/lda/si-q4.hgh")
PspHgh{Float64}(4, 0.44, [-7.33610297, 0.0, 0.0, 0.0], 1, [0.42273813, 0.48427842], [[5.90692831 -1.26189397; -1.26189397 3.25819622], [2.72701346;;]], "hgh/lda/si-q4.hgh", "Si GTH-PADE-q4 GTH-LDA-q4")

julia> irrfbz_path(create_silicon_model(psp_hgh))  # works as expected
KPath{3} (6 points, 2 paths, 8 points in paths):
 points: :U => [0.625, 0.25, 0.625]
         :W => [0.5, 0.25, 0.75]
         :K => [0.375, 0.375, 0.75]
         :Γ => [0.0, 0.0, 0.0]
         :L => [0.5, 0.5, 0.5]
         :X => [0.5, 0.0, 0.5]
  paths: [:Γ, :X, :U]
         [:K, :Γ, :L, :W, :X]
  basis: [-0.612396, 0.612396, 0.612396]
         [0.612396, -0.612396, 0.612396]
         [0.612396, 0.612396, -0.612396]

julia> psp_custom = PspHgh{Float64}(4, 0.44, [-7.33610297, 0.0, 0.0, 0.0], 1, [0.42273813, 0.48427842], [[5.90692831 -1.26189397; -1.26189397 3.25819622], [2.72701346;;]], "my_custom_identifier", "my custom name")
PspHgh{Float64}(4, 0.44, [-7.33610297, 0.0, 0.0, 0.0], 1, [0.42273813, 0.48427842], [[5.90692831 -1.26189397; -1.26189397 3.25819622], [2.72701346;;]], "my_custom_identifier", "my custom name")

julia> irrfbz_path(create_silicon_model(psp_custom))
ERROR: Could not determine pseudopotential family of 'my_custom_identifier'
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] load_psp(key::String; kwargs::@Kwargs{})
    @ DFTK ~/git/DFTK.jl/src/pseudo/load_psp.jl:18
  [3] load_psp
    @ ~/git/DFTK.jl/src/pseudo/load_psp.jl:10 [inlined]
  [4] (::DFTK.var"#905#911"{AtomsBase.Atom{3, Unitful.Quantity{…}, Unitful.Quantity{…}, Unitful.Quantity{…}}, String})()
    @ DFTK ~/git/DFTK.jl/src/external/atomsbase.jl:27
  [5] get!(default::DFTK.var"#905#911"{AtomsBase.Atom{…}, String}, h::Dict{String, ElementPsp}, key::String)
    @ Base ./dict.jl:479
  [6] (::DFTK.var"#904#910"{Dict{…}})(atom::AtomsBase.Atom{3, Unitful.Quantity{…}, Unitful.Quantity{…}, Unitful.Quantity{…}})
    @ DFTK ~/git/DFTK.jl/src/external/atomsbase.jl:25
  [7] iterate
    @ ./generator.jl:47 [inlined]
  [8] collect(itr::Base.Generator{AtomsBase.FlexibleSystem{…}, DFTK.var"#904#910"{…}})
    @ Base ./array.jl:834
  [9] map
    @ ./abstractarray.jl:3313 [inlined]
 [10] parse_system(system::AtomsBase.FlexibleSystem{3, AtomsBase.Atom{…}, Unitful.Quantity{…}})
    @ DFTK ~/git/DFTK.jl/src/external/atomsbase.jl:19
 [11] spglib_cell(system::AtomsBase.FlexibleSystem{3, AtomsBase.Atom{…}, Unitful.Quantity{…}})
    @ DFTK ~/git/DFTK.jl/src/external/spglib.jl:31
 [12] spglib_cell(model::Model{Float64, Float64}, magnetic_moments::Vector{Any})
    @ DFTK ~/git/DFTK.jl/src/external/spglib.jl:35
 [13] irrfbz_path(model::Model{Float64, Float64}, magnetic_moments::Vector{Any}; dim::Int64, space_group_number::Int64)
    @ DFTK ~/git/DFTK.jl/src/postprocess/band_structure.jl:140
 [14] irrfbz_path
    @ ~/git/DFTK.jl/src/postprocess/band_structure.jl:106 [inlined]
 [15] irrfbz_path(model::Model{Float64, Float64})
    @ DFTK ~/git/DFTK.jl/src/postprocess/band_structure.jl:106
 [16] top-level scope
    @ REPL[25]:1
Some type information was truncated. Use `show(err)` to see complete types.
mfherbst commented 8 months ago

I think the issue is that in external/spglib.jl:35 we make a. AtomsBase system and this should drop the pseudo identifier. Look at the function atomic_system(::Model) in external/atomsbase.jl. There the issue should be.