From c820c36341b609ecfd4318c2a71c35d593c4c5bb Mon Sep 17 00:00:00 2001 From: "Michael F. Herbst" Date: Thu, 19 Jan 2023 20:29:21 +0100 Subject: [PATCH] Key-based property interface (#61) --- Project.toml | 2 +- docs/src/apireference.md | 2 + docs/src/overview.md | 5 ++- docs/src/tutorial.md | 89 ++++++++++++++++++++++++++++++++++------ src/AtomsBase.jl | 2 +- src/atom.jl | 40 ++++++++---------- src/atomview.jl | 8 ++++ src/fast_system.jl | 63 +++++++++++++++++----------- src/flexible_system.jl | 42 ++++++++++++------- src/interface.jl | 25 +++++++++++ src/show.jl | 22 +++++----- test/atom.jl | 78 +++++++++++++++++++++++------------ test/fast_system.jl | 21 ++++++++++ test/interface.jl | 10 +++++ test/printing.jl | 2 +- 15 files changed, 297 insertions(+), 114 deletions(-) diff --git a/Project.toml b/Project.toml index 30f7ed9..4f731b7 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "AtomsBase" uuid = "a963bdd2-2df7-4f54-a1ee-49d51e6be12a" authors = ["JuliaMolSim community"] -version = "0.2.5" +version = "0.3.0" [deps] PeriodicTable = "7b2266bf-644c-5ea3-82d8-af4bbd25a884" diff --git a/docs/src/apireference.md b/docs/src/apireference.md index 51ad32e..d93ae1d 100644 --- a/docs/src/apireference.md +++ b/docs/src/apireference.md @@ -19,6 +19,8 @@ isinfinite n_dimensions periodicity species_type +atomkeys +hasatomkey ``` ## Species / atom properties diff --git a/docs/src/overview.md b/docs/src/overview.md index bbce24d..7fb748c 100644 --- a/docs/src/overview.md +++ b/docs/src/overview.md @@ -35,7 +35,10 @@ additional behavior depending on context. ## System state and properties The only required properties to be specified of the system is the species and implementations of standard functions accessing the properties of the species, -currently [`position`](@ref), [`velocity`](@ref), [`atomic_symbol`](@ref), [`atomic_mass`](@ref), [`atomic_number`](@ref), [`n_dimensions`](@ref), [`element`](@ref). +currently + - Geometric information: [`position`](@ref), [`velocity`](@ref), [`n_dimensions`](@ref) + - Atomic information: [`atomic_symbol`](@ref), [`atomic_mass`](@ref), [`atomic_number`](@ref), [`element`](@ref) + - Atomic and system property accessors: `getindex`, `haskey`, `getkey`, `keys`, `pairs` Based on these methods respective equivalent methods acting on an `AbstractSystem` will be automatically available, e.g. using the iteration interface of the `AbstractSystem` (see above). Most of the property accessors on the diff --git a/docs/src/tutorial.md b/docs/src/tutorial.md index 91762e3..4f7f4eb 100644 --- a/docs/src/tutorial.md +++ b/docs/src/tutorial.md @@ -80,6 +80,20 @@ using Unitful, UnitfulAtomic, AtomsBase # hide deuterium = Atom(1, atomic_symbol=:D, [0, 1, 2.]u"bohr") ```` +An equivalent dict-like interface based on `keys`, `haskey`, `getkey` and `pairs` +is also available. For example +````@example atom +keys(atom) +```` +````@example atom +atom[:atomic_symbol] +```` +````@example atom +pairs(atom) +```` +This interface seamlessly generalises to working with user-specific atomic properties +as will be discussed next. + ### Optional atomic properties Custom properties can be easily attached to an `Atom` by supplying arbitrary keyword arguments upon construction. For example to attach a pseudopotential @@ -88,11 +102,20 @@ for using the structure with [DFTK](https://dftk.org), construct the atom as using Unitful, UnitfulAtomic, AtomsBase # hide atom = Atom(:C, [0, 1, 2.]u"bohr", pseudopotential="hgh/lda/c-q4") ```` -which will make the pseudopotential identifier available as `atom.pseudopotential`. +which will make the pseudopotential identifier available as +````@example atomprop +atom[:pseudopotential] +```` +Notice that such custom properties are fully integrated with the standard atomic properties, +e.g. automatically available from the `keys`, `haskey` and `pairs` functions, e.g.: +````@example atomprop +@show haskey(atom, :pseudopotential) +pairs(atom) +```` Updating an atomic property proceeds similarly. E.g. ````@example atomprop using Unitful, UnitfulAtomic, AtomsBase # hide -newatom = Atom(;atom=atom, atomic_mass=13u"u") +newatom = Atom(atom; atomic_mass=13u"u") ```` makes a new carbon atom with all properties identical to `atom` (including custom ones), but setting the `atomic_mass` to 13 units. @@ -112,6 +135,12 @@ Property name | Unit / Type | Description `:magnetic_moments` | `Union{Float64,Vector{Float64}}` | Initial magnetic moment `:pseudopotential` | `String` | Pseudopotential or PAW keyword or `""` if Coulomb potential employed +A convenient way to iterate over all data stored in an atom offers the `pairs` function: +````@example atomprop +for (k, v) in pairs(atom) + println("$k = $v") +end +```` ## System interface and conventions Once the atoms are constructed these can be assembled into a system. @@ -119,24 +148,55 @@ For example to place a hydrogen molecule into a cubic box of `10Å` and periodic boundary conditions, use: ````@example system using Unitful, UnitfulAtomic, AtomsBase # hide -bounding_box = [[10.0, 0.0, 0.0], [0.0, 10.0, 0.0], [0.0, 0.0, 10.0]]u"Å" +box = [[10.0, 0.0, 0.0], [0.0, 10.0, 0.0], [0.0, 0.0, 10.0]]u"Å" boundary_conditions = [Periodic(), Periodic(), Periodic()] hydrogen = FlexibleSystem([Atom(:H, [0, 0, 1.]u"bohr"), Atom(:H, [0, 0, 3.]u"bohr")], - bounding_box, boundary_conditions) + box, boundary_conditions) ```` An update constructor for systems is supported as well (see [`AbstractSystem`](@ref)). For example ````@example system AbstractSystem(hydrogen; bounding_box=[[5.0, 0.0, 0.0], [0.0, 5.0, 0.0], [0.0, 0.0, 5.0]]u"Å") ```` -Note that `FlexibleSystem( ... )` would have worked as well in this example (since we are +To update the atomic composition of the system, this function supports an `atoms` (or `particles`) +keyword argument to supply the new set of atoms to be contained in the system. + +Note that in this example `FlexibleSystem( ... )` would have worked as well (since we are updating a `FlexibleSystem`). However, using the `AbstractSystem` constructor to update the system is more general as it allows for type-specific dispatching when updating other data structures implementing the `AbstractSystem` interface. -Oftentimes more convenient are the functions -[`atomic_system`](@ref), [`isolated_system`](@ref), [`periodic_system`](@ref), -which cover some standard atomic system setups. +Similar to the atoms, system objects similarly support a functional-style access to system properties +as well as a dict-style access: +````@example system +bounding_box(hydrogen) +```` +````@example system +hydrogen[:boundary_conditions] +```` +````@example system +pairs(hydrogen) +```` +Moreover atomic properties of a specific atom or all atoms can be directly queried using +the indexing notation: +````@example system +hydrogen[1, :position] # Position of first atom +```` +````@example system +hydrogen[:, :position] # All atomic symbols +```` +Finally, supported keys of atomic properties can be directly queried at the system level +using [`atomkeys`](@ref) and [`hasatomkey`](@ref). Note that these functions only apply to atomic +properties which are supported by *all* atoms of a system. In other words if a custom atomic property is only +set in a few of the contained atoms, these functions will not consider it. +````@example system +atomkeys(hydrogen) +```` + +For constructing atomic systems the functions +[`atomic_system`](@ref), [`isolated_system`](@ref), [`periodic_system`](@ref) +are oftentimes more convenient as they provide specialisations +for some standard atomic system setups. For example to setup a hydrogen system with periodic BCs, we can issue ````@example using Unitful, UnitfulAtomic, AtomsBase # hide @@ -165,11 +225,16 @@ hydrogen = isolated_system([:H => [0, 0, 1.]u"bohr", ### Optional system properties Similar to atoms, systems also support storing arbitrary data, for example -```@example +````@example sysprop using Unitful, UnitfulAtomic, AtomsBase # hide -hydrogen = isolated_system([:H => [0, 0, 1.]u"bohr", :H => [0, 0, 3.]u"bohr"]; extra_data=42) -``` -Again to simplify interoperability some optional properties are reserved, namely: +system = isolated_system([:H => [0, 0, 1.]u"bohr", :H => [0, 0, 3.]u"bohr"]; extra_data=42) +```` +Again these custom properties are fully integrated with `keys`, `haskey`, `pairs` and `getkey`. +````@example sysprop +@show keys(system) +```` +Some property names are reserved and should be considered by all libraries +supporting `AtomsBase` if possible: Property name | Unit / Type | Description :-------------- | :----------------- | :--------------------- diff --git a/src/AtomsBase.jl b/src/AtomsBase.jl index 34a3308..15eb175 100644 --- a/src/AtomsBase.jl +++ b/src/AtomsBase.jl @@ -7,8 +7,8 @@ include("interface.jl") include("properties.jl") include("show.jl") include("flexible_system.jl") -include("atom.jl") include("atomview.jl") +include("atom.jl") include("fast_system.jl") end diff --git a/src/atom.jl b/src/atom.jl index 9141e3d..3d2e4db 100644 --- a/src/atom.jl +++ b/src/atom.jl @@ -22,19 +22,20 @@ atomic_number(atom::Atom) = atom.atomic_number element(atom::Atom) = element(atomic_number(atom)) n_dimensions(::Atom{D}) where {D} = D -Base.hasproperty(at::Atom, x::Symbol) = hasfield(Atom, x) || haskey(at.data, x) -Base.getproperty(at::Atom, x::Symbol) = hasfield(Atom, x) ? getfield(at, x) : getindex(at.data, x) -function Base.propertynames(at::Atom, private::Bool=false) - if private - (fieldnames(Atom)..., keys(at.data)...) - else - (filter(!isequal(:data), fieldnames(Atom))..., keys(at.data)...) - end +Base.getindex(at::Atom, x::Symbol) = hasfield(Atom, x) ? getfield(at, x) : getindex(at.data, x) +Base.haskey(at::Atom, x::Symbol) = hasfield(Atom, x) || haskey(at.data, x) +function Base.getkey(at::Atom, x::Symbol, default) + hasfield(Atom, x) ? getfield(at, x) : getkey(at.data, x, default) +end +function Base.keys(at::Atom) + (:position, :velocity, :atomic_symbol, :atomic_number, :atomic_mass, keys(at.data)...) end +Base.pairs(at::Atom) = (k => at[k] for k in keys(at)) """ Atom(identifier::AtomId, position::AbstractVector; kwargs...) Atom(identifier::AtomId, position::AbstractVector, velocity::AbstractVector; kwargs...) + Atom(; atomic_number, position, velocity=zeros(D)u"bohr/s", kwargs...) Construct an atomic located at the cartesian coordinates `position` with (optionally) the given cartesian `velocity`. Note that `AtomId = Union{Symbol,AbstractString,Integer}`. @@ -53,18 +54,17 @@ function Atom(identifier::AtomId, atomic_number, atomic_mass, Dict(kwargs...)) end function Atom(id::AtomId, position::AbstractVector, velocity::Missing; kwargs...) - Atom(id, position; kwargs...) + Atom(id, position, zeros(length(position))u"bohr/s"; kwargs...) +end +function Atom(; atomic_symbol, position, velocity=zeros(length(position))u"bohr/s", kwargs...) + Atom(atomic_symbol, position, velocity; atomic_symbol, kwargs...) end """ Atom(atom::Atom; kwargs...) - Atom(; atom, kwargs...) Update constructor. Construct a new `Atom`, by amending the data contained -in the passed `atom` object. Note that the first version only works if `atom` is an `Atom`, -while the second version works on arbitrary species adhering to the `AtomsBase` conventions. -In library code the second version should therefore be preferred. - +in the passed `atom` object. Supported `kwargs` include `atomic_symbol`, `atomic_number`, `atomic_mass`, `charge`, `multiplicity` as well as user-specific custom properties. @@ -75,18 +75,10 @@ julia> hydrogen = Atom(:H, zeros(3)u"Å") ``` and now amend its charge and atomic mass ```julia-repl -julia> Atom(; atom, atomic_mass=1.0u"u", charge=-1.0u"e_au") +julia> Atom(atom; atomic_mass=1.0u"u", charge=-1.0u"e_au") ``` """ -function Atom(;atom, kwargs...) - extra = atom isa Atom ? atom.data : (; ) - Atom(atomic_number(atom), position(atom), velocity(atom); - atomic_symbol=atomic_symbol(atom), - atomic_number=atomic_number(atom), - atomic_mass=atomic_mass(atom), - extra..., kwargs...) -end -Atom(atom::Atom; kwargs...) = Atom(; atom, kwargs...) +Atom(atom::Atom; kwargs...) = Atom(; pairs(atom)..., kwargs...) function Base.convert(::Type{Atom}, id_pos::Pair{<:AtomId,<:AbstractVector{<:Unitful.Length}}) Atom(id_pos...) diff --git a/src/atomview.jl b/src/atomview.jl index fdfa31b..943ad07 100644 --- a/src/atomview.jl +++ b/src/atomview.jl @@ -16,3 +16,11 @@ element(atom::AtomView) = element(atomic_number(atom)) Base.show(io::IO, at::AtomView) = show_atom(io, at) Base.show(io::IO, mime::MIME"text/plain", at::AtomView) = show_atom(io, mime, at) + +Base.getindex(v::AtomView, x::Symbol) = getindex(v.system, v.index, x) +Base.haskey(v::AtomView, x::Symbol) = hasatomkey(v.system, x) +function Base.getkey(v::AtomView, x::Symbol, default) + hasatomkey(v.system, x) ? v[x] : default +end +Base.keys(v::AtomView) = atomkeys(v.system) +Base.pairs(at::AtomView) = (k => at[k] for k in keys(at)) diff --git a/src/fast_system.jl b/src/fast_system.jl index 8cf91f2..9e4d107 100644 --- a/src/fast_system.jl +++ b/src/fast_system.jl @@ -4,12 +4,12 @@ export FastSystem struct FastSystem{D, L <: Unitful.Length, M <: Unitful.Mass} <: AbstractSystem{D} - box::SVector{D, SVector{D, L}} + bounding_box::SVector{D, SVector{D, L}} boundary_conditions::SVector{D, BoundaryCondition} - positions::Vector{SVector{D, L}} - atomic_symbols::Vector{Symbol} - atomic_numbers::Vector{Int} - atomic_masses::Vector{M} + position::Vector{SVector{D, L}} + atomic_symbol::Vector{Symbol} + atomic_number::Vector{Int} + atomic_mass::Vector{M} end # Constructor to fetch the types @@ -41,22 +41,39 @@ function FastSystem(particles, box, boundary_conditions) atomic_number.(particles), atomic_mass.(particles)) end -bounding_box(sys::FastSystem) = sys.box +bounding_box(sys::FastSystem) = sys.bounding_box boundary_conditions(sys::FastSystem) = sys.boundary_conditions -Base.length(sys::FastSystem) = length(sys.positions) -Base.size(sys::FastSystem) = size(sys.positions) - -species_type(sys::FS) where {FS <: FastSystem} = AtomView{FS} -Base.getindex(sys::FastSystem, index::Int) = AtomView(sys, index) - -position(s::FastSystem) = s.positions -atomic_symbol(s::FastSystem) = s.atomic_symbols -atomic_number(s::FastSystem) = s.atomic_numbers -atomic_mass(s::FastSystem) = s.atomic_masses -velocity(s::FastSystem) = missing - -position(s::FastSystem, i) = s.positions[i] -atomic_symbol(s::FastSystem, i) = s.atomic_symbols[i] -atomic_number(s::FastSystem, i) = s.atomic_numbers[i] -atomic_mass(s::FastSystem, i) = s.atomic_masses[i] -velocity(s::FastSystem, i) = missing +Base.length(sys::FastSystem) = length(sys.position) +Base.size(sys::FastSystem) = size(sys.position) + +species_type(::FS) where {FS <: FastSystem} = AtomView{FS} +Base.getindex(sys::FastSystem, i::Integer) = AtomView(sys, i) + +position(s::FastSystem) = s.position +atomic_symbol(s::FastSystem) = s.atomic_symbol +atomic_number(s::FastSystem) = s.atomic_number +atomic_mass(s::FastSystem) = s.atomic_mass +velocity(::FastSystem) = missing + +position(s::FastSystem, i) = s.position[i] +atomic_symbol(s::FastSystem, i) = s.atomic_symbol[i] +atomic_number(s::FastSystem, i) = s.atomic_number[i] +atomic_mass(s::FastSystem, i) = s.atomic_mass[i] +velocity(::FastSystem, i) = missing + +# System property access +function Base.getindex(system::FastSystem, x::Symbol) + if x in (:bounding_box, :boundary_conditions) + getfield(system, x) + else + throw(KeyError("Key $x not found")) + end +end +Base.haskey(::FastSystem, x::Symbol) = x in (:bounding_box, :boundary_conditions) +Base.keys(::FastSystem) = (:bounding_box, :boundary_conditions) + +# Atom and atom property access +atomkeys(::FastSystem) = (:position, :atomic_symbol, :atomic_number, :atomic_mass) +hasatomkey(system::FastSystem, x::Symbol) = x in atomkeys(system) +Base.getindex(system::FastSystem, i::Integer, x::Symbol) = getfield(system, x)[i] +Base.getindex(system::FastSystem, ::Colon, x::Symbol) = getfield(system, x) diff --git a/src/flexible_system.jl b/src/flexible_system.jl index 8c73548..4f6112f 100644 --- a/src/flexible_system.jl +++ b/src/flexible_system.jl @@ -6,16 +6,33 @@ export FlexibleSystem struct FlexibleSystem{D, S, L<:Unitful.Length} <: AbstractSystem{D} particles::AbstractVector{S} - box::SVector{D, SVector{D, L}} + bounding_box::SVector{D, SVector{D, L}} boundary_conditions::SVector{D, BoundaryCondition} data::Dict{Symbol, Any} # Store arbitrary data about the atom. end -Base.hasproperty(system::FlexibleSystem, x::Symbol) = hasfield(FlexibleSystem, x) || haskey(system.data, x) -Base.getproperty(system::FlexibleSystem, x::Symbol) = hasfield(FlexibleSystem, x) ? getfield(system, x) : getindex(system.data, x) +# System property access +function Base.getindex(system::FlexibleSystem, x::Symbol) + if x in (:bounding_box, :boundary_conditions) + getfield(system, x) + else + getindex(system.data, x) + end +end +function Base.haskey(system::FlexibleSystem, x::Symbol) + x in (:bounding_box, :boundary_conditions) || haskey(system.data, x) +end +Base.keys(system::FlexibleSystem) = (:bounding_box, :boundary_conditions, keys(system.data)...) + +# Atom and atom property access +Base.getindex(system::FlexibleSystem, i::Integer) = system.particles[i] +Base.getindex(system::FlexibleSystem, i::Integer, x::Symbol) = system.particles[i][x] +Base.getindex(system::FlexibleSystem, ::Colon, x::Symbol) = [at[x] for at in system.particles] + """ - FlexibleSystem(particles, box, boundary_conditions; kwargs...) + FlexibleSystem(particles, bounding_box, boundary_conditions; kwargs...) + FlexibleSystem(particles; bounding_box, boundary_conditions, kwargs...) Construct a flexible system, a versatile data structure for atomistic systems, which puts an emphasis on flexibility rather than speed. @@ -30,7 +47,10 @@ function FlexibleSystem( if !all(length.(box) .== D) throw(ArgumentError("Box must have D vectors of length D")) end - FlexibleSystem{D, S, L}(convert.(Atom, particles), box, boundary_conditions, Dict(kwargs...)) + FlexibleSystem{D, S, L}(particles, box, boundary_conditions, Dict(kwargs...)) +end +function FlexibleSystem(particles; bounding_box, boundary_conditions, kwargs...) + FlexibleSystem(particles, bounding_box, boundary_conditions; kwargs...) end """ @@ -38,14 +58,9 @@ end Update constructor. See [`AbstractSystem`](@ref) for details. """ -function FlexibleSystem(system::AbstractSystem; - particles=nothing, atoms=nothing, - bounding_box=bounding_box(system), - boundary_conditions=boundary_conditions(system), - kwargs...) +function FlexibleSystem(system::AbstractSystem; particles=nothing, atoms=nothing, kwargs...) particles = something(particles, atoms, collect(system)) - extra = system isa FlexibleSystem ? system.data : (; ) - FlexibleSystem(particles, bounding_box, boundary_conditions; extra..., kwargs...) + FlexibleSystem(particles; pairs(system)..., kwargs...) end """ @@ -66,10 +81,9 @@ julia> AbstractSystem(system; bounding_box= ..., atoms = ... ) """ AbstractSystem(system::AbstractSystem; kwargs...) = FlexibleSystem(system; kwargs...) -bounding_box(sys::FlexibleSystem) = sys.box +bounding_box(sys::FlexibleSystem) = sys.bounding_box boundary_conditions(sys::FlexibleSystem) = sys.boundary_conditions species_type(sys::FlexibleSystem{D, S, L}) where {D, S, L} = S Base.size(sys::FlexibleSystem) = size(sys.particles) Base.length(sys::FlexibleSystem) = length(sys.particles) -Base.getindex(sys::FlexibleSystem, i::Integer) = getindex(sys.particles, i) diff --git a/src/interface.jl b/src/interface.jl index 76f34f6..43d1be2 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -5,6 +5,7 @@ export AbstractSystem export BoundaryCondition, DirichletZero, Periodic, infinite_box, isinfinite export bounding_box, boundary_conditions, periodicity, n_dimensions, species_type export position, velocity, element, atomic_mass, atomic_number, atomic_symbol +export atomkeys, hasatomkey # # Identifier for boundary conditions per dimension @@ -73,6 +74,7 @@ Base.lastindex(s::AbstractSystem) = length(s) # default to 1D indexing Base.iterate(sys::AbstractSystem, state = firstindex(sys)) = (firstindex(sys) <= state <= length(sys)) ? (@inbounds sys[state], state + 1) : nothing +Base.eltype(sys::AbstractSystem) = species_type(sys) # TODO Support similar, push, ... @@ -157,3 +159,26 @@ this may be `:D` while `atomic_number` is still `1`. """ atomic_number(sys::AbstractSystem) = atomic_number.(sys) atomic_number(sys::AbstractSystem, index) = atomic_number(sys[index]) + +""" + atomkeys(sys::AbstractSystem) + +Return the atomic properties, which are available in all atoms of the system. +""" +function atomkeys(system::AbstractSystem) + atkeys = length(system) == 0 ? () : keys(system[1]) + filter(k -> hasatomkey(system, k), atkeys) +end + +""" + hasatomkey(system::AbstractSystem, x::Symbol) + +Returns true whether the passed property available in all atoms of the passed system. +""" +hasatomkey(system::AbstractSystem, x::Symbol) = all(at -> haskey(at, x), system) + +# Defaults for system +Base.pairs(system::AbstractSystem) = (k => system[k] for k in keys(system)) +function Base.getkey(system::AbstractSystem, x::Symbol, default) + haskey(system, x) ? getindex(system, x) : default +end diff --git a/src/show.jl b/src/show.jl index 3880094..0bd66d1 100644 --- a/src/show.jl +++ b/src/show.jl @@ -51,11 +51,10 @@ function show_system(io::IO, ::MIME"text/plain", system::AbstractSystem{D}) wher end end - if system isa FlexibleSystem - for (k, v) in pairs(system.data) - extra_line = true - @printf io " %-17s : %s\n" string(k) string(v) - end + for (k, v) in pairs(system) + k in (:bounding_box, :boundary_conditions) && continue + extra_line = true + @printf io " %-17s : %s\n" string(k) string(v) end # TODO Better would be some ascii-graphical representation of the structure @@ -75,7 +74,8 @@ end function show_atom(io::IO, at) pos = [(@sprintf "%8.6g" ustrip(p)) for p in position(at)] posunit = unit(eltype(position(at))) - print(io, "Atom(", (@sprintf "%-3s" (string(atomic_symbol(at))) * ","), " [", + print(io, typeof(at).name.name, "(") + print(io, (@sprintf "%-3s" (string(atomic_symbol(at))) * ","), " [", join(pos, ", "), "]u\"$posunit\"") if ismissing(velocity(at)) || iszero(velocity(at)) print(io, ")") @@ -87,7 +87,8 @@ function show_atom(io::IO, at) end function show_atom(io::IO, ::MIME"text/plain", at) - println(io, "Atom(", atomic_symbol(at), ", atomic_number = ", atomic_number(at), + print(io, typeof(at).name.name, "(") + println(io, atomic_symbol(at), ", atomic_number = ", atomic_number(at), ", atomic_mass = ", atomic_mass(at), "):") pos = [(@sprintf "%.8g" ustrip(p)) for p in position(at)] @@ -98,9 +99,8 @@ function show_atom(io::IO, ::MIME"text/plain", at) velunit = unit(eltype(velocity(at))) @printf io " %-17s : [%s]u\"%s\"\n" "velocity" join(vel, ",") string(velunit) end - if at isa Atom - for (k, v) in pairs(at.data) - @printf io " %-17s : %s\n" string(k) string(v) - end + for (k, v) in pairs(at) + k in (:atomic_number, :atomic_mass, :atomic_symbol, :position, :velocity) && continue + @printf io " %-17s : %s\n" string(k) string(v) end end diff --git a/test/atom.jl b/test/atom.jl index 98fa173..90d976d 100644 --- a/test/atom.jl +++ b/test/atom.jl @@ -1,8 +1,8 @@ using AtomsBase using Unitful +using UnitfulAtomic using Test - @testset "atomic systems" begin @testset "Atom construction" begin at2D = Atom(:Si, zeros(2) * u"m", extradata=41) @@ -14,27 +14,32 @@ using Test @test position(at) == zeros(3) * u"m" @test velocity(at) == zeros(3) * u"bohr/s" @test element(at) == element(:Si) - @test at.atomic_symbol == :Si - @test at.atomic_number == 14 - @test hasproperty(at, :atomic_mass) - @test hasproperty(at, :atomic_symbol) - @test hasproperty(at, :atomic_number) - @test hasproperty(at, :extradata) - @test at.extradata == 42 - - @test propertynames(at) == (:position, :velocity, :atomic_symbol, - :atomic_number, :atomic_mass, :extradata) + @test atomic_symbol(at) == :Si + @test atomic_number(at) == 14 + @test at[:atomic_symbol] == :Si + @test at[:atomic_number] == 14 + @test haskey(at, :atomic_mass) + @test haskey(at, :atomic_symbol) + @test haskey(at, :atomic_number) + @test haskey(at, :extradata) + @test at[:extradata] == 42 + + @test keys(at) == (:position, :velocity, :atomic_symbol, + :atomic_number, :atomic_mass, :extradata) # Test update constructor newatom = Atom(at; extradata=43, atomic_number=15) - @test propertynames(at) == propertynames(newatom) - @test newatom.extradata == 43 - @test newatom.atomic_number == 15 - - newatom = Atom(; atom=at, extradata=43, atomic_number=15) - @test propertynames(at) == propertynames(newatom) - @test newatom.extradata == 43 - @test newatom.atomic_number == 15 + @test keys(at) == keys(newatom) + @test newatom[:extradata] == 43 + @test newatom[:atomic_number] == 15 + + newatom = Atom(at; extradata=43, atomic_number=15) + @test keys(at) == keys(newatom) + @test newatom[:extradata] == 43 + @test newatom[:atomic_number] == 15 + + newatom = Atom(:Si, ones(3)u"m", missing) + @test iszero(newatom[:velocity]) end @testset "flexible atomic systems" begin @@ -43,15 +48,36 @@ using Test dic = Dict{String, Any}("extradata_dic"=>"44") atoms = [:Si => [0.0, 1.0, 1.5]u"Å", :C => [0.0, 0.8, 1.7]u"Å", - Atom(:H, zeros(3) * u"Å", ones(3) * u"bohr/s")] + Atom(:H, zeros(3) * u"Å", ones(3) * u"bohr/s"; dat=3.0)] system = atomic_system(atoms, box, bcs, extradata=45; dic) @test length(system) == 3 @test atomic_symbol(system) == [:Si, :C, :H] @test boundary_conditions(system) == [Periodic(), DirichletZero(), DirichletZero()] @test position(system) == [[0.0, 1.0, 1.5], [0.0, 0.8, 1.7], [0.0, 0.0, 0.0]]u"Å" @test velocity(system) == [[0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [1.0, 1.0, 1.0]]u"bohr/s" - @test system.extradata == 45 - @test system.dic["extradata_dic"] == "44" + @test system[:extradata] == 45 + @test system[:dic]["extradata_dic"] == "44" + @test system[3] == atoms[3] + @test system[3, :dat] == 3.0 + @test system[2, :atomic_symbol] == :C + @test system[:, :atomic_number] == [14, 6, 1] + @test atomkeys(system) == (:position, :velocity, :atomic_symbol, + :atomic_number, :atomic_mass) + @test hasatomkey(system, :atomic_mass) + @test !hasatomkey(system, :blubber) + @test getkey(system, :blubber, :adidi) == :adidi + + @test collect(pairs(system)) == [ + :bounding_box => box, :boundary_conditions => bcs, + :extradata => 45, :dic => dic + ] + @test collect(pairs(system[1])) == [ + :position => [0.0, 1.0, 1.5]u"Å", + :velocity => zeros(3)u"Å/s", + :atomic_symbol => :Si, + :atomic_number => 14, + :atomic_mass => AtomsBase.element(:Si).atomic_mass, + ] # Test update constructor newatoms = [system[1], system[2]] @@ -75,8 +101,8 @@ using Test @test position(system) == [[0.0, 1.0, 1.5], [0.0, 0.8, 1.7], [0.0, 0.0, 0.0]]u"Å" @test velocity(system) == [[0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]u"bohr/s" @test bounding_box(system) == infinite_box(3) - @test system.extradata == 46 - @test system.dic["extradata_dic"] == "47" + @test system[:extradata] == 46 + @test system[:dic]["extradata_dic"] == "47" end @testset "periodic_system" begin @@ -91,8 +117,8 @@ using Test @test atomic_symbol(system) == [:Si, :C, :H] @test boundary_conditions(system) == [Periodic(), Periodic(), Periodic()] @test position(system) == [[0.0, -0.625, 0.0], [1.25, 0.0, 0.0], [0.0, 0.0, 0.0]]u"Å" - @test system.extradata == 48 - @test system.dic["extradata_dic"] == "49" + @test system[:extradata] == 48 + @test system[:dic]["extradata_dic"] == "49" end diff --git a/test/fast_system.jl b/test/fast_system.jl index de93d5a..e369545 100644 --- a/test/fast_system.jl +++ b/test/fast_system.jl @@ -15,8 +15,29 @@ using PeriodicTable @test atomic_mass(system) == [12.011, 12.011]u"u" @test boundary_conditions(system) == bcs @test bounding_box(system) == box + @test system[:boundary_conditions] == bcs + @test system[:bounding_box] == box @test !isinfinite(system) @test element(system[1]) == element(:C) + @test keys(system) == (:bounding_box, :boundary_conditions) + @test atomkeys(system) == (:position, :atomic_symbol, :atomic_number, :atomic_mass) + @test keys(system[1]) == (:position, :atomic_symbol, :atomic_number, :atomic_mass) + @test hasatomkey(system, :atomic_symbol) + @test system[1, :atomic_number] == 6 + @test system[:, :atomic_symbol] == [:C, :C] + @test system[2][:position] == system[2, :position] + @test system[2][:position] == [0.75, 0.75, 0.75]u"m" + @test haskey(system[1], :position) + @test !haskey(system[1], :abc) + @test getkey(system[1], :dagger, 3) == 3 + + @test collect(pairs(system)) == [(:bounding_box => box), (:boundary_conditions => bcs)] + @test collect(pairs(system[1])) == [ + :position => position(atoms[1]), + :atomic_symbol => :C, + :atomic_number => 6, + :atomic_mass => atomic_mass(atoms[1]), + ] # Test AtomView for method in (position, atomic_mass, atomic_symbol, atomic_number) diff --git a/test/interface.jl b/test/interface.jl index a450dff..7d0f6e0 100644 --- a/test/interface.jl +++ b/test/interface.jl @@ -18,6 +18,10 @@ using PeriodicTable @test atomic_number(atoms[1]) == 6 @test atomic_mass(atoms[1]) == 12.011u"u" @test element(atoms[1]) == element(:C) + @test atoms[1][:atomic_number] == 6 + @test keys(atoms[1]) == (:position, :velocity, :atomic_symbol, + :atomic_number, :atomic_mass) + @test getkey(atoms[1], :blubber, :adidi) == :adidi end @testset "System" begin @@ -43,6 +47,12 @@ using PeriodicTable @test atomic_number(flexible, 2) == 6 @test atomic_mass(flexible, 1) == 12.011u"u" + @test atomkeys(flexible) == (:position, :velocity, :atomic_symbol, + :atomic_number, :atomic_mass) + @test hasatomkey(flexible, :atomic_number) + @test flexible[1, :atomic_symbol] == :C + @test flexible[:, :atomic_symbol] == [:C, :C] + @test ismissing(velocity(fast)) @test all(position(fast) .== position(flexible)) @test all(atomic_symbol(fast) .== atomic_symbol(flexible)) diff --git a/test/printing.jl b/test/printing.jl index 2fdfe32..25b83b2 100644 --- a/test/printing.jl +++ b/test/printing.jl @@ -11,7 +11,7 @@ using Test :C => [0.125, 0.0, 0.0]] box = [[10, 0.0, 0.0], [0.0, 5, 0.0], [0.0, 0.0, 7]]u"Å" - flexible_system = periodic_system(atoms, box; fractional=true) + flexible_system = periodic_system(atoms, box; fractional=true, data=-12) @test repr(flexible_system) == """ FlexibleSystem(CSi, periodic = TTT, bounding_box = [[10.0, 0.0, 0.0], [0.0, 5.0, 0.0], [0.0, 0.0, 7.0]]u"Å")""" show(stdout, MIME("text/plain"), flexible_system)