Merge pull request #365 from SCM-NV/qmflows-0.12

BLD: Bump the minimum qmflows version >=0.12.0

.github/workflows/pythonapp.yml

@@ -167,7 +167,7 @@ jobs:
           python-version: "3.10"
 
       - name: Install linters
-        run: pip install pydocstyle pycodestyle mypy qmflows "numpy>=1.21" types-pyyaml types-pkg_resources types-setuptools
+        run: pip install pydocstyle pycodestyle mypy qmflows "numpy>=1.21" types-pyyaml types-setuptools "pyparsing>=3.0.8"
 
       - name: Python info
         run: |

conftest.py

@@ -5,6 +5,7 @@
 from collections.abc import Generator
 
 import pytest
+from nanoqm._logger import logger, stdout_handler
 
 
 @pytest.fixture(autouse=True, scope="function")
@@ -17,3 +18,12 @@ def cleunup_files() -> Generator[None, None, None]:
         os.remove("input_parameters.yml")
     if os.path.isdir("overlaps"):
         shutil.rmtree("overlaps")
+
+
+@pytest.fixture(autouse=True, scope="session")
+def prepare_logger() -> Generator[None, None, None]:
+    """Remove the logging output to stdout while running tests."""
+    assert stdout_handler in logger.handlers
+    logger.removeHandler(stdout_handler)
+    yield None
+    logger.addHandler(stdout_handler)

docs/conf.py

