crystal-structure

C++ library that helps retrieve a finite domain from the periodic atomic structure (motif) or the underlying periodic structure (lattice) of a crystal

Requiremements

You may want to have the following libraries in your cmake file as dependencies:

- [Gemmi](https://github.com/project-gemmi/gemmi)
- [Eigen3](https://gitlab.com/libeigen/eigen)

Installation

cmake -B build -S . -L -DCMAKE_INSTALL_PREFIX=install -DBUILD_CONFIG=Release
cmake --build build --target install --config Release

Usage

Periodic point cloud of a lattice

Set the cell parameters (length of a, length of b, length of c, alpha, beta, gamma)

#include <crystal/structure.h>

Crystal::Lattice lattice;
int extension = 3;
lattice.setCellParameters({10, 10, 10, 90, 90, 90});

The unit cell can be Niggli reduced with

lattice.updateParametersToReducedMode();

otherwise call the following

lattice.updateParametersToOrginalMode();

Then, extend the lattice points and get the entire point cloud in the following way

lattice.spanTheLattice(extension);
std::vector<Eigen::Vector3d> pointcloud = lattice.getPointCloud();

Periodic point cloud of a motif

Set the cell parameters and add fractional coordinates and chemical info in Crystal::Motif. The order of fractionals below must correspond to the sites order.

#include <gemmi/smcif.hpp>
#include <crystal/structure.h>

Crystal::Motif motif;
int extension = 3;
std::vector<Crystal::Chem::Atom> atoms{
    {"C1", "C", gemmi::Element("C"), Eigen::Vector3d{0.25,0.25,0.25}},
    {"C2", "C", gemmi::Element("C"), Eigen::Vector3d{0.5,0.5,0.5}},
    {"C3", "C", gemmi::Element("C"), Eigen::Vector3d{0.75,0.75,0.75}},
    {"C4", "C", gemmi::Element("C"), Eigen::Vector3d{1.0,1.0,1.0}}
};

motif.setCellParameters({10, 10, 10, 90, 90, 90});
motif.setCellMotif(atoms);

The unit cell can be Niggli reduced with

motif.updateMotifToReducedMode();

otherwise call the following

motif.updateMotifToOriginalMode();

Then, extend the motif points and get the entire point cloud

motif.repeatTheMotif(extension);
std::vector<Eigen::Vector3d> pointcloud = motif.getPointCloud();

Graph of the atomic structure

Create the atomic structure from the Sites and Bonds.

#include <gemmi/smcif.hpp>
#include <crystal/chem.h>

UndirectedGraph compound_graph;
std::map<std::string, unsigned int> mapAtomLabelToIndex;
std::vector<Crystal::Chem::Atom> atoms{
    {"C1", "C", gemmi::Element("C"), Eigen::Vector3d{0.25,0.25,0.25}},
    {"C2", "C", gemmi::Element("C"), Eigen::Vector3d{0.5,0.5,0.5}},
    {"C3", "C", gemmi::Element("C"), Eigen::Vector3d{0.75,0.75,0.75}},
    {"C4", "C", gemmi::Element("C"), Eigen::Vector3d{1.0,1.0,1.0}}
};

Retrieve the graph from it through defined Bonds

std::vector<gemmi::Restraints::Bond> bonds{
    { {0, "C1"}, {0, "C2"}, gemmi::BondType::Single, false, 1.54, 0.02 },
    { {0, "C3"}, {0, "C4"}, gemmi::BondType::Single, false, 1.54, 0.02 }
};
compound_graph = makeUndirectedGraphFromAtomicStructure(atoms, bonds, mapAtomLabelToIndex);

Note: After the last call, mapAtomLabelToIndex is populated with atom labels (e.g. C1, C2) and their indexes in the atoms array.

Get each molecule as a vector of indexes. Each index refers to the atom index in atoms vector or in motif.cell_motif if you are using the class Crystal::Motif.

std::vector<std::vector<unsigned int>> molecules_ids = getMoleculesFromUndirectedGraphCompound(compound_graph);
// molecules_ids[0] first molecule -> {0, 1}
// molecules_ids[1] second molecule -> {2, 3}

Tests

cmake -B build -S . -L -DCMAKE_INSTALL_PREFIX=install -DBUILD_TESTING=ON -DBUILD_CONFIG=Debug
cmake --build build --target install --config Debug
cd build
ctest --build-config Debug --build-target install