feat(algo): support guessing atom/bond properties from geometry

.clang-tidy

@@ -5,7 +5,7 @@
 
 # Require clang-tidy v15
 ---
-Checks: "-*,abseil-*,boost-*,bugprone-*,clang-*,concurrency-*,cppcoreguidelines-*,google-*,misc-*,modernize-*,openmp-*,performance-*,portability-*,readability-*,-bugprone-assignment-in-if-condition,-bugprone-easily-swappable-parameters,-cppcoreguidelines-avoid-do-while,-cppcoreguidelines-avoid-c-arrays,-cppcoreguidelines-init-variables,-cppcoreguidelines-narrowing-conversions,-cppcoreguidelines-pro-*,-cppcoreguidelines-rvalue-reference-param-not-moved,-google-explicit-constructor,-misc-no-recursion,-modernize-avoid-c-arrays,-modernize-return-braced-init-list,-modernize-loop-convert,-modernize-concat-nested-namespaces,-modernize-pass-by-value,-modernize-raw-string-literal,-modernize-use-trailing-return-type,-modernize-use-auto,-modernize-use-default-member-init,-modernize-use-emplace,-modernize-use-equals-delete,-modernize-use-nodiscard,-readability-isolate-declaration,-readability-identifier-length,-readability-qualified-auto,-*-magic-numbers"
+Checks: "-*,abseil-*,boost-*,bugprone-*,clang-*,concurrency-*,cppcoreguidelines-*,google-*,misc-*,modernize-*,openmp-*,performance-*,portability-*,readability-*,-bugprone-assignment-in-if-condition,-bugprone-easily-swappable-parameters,-cppcoreguidelines-avoid-do-while,-cppcoreguidelines-avoid-c-arrays,-cppcoreguidelines-init-variables,-cppcoreguidelines-narrowing-conversions,-cppcoreguidelines-pro-*,-cppcoreguidelines-rvalue-reference-param-not-moved,-cppcoreguidelines-use-default-member-init,-google-explicit-constructor,-misc-no-recursion,-modernize-avoid-c-arrays,-modernize-return-braced-init-list,-modernize-loop-convert,-modernize-concat-nested-namespaces,-modernize-pass-by-value,-modernize-raw-string-literal,-modernize-use-trailing-return-type,-modernize-use-auto,-modernize-use-default-member-init,-modernize-use-emplace,-modernize-use-equals-delete,-modernize-use-nodiscard,-readability-isolate-declaration,-readability-identifier-length,-readability-qualified-auto,-*-magic-numbers"
 FormatStyle: "file"
 # If these flags are updated, .clangd should be updated as well
 ExtraArgsBefore:
@@ -52,6 +52,10 @@ CheckOptions:
     value: google
   - key: performance-unnecessary-value-param.IncludeStyle
     value: google
+  - key: google-readability-braces-around-statements.ShortStatementLines
+    value: 3
+  - key: readability-braces-around-statements.ShortStatementLines
+    value: 3
     # The following entries are (partailly) based on the Google C++ Style Guide.
     # Namely, we use lower_case for function names, for better interoperability
     # with the Python language.

cmake/NurikitTest.cmake

@@ -10,3 +10,27 @@ if(NOT TARGET nuri_all_test)
     set_target_properties(nuri_all_test PROPERTIES EXCLUDE_FROM_ALL OFF)
   endif()
 endif()
+
+include(GoogleTest)
+
+function(nuri_add_test file)
+  get_filename_component(test_dir ${file} DIRECTORY)
+  file(RELATIVE_PATH test_dir "${PROJECT_SOURCE_DIR}/test" "${test_dir}")
+  string(REPLACE "/" "_" test_prefix ${test_dir})
+
+  get_filename_component(test_name ${file} NAME_WE)
+
+  set(target "nuri_${test_prefix}_${test_name}")
+  add_executable("${target}" "${file}")
+  target_link_libraries("${target}" PRIVATE
+    GTest::gtest GTest::gmock GTest::gtest_main
+    absl::absl_log absl::absl_check)
+
+  if(TARGET nuri_lib)
+    target_link_libraries("${target}" PRIVATE nuri_lib)
+  endif()
+
+  gtest_discover_tests("${target}"
+    WORKING_DIRECTORY "${PROJECT_SOURCE_DIR}")
+  add_dependencies(nuri_all_test "${target}")
+endfunction()

include/nuri/algo/guess.h

