➖ Remove unit library due to increasing simulation runtimes.

.gitmodules

@@ -13,6 +13,3 @@
 [submodule "libs/tinyxml2"]
 	path = libs/tinyxml2
 	url = https://github.com/leethomason/tinyxml2.git
-[submodule "libs/units"]
-	path = libs/units
-	url = https://github.com/nholthaus/units.git

docs/technology/simulation.rst

@@ -30,3 +30,15 @@ distributions of the SiDBs. Charge distribution surfaces are returned by the SiD
    :members:
 .. doxygenclass:: fiction::charge_distribution_surface< Lyt, false >
    :members:
+
+
+Physical Constants
+------------------
+
+**Header:** ``fiction/technology/physical_constants.hpp``
+
+.. doxygenvariable:: EPSILON
+.. doxygenvariable:: ELEMENTARY_CHARGE
+.. doxygenvariable:: K
+.. doxygenvariable:: POP_STABILITY_ERR
+.. doxygenvariable:: PI

docs/utils/utils.rst

@@ -124,24 +124,6 @@ Array Utils
 .. doxygenfunction:: fiction::convert_array
 .. doxygenfunction:: fiction::convert_array_of_arrays
 
-Unit Utils
------------
-
-**Header:** ``fiction/utils/unit_utils.hpp``
-
-.. doxygentypedef:: coulomb_constant_unit
-.. doxygenvariable:: POP_STABILITY_ERR
-
-There are several suffix operators exposed in this file for convenient use of units, namely
-* ``_angstrom``
-* ``_e``
-* ``_eV``
-* ``_K``
-* ``_nm``
-* ``_s``
-* ``_V``
-
-This enables the use of units as follows: e.g., ``40_K``, ``40.1_nm``.
 
 STL Extensions
 --------------

include/fiction/algorithms/path_finding/distance.hpp

@@ -92,16 +92,16 @@ template <typename Lyt, typename Dist = uint64_t>
                                                           static_cast<Dist>(std::numeric_limits<uint32_t>::max());
 }
 /**
- * Computes the distance between two SiDB cells in nanometers.
+ * Computes the distance between two SiDB cells in nanometers (unit: nm).
  *
  * @tparam Lyt SiDB cell-level layout type.
  * @tparam Dist Floating-point type for the distance.
  * @param c1 The first cell.
  * @param c2 The second cell.
- * @return The distance between the two cells in nanometers.
+ * @return The distance between the two cells in nanometers (unit: nm).
  */
 template <typename Lyt>
