diff --git a/src/AtomsBase.jl b/src/AtomsBase.jl index cf7ab22..ea2c767 100644 --- a/src/AtomsBase.jl +++ b/src/AtomsBase.jl @@ -12,6 +12,7 @@ include("flexible_system.jl") include("atomview.jl") include("atom.jl") include("fast_system.jl") +include("bonded_system.jl") function __init__() @require AtomsView="ee286e10-dd2d-4ff2-afcb-0a3cd50c8041" begin diff --git a/src/bonded_system.jl b/src/bonded_system.jl new file mode 100644 index 0000000..eca3a47 --- /dev/null +++ b/src/bonded_system.jl @@ -0,0 +1,82 @@ +# +# Implementation of AtomsBase interface in a struct-of-arrays style with bonds. +# +export BondedSystem + +struct BondedSystem{D, L <: Unitful.Length, M <: Unitful.Mass} <: AbstractSystem{D} + bounding_box::SVector{D, SVector{D, L}} + boundary_conditions::SVector{D, BoundaryCondition} + position::Vector{SVector{D, L}} + bonds::Vector{Tuple{Integer, Integer, BondOrder}} + atomic_symbol::Vector{Symbol} + atomic_number::Vector{Int} + atomic_mass::Vector{M} +end + +# Constructor to fetch the types +function BondedSystem(box, boundary_conditions, positions, bonds, atomic_symbols, atomic_numbers, atomic_masses) + BondedSystem{length(box),eltype(eltype(positions)),eltype(atomic_masses)}( + box, boundary_conditions, positions, bonds, atomic_symbols, atomic_numbers, atomic_masses + ) +end + +# Constructor to take data from another system +function BondedSystem(system::AbstractSystem) + BondedSystem(bounding_box(system), boundary_conditions(system), position(system), bonds(system), + atomic_symbol(system), atomic_number(system), atomic_mass(system)) +end + +# Convenience constructor where we don't have to preconstruct all the static stuff... +function BondedSystem(particles, box, boundary_conditions, bonds) + D = length(box) + if !all(length.(box) .== D) + throw(ArgumentError("Box must have D vectors of length D=$D.")) + end + if length(boundary_conditions) != D + throw(ArgumentError("Boundary conditions should be of length D=$D.")) + end + if !all(n_dimensions.(particles) .== D) + throw(ArgumentError("Particles must have positions of length D=$D.")) + end + BondedSystem(box, boundary_conditions, position.(particles), bonds, atomic_symbol.(particles), + atomic_number.(particles), atomic_mass.(particles)) +end + +bounding_box(sys::BondedSystem) = sys.bounding_box +boundary_conditions(sys::BondedSystem) = sys.boundary_conditions +bonds(sys::BondedSystem) = sys.bonds + +Base.length(sys::BondedSystem) = length(sys.position) +Base.size(sys::BondedSystem) = size(sys.position) + +species_type(::FS) where {FS <: BondedSystem} = AtomView{FS} +Base.getindex(sys::BondedSystem, i::Integer) = AtomView(sys, i) + +position(s::BondedSystem) = s.position +atomic_symbol(s::BondedSystem) = s.atomic_symbol +atomic_number(s::BondedSystem) = s.atomic_number +atomic_mass(s::BondedSystem) = s.atomic_mass +velocity(::BondedSystem) = missing + +# System property access +function Base.getindex(system::BondedSystem, x::Symbol) + if x === :bounding_box + bounding_box(system) + elseif x === :boundary_conditions + boundary_conditions(system) + elseif x === :bonds + bonds(system) + else + throw(KeyError(x)) + end +end +Base.haskey(::BondedSystem, x::Symbol) = x in (:bounding_box, :boundary_conditions, :bonds) +Base.keys(::BondedSystem) = (:bounding_box, :boundary_conditions, :bonds) + +# Atom and atom property access +atomkeys(::BondedSystem) = (:position, :atomic_symbol, :atomic_number, :atomic_mass) +hasatomkey(system::BondedSystem, x::Symbol) = x in atomkeys(system) +function Base.getindex(system::BondedSystem, i::Union{Integer,AbstractVector}, x::Symbol) + getfield(system, x)[i] +end +Base.getindex(system::BondedSystem, ::Colon, x::Symbol) = getfield(system, x) diff --git a/src/interface.jl b/src/interface.jl index 2d4462f..8a6739f 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -3,8 +3,10 @@ import PeriodicTable export AbstractSystem export BoundaryCondition, DirichletZero, Periodic, infinite_box, isinfinite -export bounding_box, boundary_conditions, periodicity, n_dimensions, species_type + +export bounding_box, boundary_conditions, periodicity, n_dimensions, species_type, bonds export position, velocity, element, element_symbol, atomic_mass, atomic_number, atomic_symbol + export atomkeys, hasatomkey # @@ -51,6 +53,23 @@ Return the type used to represent a species or atom. """ function species_type end +""" +The possible bond orders that can be used in the `bonds` function +""" +@enum BondOrder begin + single + double + triple +end + +""" + bonds(::AbstractSystem) + +Returns a Vector{Tuple{Integer, Integer, BondOrder}} where each entry stores the unique indices + i and j of the atoms participating in the bond and the order of the bond (e.g. [(13,21,5)]) +""" +function bonds end + """Return vector indicating whether the system is periodic along a dimension.""" periodicity(sys::AbstractSystem) = [isa(bc, Periodic) for bc in boundary_conditions(sys)]