@@ -0,0 +1,44 @@
+//
+// Project nurikit - Copyright 2023 SNU Compbio Lab.
+// SPDX-License-Identifier: Apache-2.0
+//
+
+#ifndef NURI_ALGO_GUESS_H_
+#define NURI_ALGO_GUESS_H_
+
+#include "nuri/core/molecule.h"
+
+namespace nuri {
+constexpr inline double kDefaultThreshold = 0.5;
+
+/**
+ * @brief Guess bonds, types of atoms, and number of hydrogens of a molecule.
+ * @param mut The mutator of the molecule to be guessed.
+ * @param conf The index of the conformation used for guessing.
+ * @param threshold The threshold for guessing bonds.
+ * @return true if the guessing is successful.
+ *
+ * This function assumes all connectivity information is missing, and all atom
+ * types and implicit hydrogen counts are incorrect. The information present
+ * in the molecule could be overwritten by this function.
+ *
+ * If connectivity information is already present and is correct, consider using
+ * guess_all_types().
+ */
+extern bool guess_bonds(MoleculeMutator &mut, int conf = 0,
+                        double threshold = kDefaultThreshold);
+
+/**
+ * @brief Guess types of atoms and bonds, and number of hydrogens of a molecule.
+ * @param mol The molecule to be guessed.
+ * @param conf The index of the conformation used for guessing.
+ * @return true if the guessing is successful.
+ *
+ * This function assumes all connectivity information is present and correct,
+ * and all atom/bond types and implicit hydrogen counts are incorrect. The
+ * information present in the molecule could be overwritten by this function.
+ */
+extern bool guess_all_types(Molecule &mol, int conf = 0);
+}  // namespace nuri
+
+#endif /* NURI_ALGO_GUESS_H_ */

include/nuri/algo/rings.h

@@ -0,0 +1,139 @@
+//
+// Project nurikit - Copyright 2023 SNU Compbio Lab.
+// SPDX-License-Identifier: Apache-2.0
+//
+#ifndef NURI_ALGO_RINGS_H_
+#define NURI_ALGO_RINGS_H_
+
+#include <memory>
+#include <utility>
+#include <vector>
+
+#include "nuri/core/molecule.h"
+
+namespace nuri {
+using Rings = std::vector<std::vector<int>>;
+
+/**
+ * @brief Find all elementary cycles in the molecular graph.
+ * @param mol A molecule.
+ * @return A pair of (all elementary cycles, success). If success is `false`,
+ *         the vector is in an unspecified state. This will fail if and only if
+ *         any atom is a member of more than 100 elementary cycles.
+ *
+ * This is based on the algorithm described in the following paper:
+ *    Hanser, Th. *et al.* *J. Chem. Inf. Comput. Sci.* **1996**, *36* (6),
+ *    1146-1152. DOI: [10.1021/ci960322f](https://doi.org/10.1021/ci960322f)
+ *
+ * The time complexity of this function is inherently exponential, but it is
+ * expected to run in a reasonable time (\f$\sim\mathcal{O}(V^2)\f$) for most
+ * molecules in practice.
+ */
+extern std::pair<Rings, bool> find_all_rings(const Molecule &mol);
+
+namespace internal {
+  struct FindRingsCommonData;
+}  // namespace internal
+
+/**
+ * @brief Wrapper class of the common routines of find_sssr() and
+ *        find_relevant_rings().
+ * @sa nuri::find_relevant_rings(), nuri::find_sssr()
+ *
+ * Formally, SSSR (smallest set of smallest rings) is a *minimum cycle basis*
+ * of the molecular graph. As discussed in many literatures, there is no unique
+ * SSSR for a given molecular graph (even for simple molecules such as
+ * 2-oxabicyclo[2.2.2]octane), and the SSSR is often counter-intuitive. For
+ * example, the SSSR of cubane (although unique, due to symmetry reasons)
+ * contains only five rings, which is not most chemists would expect.
+ *
+ * On the other hand, union of all SSSRs, sometimes called the *relevant
+ * rings* in the literatures, is unique for a given molecule, and is the "all
+ * smallest rings" of the molecule, chemically speaking. It is more appropriate
+ * for most applications than SSSR.
+ *
+ * We provide two functions along with this class to find the relevant rings and
+ * SSSR, respectively. If both are needed, it is recommended to construct this
+ * class first, and call find_relevant_rings() and find_sssr() member functions
+ * instead of calling the free functions directly.
+ *
+ * This is based on the algorithm described in the following paper:
+ *    Vismara, P. *Electron. J. Comb.* **1997**, *4* (1), R9.
+ *    DOI: [10.37236/1294](https://doi.org/10.37236/1294)
+ *
+ * Time complexity: theoretically \f$\mathcal{O}(\nu E^3)\f$, where \f$\nu =
+ * \mathcal{O}(E)\f$ is size of SSSR. For most molecules, however, this is
+ * \f$\mathcal{O}(V^3)\f$.
+ */
+class RingSetsFinder {
+public:
+  /**
+   * @brief Construct a new Rings Finder object.
+   * @param mol A molecule.
+   */
+  explicit RingSetsFinder(const Molecule &mol);
+
+  RingSetsFinder(const RingSetsFinder &) = delete;
+  RingSetsFinder &operator=(const RingSetsFinder &) = delete;
+  RingSetsFinder(RingSetsFinder &&) noexcept;
+  RingSetsFinder &operator=(RingSetsFinder &&) noexcept;
+
+  ~RingSetsFinder() noexcept;
+
+  /**
+   * @brief Find the relevant rings of the molecule.
+   * @return The relevant rings of the molecule.
+   * @sa nuri::find_relevant_rings()
+   */
+  Rings find_relevant_rings() const;
+
+  /**
+   * @brief Find the SSSR of the molecule.
+   * @return The smallest set of smallest rings (SSSR) of the molecule.
+   * @sa nuri::find_sssr()
+   * @note This function does not guarantee that the returned set is unique, nor
+   * that the result is reproducible even for the same molecule.
+   */
+  Rings find_sssr() const;
+
+private:
+  const Molecule *mol_;
+  std::unique_ptr<internal::FindRingsCommonData> data_;
+};
+
+/**
+ * @brief Find union of the all SSSRs in the molecular graph.
+ * @param mol A molecule.
+ * @return Union of the all SSSRs in the molecular graph.
+ * @sa find_sssr(), nuri::RingSetsFinder::find_relevant_rings()
+ *
+ * This is a convenience wrapper of the
+ * nuri::RingSetsFinder::find_relevant_rings() member function.
+ *
+ * @note If both relevant rings and SSSR are needed, it is recommended to use
+ * the nuri::RingSetsFinder class instead of the free functions.
+ */
+inline Rings find_relevant_rings(const Molecule &mol) {
+  return RingSetsFinder(mol).find_relevant_rings();
+}
+
+/**
+ * @brief Find a smallest set of smallest rings (SSSR) of the molecular graph.
+ * @param mol A molecule.
+ * @return *A* smallest set of smallest rings (SSSR) of the molecular graph.
+ * @sa find_relevant_rings(), nuri::RingSetsFinder::find_sssr()
+ * @note This function does not guarantee that the returned set is unique, nor
+ *       that the result is reproducible even for the same molecule.
+ *
+ * This is a convenience wrapper of the nuri::RingSetsFinder::find_sssr() member
+ * function.
+ *
+ * @note If both relevant rings and SSSR are needed, it is recommended to use
+ * the nuri::RingSetsFinder class instead of the free functions.
+ */
+inline Rings find_sssr(const Molecule &mol) {
+  return RingSetsFinder(mol).find_sssr();
+}
+}  // namespace nuri
+
+#endif /* NURI_ALGO_RINGS_H_ */