-[[nodiscard]] constexpr units::length::nanometer_t
+[[nodiscard]] constexpr double
 sidb_nanometer_distance([[maybe_unused]] const Lyt& lyt, const coordinate<Lyt>& source, const coordinate<Lyt>& target,
                         const sidb_simulation_parameters& sp = sidb_simulation_parameters{}) noexcept
 {
@@ -115,7 +115,7 @@ sidb_nanometer_distance([[maybe_unused]] const Lyt& lyt, const coordinate<Lyt>&
     const auto x = pos_c1.first - pos_c2.first;
     const auto y = pos_c1.second - pos_c2.second;
 
-    return units::math::hypot(x, y);
+    return std::hypot(x, y);
 }
 
 // NOLINTBEGIN(*-special-member-functions): virtual destructor is prudent

include/fiction/algorithms/simulation/sidb/calculate_energy_and_state_type.hpp

@@ -22,7 +22,7 @@ namespace fiction
  *  Data type to collect electrostatic potential energies (in eV) of charge distributions with corresponding state types
  * (i.e., true = transparent, false = erroneous).
  */
-using sidb_energy_and_state_type = std::vector<std::pair<units::energy::electron_volt_t, bool>>;
+using sidb_energy_and_state_type = std::vector<std::pair<double, bool>>;
 
 /**
  * This function takes in an SiDB energy distribution. For each charge distribution, the state type is determined (i.e.

include/fiction/algorithms/simulation/sidb/critical_temperature.hpp

@@ -19,12 +19,10 @@
 #include "fiction/utils/hash.hpp"
 #include "fiction/utils/math_utils.hpp"
 #include "fiction/utils/truth_table_utils.hpp"
-#include "fiction/utils/units_utils.hpp"
 
 #include <fmt/format.h>
 #include <kitty/bit_operations.hpp>
 #include <kitty/dynamic_truth_table.hpp>
-#include <units.h>
 
 #include <algorithm>
 #include <cassert>
@@ -97,9 +95,9 @@ struct critical_temperature_params
      */
     double confidence_level{0.99};
     /**
-     * Simulation stops at max_temperature (~ 126 °C by default).
+     * Simulation stops at max_temperature (~ 126 °C by default) (unit: K).
      */
-    units::temperature::kelvin_t max_temperature{400_K};
+    double max_temperature{400};
     /**
      * Truth table of the given gate (if layout is simulated in `gate-based` mode).
      */
@@ -123,32 +121,31 @@ struct critical_temperature_stats
      */
     std::string algorithm_name{};
     /**
-     * Critical Temperature of the given layout.
+     * Critical Temperature of the given layout (unit: K).
      */
-    units::temperature::kelvin_t critical_temperature{0_K};
+    double critical_temperature{0};
     /**
      * Number of physically valid charge configurations.
      */
     uint64_t num_valid_lyt{};
     /**
-     * Energy difference between the ground state and the first (erroneous) excited state.
+     * Energy difference between the ground state and the first (erroneous) excited state (unit: eV).
      */
-    units::energy::electron_volt_t energy_between_ground_state_and_first_erroneous =
-        units::energy::electron_volt_t(std::numeric_limits<double>::infinity());
+    double energy_between_ground_state_and_first_erroneous = std::numeric_limits<double>::infinity();
     /**
      * Prints the simulation results to the given output stream.
      *
      * @param out Output stream.
      */
     void report(std::ostream& out = std::cout) const
     {
-        out << fmt::format("Critical Temperature  = {:.2f} K\n", critical_temperature.value());
+        out << fmt::format("Critical Temperature  = {:.2f} K\n", critical_temperature);
 
         if (num_valid_lyt != 0)
         {
             out << fmt::format("'# of physically valid charge configurations': {} | Energy between ground state and "
                                "first erroneous: {}\n",
-                               num_valid_lyt, energy_between_ground_state_and_first_erroneous.value());
+                               num_valid_lyt, energy_between_ground_state_and_first_erroneous);
         }
         else
         {
@@ -287,8 +284,8 @@ class critical_temperature_impl
 
             else
             {
-                temperature_stats.critical_temperature = 0.0_K;  // If no ground state fulfills the logic, the Critical
-                                                                 // Temperature is zero. May be worth it to change µ_.
+                temperature_stats.critical_temperature = 0.0;  // If no ground state fulfills the logic, the Critical
+                                                               // Temperature is zero. May be worth it to change µ_.
             }
         }
 
@@ -329,10 +326,10 @@ class critical_temperature_impl
                 (first_excited_state_energy - ground_state_energy) * 1000;
         }
 
-        std::vector<units::temperature::kelvin_t> temp_values{};
-        temp_values.reserve(static_cast<uint64_t>(parameter.max_temperature.value() * 100));
+        std::vector<double> temp_values{};  // unit: K
+        temp_values.reserve(static_cast<uint64_t>(parameter.max_temperature * 100));
 
-        for (uint64_t i = 1; i <= static_cast<uint64_t>(parameter.max_temperature.value() * 100); i++)
+        for (uint64_t i = 1; i <= static_cast<uint64_t>(parameter.max_temperature * 100); i++)
         {
             temp_values.emplace_back(static_cast<double>(i) / 100.0);
         }
@@ -349,7 +346,7 @@ class critical_temperature_impl
                 break;
             }
 