@@ -210,7 +210,7 @@
 autodoc_mock_imports = [
     'h5py',
     'qmflows',
-    'compute_integrals',
+    'nanoqm.compute_integrals',
     'mendeleev',
     'scm',
     'noodles',

libint/compute_integrals.cc

@@ -82,20 +82,33 @@ template <typename Lambda> void parallel_do(Lambda &lambda) {
  * \brief Set the number of thread to use
  */
 void set_nthread() {
-
+#if defined(_OPENMP)
   using libint2::nthreads;
   nthreads = std::thread::hardware_concurrency();
-
-#if defined(_OPENMP)
   omp_set_num_threads(nthreads);
 #endif
-  std::cout << "Will scale over " << nthreads
-#if defined(_OPENMP)
-            << " OpenMP"
-#else
-            << " C++11"
-#endif
-            << " threads" << std::endl;
+}
+
+/**
+ * \brief Get the number of threads and the thread type as a 2-tuple
+ */
+int get_thread_count() {
+  using libint2::nthreads;
+  return std::thread::hardware_concurrency();
+}
+
+/**
+ * \brief Get the type of threads
+ */
+string get_thread_type() {
+  string ret;
+
+  #if defined(_OPENMP)
+    ret = "OpenMP";
+  #else
+    ret = "C++11";
+  #endif
+  return ret;
 }
 
 /**
@@ -287,13 +300,9 @@ std::vector<Matrix> compute_multipoles(
 CP2K_Basis_Atom read_basis_from_hdf5(const string &path_file,
                                      const string &symbol,
                                      const string &basis) {
-  std::vector<std::vector<double>> coefficients;
-  std::vector<double> exponents;
-  std::vector<int64_t> format;
-
-  libint2::svector<double> small_exp;
+  libint2::svector<libint2::svector<double>> small_exp;
   libint2::svector<libint2::svector<double>> small_coef;
-  libint2::svector<CP2K_Contractions> small_fmt;
+  libint2::svector<libint2::svector<CP2K_Contractions>> small_fmt;
 
   try {
     // Open an existing HDF5 File
@@ -310,6 +319,10 @@ CP2K_Basis_Atom read_basis_from_hdf5(const string &path_file,
     // only a single set of exponents, but there are exception such
     // as BASIS_ADMM_MOLOPT
     for (const auto &name : dset_names) {
+      std::vector<std::vector<double>> coefficients;
+      std::vector<double> exponents;
+      std::vector<int64_t> format;
+
       const string path_coefficients = root + "/" + name + "/coefficients";
       const string path_exponents = root + "/" + name + "/exponents";
 
@@ -323,14 +336,16 @@ CP2K_Basis_Atom read_basis_from_hdf5(const string &path_file,
       dataset_es.read(exponents);
       attr.read(format);
 
-      // Move data to small vectors and keep extending them as iteration
-      // over `dset_names` continues
-      std::move(exponents.begin(), exponents.end(), std::back_inserter(small_exp));
+      // Move data to small vectors and keep appending or extending them as
+      // iteration over `dset_names` continues
+      libint2::svector<double> small_exp_1d;
+      std::move(exponents.begin(), exponents.end(), std::back_inserter(small_exp_1d));
+      small_exp.push_back(small_exp_1d);
 
       for (const auto &v : coefficients) {
-        libint2::svector<double> small;
-        std::move(v.begin(), v.end(), std::back_inserter(small));
-        small_coef.push_back(small);
+        libint2::svector<double> small_coef_1d;
+        std::move(v.begin(), v.end(), std::back_inserter(small_coef_1d));
+        small_coef.push_back(small_coef_1d);
       }
 
       // The CP2K basis format is defined by a vector of integers, for each atom.
@@ -349,15 +364,12 @@ CP2K_Basis_Atom read_basis_from_hdf5(const string &path_file,
       // Note: Elements 4 and onwards define the number of contracted for each
       // angular momentum quantum number (all prior elements are disgarded).
       int l, i;
+      libint2::svector<CP2K_Contractions> small_fmt_1d;
       for (i=4, l=format[1]; i != static_cast<int>(format.size()); i++, l++) {
         int count = format[i];
-        small_fmt.push_back(CP2K_Contractions{l, count});
+        small_fmt_1d.push_back({l, count});
       }
-
-      // Clear the temp vectors for the next iteration cycle
-      coefficients.clear();
-      exponents.clear();
-      format.clear();
+      small_fmt.push_back(small_fmt_1d);
     }
   } catch (HighFive::Exception &err) {
     // catch and print any HDF5 error
@@ -405,18 +417,23 @@ create_map_symbols_basis(const string &path_hdf5,
  */
 libint2::svector<Shell> create_shells_for_atom(const CP2K_Basis_Atom &data,
                                                const Atom &atom) {
-  libint2::svector<CP2K_Contractions> basis_format = data.basis_format;
   libint2::svector<Shell> shells;
+  libint2::svector<double> exponents;
+  libint2::svector<CP2K_Contractions> basis_format;
 
   int acc = 0;
-  for (auto contractions : basis_format) {
-    for (int i = 0; i < contractions.count; i++) {
-      shells.push_back({
-        data.exponents,
-        {{contractions.l, true, data.coefficients[acc]}},  // compute integrals in sphericals
-        {{atom.x, atom.y, atom.z}}  // Atomic Coordinates
-      });
-      acc += 1;
+  for (int i = 0; i != static_cast<int>(data.exponents.size()); i++) {
+    exponents = data.exponents[i];
+    basis_format = data.basis_format[i];
+    for (auto contractions : basis_format) {
+      for (int j = 0; j < contractions.count; j++) {
+        shells.push_back({
+          exponents,
+          {{contractions.l, true, data.coefficients[acc]}},  // compute integrals in sphericals
+          {{atom.x, atom.y, atom.z}}  // Atomic Coordinates
+        });
+        acc += 1;
+      }
     }
   }
   return shells;
@@ -551,4 +568,10 @@ PYBIND11_MODULE(compute_integrals, m) {
 
   m.def("compute_integrals_multipole", &compute_integrals_multipole,
         py::return_value_policy::reference_internal);
+
+  m.def("get_thread_count", &get_thread_count,
+        py::return_value_policy::reference_internal);
+
+  m.def("get_thread_type", &get_thread_type,
+        py::return_value_policy::reference_internal);
 }

libint/include/namd.hpp

@@ -57,8 +57,8 @@ struct CP2K_Basis_Atom {
   // Contains the basis specificationf for a given atom
   std::string symbol;
   libint2::svector<libint2::svector<double>> coefficients;
-  libint2::svector<double> exponents;
-  libint2::svector<CP2K_Contractions> basis_format;
+  libint2::svector<libint2::svector<double>> exponents;
+  libint2::svector<libint2::svector<CP2K_Contractions>> basis_format;
 };
 
 // Map from atomic_number to symbol

nanoqm/__init__.py

@@ -3,6 +3,8 @@
 chemical (QM) packages."""
 
 from ._version import __version__ as __version__
+from ._version_info import version_info as version_info
+from ._logger import logger as logger
 
 from .analysis import (
     autocorrelate, dephasing, convolute, func_conv, gauss_function,

nanoqm/_logger.py

@@ -0,0 +1,79 @@
+"""The Nano-QMFlows logger."""
+
+from __future__ import annotations
+
+import os
+import sys
+import types
+import logging
+import contextlib
+from typing import ClassVar
+
+from qmflows.type_hints import PathLike
+
+__all__ = ["logger", "stdout_handler", "EnableFileHandler"]
+
+#: The Nano-QMFlows logger.
+logger = logging.getLogger("nanoqm")
+logger.setLevel(logging.DEBUG)
+
+qmflows_logger = logging.getLogger("qmflows")
+noodles_logger = logging.getLogger("noodles")
+noodles_logger.setLevel(logging.WARNING)
+
+#: The Nano-QMFlows stdout handler.
+stdout_handler = logging.StreamHandler(stream=sys.stdout)
+stdout_handler.setLevel(logging.DEBUG)
+stdout_handler.setFormatter(logging.Formatter(
+    fmt='[%(asctime)s] %(levelname)s: %(message)s',
+    datefmt='%H:%M:%S',
+))
+logger.addHandler(stdout_handler)
+
+
+class EnableFileHandler(contextlib.ContextDecorator):
+    """Add a file handler to the noodles, qmflows and nanoqm loggers.
+
+    Attributes
+    ----------
+    handler : logging.FileHandler
+        The relevant titular handler.
+
+    """
+
+    __slots__ = ("handler",)
+
+    LOGGERS: ClassVar = (logger, qmflows_logger, noodles_logger)
+
+    def __init__(self, path: PathLike) -> None:
+        """Initialize the context manager.
+
+        Parameters
+        ----------
+        path : path-like object
+            Path to the log file.
+
+        """
+        self.handler = logging.FileHandler(os.fsdecode(path))
+        self.handler.setLevel(logging.DEBUG)
+        self.handler.setFormatter(logging.Formatter(
+            fmt='%(asctime)s---%(levelname)s\n%(message)s\n',
+            datefmt='%H:%M:%S',
+        ))
+
+    def __enter__(self) -> None:
+        """Add the file handler."""
+        for logger in self.LOGGERS:
+            if self.handler not in logger.handlers:
+                logger.addHandler(self.handler)
+
+    def __exit__(
+        self,
+        exc_type: type[BaseException] | None,
+        exc_value: BaseException | None,
+        exc_traceback: types.TracebackType | None,
+    ) -> None:
+        """Remove the file handler."""
+        for logger in self.LOGGERS:
+            if self.handler in logger.handlers:
+                logger.removeHandler(self.handler)

nanoqm/_version_info.py

@@ -0,0 +1,11 @@
+"""The Nano-QMFlows version tuple."""
+
+from nanoutils import VersionInfo
+from packaging.version import Version
+
+from ._version import __version__
+
+__all__ = ["version_info"]
+
+VERSION = Version(__version__)
+version_info = VersionInfo._make(VERSION.release[:3])

nanoqm/common.py

@@ -32,11 +32,11 @@
 
 import os
 import json
-from itertools import chain, repeat
+from itertools import chain
+from pathlib import Path
 from typing import (Any, Dict, Iterable, List, Mapping, NamedTuple, Tuple,
                     Sequence, overload, TypeVar, TYPE_CHECKING, Iterator)
 
-import pkg_resources as pkg
 import h5py
 import mendeleev
 import numpy as np
@@ -46,15 +46,16 @@
 from scm.plams import Atom, Molecule
 
 from qmflows.yaml_utils import UniqueSafeLoader
+from . import __path__ as nanoqm_path
 
 if TYPE_CHECKING:
     import numpy.typing as npt
 
 _T = TypeVar("_T")
 
-_path_valence_electrons = pkg.resource_filename(
-    "nanoqm", "basis/valence_electrons.json")
-_path_aux_fit = pkg.resource_filename("nanoqm", "basis/aux_fit.json")
+
+_path_valence_electrons = Path(nanoqm_path[0]) / "basis" / "valence_electrons.json"
+_path_aux_fit = Path(nanoqm_path[0]) / "basis" / "aux_fit.json"
 
 with open(_path_valence_electrons, 'r') as f1, open(_path_aux_fit, 'r') as f2:
     valence_electrons: "dict[str, int]" = json.load(f1)

libint/compute_integrals.pyi → nanoqm/compute_integrals.pyi

@@ -14,3 +14,7 @@ def compute_integrals_multipole(
     __basis_name: str,
     __multipole: str,
 ) -> npt.NDArray[np.float64]: ...
+
+def get_thread_count() -> int: ...
+
+def get_thread_type() -> str: ...

nanoqm/integrals/multipole_matrices.py

@@ -18,26 +18,24 @@
 
 from __future__ import annotations
 
-import logging
 import os
 import uuid
 from os.path import join
 from pathlib import Path
 from typing import List, Union, TYPE_CHECKING
 
 import numpy as np
-from compute_integrals import compute_integrals_multipole
 from qmflows.common import AtomXYZ
 
+from .. import logger
 from ..common import (DictConfig, is_data_in_hdf5, retrieve_hdf5_data,
                       store_arrays_in_hdf5, tuplesXYZ_to_plams)
+from ..compute_integrals import compute_integrals_multipole, get_thread_count, get_thread_type
 
 if TYPE_CHECKING:
     from numpy.typing import NDArray
     from numpy import float64 as f8
 
-logger = logging.getLogger(__name__)
-
 
 def get_multipole_matrix(config: DictConfig, inp: DictConfig, multipole: str) -> NDArray[f8]:
     """Retrieve the `multipole` number `i` from the trajectory. Otherwise compute it.
@@ -114,6 +112,9 @@ def compute_matrix_multipole(
 
     # name of the basis set
     basis_name = config["cp2k_general_settings"]["basis"]
+    thread_count = get_thread_count()
+    thread_type = get_thread_type()
+    logger.info(f"Will scale over {thread_count} {thread_type} threads")
 
     if multipole == 'overlap':
         matrix_multipole = compute_integrals_multipole(
@@ -129,7 +130,6 @@ def compute_matrix_multipole(
 
     elif multipole == 'quadrupole':
         # The tensor contains the overlap + {xx, xy, xz, yy, yz, zz} quadrupole matrices
-        print("super_matrix: ", path, path_hdf5, basis_name, multipole)
         super_matrix = compute_integrals_multipole(
             path, path_hdf5, basis_name, multipole)
         dim = super_matrix.shape[1]