Dev/atomsbase integration (#18)

* set up the integration for AtomsBase

* add atomsbase functions

* example update before I formalize it

* fix dead link

* add atomsbase integration example doc

* fix dependency

* remove unnecessary file

* update tests

* update more tests

* add docs to ccread func

* add docs to the rest of the functions

* update docs

* Update tests and atomsbase funcs

* update as per comment

Author: vcanogil

Project.toml

@@ -4,14 +4,15 @@ authors = ["cclib development team"]
 version = "0.2.0"
 
 [deps]
+AtomsBase = "a963bdd2-2df7-4f54-a1ee-49d51e6be12a"
 CondaPkg = "992eb4ea-22a4-4c89-a5bb-47a3300528ab"
 Logging = "56ddb016-857b-54e1-b83d-db4d58db5568"
-PyCall = "438e738f-606a-5dbb-bf0a-cddfbfd45ab0"
 PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d"
+Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
+UnitfulAtomic = "a7773ee8-282e-5fa2-be4e-bd808c38a91a"
 
 [compat]
 CondaPkg = "0.2"
-PyCall = "1"
 PythonCall = "0.9"
 julia = "1"
 

docs/make.jl

@@ -16,6 +16,7 @@ makedocs(;
     pages=[
         "Home" => "index.md",
         "How to parse files" => "io.md",
+        "AtomsBase Integration" => "atombase.md"
     ],
 )
 

docs/src/atombase.md

@@ -0,0 +1,74 @@
+# AtomsBase Integration
+
+Cclib.jl allows loading atom information into [AtomsBase.jl](https://github.com/JuliaMolSim/AtomsBase.jl) objects using `get_atom_objects` function:
+
+```Julia
+# Input file can be found in the in the repo under "test" folder
+julia> using Cclib
+julia> using AtomsBase
+julia> mol = ccread("./Trp_polar.fchk")
+julia> atoms = get_atom_objects(mol)
+julia> atoms[1]
+
+Atom(N, atomic_number = 7, atomic_mass = 14.007 u):
+    position          : [-0.069982688,0.33219872,0.28212832]u"Å"
+```
+
+Like before, all the stored information for each atom is still accessible:
+
+```Julia
+julia> keys(atoms[1])
+(:position, :velocity, :atomic_symbol, :atomic_number, :atomic_mass)
+```
+# Interoperability With Other Tools
+
+We can use data loaded with Cclib.jl to perform calculations using other libraries that use AtomsBase.jl, such as [InteratomicPotentials.jl](https://github.com/cesmix-mit/InteratomicPotentials.jl)
+
+First, we load the data and define the system:
+```Julia
+julia> using Cclib
+julia> using AtomsBase
+julia> using Unitful
+julia> using UnitfulAtomic
+julia> using InteratomicPotentials
+
+julia> atoms = get_atom_objects("./Trp_polar.fchk")
+julia> boundaryconditions = [Periodic(), Periodic(), Periodic()]
+julia> box = [[10.0, 0.0, 0.0], [0.0, 10.0, 0.0], [0.0, 0.0, 10.0]]u"Å"
+julia> system = FlexibleSystem(atoms, box, boundaryconditions)
+
+FlexibleSystem(C₁₁H₁₂N₂O₂, periodic = TTT):
+    bounding_box      : [      10        0        0;
+                                0       10        0;
+                                0        0       10]u"Å"
+
+          .------------------------.
+         /|  H                     |
+        / |                        |
+       /  |                        |
+      /   |    H  C      C  H      |
+     /    H                        |
+    *  C  C     C  C H C    H      |
+    |   H |  C    HN    H          |
+    |     |                        |
+    |     |   C  C                 |
+    |     |    C                   |
+    |     |      O                 |
+    |     .----------------------H-.
+    |    /                        /
+    |   /                        /
+    | H/O                       /
+    | /          H             /
+    |/                       N/
+    *------------------------*
+```
+We can now perform calculate various properties of the system:
+```Julia
+julia> ϵ = 1.0 * 1u"eV"
+julia> σ = 0.25 * 1u"Å"
+julia> rcutoff  = 2.25 * 1u"Å"
+julia> lj = LennardJones(ϵ, σ, rcutoff, [:N, :C, :O, :H])
+julia> potential_energy(system, lj).
+-0.00039418327614247183 Eₕ
+```
+Refer to InteratomicPotentials.jl documentation for more details.

docs/src/io.md

@@ -2,11 +2,12 @@
 
 This page outlines how to access information store in computational chemistry output files.
 
-### Supported formats
-Properties that can be parsed and supported file formats can be found [here](https://cclib.github.io/data.html").
+# Supported formats
+Properties that can be parsed and supported file formats can be found [here](https://cclib.github.io/data.html#details-of-current-implementation").
 
-### How to read files
+# How to read files
 ```Julia
+# Input file can be found in the in the repo under "test" folder
 julia> using Cclib
 julia> mol = ccread("uracil_two.xyz")
 ```

src/Cclib.jl

@@ -1,10 +1,14 @@
 module Cclib
 using PythonCall
 using CondaPkg
+using AtomsBase
+using Unitful
+using UnitfulAtomic
 using Logging
 
-include("functions.jl")
 include("config.jl")
+include("functions.jl")
+include("ab_integration.jl")
 
 const cclib = Ref{Py}()
 

src/ab_integration.jl

@@ -0,0 +1,41 @@
+#
+# Functions for AtomsBase integration https://github.com/JuliaMolSim/AtomsBase.jl/tree/master
+#
+export get_atom_objects
+
+"""
+    get_atom_objects(mol::Dict)
+
+Load atom data contained in a cclib dictionary into a list
+of AtomsBase atom objects.
+
+# Arguments
+- `mol::Dict`: Dictionary that was returned by `ccread` function
+
+# Returns
+A list of AtomsBase atom objects
+"""
+function get_atom_objects(mol::Dict)
+    atoms = [
+             Atom(mol["atomnos"][i],
+             vec(mol["atomcoords"][:, i, :])u"Å")
+             for i in 1:mol["natom"]
+            ]
+    return atoms
+end
+
+"""
+    get_atom_objects(mol::String)
+
+Load atom data contained in a cclib-supported file format
+into a list of AtomsBase atom objects.
+
+# Arguments
+- `mol::String`: File supported by `ccread` function
+
+# Returns
+A list of AtomsBase atom objects
+"""
+function get_atom_objects(mol::String)
+    return get_atom_objects(ccread(mol))
+end

src/functions.jl

@@ -1,3 +1,7 @@
+#
+# General functions for parsing chemical outputs
+#
+
 export ccread
 
 function pyccread(file)
@@ -16,7 +20,31 @@ function get_data(file)
     return datadict
 end
 
-function ccread(file)
+"""
+    ccread(file::String)
+
+Read in the file of supported format and store the data
+it contains in a dictionary.
+Properties read in are accessible as dictionary keys.
+
+# Arguments
+- `file::String`: Cclib-supported file format
+
+# Returns
+A dictionary containing stored information
+
+# Example
+```Julia
+julia> mol = ccread("uracil_two.xyz")
+julia> keys(mol)
+KeySet for a Dict{String, Any} with 4 entries. Keys:
+  "atomcoords"
+  "natom"
+  "atomnos"
+  "metadata"
+```
+"""
+function ccread(file::String)
     try
         if !isfile(file)
             @error "$(file) is not file"

test/runtests.jl

@@ -1,21 +1,33 @@
 using Cclib
 using PythonCall
 using Test
+using AtomsBase
+
 
 
 @testset "Cclib.jl" begin
+    test_file = "./data/uracil_two.xyz"
+    test_ccdata = ccread(test_file)
 
     # Check that input is read correctly
-    @test ccread("./data/uracil_two.xyz")["natom"] == 12
-    @test ccread("./data/uracil_two.xyz")["atomnos"] == [7, 6, 1, 8, 7, 6, 8, 6, 6, 1, 1, 1]
-    @test ccread("./data/uracil_two.xyz")["atomcoords"] == [5.435 4.091 6.109 3.647 3.253 3.636 2.79 5.065 5.906 2.243 5.448 7.003;
-                                                                 5.42432 4.08846 6.05396 3.68445 3.24195 3.62987 2.83542 5.0864 5.91213 2.2534 5.42415 6.9915;;;
-                                                                 3.248 3.266 3.347 3.394 3.118 2.955 2.842 2.953 3.101 3.138 2.833 3.119;
-                                                                 3.24732 3.26356 3.354 3.39863 3.11963 2.96568 2.83961 2.95832 3.09904 3.12801 2.83814 3.10206;;;
-                                                                 -0.916 -0.635 -0.144 0.495 -1.717 -3.039 -3.923 -3.239 -2.202 -1.514 -4.258 -2.345;
-                                                                 -0.9432 -0.64127 -0.15563 0.50851 -1.70419 -3.0054 -3.9297 -3.25044 -2.21101 -1.51035 -4.27067 -2.32366]
+    @test test_ccdata["natom"] == 12
+    @test test_ccdata["atomnos"] == [7, 6, 1, 8, 7, 6, 8, 6, 6, 1, 1, 1]
+    @test test_ccdata["atomcoords"] == [5.435 4.091 6.109 3.647 3.253 3.636 2.79 5.065 5.906 2.243 5.448 7.003;
+                                        5.42432 4.08846 6.05396 3.68445 3.24195 3.62987 2.83542 5.0864 5.91213 2.2534 5.42415 6.9915;;;
+                                        3.248 3.266 3.347 3.394 3.118 2.955 2.842 2.953 3.101 3.138 2.833 3.119;
+                                        3.24732 3.26356 3.354 3.39863 3.11963 2.96568 2.83961 2.95832 3.09904 3.12801 2.83814 3.10206;;;
+                                        -0.916 -0.635 -0.144 0.495 -1.717 -3.039 -3.923 -3.239 -2.202 -1.514 -4.258 -2.345;
+                                        -0.9432 -0.64127 -0.15563 0.50851 -1.70419 -3.0054 -3.9297 -3.25044 -2.21101 -1.51035 -4.27067 -2.32366]
 
     # Check that it handles unsupported files correctly
     @test isnothing(ccread("./data/invalid_file.txt"))
     @test isnothing(ccread("./data/"))
+
+    # Check AtomsBase Integration
+    test_atom_objects = get_atom_objects(test_file)
+    @test atomic_number(test_atom_objects[1]) == 7
+    @test atomic_symbol(test_atom_objects[1]) == :N
+    @test atomic_number(test_atom_objects[6]) == 6
+    @test atomic_symbol(test_atom_objects[6]) == :C
+
 end