-            if (units::math::abs(temp - parameter.max_temperature) < 0.001_K)
+            if (std::abs(temp - parameter.max_temperature) < 0.001)
             {
                 // Maximal temperature is stored as the Critical Temperature.
                 temperature_stats.critical_temperature = parameter.max_temperature;
@@ -366,11 +363,11 @@ class critical_temperature_impl
      *
      * @param energy_and_state_type All energies of all physically valid charge distributions with the corresponding
      * state type (i.e. transparent, erroneous).
-     * @param min_energy Minimal energy of all physically valid charge distributions of a given layout.
+     * @param min_energy Minimal energy of all physically valid charge distributions of a given layout (unit: eV).
      * @return State type (i.e. transparent, erroneous) of the ground state is returned.
      */
-    bool energy_between_ground_state_and_first_erroneous(const sidb_energy_and_state_type&    energy_and_state_type,
-                                                         const units::energy::electron_volt_t min_energy)
+    bool energy_between_ground_state_and_first_erroneous(const sidb_energy_and_state_type& energy_and_state_type,
+                                                         const double                      min_energy)
     {
         bool ground_state_is_transparent = false;
         for (const auto& [energy, state_type] : energy_and_state_type)
@@ -401,10 +398,10 @@ class critical_temperature_impl
     void determine_critical_temperature(const sidb_energy_and_state_type& energy_state_type)
     {
         // Vector with temperature values from 0.01 to max_temperature * 100 K in 0.01 K steps is generated.
-        std::vector<units::temperature::kelvin_t> temp_values{};
-        temp_values.reserve(static_cast<uint64_t>(parameter.max_temperature.value() * 100));
+        std::vector<double> temp_values{};
+        temp_values.reserve(static_cast<uint64_t>(parameter.max_temperature * 100));
 
-        for (uint64_t i = 1; i <= static_cast<uint64_t>(parameter.max_temperature.value() * 100); i++)
+        for (uint64_t i = 1; i <= static_cast<uint64_t>(parameter.max_temperature * 100); i++)
         {
             temp_values.emplace_back(static_cast<double>(i) / 100.0);
         }
@@ -420,7 +417,7 @@ class critical_temperature_impl
                 break;
             }
 
-            if (units::math::abs(temp - parameter.max_temperature) < 0.001_K)
+            if (std::abs(temp - parameter.max_temperature) < 0.001)
             {
                 // Maximal temperature is stored as Critical Temperature.
                 temperature_stats.critical_temperature = parameter.max_temperature;

include/fiction/algorithms/simulation/sidb/energy_distribution.hpp

@@ -7,7 +7,6 @@
 
 #include "fiction/technology/charge_distribution_surface.hpp"
 #include "fiction/utils/math_utils.hpp"
-#include "fiction/utils/units_utils.hpp"
 
 #include <cmath>
 #include <cstdint>
@@ -18,10 +17,10 @@ namespace fiction
 {
 
 /**
- *  Data type to collect electrostatic potential energies of charge distributions with corresponding degeneracy (i.e.
- * how often a certain energy value occurs).
+ * Data type to collect electrostatic potential energies (unit: eV) of charge distributions with corresponding
+ * degeneracy (i.e. how often a certain energy value occurs).
  */
-using sidb_energy_distribution = std::map<units::energy::electron_volt_t, uint64_t>;
+using sidb_energy_distribution = std::map<double, uint64_t>;  // unit: (eV, unitless)
 
 /**
  * This function takes in a vector of charge_distribution_surface objects and returns a map containing the system energy
@@ -39,7 +38,7 @@ energy_distribution(const std::vector<charge_distribution_surface<Lyt>>& input_v
     static_assert(is_cell_level_layout_v<Lyt>, "Lyt is not a cell-level layout");
     static_assert(has_sidb_technology_v<Lyt>, "Lyt is not an SiDB layout");
 
-    std::map<units::energy::electron_volt_t, uint64_t> distribution{};
+    std::map<double, uint64_t> distribution{};  // unit: (eV, unitless)
 
     for (const auto& lyt : input_vec)
     {

include/fiction/algorithms/simulation/sidb/is_ground_state.hpp

@@ -12,8 +12,6 @@
 #include "fiction/traits.hpp"
 #include "fiction/utils/math_utils.hpp"
 
-#include <units.h>
-
 #include <cmath>
 
 namespace fiction
@@ -44,7 +42,7 @@ template <typename Lyt>
     const auto min_energy_exact  = minimum_energy(exhaustive_results.charge_distributions);
     const auto min_energy_new_ap = minimum_energy(quicksim_results.charge_distributions);
 
-    return round_to_n_decimal_places(units::math::abs(min_energy_exact - min_energy_new_ap), 6).value() == 0;
+    return round_to_n_decimal_places(std::abs(min_energy_exact - min_energy_new_ap), 6) == 0;
 }
 
 }  // namespace fiction

include/fiction/algorithms/simulation/sidb/minimum_energy.hpp

@@ -19,18 +19,16 @@ namespace fiction
  *
  * @tparam Lyt Cell-level layout type.
  * @param charge_lyts Vector of charge_distribution_surface objects.
- * @return Value of the minimum energy found in the input vector.
+ * @return Value of the minimum energy found in the input vector (unit: eV).
  */
 template <typename Lyt>
-[[nodiscard]] units::energy::electron_volt_t
-minimum_energy(const std::vector<charge_distribution_surface<Lyt>>& charge_lyts) noexcept
+[[nodiscard]] double minimum_energy(const std::vector<charge_distribution_surface<Lyt>>& charge_lyts) noexcept
 {
     static_assert(is_cell_level_layout_v<Lyt>, "Lyt is not a cell-level layout");
     static_assert(has_sidb_technology_v<Lyt>, "Lyt is not an SiDB layout");
 
-    return units::energy::electron_volt_t(
-        std::accumulate(charge_lyts.cbegin(), charge_lyts.cend(), std::numeric_limits<double>::max(),
-                        [](const double a, const auto& lyt) { return std::min(a, lyt.get_system_energy().value()); }));
+    return std::accumulate(charge_lyts.cbegin(), charge_lyts.cend(), std::numeric_limits<double>::max(),
+                           [](const double a, const auto& lyt) { return std::min(a, lyt.get_system_energy()); });
 }
 
 }  // namespace fiction

include/fiction/algorithms/simulation/sidb/occupation_probability_of_excited_states.hpp

@@ -7,7 +7,6 @@
 
 #include "fiction/algorithms/simulation/sidb/calculate_energy_and_state_type.hpp"
 #include "fiction/utils/math_utils.hpp"
-#include "fiction/utils/units_utils.hpp"
 
 #include <algorithm>
 #include <cassert>
@@ -27,20 +26,20 @@ namespace fiction
  *
  * @param energy_and_state_type This contains the energies of all possible charge distributions together with the
  * information if the charge distribution (state) is transparent or erroneous.
- * @param temperature System temperature to assume.
+ * @param temperature System temperature to assume (unit: K).
  * @return The occupation probability of all erroneous states is returned.
  */
-[[nodiscard]] inline double occupation_probability_gate_based(const sidb_energy_and_state_type&   energy_and_state_type,
-                                                              const units::temperature::kelvin_t& temperature) noexcept
+[[nodiscard]] inline double occupation_probability_gate_based(const sidb_energy_and_state_type& energy_and_state_type,
+                                                              const double                      temperature) noexcept
 {
-    assert((temperature > 0.0_K) && "temperature should be slightly above 0 K");
+    assert((temperature > 0.0) && "temperature should be slightly above 0 K");
 
     if (energy_and_state_type.empty())
     {
         return 0.0;
     }
 
-    auto min_energy = units::energy::electron_volt_t(std::numeric_limits<double>::infinity());
+    auto min_energy = std::numeric_limits<double>::infinity();  // unit: eV
 
     // Determine the minimal energy.
     const auto [energy, state_type] = *std::min_element(energy_and_state_type.cbegin(), energy_and_state_type.cend(),
@@ -51,7 +50,7 @@ namespace fiction
     const double partition_function =
         std::accumulate(energy_and_state_type.cbegin(), energy_and_state_type.cend(), 0.0,
                         [&](const double sum, const auto& it)
-                        { return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature).value()); });
+                        { return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature)); });
 
     // All Boltzmann factors of the erroneous states are summed.
     double p = 0;