include/nuri/core/element.h

@@ -132,9 +132,44 @@ constexpr inline bool operator!=(const Isotope &lhs,
  *
  *  - Retrieved from the BODR (Blue Obelisk Data Repository):
  *    https://github.com/BlueObelisk/bodr/blob/29ce17071c71b2d4d5ee81a2a28f0407331f1624/bodr/elements/elements.xml
+ *
+ * @subsection type-state About the type and standard state data
+ *
+ * The element type and standard state data were obtained from the PubChem
+ * periodic table.
+ *
+ * The element type represents classification of the element into one of the
+ * following categories. The categories are mutually exclusive, and elements
+ * that do not fall into any of the categories are classified as metals.
+ *
+ *  - Unknown (only dummy atoms),
+ *  - Nonmetal,
+ *  - Metalloid
+ *
+ * The standard state represents the state of the element at 298.15 K and 1 atm,
+ * and is one of gas, liquid, or solid. Unknown is only used for dummy atoms.
+ *
+ * @subsubsection type-state-ref References
+ *
+ *  - National Center for Biotechnology Information. Periodic Table of Elements.
+ *    https://pubchem.ncbi.nlm.nih.gov/periodic-table. (Accessed 2023-12-05)
  */
 class Element {
 public:
+  enum class Type : std::uint8_t {
+    kUnknown,
+    kMetal,
+    kMetalloid,
+    kNonmetal,
+  };
+
+  enum class State : std::uint8_t {
+    kUnknown,
+    kSolid,
+    kLiquid,
+    kGas,
+  };
+
   Element() = delete;
   ~Element() noexcept = default;
 
@@ -198,6 +233,18 @@ class Element {
     return internal::check_flag(flags_, ElementFlags::kActinide);
   }
 
+  /**
+   * @brief Get the type of the element.
+   * @return The type of the element.
+   */
+  constexpr Type type() const noexcept { return type_; }
+
+  /**
+   * @brief Get the standard state of the element.
+   * @return The standard state of the element.
+   */
+  constexpr State state() const noexcept { return state_; }
+
   /**
    * @brief Get the IUPAC Symbol of the atom.
    * @return The IUPAC Symbol.
@@ -299,6 +346,8 @@ class Element {
   std::int16_t period_;
   std::int16_t group_;
   ElementFlags flags_;
+  Type type_;
+  State state_;
   std::string_view symbol_;
   std::string_view name_;
   double atomic_weight_;
@@ -447,6 +496,9 @@ class PeriodicTable final {
   absl::flat_hash_map<std::string_view, const Element *> symbol_to_element_;
   absl::flat_hash_map<std::string_view, const Element *> name_to_element_;
 };
+
+// NOLINTNEXTLINE(readability-identifier-naming)
+static const PeriodicTable &kPt = PeriodicTable::get();
 }  // namespace nuri
 
 #endif /* NURI_CORE_ELEMENT_H_ */