From 0331a85c4f0acccf0e7da87d2c4c8d5905c4ce74 Mon Sep 17 00:00:00 2001 From: "Michael F. Herbst" Date: Thu, 15 Dec 2022 22:16:23 +0100 Subject: [PATCH] More detailed printing --- Project.toml | 1 + src/AtomsBase.jl | 4 +- src/atom.jl | 19 +++----- src/atomview.jl | 4 +- src/element.jl | 7 +++ src/fast_system.jl | 5 -- src/flexible_system.jl | 5 -- src/interface.jl | 16 ++++++- src/properties.jl | 18 ------- src/show.jl | 104 +++++++++++++++++++++++++++++++++++++++++ test/fast_system.jl | 2 +- test/interface.jl | 2 +- test/printing.jl | 12 +++-- 13 files changed, 150 insertions(+), 49 deletions(-) create mode 100644 src/element.jl create mode 100644 src/show.jl diff --git a/Project.toml b/Project.toml index 76d0966..f179757 100644 --- a/Project.toml +++ b/Project.toml @@ -5,6 +5,7 @@ version = "0.2.3" [deps] PeriodicTable = "7b2266bf-644c-5ea3-82d8-af4bbd25a884" +Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" UnitfulAtomic = "a7773ee8-282e-5fa2-be4e-bd808c38a91a" diff --git a/src/AtomsBase.jl b/src/AtomsBase.jl index a8d5918..f7c057f 100644 --- a/src/AtomsBase.jl +++ b/src/AtomsBase.jl @@ -1,11 +1,13 @@ module AtomsBase using Unitful using UnitfulAtomic -import PeriodicTable: elements +using Printf using StaticArrays include("interface.jl") +include("element.jl") include("properties.jl") +include("show.jl") include("flexible_system.jl") include("atom.jl") include("atomview.jl") diff --git a/src/atom.jl b/src/atom.jl index 8d8852c..26ff2e7 100644 --- a/src/atom.jl +++ b/src/atom.jl @@ -3,9 +3,6 @@ # export Atom, atomic_system, periodic_system, isolated_system -# Valid types for atom identifiers -const AtomId = Union{Symbol,AbstractString,Integer} - struct Atom{D, L<:Unitful.Length, V<:Unitful.Velocity, M<:Unitful.Mass} position::SVector{D, L} velocity::SVector{D, V} @@ -19,8 +16,7 @@ position(atom::Atom) = atom.position atomic_mass(atom::Atom) = atom.atomic_mass atomic_symbol(atom::Atom) = atom.atomic_symbol atomic_number(atom::Atom) = atom.atomic_number -element(atom::Atom) = elements[atomic_symbol(atom)] -n_dimensions(atom::Atom{D}) where {D} = D +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) @@ -35,9 +31,9 @@ end function Atom(identifier::AtomId, position::AbstractVector{L}, velocity::AbstractVector{V}=zeros(length(position))u"bohr/s"; - atomic_symbol=Symbol(elements[identifier].symbol), - atomic_number=elements[identifier].number, - atomic_mass::M=elements[identifier].atomic_mass, + atomic_symbol=Symbol(element(identifier).symbol), + atomic_number=element(identifier).number, + atomic_mass::M=element(identifier).atomic_mass, kwargs...) where {L <: Unitful.Length, V <: Unitful.Velocity, M <: Unitful.Mass} Atom{length(position), L, V, M}(position, velocity, atomic_symbol, atomic_number, atomic_mass, Dict(kwargs...)) @@ -58,10 +54,9 @@ function Base.convert(::Type{Atom}, id_pos::Pair{<:AtomId,<:AbstractVector{<:Uni Atom(id_pos...) end -function Base.show(io::IO, at::Atom{D, L}) where {D, L} - pos = ustrip.(at.position) - print(io, "Atom($(at.atomic_symbol), [", join(pos, ", "), "]u\"$(unit(L))\")") -end +Base.show(io::IO, at::Atom) = show_atom(io, at) +Base.show(io::IO, mime::MIME"text/plain", at::Atom) = show_atom(io, mime, at) + # # Special high-level functions to construct atomic systems diff --git a/src/atomview.jl b/src/atomview.jl index bd937bf..7208b61 100644 --- a/src/atomview.jl +++ b/src/atomview.jl @@ -11,5 +11,7 @@ position(v::AtomView) = position(v.system, v.index) atomic_mass(v::AtomView) = atomic_mass(v.system, v.index) atomic_symbol(v::AtomView) = atomic_symbol(v.system, v.index) atomic_number(v::AtomView) = atomic_number(v.system, v.index) -element(v::AtomView) = elements[atomic_symbol(v)] n_dimensions(v::AtomView) = n_dimensions(v.system) + +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) diff --git a/src/element.jl b/src/element.jl new file mode 100644 index 0000000..6a72d78 --- /dev/null +++ b/src/element.jl @@ -0,0 +1,7 @@ +import PeriodicTable + +# Valid types for atom identifiers +const AtomId = Union{Symbol,AbstractString,Integer} + +# Lookup element in PeriodicTable +element(id::AtomId) = PeriodicTable.elements[id] diff --git a/src/fast_system.jl b/src/fast_system.jl index 6120826..8cf91f2 100644 --- a/src/fast_system.jl +++ b/src/fast_system.jl @@ -41,11 +41,6 @@ function FastSystem(particles, box, boundary_conditions) atomic_number.(particles), atomic_mass.(particles)) end -function Base.show(io::IO, system::FastSystem) - print(io, "FastSystem") - show_system(io, system) -end - bounding_box(sys::FastSystem) = sys.box boundary_conditions(sys::FastSystem) = sys.boundary_conditions Base.length(sys::FastSystem) = length(sys.positions) diff --git a/src/flexible_system.jl b/src/flexible_system.jl index d80d355..fc523e8 100644 --- a/src/flexible_system.jl +++ b/src/flexible_system.jl @@ -40,11 +40,6 @@ function FlexibleSystem(system::AbstractSystem; end FlexibleSystem(;system::FlexibleSystem, kwargs...) = FlexibleSystem(system; kwargs...) -function Base.show(io::IO, system::FlexibleSystem) - print(io, "FlexibleSystem") - show_system(io, system) -end - bounding_box(sys::FlexibleSystem) = sys.box boundary_conditions(sys::FlexibleSystem) = sys.boundary_conditions species_type(sys::FlexibleSystem{D, S, L}) where {D, S, L} = S diff --git a/src/interface.jl b/src/interface.jl index cdf3c08..81d8c50 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -1,7 +1,7 @@ import Base.position export AbstractSystem -export BoundaryCondition, DirichletZero, Periodic, infinite_box +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 @@ -52,6 +52,10 @@ function species_type end """Return vector indicating whether the system is periodic along a dimension.""" periodicity(sys::AbstractSystem) = [isa(bc, Periodic) for bc in boundary_conditions(sys)] +"""Returns true if the given system is infinite""" +isinfinite(sys::AbstractSystem{D}) where {D} = bounding_box(sys) == infinite_box(D) + + """ n_dimensions(::AbstractSystem) n_dimensions(atom) @@ -126,6 +130,11 @@ atomic_mass(sys::AbstractSystem, index) = atomic_mass(sys[index]) Vector of atomic symbols in the system `sys` or the atomic symbol of a particular `species` / the `i`th species in `sys`. + +The intention is that [`atomic_number`](@ref) carries the meaning +of identifying the type of a `species` (e.g. the element for the case of an atom), whereas +[`atomic_symbol`](@ref) may return a more unique identifier. For example for a deuterium atom +this may be `:D` while `atomic_number` is still `1`. """ atomic_symbol(sys::AbstractSystem) = atomic_symbol.(sys) atomic_symbol(sys::AbstractSystem, index) = atomic_symbol(sys[index]) @@ -138,6 +147,11 @@ atomic_symbol(sys::AbstractSystem, index) = atomic_symbol(sys[index]) Vector of atomic numbers in the system `sys` or the atomic number of a particular `species` / the `i`th species in `sys`. + +The intention is that [`atomic_number`](@ref) carries the meaning +of identifying the type of a `species` (e.g. the element for the case of an atom), whereas +[`atomic_symbol`](@ref) may return a more unique identifier. For example for a deuterium atom +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]) diff --git a/src/properties.jl b/src/properties.jl index 5e859cd..d35d539 100644 --- a/src/properties.jl +++ b/src/properties.jl @@ -1,6 +1,5 @@ export chemical_formula - """ Returns the chemical formula of an AbstractSystem as a string. """ @@ -18,20 +17,3 @@ function chemical_formula(symbols::AbstractVector{Symbol}) join(sort(parts)) end chemical_formula(system) = chemical_formula(atomic_symbol(system)) - - -function show_system(io::IO, system::AbstractSystem{D}) where {D} - print(io, "($(chemical_formula(system)), ") - bc = boundary_conditions(system) - if all(isequal(bc[1]), bc) - print(io, typeof(bc[1]), ", ") - end - box = bounding_box(system) - if box != infinite_box(D) - box_str = ["[" * join(ustrip.(bvector), ", ") * "]" for bvector in box] - print(io, "box=[", join(box_str, ", "), "]u\"$(unit(box[1][1]))\"") - else - print(io, "box=infinite") - end - print(io, ")") -end diff --git a/src/show.jl b/src/show.jl new file mode 100644 index 0000000..004c4a6 --- /dev/null +++ b/src/show.jl @@ -0,0 +1,104 @@ + +""" +Suggested function to print AbstractSystem objects to screen +""" +function show_system(io::IO, system::AbstractSystem{D}) where {D} + bc = boundary_conditions(system) + + print(io, typeof(system).name.name, "($(chemical_formula(system))") + if isinfinite(system) + print(io, ", infinite") + else + perstr = [p ? "T" : "F" for p in periodicity(system)] + print(io, ", periodic = ", join(perstr, "")) + end + + if !isinfinite(system) + box_str = ["[" * join(ustrip.(bvector), ", ") * "]" + for bvector in bounding_box(system)] + bunit = unit(eltype(first(bounding_box(system)))) + print(io, ", bounding_box = [", join(box_str, ", "), "]u\"$bunit\"") + end + print(io, ")") +end +function show_system(io::IO, ::MIME"text/plain", system::AbstractSystem{D}) where {D} + bc = boundary_conditions(system) + box = bounding_box(system) + print(io, typeof(system).name.name, "($(chemical_formula(system))") + if isinfinite(system) + print(io, ", infinite") + else + perstr = [p ? "T" : "F" for p in periodicity(system)] + print(io, ", periodic = ", join(perstr, "")) + end + println(io, "):") + + extra_line = false + if !isinfinite(system) + extra_line = true + box = bounding_box(system) + bunit = unit(eltype(first(bounding_box(system)))) + for (i, bvector) in enumerate(box) + if i == 1 + @printf io " %-17s : [" "bounding_box" + else + print(io, " "^25) + end + boxstr = [(@sprintf "%8.6g" ustrip(b)) for b in bvector] + print(io, join(boxstr, " ")) + println(io, i == D ? "]u\"$bunit\"" : ";") + 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 + end + + # TODO Better would be some ascii-graphical representation of the structure + if length(system) < 20 + extra_line && println(io) + for atom in system + println(io, " ", atom) + end + end +end + +Base.show(io::IO, system::AbstractSystem) = show_system(io, system) +function Base.show(io::IO, mime::MIME"text/plain", system::AbstractSystem) + show_system(io, mime, system) +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($(atomic_symbol(at)), [", join(pos, ", "), "]u\"$posunit\"") + if ismissing(velocity(at)) || iszero(velocity(at)) + print(io, ")") + else + vel = [(@sprintf "%8.6g" ustrip(p)) for p in velocity(at)] + velunit = unit(eltype(velocity(at))) + print(io, ", [", join(vel, ", "), "]u\"$velunit\")") + end +end + +function show_atom(io::IO, ::MIME"text/plain", at) + println(io, "Atom(", atomic_symbol(at), ", atomic_number = ", atomic_number(at), + ", atomic_mass = ", atomic_mass(at), "):") + + pos = [(@sprintf "%.8g" ustrip(p)) for p in position(at)] + posunit = unit(eltype(position(at))) + @printf io " %-17s : [%s]u\"%s\"\n" "position" join(pos, ",") string(posunit) + if !ismissing(velocity(at)) && !iszero(velocity(at)) + vel = [(@sprintf "%.8g" ustrip(p)) for p in velocity(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 + end +end diff --git a/test/fast_system.jl b/test/fast_system.jl index 9736d38..61e7f1c 100644 --- a/test/fast_system.jl +++ b/test/fast_system.jl @@ -15,6 +15,7 @@ using PeriodicTable @test atomic_mass(system) == [12.011, 12.011]u"u" @test boundary_conditions(system) == bcs @test bounding_box(system) == box + @test !isinfinite(system) # Test AtomView @@ -24,5 +25,4 @@ using PeriodicTable end @test ismissing(velocity(system[1])) @test n_dimensions(system[1]) == n_dimensions(system) - @test element(system[1]) == elements[:C] end diff --git a/test/interface.jl b/test/interface.jl index 269a437..ca3c420 100644 --- a/test/interface.jl +++ b/test/interface.jl @@ -14,7 +14,6 @@ using PeriodicTable @testset "Atoms" begin @test position(atoms[1]) == [0.25, 0.25, 0.25]u"m" @test velocity(atoms[1]) == [0.0, 0.0, 0.0]u"bohr/s" - @test element(atoms[1]) == PeriodicTable.elements[:C] @test atomic_symbol(atoms[1]) == :C @test atomic_number(atoms[1]) == 6 @test atomic_mass(atoms[1]) == 12.011u"u" @@ -29,6 +28,7 @@ using PeriodicTable @test bounding_box(flexible) == [[1, 0, 0], [0, 1, 0], [0, 0, 1]]u"m" @test boundary_conditions(flexible) == [Periodic(), Periodic(), DirichletZero()] @test periodicity(flexible) == [1, 1, 0] + @test !isinfinite(flexible) @test n_dimensions(flexible) == 3 @test position(flexible) == [[0.25, 0.25, 0.25], [0.75, 0.75, 0.75]]u"m" @test position(flexible, 1) == [0.25, 0.25, 0.25]u"m" diff --git a/test/printing.jl b/test/printing.jl index 68b53c5..2fdfe32 100644 --- a/test/printing.jl +++ b/test/printing.jl @@ -4,17 +4,21 @@ using Test @testset "Printing atomic systems" begin at = Atom(:Si, zeros(3) * u"m", extradata=42) - @test repr(at) == "Atom(Si, [0.0, 0.0, 0.0]u\"m\")" - + @test repr(at) == "Atom(Si, [ 0, 0, 0]u\"m\")" + show(stdout, MIME("text/plain"), at) + atoms = [:Si => [0.0, -0.125, 0.0], :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) @test repr(flexible_system) == """ - FlexibleSystem(CSi, Periodic, box=[[10.0, 0.0, 0.0], [0.0, 5.0, 0.0], [0.0, 0.0, 7.0]]u"Å")""" + 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) fast_system = FastSystem(flexible_system) @test repr(fast_system) == """ - FastSystem(CSi, Periodic, box=[[10.0, 0.0, 0.0], [0.0, 5.0, 0.0], [0.0, 0.0, 7.0]]u"Å")""" + FastSystem(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"), fast_system) + show(stdout, MIME("text/plain"), fast_system[1]) end