@@ -61,7 +60,7 @@ namespace fiction
     {
         if (!state_transparent_erroneous)
         {
-            p += std::exp(-((energies - min_energy) * 12'000 / temperature).value());
+            p += std::exp(-((energies - min_energy) * 12'000 / temperature));
         }
     }
 
@@ -72,31 +71,30 @@ namespace fiction
  * This function computes the occupation probability of excited states (charge distributions with energy higher than the
  * ground state) at a given temperature.
  *
- * @param energy_distribution This contains the energies of all possible charge distributions with the degeneracy.
- * @param temperature System temperature to assume.
+ * @param energy_distribution This contains the energies in eV of all possible charge distributions with the degeneracy.
+ * @param temperature System temperature to assume (unit: K).
  * @return The total occupation probability of all excited states is returned.
  */
-[[nodiscard]] inline double
-occupation_probability_non_gate_based(const sidb_energy_distribution&     energy_distribution,
-                                      const units::temperature::kelvin_t& temperature) noexcept
+[[nodiscard]] inline double occupation_probability_non_gate_based(const sidb_energy_distribution& energy_distribution,
+                                                                  const double                    temperature) noexcept
 {
-    assert((temperature > 0.0_K) && "Temperature should be slightly above 0 K");
+    assert((temperature > 0.0) && "Temperature should be slightly above 0 K");
 
     if (energy_distribution.empty())
     {
         return 0.0;
     }
 
-    auto min_energy = units::energy::electron_volt_t(std::numeric_limits<double>::infinity());
+    auto min_energy = std::numeric_limits<double>::infinity();
 
     const auto& [energy, degeneracy] = *(energy_distribution.begin());
-    min_energy                       = energy;
+    min_energy                       = energy;  // unit: eV
 
     // The partition function is obtained by summing up all the Boltzmann factors.
     const double partition_function =
         std::accumulate(energy_distribution.cbegin(), energy_distribution.cend(), 0.0,
                         [&](const double sum, const auto& it)
-                        { return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature).value()); });
+                        { return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature)); });
 
     // All Boltzmann factors of the excited states are summed.
     const double p =
@@ -107,7 +105,7 @@ occupation_probability_non_gate_based(const sidb_energy_distribution&     energy
                             // possible rounding errors and for comparability with the min_energy.
                             if (round_to_n_decimal_places(it.first, 6) != round_to_n_decimal_places(min_energy, 6))
                             {
-                                return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature).value());
+                                return sum + std::exp(-((it.first - min_energy) * 12'000 / temperature));
                             }
                             return sum;
                         });