clease - npm Package Compare versions

Comparing version 1.1.0 to 1.2.0

cxx/src/ce_updater.cpp

+#include "ce_updater.hpp"
+
+#include <algorithm>
+#include <cassert>
+#include <iostream>
+#include <iterator>
+#include <sstream>
+#include <stdexcept>
+
+#include "additional_tools.hpp"
+#include "atomic_numbers.hpp"
+#include "basis_function.hpp"
+#include "cluster_name.hpp"
+#include "config.hpp"
+
+using namespace std;
+
+CEUpdater::CEUpdater(){};
+CEUpdater::~CEUpdater(){};
+
+void CEUpdater::init(PyObject *py_atoms, PyObject *settings, PyObject *corrFunc, PyObject *pyeci,
+                     PyObject *cluster_list) {
+    atoms = py_atoms;
+    if (settings == nullptr) {
+        throw invalid_argument("Settings object is nullptr!");
+    }
+
+#ifdef PRINT_DEBUG
+    cout << "Getting symbols from settings object\n";
+#endif
+    PyObject *py_ignore_bck = get_attr(settings, "ignore_background_atoms");
+    ignore_background_indices = PyObject_IsTrue(py_ignore_bck);
+    Py_DECREF(py_ignore_bck);
+
+    unsigned int n_atoms = PyObject_Length(atoms);
+    if (n_atoms < 0) {
+        throw invalid_argument("Could not retrieve the length of the atoms object!");
+    }
+
+    // Read the atomic symbols
+    std::vector<std::string> symbols = get_symbols_from_atoms(py_atoms);
+
+    trans_symm_group.resize(n_atoms);
+    set<string> unique_symbols;
+
+    // Extract unique symbols from settings
+    PyObject *py_unique_symb = get_attr(settings, "unique_elements");
+    for (unsigned int i = 0; i < list_size(py_unique_symb); i++) {
+        unique_symbols.insert(py2string(PyList_GetItem(py_unique_symb, i)));
+    }
+    Py_DECREF(py_unique_symb);
+
+    insert_in_set(symbols, unique_symbols);
+    symbols_with_id = std::make_unique<Symbols>(symbols, unique_symbols);
+
+    // Build read the translational sites
+    PyObject *py_trans_symm_group = get_attr(settings, "index_by_sublattice");
+
+    if (py_trans_symm_group == nullptr) {
+        throw invalid_argument("index_by_sublattice is nullptr!");
+    }
+
+#ifdef PRINT_DEBUG
+    cout << "Reading background indices\n";
+#endif
+    // Read the backgound indices from settings
+    PyObject *bkg_indx = get_attr(settings, "background_indices");
+    read_background_indices(bkg_indx);
+    Py_DECREF(bkg_indx);
+
+    count_non_bkg_sites();
+    build_trans_symm_group(py_trans_symm_group);
+    Py_DECREF(py_trans_symm_group);
+
+#ifdef PRINT_DEBUG
+    cout << "Getting cluster names from atoms object\n";
+#endif
+
+    // Read cluster names
+    create_cname_with_dec(corrFunc);
+
+#ifdef PRINT_DEBUG
+    cout << "Cluster names with decoration number created...\n";
+#endif
+
+    PyObject *py_num_elements = get_attr(settings, "num_unique_elements");
+
+    if (py_num_elements == nullptr) {
+        throw invalid_argument("num_unique_elements is nullptr!");
+    }
+    int num_bfs = py2int(py_num_elements) - 1;
+    Py_DECREF(py_num_elements);
+
+    if (cluster_list == nullptr) {
+        throw invalid_argument("cluster_list is nullptr!");
+    }
+    // unsigned int num_trans_symm = list_size(cluster_info);
+    unsigned int num_clusters = PySequence_Size(cluster_list);
+
+#ifdef PRINT_DEBUG
+    cout << "Parsing cluster list...\n";
+#endif
+
+    PyObject *py_no_si = get_attr(cluster_list, "assume_no_self_interactions");
+    assume_no_self_interactions = PyObject_IsTrue(py_no_si);
+    Py_DECREF(py_no_si);
+
+#ifdef PRINT_DEBUG
+    std::cout << "Assuming no self-interaction?: " << assume_no_self_interactions << std::endl;
+#endif
+
+    for (unsigned int i = 0; i < num_clusters; i++) {
+        PyObject *py_cluster = PySequence_GetItem(cluster_list, i);
+
+        Cluster new_clst(py_cluster);
+        PyObject *py_cluster_name = get_attr(py_cluster, "name");
+        string cluster_name = py2string(py_cluster_name);
+        Py_DECREF(py_cluster_name);
+        Py_DECREF(py_cluster);
+
+        new_clst.construct_equivalent_deco(num_bfs);
+        clusters.append(new_clst);
+
+        int norm_fact = new_clst.get().size() * trans_symm_group_count[new_clst.symm_group];
+        if (normalisation_factor.find(cluster_name) == normalisation_factor.end()) {
+            normalisation_factor[cluster_name] = norm_fact;
+        } else {
+            normalisation_factor[cluster_name] += norm_fact;
+        }
+    }
+#ifdef PRINT_DEBUG
+    cout << "Finished reading cluster_info\n";
+#endif
+
+#ifdef PRINT_DEBUG
+    cout << "Reading basis functions from settings object\n";
+#endif
+
+    PyObject *bfs = get_attr(settings, "basis_functions");
+    if (bfs == NULL) {
+        status = Status_t::INIT_FAILED;
+        return;
+    }
+
+    // Reading basis functions from python object
+    PyObject *key;
+    PyObject *value;
+    unsigned int n_bfs = list_size(bfs);
+    bf_list basis_func_raw;
+    for (unsigned int i = 0; i < n_bfs; i++) {
+        Py_ssize_t pos = 0;
+        map<string, double> new_entry;
+        PyObject *bf_dict = PyList_GetItem(bfs, i);
+        while (PyDict_Next(bf_dict, &pos, &key, &value)) {
+            new_entry[py2string(key)] = PyFloat_AS_DOUBLE(value);
+        }
+        basis_func_raw.push_back(new_entry);
+    }
+
+    this->basis_functions = std::make_unique<BasisFunction>(basis_func_raw, *symbols_with_id);
+
+#ifdef PRINT_DEBUG
+    cout << "Reading translation matrix from settings\n";
+#endif
+    // Retrieve the TransMatrix object
+    PyObject *trans_mat_obj = get_attr(settings, "trans_matrix");
+    if (trans_mat_obj == NULL) {
+        status = Status_t::INIT_FAILED;
+        return;
+    }
+
+    // Get the internal trans_matrix object from within the TransMatrix object
+    PyObject *trans_mat_orig = get_attr(trans_mat_obj, "trans_matrix");
+
+    read_trans_matrix(trans_mat_orig);
+    Py_DECREF(trans_mat_obj);
+    Py_DECREF(trans_mat_orig);
+
+    // Read the ECIs, and parse the names.
+    this->set_eci(pyeci);
+#ifdef PRINT_DEBUG
+    cout << "Parsing correlation function\n";
+#endif
+
+    vector<string> flattened_cnames;
+    flattened_cf_names(flattened_cnames);
+
+    history = std::make_unique<CFHistoryTracker>(eci.get_names());
+    history->insert(corrFunc, nullptr);
+
+    // Store the singlets names
+    for (unsigned int i = 0; i < flattened_cnames.size(); i++) {
+        std::string name = flattened_cnames[i];
+        // Fetch the pre-parsed version of the name.
+        const ParsedName parsed = this->m_parsed_names[i];
+        if (parsed.size == 1) {
+            singlets.push_back(name);
+        }
+    }
+
+    status = Status_t::READY;
+    clear_history();
+#ifdef PRINT_DEBUG
+    cout << "CEUpdater initialized sucessfully!\n";
+#endif
+
+    // Verify that the ECIs given corresponds to a correlation function
+    if (!all_eci_corresponds_to_cf()) {
+        throw invalid_argument("All ECIs does not correspond to a correlation function!");
+    }
+}
+
+void CEUpdater::parse_eci_names() {
+    std::size_t num = eci.size();
+    this->m_parsed_names.clear();
+    this->m_parsed_names.reserve(num);
+
+    for (unsigned int i = 0; i < num; i++) {
+        std::string name = eci.name(i);
+        ClusterName c_name = ClusterName(name);
+        ParsedName parsed = c_name.get_parsed();
+        this->m_parsed_names.emplace_back(parsed);
+    }
+}
+
+double CEUpdater::get_energy() {
+    double energy = 0.0;
+    cf &corr_func = history->get_current();
+    energy = eci.dot(corr_func);
+    return energy * symbols_with_id->size();
+}
+
+double CEUpdater::spin_product_one_atom(int ref_indx, const Cluster &cluster,
+                                        const vector<int> &dec, int ref_id) const {
+    double sp = 0.0;
+
+    // Note: No duplication factor is included here, since this method is used for calculating the
+    // CF from scratch (which will account for self-interactions with no issue), and not updating
+    // (which must account for self-interactions via the duplication factor).
+    // i.e. spin_product_one_atom_delta must accounts for the self-interaction.
+
+    const vector<vector<int>> &indx_list = cluster.get();
+    unsigned int num_indx = indx_list.size();
+    unsigned int n_memb = indx_list[0].size();
+
+    // Cache the relevant row from the trans matrix.
+    int *trans_matrix_row = trans_matrix.get_row(ref_indx);
+
+    for (unsigned int i = 0; i < num_indx; i++) {
+        double sp_temp = 1.0;
+
+        // Use pointer arithmetics in the inner most loop
+        const int *indices = &indx_list[i][0];
+
+        for (unsigned int j = 0; j < n_memb; j++) {
+            int trans_index = trans_matrix.lookup_in_row(trans_matrix_row, indices[j]);
+            int id = (trans_index == ref_indx) ? ref_id : symbols_with_id->id(trans_index);
+            sp_temp *= basis_functions->get(dec[j], id);
+        }
+        sp += sp_temp;
+    }
+    return sp;
+}
+
+int CEUpdater::get_original_index(int ref_indx) const {
+    int *trans_matrix_row = trans_matrix.get_row(ref_indx);
+    int *allowed_lu = trans_matrix.get_allowed_lookup_values();
+    for (unsigned int j = 0; j < trans_matrix.get_num_non_zero(); j++) {
+        int col = allowed_lu[j];
+        int indx = trans_matrix.lookup_in_row(trans_matrix_row, col);
+        if (indx == ref_indx) {
+            return col;
+        }
+    }
+    std::stringstream err;
+    err << "Did not locate original index for ref index: " << ref_indx;
+    throw std::runtime_error(err.str());
+}
+
+double CEUpdater::spin_product_one_atom_delta_no_si(const SpinProductCache &sp_cache,
+                                                    const Cluster &cluster,
+                                                    const deco_t &deco) const {
+    /* Note: This function assumes no self-interaction within a cluster. */
+
+    // Figure out how many times we need to iterate
+    unsigned int num_indx = cluster.get_num_figures();  // Outer loop count
+    // Assume 1 ref site in a figure, so we iterate 1 less
+    unsigned int n_non_ref = cluster.get_size() - 1;  // Inner loop count
+
+    int *tm_row = sp_cache.trans_matrix_row;
+
+    /* There are n_non_ref sites per each ref site, so the non_ref_site_ptr
+    iterates faster than the ref_site_ptr. Figures are placed contiguously
+    in a 1D array.*/
+    const ClusterSite *non_ref_site_ptr = &cluster.get_non_ref_sites()[0];
+    const int *ref_site_ptr = &cluster.get_ref_cluster_sites()[0];
+
+    // Keep track of the change in spin-product
+    double sp_delta = 0.0;
+
+    // Iterate each figure in the cluster. 1 reference site is assumed per cluster
+    for (unsigned int i = 0; i < num_indx; i++, ++ref_site_ptr) {
+        /* Calculate the spin product for both new and the old (ref)
+        The constant term to the spin product from the sites which didn't change.*/
+        const int dec_ref = deco[*ref_site_ptr];
+        double new_bf = basis_functions->get(dec_ref, sp_cache.new_symb_id);
+        double old_bf = basis_functions->get(dec_ref, sp_cache.old_symb_id);
+        double sp_change = new_bf - old_bf;
+
+        /* Iterate the remaining non-reference sites, as we already took care of the reference
+        site (assuming no self-interaction)*/
+        for (unsigned int j = 0; j < n_non_ref; j++, ++non_ref_site_ptr) {
+            const ClusterSite &site = *non_ref_site_ptr;
+            const int dec_j = deco[site.cluster_index];
+
+            const int trans_index = trans_matrix.lookup_in_row(tm_row, site.lattice_index);
+            sp_change *= basis_functions->get(dec_j, symbols_with_id->id(trans_index));
+        }
+
+        sp_delta += sp_change;
+    }
+    return sp_delta;
+}
+
+double CEUpdater::spin_product_one_atom_delta(const SpinProductCache &sp_cache,
+                                              const Cluster &cluster, const deco_t &deco) const {
+    // Keep track of the change in spin-product
+    double sp_delta = 0.0;
+
+    // List of figures in the cluster
+    const vector<vector<int>> &indx_list = cluster.get();
+    // Account for the self-interaction, in case we updated 2 sites with 1 change
+    const std::vector<double> &dup_factors = cluster.get_duplication_factors();
+    unsigned int num_indx = indx_list.size();
+    unsigned int n_memb = indx_list[0].size();
+
+    int *tm_row = sp_cache.trans_matrix_row;
+
+    // Iterate each site in the cluster
+    for (unsigned int i = 0; i < num_indx; i++) {
+        // Use pointer arithmetics in the inner most loop
+        const int *indices = &indx_list[i][0];
+
+        // Calculate the spin product for both new and the old (ref)
+        double sp_temp_new = 1.0, sp_temp_ref = 1.0;
+        // The constant term to the spin product from the sites which didn't change.
+        double sp_const = dup_factors[i];
+
+        for (unsigned int j = 0; j < n_memb; j++) {
+            const int site_index = indices[j];
+            const int dec_j = deco[j];
+            if (site_index == sp_cache.original_index) {
+                // This site is the reference index.
+                // Look up the BF values directly for the new and old symbols
+                sp_temp_new *= basis_functions->get(dec_j, sp_cache.new_symb_id);
+                sp_temp_ref *= basis_functions->get(dec_j, sp_cache.old_symb_id);
+            } else {
+                // Look up the symbol of the non-reference site, which hasn't changed.
+                const int trans_index = trans_matrix.lookup_in_row(tm_row, site_index);
+                sp_const *= basis_functions->get(dec_j, symbols_with_id->id(trans_index));
+            }
+        }
+        // The change in spin-product is the difference in SP between the site(s) which
+        // changed, multiplied by the constant SP from the other un-changed sites (since
+        // these contribute equally before and after the change).
+        sp_delta += (sp_temp_new - sp_temp_ref) * sp_const;
+    }
+    return sp_delta;
+}
+
+void CEUpdater::update_cf(PyObject *single_change) {
+    SymbolChange symb_change = SymbolChange(single_change);
+    update_cf(symb_change);
+}
+
+void CEUpdater::py_changes2_symb_changes(PyObject *all_changes,
+                                         vector<SymbolChange> &symb_changes) {
+    unsigned int size = list_size(all_changes);
+    for (unsigned int i = 0; i < size; i++) {
+        SymbolChange symb_change = SymbolChange(PyList_GetItem(all_changes, i));
+        symb_changes.push_back(symb_change);
+    }
+}
+
+SpinProductCache CEUpdater::build_sp_cache(const SymbolChange &symb_change) const {
+    int ref_indx = symb_change.indx;
+    // Look up the untranslated index of the reference index.
+    int orig_indx = this->get_original_index(ref_indx);
+    // Cache the relevant row from the trans matrix.
+    int *trans_matrix_row = this->trans_matrix.get_row(ref_indx);
+
+    unsigned int old_symb_id = symbols_with_id->get_symbol_id(symb_change.old_symb);
+    unsigned int new_symb_id = symbols_with_id->get_symbol_id(symb_change.new_symb);
+
+    SpinProductCache sp_cache = {ref_indx, orig_indx, trans_matrix_row, new_symb_id, old_symb_id};
+    return sp_cache;
+}
+
+cf &CEUpdater::get_next_cf(SymbolChange &symb_change) {
+    SymbolChange *symb_change_track;
+    cf *next_cf_ptr = nullptr;
+    history->get_next(&next_cf_ptr, &symb_change_track);
+    cf &next_cf = *next_cf_ptr;
+
+    symb_change_track->indx = symb_change.indx;
+    symb_change_track->old_symb = symb_change.old_symb;
+    symb_change_track->new_symb = symb_change.new_symb;
+    symb_change_track->track_indx = symb_change.track_indx;
+    return next_cf;
+}
+
+void CEUpdater::update_cf(SymbolChange &symb_change) {
+    if (symb_change.old_symb == symb_change.new_symb) {
+        return;
+    }
+
+    if (is_background_index[symb_change.indx]) {
+        throw runtime_error("Attempting to move a background atom!");
+    }
+
+    cf &current_cf = history->get_current();
+    cf &next_cf = get_next_cf(symb_change);
+
+    symbols_with_id->set_symbol(symb_change.indx, symb_change.new_symb);
+
+    /* The following prepares a range of properties which will be
+     useful for all of the clusters, so we don't compute more
+     than we have to inside the main ECI loop */
+    SpinProductCache sp_cache = this->build_sp_cache(symb_change);
+
+    if (atoms != nullptr) {
+        set_symbol_in_atoms(atoms, symb_change.indx, symb_change.new_symb);
+    }
+
+    int symm = trans_symm_group[symb_change.indx];
+    const std::vector<ClusterCache> &clusters_cache = m_cluster_by_symm[symm];
+
+    // Loop over all ECIs
+    // As work load for different clusters are different due to a different
+    // multiplicity factor, we need to apply a dynamic schedule
+#ifdef HAS_OMP
+    bool is_par = this->cf_update_num_threads > 1;
+#pragma omp parallel for if (is_par) num_threads(this->cf_update_num_threads) schedule(dynamic)
+#endif
+    for (unsigned int i = 0; i < eci.size(); i++) {
+        // The pre-parsed version of the cluster name.
+        const ParsedName &parsed = this->m_parsed_names[i];
+
+        // 0-body
+        if (parsed.size == 0) {
+            next_cf[i] = current_cf[i];
+            continue;
+        }
+
+        // Singlet
+        if (parsed.size == 1) {
+            unsigned int dec = parsed.dec_num;
+            double new_bf = basis_functions->get(dec, sp_cache.new_symb_id);
+            double old_bf = basis_functions->get(dec, sp_cache.old_symb_id);
+            next_cf[i] = current_cf[i] + (new_bf - old_bf) / num_non_bkg_sites;
+            continue;
+        }
+
+        // n-body
+        const ClusterCache &cluster_cache = clusters_cache[i];
+        const Cluster *cluster_ptr = cluster_cache.cluster_ptr;
+
+        if (cluster_ptr == nullptr) {
+            // This cluster was not present in the symmetry group.
+            next_cf[i] = current_cf[i];
+            continue;
+        }
+        // The cluster is in the symmetry group, so calculate the spin product
+        // change for this cluster.
+        const Cluster &cluster = *cluster_ptr;
+        const equiv_deco_t &equiv_deco = *cluster_cache.equiv_deco_ptr;
+
+        double delta_sp = 0.0;
+        // Calculate the change (delta) in spin product
+        for (const deco_t &deco : equiv_deco) {
+            if (this->assume_no_self_interactions) {
+                // Faster version for large cells with no self interaction
+                delta_sp += spin_product_one_atom_delta_no_si(sp_cache, cluster, deco);
+            } else {
+                // Safe fall-back version
+                delta_sp += spin_product_one_atom_delta(sp_cache, cluster, deco);
+            }
+        }
+
+        delta_sp *= cluster_cache.normalization;
+        next_cf[i] = current_cf[i] + delta_sp;
+    }
+}
+
+void CEUpdater::undo_changes() {
+    unsigned int buf_size = history->history_size();
+    undo_changes(buf_size - 1);
+}
+
+void CEUpdater::undo_changes(int num_steps) {
+    int buf_size = history->history_size();
+
+    if (num_steps > buf_size - 1) {
+        throw invalid_argument("Can't reset history beyond the buffer size!");
+    }
+
+    SymbolChange *last_changes;
+    for (int i = 0; i < num_steps; i++) {
+        history->pop(&last_changes);
+        symbols_with_id->set_symbol(last_changes->indx, last_changes->old_symb);
+
+        if (atoms != nullptr) {
+            set_symbol_in_atoms(atoms, last_changes->indx, last_changes->old_symb);
+        }
+    }
+}
+
+double CEUpdater::calculate(PyObject *system_changes) {
+    unsigned int size = list_size(system_changes);
+    if (size == 0) {
+        return get_energy();
+    } else if (size == 1) {
+        for (unsigned int i = 0; i < size; i++) {
+            update_cf(PyList_GetItem(system_changes, i));
+        }
+        return get_energy();
+    }
+
+    if (size % 2 == 0) {
+        bool sequence_arbitrary_moves = false;
+        vector<swap_move> sequence;
+        for (unsigned int i = 0; i < size / 2; i++) {
+            swap_move changes;
+            changes[0] = SymbolChange(PyList_GetItem(system_changes, 2 * i));
+            changes[1] = SymbolChange(PyList_GetItem(system_changes, 2 * i + 1));
+
+            if (!is_swap_move(changes)) {
+                sequence_arbitrary_moves = true;
+                break;
+            }
+            sequence.push_back(changes);
+        }
+
+        if (!sequence_arbitrary_moves) {
+            return calculate(sequence);
+        }
+    }
+
+    // Last option is that this is a sequence of arbitrary moves
+    vector<SymbolChange> changes(size);
+    for (unsigned int i = 0; i < size; i++) {
+        changes[i] = SymbolChange(PyList_GetItem(system_changes, i));
+    }
+    return calculate(changes);
+}
+
+double CEUpdater::calculate(swap_move &system_changes) {
+    if (symbols_with_id->id(system_changes[0].indx) ==
+        symbols_with_id->id(system_changes[1].indx)) {
+        cout << system_changes[0] << endl;
+        cout << system_changes[1] << endl;
+        throw runtime_error(
+            "This version of the calculate function assumes that the provided update is swapping "
+            "two atoms\n");
+    }
+
+    if (symbols_with_id->get_symbol(system_changes[0].indx) != system_changes[0].old_symb) {
+        throw runtime_error("The atom position tracker does not match the current state\n");
+    } else if (symbols_with_id->get_symbol(system_changes[1].indx) != system_changes[1].old_symb) {
+        throw runtime_error("The atom position tracker does not match the current state\n");
+    }
+
+    // Update correlation function
+    update_cf(system_changes[0]);
+    update_cf(system_changes[1]);
+
+    return get_energy();
+}
+
+void CEUpdater::clear_history() {
+    history->clear();
+}
+
+void CEUpdater::flattened_cf_names(vector<string> &flattened) {
+    flattened = eci.get_names();
+
+    // Sort the cluster names for consistency
+    sort(flattened.begin(), flattened.end());
+}
+
+PyObject *CEUpdater::get_cf() {
+    PyObject *cf_dict = PyDict_New();
+    cf &corrfunc = history->get_current();
+
+    for (unsigned int i = 0; i < corrfunc.size(); i++) {
+        PyObject *pyvalue = PyFloat_FromDouble(corrfunc[i]);
+        PyDict_SetItemString(cf_dict, corrfunc.name(i).c_str(), pyvalue);
+        Py_DECREF(pyvalue);
+    }
+    return cf_dict;
+}
+
+void CEUpdater::set_symbols(const vector<string> &new_symbs) {
+    if (new_symbs.size() != symbols_with_id->size()) {
+        throw runtime_error(
+            "The number of atoms in the updater cannot be changed via the set_symbols function\n");
+    }
+    symbols_with_id->set_symbols(new_symbs);
+}
+
+void CEUpdater::cluster_by_symm_group() {
+    m_cluster_by_symm.clear();
+
+    // Find unique symmetry groups
+    std::set<int> unique;
+    insert_in_set(this->trans_symm_group, unique);
+
+    for (const int symm : unique) {
+        if (symm == -1) {
+            // Background symmetry group
+            continue;
+        }
+        // 1 ClusterCache per ECI value
+        std::vector<ClusterCache> cluster_cache;
+        cluster_cache.reserve(this->m_parsed_names.size());
+        for (const ParsedName &parsed : this->m_parsed_names) {
+            ClusterCache cache;
+            if (parsed.size == 0 || parsed.size == 1 ||
+                !clusters.is_in_symm_group(parsed.prefix, symm)) {
+                /* Either 0- or 1-body cluster, or cluster is not in this
+                symmetry group. */
+                cluster_cache.push_back(cache);
+                continue;
+            }
+            Cluster *cluster_ptr = clusters.get_ptr(parsed.prefix, symm);
+            equiv_deco_t *equiv_ptr = cluster_ptr->get_equiv_deco_ptr(parsed.dec_str);
+            // Calculate the normalization of the resulting cluster functions
+            double normalization = static_cast<double>(cluster_ptr->size);
+            normalization /= equiv_ptr->size();
+            normalization /= normalisation_factor.at(parsed.prefix);
+
+            // Populate the new cache object
+            cache.cluster_ptr = cluster_ptr;
+            cache.equiv_deco_ptr = equiv_ptr;
+            cache.normalization = normalization;
+
+            cluster_cache.push_back(cache);
+        }
+        m_cluster_by_symm.insert({symm, cluster_cache});
+    }
+}
+void CEUpdater::set_eci(PyObject *pyeci) {
+    PyObject *key;
+    PyObject *value;
+
+    // Read the ECIs
+    Py_ssize_t pos = 0;
+    std::map<std::string, double> temp_eci;
+
+    while (PyDict_Next(pyeci, &pos, &key, &value)) {
+        temp_eci[py2string(key)] = PyFloat_AS_DOUBLE(value);
+    }
+    this->eci.init(temp_eci);
+    // Pre-parse the names of the clusters.
+    this->parse_eci_names();
+
+    // If status is not READY, then we're still initializing, and CF's are missing.
+    if (this->status == Status_t::READY && !all_eci_corresponds_to_cf()) {
+        throw invalid_argument("All ECIs has to correspond to a correlation function!");
+    }
+    // Update the cluster pointers to match the order with the ECI's.
+    cluster_by_symm_group();
+}
+
+bool CEUpdater::all_decoration_nums_equal(const vector<int> &dec_nums) const {
+    for (unsigned int i = 1; i < dec_nums.size(); i++) {
+        if (dec_nums[i] != dec_nums[0]) {
+            return false;
+        }
+    }
+    return true;
+}
+
+std::vector<double> CEUpdater::get_singlets() const {
+    size_t n = singlets.size();
+    std::vector<double> singlets_out(n, 0.0);
+
+    cf &cfs = history->get_current();
+    for (unsigned int i = 0; i < n; i++) {
+        singlets_out[i] = cfs[singlets[i]];
+    }
+    return singlets_out;
+}
+
+void CEUpdater::create_cname_with_dec(PyObject *cf) {
+    if (!PyDict_Check(cf)) {
+        throw invalid_argument("Correlation functons has to be dictionary!");
+    }
+    Py_ssize_t pos = 0;
+    PyObject *key;
+    PyObject *value;
+    while (PyDict_Next(cf, &pos, &key, &value)) {
+        string new_key = py2string(key);
+        ClusterName c_name = ClusterName(new_key);
+        unsigned int size = c_name.get_size();
+#ifdef PRINT_DEBUG
+        cout << "Read CF: " << new_key << endl;
+#endif
+        if ((size == 0 || size == 1)) {
+            cname_with_dec[new_key] = new_key;
+        } else {
+            std::string prefix = c_name.get_prefix();
+            cname_with_dec[prefix] = new_key;
+        }
+    }
+}
+
+void CEUpdater::build_trans_symm_group(PyObject *py_trans_symm_group) {
+    // Fill the symmetry group array with -1 indicating an invalid value
+    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
+        trans_symm_group[i] = -1;
+    }
+
+    unsigned int py_list_size = list_size(py_trans_symm_group);
+    for (unsigned int i = 0; i < py_list_size; i++) {
+        PyObject *sublist = PyList_GetItem(py_trans_symm_group, i);
+        unsigned int n_sites = list_size(sublist);
+        for (unsigned int j = 0; j < n_sites; j++) {
+            int indx = py2int(PyList_GetItem(sublist, j));
+            if (trans_symm_group[indx] != -1) {
+                throw runtime_error(
+                    "One site appears to be present in more than one translation symmetry "
+                    "group!");
+            }
+            trans_symm_group[indx] = i;
+        }
+    }
+
+    // Check that all sites belongs to one translational symmetry group
+    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
+        if ((trans_symm_group[i] == -1) && !is_background_index[i]) {
+            stringstream msg;
+            msg << "Site " << i << " has not been assigned to any translational symmetry group!";
+            throw runtime_error(msg.str());
+        }
+    }
+
+    // Count the number of atoms in each symmetry group
+    trans_symm_group_count.resize(py_list_size);
+    fill(trans_symm_group_count.begin(), trans_symm_group_count.end(), 0);
+    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
+        if (trans_symm_group[i] >= 0) {
+            trans_symm_group_count[trans_symm_group[i]] += 1;
+        }
+    }
+}
+
+bool CEUpdater::all_eci_corresponds_to_cf() {
+    cf &corrfunc = history->get_current();
+    return eci.names_are_equal(corrfunc);
+}
+
+double CEUpdater::calculate(vector<swap_move> &sequence) {
+    if (sequence.size() >= history->max_history / 2) {
+        throw invalid_argument(
+            "The length of sequence of swap move exceeds the buffer size for the history "
+            "tracker");
+    }
+
+    for (unsigned int i = 0; i < sequence.size(); i++) {
+        calculate(sequence[i]);
+    }
+    return get_energy();
+}
+
+double CEUpdater::calculate(vector<SymbolChange> &sequence) {
+    for (auto &change : sequence) {
+        update_cf(change);
+    }
+    return get_energy();
+}
+
+void CEUpdater::read_trans_matrix(PyObject *py_trans_mat) {
+    bool is_list = PyList_Check(py_trans_mat);
+
+#ifdef PRINT_DEBUG
+    cout << "read_trans_matrix: Extracting unique indices" << endl;
+#endif
+
+    set<int> unique_indx;
+    clusters.unique_indices(unique_indx);
+    vector<int> unique_indx_vec;
+    set2vector(unique_indx, unique_indx_vec);
+
+    // Compute the max index that is ever going to be checked
+    unsigned int max_indx = clusters.max_index();
+    if (is_list) {
+        unsigned int size = list_size(py_trans_mat);
+
+#ifdef PRINT_DEBUG
+        cout << "read_trans_matrix: Updating size of trans_matrix" << endl;
+#endif
+
+        trans_matrix.set_size(size, unique_indx_vec.size(), max_indx);
+
+#ifdef PRINT_DEBUG
+        cout << "read_trans_matrix: Setting lookup values in trans_matrix" << endl;
+#endif
+        trans_matrix.set_lookup_values(unique_indx_vec);
+
+#ifdef PRINT_DEBUG
+        cout << "read_trans_matrix: Reading translation matrix from list of dictionaries" << endl;
+#endif
+        unsigned int n_elements_insterted = 0;
+        for (unsigned int i = 0; i < size; i++) {
+            // Background atoms are ignored (and should never be accessed)
+            if (is_background_index[i] && ignore_background_indices) {
+                continue;
+            }
+
+            PyObject *dict = PyList_GetItem(py_trans_mat, i);
+            for (unsigned int j = 0; j < unique_indx_vec.size(); j++) {
+                int col = unique_indx_vec[j];
+                PyObject *value = PyDict_GetItem(dict, int2py(col));
+
+                if (value == NULL) {
+                    stringstream ss;
+                    ss << "Requested value " << col << " is not a key in the dictionary!";
+                    throw invalid_argument(ss.str());
+                }
+                trans_matrix(i, col) = py2int(value);
+                n_elements_insterted++;
+            }
+        }
+#ifdef PRINT_DEBUG
+        cout << "Inserted " << n_elements_insterted << " into the translation matrix\n";
+#endif
+    } else {
+        // Safety, this should no longer be possible to hit.
+        throw std::runtime_error("trans_matrix must be a list");
+    }
+}
+
+void CEUpdater::sort_indices(int indices[], const vector<int> &order, unsigned int n_indices) {
+    // This function is called many times
+    // profiling (with YEP) revealed that
+    // [] operator of the vector used quite a bit of time
+    // Therefore we here use raw C-arrays or pointer arithmetics
+    int sorted[4];
+    const int *ptr = &order[0];
+    for (unsigned int i = 0; i < n_indices; i++) {
+        sorted[i] = indices[*(ptr + i)];
+    }
+    memcpy(indices, sorted, n_indices * sizeof(int));
+}
+
+bool CEUpdater::is_swap_move(const swap_move &move) const {
+    return (move[0].old_symb == move[1].new_symb) && (move[1].new_symb == move[1].old_symb);
+}
+
+void CEUpdater::read_background_indices(PyObject *bkg_indices) {
+    // Fill array with false
+    is_background_index.resize(symbols_with_id->size());
+    fill(is_background_index.begin(), is_background_index.end(), false);
+
+    // Set to true if index is in bkg_indices
+    int size = list_size(bkg_indices);
+    for (int i = 0; i < size; i++) {
+        PyObject *py_indx = PyList_GetItem(bkg_indices, i);
+        int indx = py2int(py_indx);
+        is_background_index[indx] = true;
+    }
+}
+
+void CEUpdater::count_non_bkg_sites() {
+    // Count and store the number of non-background sites
+    num_non_bkg_sites = 0;
+    for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
+        if (!is_background_index[atom_no] || !ignore_background_indices) {
+            num_non_bkg_sites += 1;
+        }
+    }
+}
+
+void CEUpdater::get_changes(const std::vector<std::string> &new_symbols,
+                            std::vector<unsigned int> &changed_sites) const {
+    if (new_symbols.size() != symbols_with_id->size()) {
+        throw invalid_argument("Size of passed atoms does not match!");
+    }
+
+    for (unsigned int i = 0; i < new_symbols.size(); i++) {
+        unsigned int symb_id = symbols_with_id->id(i);
+        if (symbols_with_id->get_symbol_id(new_symbols[i]) != symb_id) {
+            changed_sites.push_back(i);
+        }
+    }
+}
+
+void CEUpdater::calculate_cf_from_scratch(const vector<string> &cf_names, map<string, double> &cf) {
+    cf.clear();
+
+    // Initialise all cluster names
+    for (const string &name : cf_names) {
+        cf[name] = 0.0;
+    }
+
+    // Loop over all clusters
+    for (const string &name : cf_names) {
+        ClusterName c_name = ClusterName(name);
+        unsigned int cluster_size = c_name.get_size();
+
+        // Handle empty cluster
+        if (cluster_size == 0) {
+            cf[name] = 1.0;
+            continue;
+        }
+
+        // Handle singlet cluster
+        if (cluster_size == 1) {
+            unsigned int dec = c_name.get_dec_num();
+            double new_value = 0.0;
+            // Normalise with respect to the actual number of atoms included
+            for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
+                if (!is_background_index[atom_no] || !ignore_background_indices) {
+                    new_value += basis_functions->get(dec, symbols_with_id->id(atom_no));
+                }
+            }
+            cf[name] = new_value / num_non_bkg_sites;
+            continue;
+        }
+
+        // Handle the rest of the clusters
+        std::string prefix, dec_str;
+        c_name.get_prefix_and_dec_str(prefix, dec_str);
+
+        double sp = 0.0;
+        double count = 0;
+        for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
+            int symm = trans_symm_group[atom_no];
+
+            if ((!clusters.is_in_symm_group(prefix, symm)) ||
+                (is_background_index[atom_no] && ignore_background_indices)) {
+                continue;
+            }
+
+            const Cluster &cluster = clusters.get(prefix, symm);
+
+            const equiv_deco_t &equiv_deco = cluster.get_equiv_deco(dec_str);
+            unsigned int ref_id = symbols_with_id->id(atom_no);
+            double sp_temp = 0.0;
+            for (const vector<int> &deco : equiv_deco) {
+                sp_temp += spin_product_one_atom(atom_no, cluster, deco, ref_id);
+            }
+            sp += sp_temp / equiv_deco.size();
+            count += cluster.get().size();
+        }
+
+        if (count == 0) {
+            cf[name] = 0.0;
+        } else {
+            cf[name] = sp / count;
+        }
+    }
+
+    history->get_current().init(cf);
+}
+
+void CEUpdater::set_atoms(PyObject *py_atoms) {
+    unsigned int num_atoms = PySequence_Length(py_atoms);
+
+    if (num_atoms != symbols_with_id->size()) {
+        throw invalid_argument("Length of passed atoms object is different from current");
+    }
+
+    std::vector<std::string> symbols = get_symbols_from_atoms(py_atoms);
+
+    this->atoms = py_atoms;
+    symbols_with_id->set_symbols(symbols);
+}

clease.egg-info/PKG-INFO

@@ -1,4 +0,4 @@
-Metadata-Version: 2.1
+Metadata-Version: 2.4
 Name: clease
-Version: 1.1.0
+Version: 1.2.0
 Summary: CLuster Expansion in Atomistic Simulation Environment
@@ -14,13 +14,10 @@ Home-page: https://gitlab.com/computationalmaterials/clease/
 Classifier: License :: OSI Approved :: Mozilla Public License 2.0 (MPL 2.0)
-Classifier: Programming Language :: Python :: 3.7
-Classifier: Programming Language :: Python :: 3.8
-Classifier: Programming Language :: Python :: 3.9
-Classifier: Programming Language :: Python :: 3.10
+Classifier: Programming Language :: Python :: 3
 Classifier: Topic :: Scientific/Engineering :: Physics
 Classifier: Topic :: Scientific/Engineering :: Chemistry
-Requires-Python: >=3.7
+Requires-Python: >=3.10
 Description-Content-Type: text/markdown
 License-File: LICENSE.md
 Requires-Dist: ase>=3.22
-Requires-Dist: numpy<2
+Requires-Dist: numpy
 Requires-Dist: cython
@@ -52,5 +49,4 @@ Requires-Dist: matplotlib
 Requires-Dist: twine; extra == "dev"
-Requires-Dist: black>=22.1.0; extra == "dev"
-Requires-Dist: clang-format>=14.0.3; extra == "dev"
-Requires-Dist: ruff; extra == "dev"
+Requires-Dist: clang-format==21.1.5; extra == "dev"
+Requires-Dist: ruff==0.14.5; extra == "dev"
 Requires-Dist: pyclean>=2.0.0; extra == "dev"
@@ -62,21 +58,22 @@ Requires-Dist: pytest-cov; extra == "dev"
 Provides-Extra: all
-Requires-Dist: black>=22.1.0; extra == "all"
+Requires-Dist: pyclean>=2.0.0; extra == "all"
+Requires-Dist: pre-commit; extra == "all"
+Requires-Dist: ipython; extra == "all"
+Requires-Dist: sphinx; extra == "all"
+Requires-Dist: pytest; extra == "all"
+Requires-Dist: twine; extra == "all"
+Requires-Dist: cython; extra == "all"
+Requires-Dist: ruff==0.14.5; extra == "all"
+Requires-Dist: pytest-benchmark[histogram]>=3.4.1; extra == "all"
 Requires-Dist: mock; extra == "all"
-Requires-Dist: twine; extra == "all"
-Requires-Dist: ruff; extra == "all"
+Requires-Dist: clang-format==21.1.5; extra == "all"
 Requires-Dist: sphinx_rtd_theme; extra == "all"
-Requires-Dist: pre-commit; extra == "all"
-Requires-Dist: pyclean>=2.0.0; extra == "all"
 Requires-Dist: pip; extra == "all"
-Requires-Dist: pytest-cov; extra == "all"
 Requires-Dist: clease-gui; extra == "all"
+Requires-Dist: tox>=4; extra == "all"
+Requires-Dist: pytest-mock; extra == "all"
 Requires-Dist: build; extra == "all"
-Requires-Dist: sphinx; extra == "all"
-Requires-Dist: pytest-mock; extra == "all"
-Requires-Dist: pytest; extra == "all"
-Requires-Dist: cython; extra == "all"
-Requires-Dist: ipython; extra == "all"
-Requires-Dist: clang-format>=14.0.3; extra == "all"
-Requires-Dist: tox>=4; extra == "all"
-Requires-Dist: pytest-benchmark[histogram]>=3.4.1; extra == "all"
+Requires-Dist: pytest-cov; extra == "all"
+Dynamic: license-file
+Dynamic: provides-extra
 
@@ -94,2 +91,10 @@ # CLEASE
 
+# Documentation
+
+   The updated version of the CLEASE documentation is now available at:
+
+   https://www.clease.org
+
+   Please visit the new site for the latest documentation, tutorials, and resources.
+
 # Installation
@@ -96,0 +101,0 @@

clease.egg-info/requires.txt

 ase>=3.22
-numpy<2
+numpy
 cython
@@ -16,21 +16,20 @@ matplotlib
 [all]
-black>=22.1.0
+pyclean>=2.0.0
+pre-commit
+ipython
+sphinx
+pytest
+twine
+cython
+ruff==0.14.5
+pytest-benchmark[histogram]>=3.4.1
 mock
-twine
-ruff
+clang-format==21.1.5
 sphinx_rtd_theme
-pre-commit
-pyclean>=2.0.0
 pip
-pytest-cov
 clease-gui
+tox>=4
+pytest-mock
 build
-sphinx
-pytest-mock
-pytest
-cython
-ipython
-clang-format>=14.0.3
-tox>=4
-pytest-benchmark[histogram]>=3.4.1
+pytest-cov
 
@@ -43,5 +42,4 @@ [dev]
 twine
-black>=22.1.0
-clang-format>=14.0.3
-ruff
+clang-format==21.1.5
+ruff==0.14.5
 pyclean>=2.0.0
@@ -48,0 +46,0 @@ pytest-cov

clease.egg-info/SOURCES.txt

@@ -153,2 +153,3 @@ LICENSE.md
 cxx/src/basis_function.cpp
+cxx/src/ce_updater.cpp
 cxx/src/cf_history_tracker.cpp
@@ -162,3 +163,2 @@ cxx/src/cluster.cpp
 cxx/src/symbols_with_numbers.cpp
-cxx/src/with_numpy/ce_updater.cpp
 doc/source/acknowledgements.rst
@@ -165,0 +165,0 @@ doc/source/api_doc.rst

clease/_version.txt

@@ -1,1 +0,1 @@
-1.1.0
+1.2.0

clease/basis_function.py

 """Module for setting up pseudospins and basis functions."""
+
 from abc import ABC, abstractmethod
+from collections.abc import Sequence
 import math
-from typing import Dict, List, Optional, Sequence
 
@@ -38,3 +39,3 @@ import numpy as np
     @property
-    def unique_elements(self) -> List[str]:
+    def unique_elements(self) -> list[str]:
         return self._unique_elements
@@ -51,7 +52,7 @@
     @property
-    def spin_dict(self) -> Dict[str, int]:
+    def spin_dict(self) -> dict[str, int]:
         return self.get_spin_dict()
 
     @property
-    def basis_functions(self) -> List[Dict[str, float]]:
+    def basis_functions(self) -> list[dict[str, float]]:
         """Property access to :meth:`get_basis_functions`."""
@@ -89,3 +90,3 @@ return self.get_basis_functions()
 
-    def get_spin_dict(self) -> Dict[str, int]:
+    def get_spin_dict(self) -> dict[str, int]:
         """Define pseudospins for all consistuting elements."""
@@ -99,3 +100,3 @@ gram_schmidt = GramSchmidtMonimial(self.num_unique_elements)
 
-    def get_basis_functions(self) -> List[Dict[str, float]]:
+    def get_basis_functions(self) -> list[dict[str, float]]:
         """Create basis functions to guarantee the orthonormality."""
@@ -118,3 +119,3 @@ gram_schmidt = GramSchmidtMonimial(self.num_unique_elements)
 
-    def get_spin_dict(self) -> Dict[str, int]:
+    def get_spin_dict(self) -> dict[str, int]:
         """Define pseudospins for all consistuting elements."""
@@ -127,3 +128,3 @@ spin_values = list(range(self.num_unique_elements))
 
-    def get_basis_functions(self) -> List[Dict[str, float]]:
+    def get_basis_functions(self) -> list[dict[str, float]]:
         """Create basis functions to guarantee the orthonormality."""
@@ -175,3 +176,3 @@ alpha = list(range(1, self.num_unique_elements))
 
-    def __init__(self, unique_elements: List[str], redundant_element: Optional[str] = "auto"):
+    def __init__(self, unique_elements: list[str], redundant_element: str | None = "auto"):
         super().__init__(unique_elements)
@@ -183,3 +184,3 @@ if redundant_element == "auto":
 
-    def get_spin_dict(self) -> Dict[str, int]:
+    def get_spin_dict(self) -> dict[str, int]:
         """Define pseudospins for all consistuting elements."""
@@ -192,3 +193,3 @@ spin_values = list(range(self.num_unique_elements))
 
-    def get_basis_functions(self) -> List[Dict[str, float]]:
+    def get_basis_functions(self) -> list[dict[str, float]]:
         """Create orthonormal basis functions.
@@ -195,0 +196,0 @@

clease/calculator/clease_vol_dep.py

@@ -1,2 +0,2 @@
-from typing import Dict, Iterable, List, Optional, Sequence, Set, Union
+from collections.abc import Iterable, Sequence
 
@@ -83,5 +83,5 @@ from ase import Atoms
         settings: ClusterExpansionSettings,
-        eci: Dict[str, float],
-        vol_coeff: Dict[str, float],
-        init_cf: Optional[Dict[str, float]] = None,
+        eci: dict[str, float],
+        vol_coeff: dict[str, float],
+        init_cf: dict[str, float] | None = None,
     ):
@@ -92,5 +92,3 @@ if not eci_format_ok(eci.keys()):
         if not vol_coeff_format_ok(vol_coeff.keys()):
-            raise ValueError(
-                "Invalid format of volume coefficient names. " f"Got\n{vol_coeff.keys()}"
-            )
+            raise ValueError(f"Invalid format of volume coefficient names. Got\n{vol_coeff.keys()}")
 
@@ -116,3 +114,3 @@ cf_to_track = unique_eci_names_no_vol(eci.keys())
 
-    def get_volume(self, cf: Optional[Dict[str, float]] = None) -> float:
+    def get_volume(self, cf: dict[str, float] | None = None) -> float:
         """
@@ -133,3 +131,3 @@ Returns the volume per atom
 
-    def get_pressure(self, cf: Optional[Dict[str, float]] = None) -> float:
+    def get_pressure(self, cf: dict[str, float] | None = None) -> float:
         """
@@ -155,3 +153,3 @@ Return the pressure
 
-    def get_bulk_modulus(self, cf: Optional[Dict[str, float]] = None) -> float:
+    def get_bulk_modulus(self, cf: dict[str, float] | None = None) -> float:
         """
@@ -173,3 +171,3 @@ Return the bulk modulus of the current atoms object
 
-    def _d2EdV2(self, cf: Dict[str, float], vol: float) -> float:
+    def _d2EdV2(self, cf: dict[str, float], vol: float) -> float:
         """
@@ -188,3 +186,3 @@ Return the double derivative of the energy with respect
 
-    def _d3EdV3(self, cf: Dict[str, float], vol: float) -> float:
+    def _d3EdV3(self, cf: dict[str, float], vol: float) -> float:
         """
@@ -208,3 +206,3 @@ Return the third derivative of the energy with respect to
 
-    def get_dBdP(self, cf: Optional[Dict[str, float]] = None) -> float:
+    def get_dBdP(self, cf: dict[str, float] | None = None) -> float:
         """
@@ -243,5 +241,5 @@ Return the pressure derivative of the bulk modulus of the
         self,
-        atoms: Optional[Atoms] = None,
-        properties: Optional[List[str]] = None,
-        system_changes: Union[Sequence[SystemChange], None] = None,
+        atoms: Atoms | None = None,
+        properties: list[str] | None = None,
+        system_changes: Sequence[SystemChange] | None = None,
     ) -> float:
@@ -272,3 +270,3 @@ """Calculate the energy of the passed Atoms object.
 
-def unique_eci_names_no_vol(eci_names: Iterable[str]) -> Set[str]:
+def unique_eci_names_no_vol(eci_names: Iterable[str]) -> set[str]:
     """
@@ -275,0 +273,0 @@ Return a set with the unique ECI names without any volume tag

clease/calculator/clease.py

 """Calculator for Cluster Expansion."""
+
+from collections.abc import Sequence
 import contextlib
 import sys
-from typing import Any, Dict, List, Optional, Sequence, TextIO, Union
+from typing import Any, TextIO
 
@@ -52,5 +54,5 @@ from ase import Atoms
         settings: ClusterExpansionSettings,
-        eci: Dict[str, float],
-        init_cf: Optional[Dict[str, float]] = None,
-        logfile: Union[TextIO, str, None] = None,
+        eci: dict[str, float],
+        init_cf: dict[str, float] | None = None,
+        logfile: TextIO | str | None = None,
     ) -> None:
@@ -146,3 +148,3 @@ if not isinstance(settings, ClusterExpansionSettings):
 
-    def calculate_cf_from_scratch(self) -> Dict[str, float]:
+    def calculate_cf_from_scratch(self) -> dict[str, float]:
         """Calculate correlation functions from scratch."""
@@ -162,5 +164,5 @@ self.require_updater()
         self,
-        atoms: Optional[Atoms] = None,
-        properties: Optional[List[str]] = None,
-        system_changes: Optional[SystemChanges] = None,
+        atoms: Atoms | None = None,
+        properties: list[str] | None = None,
+        system_changes: SystemChanges | None = None,
     ) -> float:
@@ -202,5 +204,3 @@ """Calculate the energy of the passed Atoms object.
 
-    def get_property(
-        self, name: str, atoms: Optional[Atoms] = None, allow_calculation: bool = True
-    ):
+    def get_property(self, name: str, atoms: Atoms | None = None, allow_calculation: bool = True):
         """Get a property from the calculator.
@@ -216,3 +216,3 @@
 
-    def get_potential_energy(self, atoms: Optional[Atoms] = None) -> float:
+    def get_potential_energy(self, atoms: Atoms | None = None) -> float:
         """Calculate the energy from scratch with an atoms object"""
@@ -231,5 +231,3 @@ # self.set_atoms(atoms)
 
-    def calculation_required(
-        self, atoms: Atoms, properties: Optional[Sequence[str]] = None
-    ) -> bool:
+    def calculation_required(self, atoms: Atoms, properties: Sequence[str] | None = None) -> bool:
         """Check whether a calculation is required for a given atoms object.
@@ -247,3 +245,3 @@ The ``properties`` argument only exists for compatibility reasons, and has no effect.
 
-    def check_state(self, atoms: Atoms) -> List[str]:
+    def check_state(self, atoms: Atoms) -> list[str]:
         """Method for checking if energy needs calculation.
@@ -266,3 +264,3 @@ Primarily for ASE compatibility.
     @property
-    def indices_of_changed_atoms(self) -> List[int]:
+    def indices_of_changed_atoms(self) -> list[int]:
         """Return the indices of atoms that have been changed."""
@@ -277,3 +275,3 @@ changed = self.get_changed_sites(self.atoms)
 
-    def get_changed_sites(self, atoms: Atoms) -> List[int]:
+    def get_changed_sites(self, atoms: Atoms) -> list[int]:
         """Return the list of indices which differ from the internal ones."""
@@ -283,7 +281,7 @@ self.require_updater()
 
-    def get_cf(self) -> Dict[str, float]:
+    def get_cf(self) -> dict[str, float]:
         """Return the correlation functions as a dict"""
         return self.updater.get_cf()
 
-    def update_cf(self, system_changes: Optional[SystemChanges] = None) -> None:
+    def update_cf(self, system_changes: SystemChanges | None = None) -> None:
         """Update correlation function based on the reference value.
@@ -315,3 +313,3 @@
     @property
-    def cf(self) -> List[float]:
+    def cf(self) -> list[float]:
         temp_cf = self.updater.get_cf()
@@ -327,3 +325,3 @@ return [temp_cf[x] for x in self.cf_names]
 
-    def update_eci(self, eci: Dict[str, float]) -> None:
+    def update_eci(self, eci: dict[str, float]) -> None:
         """Update the ECI values.
@@ -338,3 +336,4 @@
     def get_singlets(self) -> np.ndarray:
-        return self.updater.get_singlets()
+        singlets = self.updater.get_singlets()
+        return np.array(singlets, dtype=float)
 
@@ -431,3 +430,3 @@ def get_energy(self) -> float:
     @property
-    def parameters(self) -> Dict[str, Any]:
+    def parameters(self) -> dict[str, Any]:
         """Return a dictionary with relevant parameters."""
@@ -443,3 +442,3 @@ return {"eci": self.eci}
 
-def _check_properties(properties: Optional[List[str]], implemented_properties: List[str]) -> None:
+def _check_properties(properties: list[str] | None, implemented_properties: list[str]) -> None:
     """Check whether the passed properties is supported. If it is None, nothing is checked.
@@ -446,0 +445,0 @@ Raises PropertyNotImplementedError upon finding a bad property.

clease/calculator/util.py

@@ -1,3 +0,1 @@
-from typing import Dict, Optional
-
 from ase import Atoms
@@ -14,4 +12,4 @@
     atoms: Atoms,
-    eci: Optional[Dict[str, float]] = None,
-    num_threads: Optional[int] = None,
+    eci: dict[str, float] | None = None,
+    num_threads: int | None = None,
 ) -> Atoms:
@@ -46,3 +44,3 @@ """Utility function for an efficient initialization of large cells. Will set the atoms
 
-def get_ce_energy(settings: ClusterExpansionSettings, atoms: Atoms, eci: Dict[str, float]) -> float:
+def get_ce_energy(settings: ClusterExpansionSettings, atoms: Atoms, eci: dict[str, float]) -> float:
     """Get energy of the ASE Atoms object based on given ECI values.
@@ -49,0 +47,0 @@

clease/cluster_coverage.py

@@ -0,4 +1,4 @@
+from collections.abc import Sequence
 from copy import deepcopy
 import sys
-from typing import Dict, List, Optional, Sequence, Tuple
 
@@ -56,3 +56,3 @@ from ase import Atoms
         settings: ClusterExpansionSettings,
-        select_cond: Optional[Sequence[Tuple[str, str, str]]] = None,
+        select_cond: Sequence[tuple[str, str, str]] | None = None,
     ):
@@ -67,3 +67,3 @@ # Make copy such that we don't alter the settings object
 
-    def _unique_templates(self) -> List[Atoms]:
+    def _unique_templates(self) -> list[Atoms]:
         unique_templates = []
@@ -78,3 +78,3 @@
 
-    def coverage(self, template: Atoms) -> Dict[str, float]:
+    def coverage(self, template: Atoms) -> dict[str, float]:
         """
@@ -101,3 +101,3 @@ Return the cluster coverage for the passed atoms object.
 
-    def max_coverage(self) -> Dict[str, float]:
+    def max_coverage(self) -> dict[str, float]:
         """
@@ -115,3 +115,3 @@ Return the maximum cluster coverage among all templates in the database.
 
-    def print_report(self, coverage: Optional[Dict[str, float]] = None, file=sys.stdout) -> None:
+    def print_report(self, coverage: dict[str, float] | None = None, file=sys.stdout) -> None:
         """
@@ -155,3 +155,3 @@ Prints a nicely formatted report of coverage.
 
-def _grade(coverage: Dict[str, float]) -> str:
+def _grade(coverage: dict[str, float]) -> str:
     grade = "Very good (coverage > 0.75)"
@@ -158,0 +158,0 @@ for v in coverage.values():

clease/cluster/cluster_fingerprint.py

@@ -0,3 +1,4 @@
+from collections.abc import Iterable
 from functools import total_ordering
-from typing import Any, Iterable
+from typing import Any
 
@@ -4,0 +5,0 @@ import attr

clease/cluster/cluster_generator.py

@@ -0,4 +1,4 @@
+from collections.abc import Iterable, Iterator, Sequence
 from itertools import product
 from math import sqrt
-from typing import Dict, Iterable, Iterator, List, Optional, Sequence, Set, Tuple, Union
 
@@ -52,3 +52,3 @@ from ase import Atoms
 
-    def _as_euc(self, x: Union[np.ndarray, FourVector]) -> np.ndarray:
+    def _as_euc(self, x: np.ndarray | FourVector) -> np.ndarray:
         """Helper function to translate vector to NumPy array format"""
@@ -62,3 +62,3 @@ if isinstance(x, FourVector):
         self, cartesian: np.ndarray, sublattices: Sequence[int]
-    ) -> List[FourVector]:
+    ) -> list[FourVector]:
         """Translate many positions into FourVector's"""
@@ -90,3 +90,3 @@
 
-    def to_four_vector(self, cartesian: np.ndarray, sublattice: Optional[int] = None) -> FourVector:
+    def to_four_vector(self, cartesian: np.ndarray, sublattice: int | None = None) -> FourVector:
         """Translate a position in Cartesian coordinates to its FourVector"""
@@ -128,4 +128,4 @@ if cartesian.ndim != 1:
         self,
-        x1: Union[np.ndarray, FourVector],
-        x2: Union[np.ndarray, FourVector],
+        x1: np.ndarray | FourVector,
+        x2: np.ndarray | FourVector,
     ) -> float:
@@ -204,3 +204,3 @@ """
         self, cutoff: float, lattice: int
-    ) -> Dict[FourVector, Set[FourVector]]:
+    ) -> dict[FourVector, set[FourVector]]:
         """Prepare all sites which are within the cutoff sphere. Note, this only prepares sites
@@ -255,3 +255,3 @@ which are pair-wise within the cutoff, and does not consider the distance to the
         self, size: int, cutoff: float, ref_lattice: int
-    ) -> Tuple[List[List[Figure]], List[ClusterFingerprint]]:
+    ) -> tuple[list[list[Figure]], list[ClusterFingerprint]]:
         """Generate all possible figures of a given size, are within a given cutoff radius
@@ -267,7 +267,7 @@ (from the center of mass of the figure), and from a reference lattice.
         Returns:
-            List[List[Figure]], List[ClusterFingerprint]: The collection of figures and
+            list[list[Figure]], list[ClusterFingerprint]: The collection of figures and
                 their corresponding fingerprints.
         """
-        clusters: List[List[Figure]] = []
-        all_fps: List[ClusterFingerprint] = []
+        clusters: list[list[Figure]] = []
+        all_fps: list[ClusterFingerprint] = []
 
@@ -325,5 +325,5 @@ for new_figure in self.figure_iterator(size, cutoff, ref_lattice):
         order = np.lexsort(dists.T)[::-1]
-        return Figure((fvs[i] for i in order))
+        return Figure(fvs[i] for i in order)
 
-    def to_atom_index(self, cluster: Cluster, lut: Dict[FourVector, int]) -> List[List[int]]:
+    def to_atom_index(self, cluster: Cluster, lut: dict[FourVector, int]) -> list[list[int]]:
         """
@@ -340,3 +340,3 @@ Convert the integer vector representation to an atomic index
 
-    def equivalent_sites(self, figure: Figure) -> List[List[int]]:
+    def equivalent_sites(self, figure: Figure) -> list[list[int]]:
         """Find the equivalent sites of a figure."""
@@ -352,3 +352,3 @@ dists = self._get_internal_distances(figure, sort=True)
         # Merge pairs into groups
-        merged: List[Set[int]] = []
+        merged: list[set[int]] = []
         for equiv in equiv_sites:
@@ -380,3 +380,3 @@ found_group = False
         # Dictionary mapping a ref lattice to a list of four-vectors
-        self.pre_calc: Dict[int, List[FourVector]] = {}
+        self.pre_calc: dict[int, list[FourVector]] = {}
 
@@ -391,3 +391,3 @@ def must_generate(self, cutoff: float, ref_lattice: int) -> bool:
 
-    def get(self, cutoff: float, ref_lattice: int) -> List[FourVector]:
+    def get(self, cutoff: float, ref_lattice: int) -> list[FourVector]:
         """
@@ -404,3 +404,3 @@ Return sites within the cutoff
 def site_iterator(
-    within_cutoff: Dict[FourVector, Set[FourVector]], size: int, ref_lattice: int
+    within_cutoff: dict[FourVector, set[FourVector]], size: int, ref_lattice: int
 ) -> Iterator[Figure]:
@@ -407,0 +407,0 @@ """

clease/cluster/cluster_list.py

@@ -6,3 +6,3 @@ from __future__ import annotations
 import logging
-from typing import Any, Dict, List
+from typing import Any
 
@@ -30,9 +30,9 @@ from ase import Atoms
         # Format of the names cache: {num_bf: names}
-        self._all_cf_name_cache: Dict[int, List[str]] = {}
+        self._all_cf_name_cache: dict[int, list[str]] = {}
 
-    def todict(self) -> Dict[str, Any]:
+    def todict(self) -> dict[str, Any]:
         return {"clusters": self.clusters}
 
     @classmethod
-    def from_dict(cls, dct: Dict[str, Any]) -> ClusterList:
+    def from_dict(cls, dct: dict[str, Any]) -> ClusterList:
         cluster_lst = cls()
@@ -45,3 +45,3 @@ clusters = dct["clusters"]
     @property
-    def clusters(self) -> List[Cluster]:
+    def clusters(self) -> list[Cluster]:
         return self._clusters
@@ -63,7 +63,7 @@
     @property
-    def names(self) -> List[str]:
+    def names(self) -> list[str]:
         """Get all names in the cluster list"""
         return [cluster.name for cluster in self.clusters]
 
-    def get_by_name(self, name) -> List[Cluster]:
+    def get_by_name(self, name) -> list[Cluster]:
         return [cluster for cluster in self.clusters if cluster.name == name]
@@ -79,3 +79,3 @@
 
-    def get_by_size(self, size) -> List[Cluster]:
+    def get_by_size(self, size) -> list[Cluster]:
         # Return all clusters with a given size
@@ -120,3 +120,3 @@ return [c for c in self.clusters if c.size == size]
 
-    def get_all_cf_names(self, num_bf: int) -> List[str]:
+    def get_all_cf_names(self, num_bf: int) -> list[str]:
         """
@@ -138,3 +138,3 @@ Return a list of all correlation function names
 
-    def _build_all_cf_names(self, num_bf: int) -> List[str]:
+    def _build_all_cf_names(self, num_bf: int) -> list[str]:
         if not isinstance(num_bf, int):
@@ -199,3 +199,3 @@ raise TypeError(f"Number of basis functions must be integer, got {num_bf}")
 
-    def tolist(self) -> List[Cluster]:
+    def tolist(self) -> list[Cluster]:
         """Returns a copy of the ClusterList as a regular list."""
@@ -224,3 +224,3 @@ return list(self.clusters)
 
-    def multiplicity_factors(self, num_sites_per_group: List[int]) -> Dict[str, float]:
+    def multiplicity_factors(self, num_sites_per_group: list[int]) -> dict[str, float]:
         mult_factors = {}
@@ -245,5 +245,5 @@ norm = {}
 
-    def get_figures(self, generator: ClusterGenerator) -> List[Atoms]:
+    def get_figures(self, generator: ClusterGenerator) -> list[Atoms]:
         """Get the figures (in their ASE Atoms object representation)"""
-        figures: List[Atoms] = []
+        figures: list[Atoms] = []
         self.sort()
@@ -250,0 +250,0 @@ # We want to skip c0 and c1 anyways

clease/cluster/cluster_manager.py

@@ -0,1 +1,2 @@
+from collections.abc import Callable, Iterator, Sequence
 from copy import deepcopy
@@ -5,3 +6,2 @@ import functools
 import logging
-from typing import Callable, Dict, Iterator, List, Optional, Sequence, Set, Tuple
 
@@ -35,3 +35,3 @@ import ase
 
-    def __init__(self, prim_cell: ase.Atoms, background_syms: Optional[Set[str]] = None):
+    def __init__(self, prim_cell: ase.Atoms, background_syms: set[str] | None = None):
         self._background_syms = background_syms or set()
@@ -56,3 +56,3 @@
     @property
-    def background_syms(self) -> Set[str]:
+    def background_syms(self) -> set[str]:
         """The symbols which are considered background."""
@@ -118,4 +118,4 @@ return self._background_syms
         num_lattices = range(len(self.prim))
-        all_fps: List[ClusterFingerprint] = []
-        all_figures: List[List[Figure]] = []
+        all_fps: list[ClusterFingerprint] = []
+        all_figures: list[list[Figure]] = []
         lattices = []
@@ -233,3 +233,3 @@
 
-    def unique_four_vectors(self) -> Set[FourVector]:
+    def unique_four_vectors(self) -> set[FourVector]:
         """
@@ -250,3 +250,3 @@ Return a list with all unique 4-vectors which are
 
-    def get_figures(self) -> List[ase.Atoms]:
+    def get_figures(self) -> list[ase.Atoms]:
         """
@@ -257,3 +257,3 @@ Return a list of atoms object representing the clusters
 
-    def make_all_four_vectors(self, atoms: ase.Atoms) -> Tuple[np.ndarray, List[FourVector]]:
+    def make_all_four_vectors(self, atoms: ase.Atoms) -> tuple[np.ndarray, list[FourVector]]:
         """Construct the FourVector to every site which is not a background site.
@@ -273,3 +273,3 @@
 
-    def create_four_vector_lut(self, template: ase.Atoms) -> Dict[FourVector, int]:
+    def create_four_vector_lut(self, template: ase.Atoms) -> dict[FourVector, int]:
         """
@@ -290,3 +290,3 @@ Construct a lookup table (LUT) for the index in template given the
         self, template: ase.Atoms, unique: Sequence[FourVector]
-    ) -> Dict[FourVector, int]:
+    ) -> dict[FourVector, int]:
         """Translate a set of unique FourVectors into their corresponding index
@@ -307,3 +307,3 @@ in a template atoms object."""
 
-    def translation_matrix(self, template: ase.Atoms) -> List[Dict[int, int]]:
+    def translation_matrix(self, template: ase.Atoms) -> list[dict[int, int]]:
         """
@@ -336,3 +336,3 @@ Construct the translation matrix.
 
-        def _make_site_mapping(index: int) -> Dict[int, int]:
+        def _make_site_mapping(index: int) -> dict[int, int]:
             """Helper function to calculate the translation mapping for each
@@ -399,3 +399,3 @@ atomic site."""
         unique_xyz: np.ndarray,
-        sublattices: List[int],
+        sublattices: list[int],
         rep_arr: np.ndarray,
@@ -415,3 +415,3 @@ ) -> Iterator[FourVector]:
         cell: np.ndarray,
-        cell_T_inv: Optional[np.ndarray] = None,
+        cell_T_inv: np.ndarray | None = None,
     ) -> Iterator[FourVector]:
@@ -436,4 +436,4 @@ """Generalized FourVector wrapping function."""
         max_cluster_dia: Sequence[float],
-        index_by_sublattice: List[List[int]],
-    ) -> Tuple[ClusterList, TransMatrix]:
+        index_by_sublattice: list[list[int]],
+    ) -> tuple[ClusterList, TransMatrix]:
         """Create a ClusterList and a TransMatrix object, and calculate the norm factors."""
@@ -458,3 +458,3 @@ # Ensure that we have built the clusters for the cutoff
 def _set_norm_factors(
-    index_by_sublattice: List[List[int]],
+    index_by_sublattice: list[list[int]],
     trans_matrix: TransMatrix,
@@ -461,0 +461,0 @@ cluster_list: ClusterList,

clease/cluster/cluster.py

@@ -0,3 +1,4 @@
+from collections.abc import Sequence
 from functools import total_ordering
-from typing import Any, Dict, List, Sequence
+from typing import Any
 
@@ -34,3 +35,3 @@ from ase import Atoms
 
-    info: Dict[str, Any] = attr.field(default=attr.Factory(dict))
+    info: dict[str, Any] = attr.field(default=attr.Factory(dict))
     # "indices" are the integer index representation of the Figures.
@@ -48,3 +49,3 @@ # therefore, "indices" and "ref_indx" depend on the currently active template,
                 raise TypeError(
-                    f"All values must Figure type, got {value} " f"of type {type(v)} in index {ii}."
+                    f"All values must Figure type, got {value} of type {type(v)} in index {ii}."
                 )
@@ -94,3 +95,3 @@
 
-    def _order_equiv_sites(self, figure: Sequence[int]) -> List[int]:
+    def _order_equiv_sites(self, figure: Sequence[int]) -> list[int]:
         """Sort equivalent sites of a figure in index representation."""
@@ -106,3 +107,3 @@ figure_cpy = list(figure)
     @property
-    def num_fig_occurences(self) -> Dict[str, int]:
+    def num_fig_occurences(self) -> dict[str, int]:
         """Number of ocurrences for each figures."""
@@ -120,3 +121,3 @@ occ_count = {}
         target_figure: Sequence[int],
-        trans_matrix: List[Dict[int, int]],
+        trans_matrix: list[dict[int, int]],
     ):
@@ -146,7 +147,6 @@ """Find figures that correspond to another reference index.
         raise RuntimeError(
-            f"There are no matching figure for ref_indx: "
-            f"{ref_indx} and figure: {target_figure}!"
+            f"There are no matching figure for ref_indx: {ref_indx} and figure: {target_figure}!"
         )
 
-    def get_all_figure_keys(self) -> List[str]:
+    def get_all_figure_keys(self) -> list[str]:
         return [self.get_figure_key(fig) for fig in self.indices]
@@ -153,0 +153,0 @@

clease/convexhull.py

@@ -1,3 +0,1 @@
-from typing import Optional
-
 from ase.db import connect
@@ -48,3 +46,3 @@ import matplotlib.pyplot as plt
         conc_scale=1.0,
-        conc_ranges: Optional[dict] = None,
+        conc_ranges: dict | None = None,
     ):
@@ -104,3 +102,3 @@ if conc_ranges is None:
                 v[f"{k2}_conc"] = 0.0
-            v["energy"] = np.infty
+            v["energy"] = np.inf
 
@@ -107,0 +105,0 @@ # Iterate through all rows, update

clease/corr_func.py

 """Module for calculating correlation functions."""
+
+from collections.abc import Iterator
 import logging
-from typing import Any, Dict, Iterator, Tuple
+from typing import Any
 
@@ -16,3 +18,3 @@ from ase.atoms import Atoms
 # Type alias for a Correlation function
-CF_T = Dict[str, float]
+CF_T = dict[str, float]
 
@@ -148,3 +150,3 @@
 
-    def iter_reconfigure_db_entries(self, select_cond=None) -> Iterator[Tuple[int, int, int]]:
+    def iter_reconfigure_db_entries(self, select_cond=None) -> Iterator[tuple[int, int, int]]:
         """Iterator which reconfigures the correlation function values in the DB,
@@ -156,3 +158,3 @@ which yields after each reconfiguration and reports on the progress.
         Yields:
-            Tuple[int, int, int]: (row_id, count, total) A tuple containing the ID
+            tuple[int, int, int]: (row_id, count, total) A tuple containing the ID
             of the row which was just reconfigured, current
@@ -159,0 +161,0 @@ count which has been reconfigured, as well as the total number of

clease/data_manager.py

 from abc import ABC, abstractmethod
 from collections import defaultdict
+from collections.abc import Sequence
 from itertools import combinations_with_replacement as cwr
 import logging
 import sqlite3
-from typing import Dict, List, Optional, Sequence, Set, Tuple
 
@@ -81,3 +81,3 @@ from ase.db import connect
     @abstractmethod
-    def get_data(self, select_cond: List[tuple]) -> Tuple[np.ndarray, np.ndarray]:
+    def get_data(self, select_cond: list[tuple]) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -88,4 +88,4 @@ Return the design matrix X and the target data y
     def _get_data_from_getters(
-        self, select_cond: List[tuple], feature_getter: FeatureGetter, target_getter: TargetGetter
-    ) -> Tuple[np.ndarray, np.ndarray]:
+        self, select_cond: list[tuple], feature_getter: FeatureGetter, target_getter: TargetGetter
+    ) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -178,3 +178,3 @@ Return the design matrix X and the target data y
 
-    def get_matching_names(self, pattern: str) -> List[str]:
+    def get_matching_names(self, pattern: str) -> list[str]:
         """
@@ -192,3 +192,3 @@ Get names that matches pattern
 
-    def get_cols(self, names: List[str]) -> np.ndarray:
+    def get_cols(self, names: list[str]) -> np.ndarray:
         """
@@ -203,3 +203,3 @@ Get all columns corresponding to the names
 
-    def groups(self) -> List[int]:
+    def groups(self) -> list[int]:
         """
@@ -235,3 +235,3 @@ Returns the group of each item in the X matrix. In the top-level
         tab_name: str,
-        cf_names: Optional[List[str]] = None,
+        cf_names: list[str] | None = None,
         order: int = 1,
@@ -244,3 +244,3 @@ ) -> None:
 
-    def get_data(self, select_cond: List[tuple]) -> Tuple[np.ndarray, np.ndarray]:
+    def get_data(self, select_cond: list[tuple]) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -284,3 +284,3 @@ Return X and y, where X is the design matrix containing correlation
         tab_name: str,
-        cf_names: Optional[List[str]] = None,
+        cf_names: list[str] | None = None,
         order: int = 1,
@@ -294,3 +294,3 @@ ) -> None:
 
-    def get_data(self, select_cond: List[tuple]) -> Tuple[np.ndarray, np.ndarray]:
+    def get_data(self, select_cond: list[tuple]) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -328,3 +328,3 @@ Return X and y, where X is the design matrix containing correlation
     def __init__(
-        self, db_name: str, tab_name: str, cf_names: Optional[List[str]] = None, order: int = 1
+        self, db_name: str, tab_name: str, cf_names: list[str] | None = None, order: int = 1
     ) -> None:
@@ -350,3 +350,3 @@ super().__init__(db_name)
     @staticmethod
-    def _is_matrix_representable(id_cf_names: Dict[int, List[str]]) -> bool:
+    def _is_matrix_representable(id_cf_names: dict[int, list[str]]) -> bool:
         """
@@ -356,3 +356,3 @@ Check if the extracted correlation functions can be represented as a
         """
-        reference = sorted(list(id_cf_names.values())[0])
+        reference = sorted(next(iter(id_cf_names.values())))
 
@@ -370,3 +370,3 @@ # Check the names of the correlation functions is equal for
     @staticmethod
-    def _minimum_common_cf_set(id_cf_names: Dict[int, List[str]]) -> Set[str]:
+    def _minimum_common_cf_set(id_cf_names: dict[int, list[str]]) -> set[str]:
         """
@@ -376,3 +376,3 @@ Returns the minimum set of correlation functions that exists for all
         """
-        common_cf = set(list(id_cf_names.values())[0])
+        common_cf = set(next(iter(id_cf_names.values())))
         for v in id_cf_names.values():
@@ -665,3 +665,3 @@ common_cf = common_cf.intersection(v)
         tab_name: str,
-        cf_names: Optional[List[str]] = None,
+        cf_names: list[str] | None = None,
         order: int = 1,
@@ -674,3 +674,3 @@ ) -> None:
 
-    def get_data(self, select_cond: List[tuple]) -> Tuple[np.ndarray, np.ndarray]:
+    def get_data(self, select_cond: list[tuple]) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -693,3 +693,3 @@ Return X and y, where X is the design matrix containing correlation
 
-def extract_num_atoms(cur: sqlite3.Cursor, ids: Set[int]) -> Dict[int, int]:
+def extract_num_atoms(cur: sqlite3.Cursor, ids: set[int]) -> dict[int, int]:
     """
@@ -744,5 +744,5 @@ Extract the number of atoms for all ids
         tab_name: str,
-        cf_names: List[str],
-        order: Optional[int] = 0,
-        properties: Tuple[str] = ("energy", "pressure"),
+        cf_names: list[str],
+        order: int | None = 0,
+        properties: tuple[str] = ("energy", "pressure"),
         cf_order: int = 1,
@@ -758,3 +758,3 @@ ) -> None:
 
-    def build(self, ids: List[int]) -> np.ndarray:
+    def build(self, ids: list[int]) -> np.ndarray:
         """
@@ -877,3 +877,3 @@ Construct the design matrix and the target value required to fit a
 
-    def _extract_key(self, ids: Set[int], key: str) -> Dict[int, float]:
+    def _extract_key(self, ids: set[int], key: str) -> dict[int, float]:
         """
@@ -900,3 +900,3 @@ Extract a key from the database for the ids in the passed set that
 
-    def get_data(self, select_cond: List[tuple]) -> Tuple[np.ndarray, np.ndarray]:
+    def get_data(self, select_cond: list[tuple]) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -914,3 +914,3 @@ Return the design matrix and the target values for the entries
 
-    def groups(self) -> List[int]:
+    def groups(self) -> list[int]:
         """
@@ -917,0 +917,0 @@ Return the group of each rows.

clease/data_normalizer.py

@@ -1,3 +0,1 @@
-from typing import Tuple
-
 import numpy as np
@@ -24,3 +22,3 @@
 
-    def normalize(self, X: np.ndarray, y: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
+    def normalize(self, X: np.ndarray, y: np.ndarray) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -27,0 +25,0 @@ Normalizes the each column of X to zero mean and unit variance. y is

clease/datastructures/figure.py

 """This module defines the "Figure" class, which is a collection of FourVector objects."""
-from typing import Any, Iterable, Optional, Tuple
 
+from collections.abc import Iterable
+from typing import Any
+
 import ase
@@ -45,3 +47,3 @@ import attr
 
-    components: Tuple[FourVector] = attr.field(
+    components: tuple[FourVector] = attr.field(
         converter=_convert_figure, validator=attr.validators.instance_of(tuple)
@@ -57,8 +59,7 @@ )
                 raise TypeError(
-                    f"All values must FourVector type, got {value} "
-                    f"of type {type(v)} in index {ii}."
+                    f"All values must FourVector type, got {value} of type {type(v)} in index {ii}."
                 )
 
     def to_cartesian(
-        self, prim: ase.Atoms, transposed_cell: Optional[np.ndarray] = None
+        self, prim: ase.Atoms, transposed_cell: np.ndarray | None = None
     ) -> np.ndarray:
@@ -84,3 +85,3 @@ """Get the Figure in terms of the cartesian coordinates, as defined
 
-    def get_diameter(self, prim: ase.Atoms, transposed_cell: Optional[np.ndarray] = None) -> float:
+    def get_diameter(self, prim: ase.Atoms, transposed_cell: np.ndarray | None = None) -> float:
         """Calculate the diameter of the figure, as the maximum distance to the
@@ -87,0 +88,0 @@ geometric center of the figure in cartesian coordinates.

clease/datastructures/four_vector.py

@@ -6,3 +6,3 @@ from __future__ import annotations
 from itertools import product
-from typing import List, NamedTuple, Optional, Tuple
+from typing import NamedTuple
 
@@ -34,3 +34,3 @@ from ase import Atoms
 
-    def to_cartesian(self, prim: Atoms, transposed_cell: Optional[np.ndarray] = None) -> np.ndarray:
+    def to_cartesian(self, prim: Atoms, transposed_cell: np.ndarray | None = None) -> np.ndarray:
         """Convert the four vector into cartesian coordinates
@@ -67,3 +67,3 @@
 
-    def to_tuple(self) -> Tuple[int]:
+    def to_tuple(self) -> tuple[int]:
         """Get the tuple representation of the four-vector"""
@@ -136,3 +136,3 @@ return attr.astuple(self)
 
-def _make_grid(bbox: _Box, prim: Atoms) -> Tuple[List[FourVector], np.ndarray]:
+def _make_grid(bbox: _Box, prim: Atoms) -> tuple[list[FourVector], np.ndarray]:
     """
@@ -165,3 +165,3 @@ Make a grid bounded by bbox using repetitions of the atomic positions in prim
 
-def construct_four_vectors(prim: Atoms, supercell: Atoms) -> List[FourVector]:
+def construct_four_vectors(prim: Atoms, supercell: Atoms) -> list[FourVector]:
     """
@@ -168,0 +168,0 @@ Calculate the all four vector corresponding to the passed Atoms object. A four-vector u is

clease/datastructures/mc_step.py

@@ -1,2 +0,2 @@
-from typing import Any, Dict
+from typing import Any
 
@@ -38,3 +38,3 @@ import attr
     # Is not included in equality comparison.
-    other: Dict[str, Any] = attr.field(default=attr.Factory(dict), eq=False)
+    other: dict[str, Any] = attr.field(default=attr.Factory(dict), eq=False)
 
@@ -41,0 +41,0 @@ @property

clease/datastructures/system_changes.py

@@ -0,3 +1,3 @@
+from collections.abc import Sequence
 import typing
-from typing import Sequence, Tuple
 
@@ -54,3 +54,3 @@ import ase
 
-    def astuple(self) -> Tuple[int, str, str, str]:
+    def astuple(self) -> tuple[int, str, str, str]:
         """Return the tuple representation of the SystemChange."""
@@ -57,0 +57,0 @@ return attr.astuple(self, recurse=False)

clease/datastructures/trans_matrix.py

 from __future__ import annotations
 
-from typing import Dict, List
-
 import attr
@@ -19,3 +17,3 @@ import numpy as np
 
-    trans_matrix: List[Dict[int, int]] = attr.field()
+    trans_matrix: list[dict[int, int]] = attr.field()
 
@@ -22,0 +20,0 @@ def key_array(self) -> np.ndarray:

clease/db_util.py

 from collections import defaultdict
+from collections.abc import Sequence
 from datetime import datetime
 import logging
-from typing import Dict, List, Sequence, Set, Union
 
@@ -51,3 +51,3 @@ import ase
 # Type aliasing
-_TABLE = Dict[str, Union[int, float, str, np.integer, np.floating]]  # Alias for a table
+_TABLE = dict[str, int | float | str | np.integer | np.floating]  # Alias for a table
 _ROW = ase.db.row.AtomsRow
@@ -239,6 +239,4 @@ _CONNECTION = ase.db.core.Database
         raise ValueError(
-            (
-                f"Number of table names should match number of tables. Got {len(table_names)} "
-                f"table names, but {len(tables)} tables."
-            )
+            f"Number of table names should match number of tables. Got {len(table_names)} "
+            f"table names, but {len(tables)} tables."
         )
@@ -287,3 +285,3 @@
 
-def get_all_cf_names(db_name: str, tab_name: str) -> Set[str]:
+def get_all_cf_names(db_name: str, tab_name: str) -> set[str]:
     """
@@ -309,3 +307,3 @@ Return a list with all correlation function names
 
-def get_all_cf(db_name: str, tab_name: str, db_id: int) -> Dict[str, float]:
+def get_all_cf(db_name: str, tab_name: str, db_id: int) -> dict[str, float]:
     """
@@ -333,3 +331,3 @@ Return all correlation functions associated with an entry in the database
 
-def get_cf_tables(db_name: str) -> List[str]:
+def get_cf_tables(db_name: str) -> list[str]:
     """
@@ -336,0 +334,0 @@ Return a list with table names that contain correlation functions

clease/evaluate.py

 """Module that fits ECIs to energy data."""
+
 from collections import Counter, defaultdict
+from collections.abc import Sequence
 import json
@@ -8,3 +10,2 @@ import logging as lg
 import sys
-from typing import Dict, List, Optional, Sequence
 
@@ -107,3 +108,3 @@ from ase.db import connect
         num_repetitions=1,
-        normalization_symbols: Optional[Sequence[str]] = None,
+        normalization_symbols: Sequence[str] | None = None,
     ):
@@ -179,3 +180,3 @@ """Initialize the Evaluate class."""
 
-    def set_normalization(self, normalization_symbols: Optional[Sequence[str]] = None) -> None:
+    def set_normalization(self, normalization_symbols: Sequence[str] | None = None) -> None:
         """Set the energy normalization factor, e.g. to normalize the final energy reports
@@ -345,3 +346,3 @@ in energy per metal atom, rather than energy per atom (i.e. every atom).
 
-    def get_eci_dict(self, cutoff_tol: float = 1e-14) -> Dict[str, float]:
+    def get_eci_dict(self, cutoff_tol: float = 1e-14) -> dict[str, float]:
         """Determine cluster names and their corresponding ECI value and return
@@ -355,3 +356,3 @@ them in a dictionary format.
         Returns:
-            Dict[str, float]: Dictionary with the CF names and the corresponding
+            dict[str, float]: Dictionary with the CF names and the corresponding
                 ECI value.
@@ -369,3 +370,3 @@ """
 
-    def load_eci_dict(self, eci_dict: Dict[str, float]) -> None:
+    def load_eci_dict(self, eci_dict: dict[str, float]) -> None:
         """Load the ECI's from a dictionary. Any ECI's which are missing
@@ -392,3 +393,3 @@ from the internal cf_names list are assumed to be 0.
         full_fname = add_file_extension(fname, ".json")
-        with open(full_fname, "r") as fd:
+        with open(full_fname) as fd:
             self.load_eci_dict(json.load(fd))
@@ -485,3 +486,3 @@
             0.01,
-            cv_name + f" = {cv:.3f} meV/atom \n" f"RMSE = {self.rmse() * 1000.0:.3f} meV/atom",
+            cv_name + f" = {cv:.3f} meV/atom \nRMSE = {self.rmse() * 1000.0:.3f} meV/atom",
             verticalalignment="bottom",
@@ -669,7 +670,7 @@ horizontalalignment="right",
                 raise TypeError(
-                    "Each entry in fitting_schemes should be an " "instance of LinearRegression"
+                    "Each entry in fitting_schemes should be an instance of LinearRegression"
                 )
             elif not scheme.is_scalar():
                 raise TypeError(
-                    "plot_CV only supports the fitting schemes " "with a scalar paramater."
+                    "plot_CV only supports the fitting schemes with a scalar paramater."
                 )
@@ -684,7 +685,7 @@
                 raise TypeError(
-                    "Each entry in fitting_schemes should be an " "instance of LinearRegression"
+                    "Each entry in fitting_schemes should be an instance of LinearRegression"
                 )
             elif not scheme.is_scalar():
                 raise TypeError(
-                    "plot_CV only supports the fitting schemes " "with a scalar paramater."
+                    "plot_CV only supports the fitting schemes with a scalar paramater."
                 )
@@ -728,3 +729,3 @@
 
-    def cv_for_alpha(self, alphas: List[float]) -> None:
+    def cv_for_alpha(self, alphas: list[float]) -> None:
         """
@@ -840,5 +841,3 @@ Calculate the CV scores for alphas using the fitting scheme
             0.01,
-            f"min. CV score:\n"
-            f"alpha = {min_alpha:.10f} \n"
-            f"CV = {min_cv * 1000.0:.3f} meV/atom",
+            f"min. CV score:\nalpha = {min_alpha:.10f} \nCV = {min_cv * 1000.0:.3f} meV/atom",
             verticalalignment="bottom",
@@ -1124,4 +1123,5 @@ horizontalalignment="left",
                 raise ValueError(
-                    "Inconsistent length of max_dia, size: {}, expected at most: {}".format(
-                        size, max_size_expected
+                    (
+                        f"Inconsistent length of max_dia, size: {size}, ",
+                        f"expected at most: {max_size_expected}",
                     )
@@ -1143,3 +1143,3 @@ )
 
-    def generalization_error(self, validation_id: List[int]):
+    def generalization_error(self, validation_id: list[int]):
         """
@@ -1223,3 +1223,3 @@ Estimate the generalization error to new datapoints
 
-    def get_eci_by_size(self) -> Dict[str, Dict[str, list]]:
+    def get_eci_by_size(self) -> dict[str, dict[str, list]]:
         """
@@ -1226,0 +1226,0 @@ Classify distance, eci and cf_name according to cluster body size

clease/geometry.py

 """Module for tools pertaining to geometry of atoms and cells."""
+
 from ase import Atoms
@@ -3,0 +4,0 @@ import numpy as np

clease/gramSchmidthMonomials.py

@@ -1,3 +0,1 @@
-from typing import Dict, List
-
 import numpy as np
@@ -49,3 +47,3 @@
 
-    def spin_dict(self, symbols: List[str]) -> Dict[str, int]:
+    def spin_dict(self, symbols: list[str]) -> dict[str, int]:
         """
@@ -57,3 +55,3 @@ Return a dictionary with the spin variable for each of the
 
-    def basis_functions(self, symbols: List[str]) -> List[Dict[str, float]]:
+    def basis_functions(self, symbols: list[str]) -> list[dict[str, float]]:
         """
@@ -142,3 +140,3 @@ Construct the spin dictionary from a list of symbols
 
-    def norm(self, bf: List[float]) -> float:
+    def norm(self, bf: list[float]) -> float:
         """
@@ -154,3 +152,3 @@ Computes the norm of a basis function
 
-    def eval(self, bf: List[float], valueIndex: int) -> float:
+    def eval(self, bf: list[float], valueIndex: int) -> float:
         """
@@ -157,0 +155,0 @@ Evaluates the basis function

clease/interactive_plot.py

@@ -0,3 +1,3 @@
+from collections.abc import Callable, Sequence
 from tkinter import TclError
-from typing import Callable, List, Sequence, Union
 
@@ -64,3 +64,3 @@ from ase.db import connect
 
-    def __init__(self, fig, annotated_axes: Union[AnnotatedAx, Sequence[AnnotatedAx]]):
+    def __init__(self, fig, annotated_axes: AnnotatedAx | Sequence[AnnotatedAx]):
         self.fig = fig
@@ -89,3 +89,3 @@ self.annotated_axes = annotated_axes
 
-    def get_mpl_events(self) -> List[MPLEvent]:
+    def get_mpl_events(self) -> list[MPLEvent]:
         event = MPLEvent("motion_notify_event", self.hover)
@@ -196,3 +196,3 @@ return [event]
 class ShowStructureOnClick(InteractivePlot):
-    def __init__(self, fig, axes: Union[AnnotatedAx, Sequence[AnnotatedAx]], db_name: str):
+    def __init__(self, fig, axes: AnnotatedAx | Sequence[AnnotatedAx], db_name: str):
         self.db_name = db_name
@@ -205,3 +205,3 @@ self.active_images = Images()
 
-    def get_mpl_events(self) -> List[MPLEvent]:
+    def get_mpl_events(self) -> list[MPLEvent]:
         event = MPLEvent("button_press_event", self.on_click)
@@ -208,0 +208,0 @@ events = super().get_mpl_events()

clease/jsonio.py

 import json
-from typing import Any, Dict
+from typing import Any
 
@@ -166,3 +166,3 @@ from ase.io import jsonio as aseio
 
-    def todict(self) -> Dict[str, Any]:
+    def todict(self) -> dict[str, Any]:
         """Convert into a dictionary representation."""
@@ -173,4 +173,4 @@ # Disable recursive, since the json module takes care of that
     @classmethod
-    def from_dict(cls, dct: Dict[str, Any]):
+    def from_dict(cls, dct: dict[str, Any]):
         """Load an instance of the class from a dictionary."""
         return cls(**dct)

clease/logging_utils.py

@@ -0,5 +1,6 @@
+from collections.abc import Iterator
 from contextlib import contextmanager
 import logging
 import sys
-from typing import Iterator, Optional, TextIO
+from typing import TextIO
 
@@ -15,3 +16,3 @@ __all__ = ("log_stream", "log_stream_context", "get_root_clease_logger")
 def log_stream(
-    level: Optional[int] = None, stream: TextIO = sys.stdout, fmt: Optional[str] = None
+    level: int | None = None, stream: TextIO = sys.stdout, fmt: str | None = None
 ) -> logging.StreamHandler:
@@ -54,3 +55,3 @@ """Helper function to enable CLEASE logging to a stream. Default stream is stdout.
 def log_stream_context(
-    level: Optional[int] = None, stream: TextIO = sys.stdout, fmt: Optional[str] = None
+    level: int | None = None, stream: TextIO = sys.stdout, fmt: str | None = None
 ) -> Iterator[logging.StreamHandler]:
@@ -98,3 +99,3 @@ """Context which temporarily adds a stream handler to the root CLEASE logger.
 def _make_stream_handler(
-    level: Optional[int], stream: TextIO, fmt: Optional[str] = None
+    level: int | None, stream: TextIO, fmt: str | None = None
 ) -> logging.StreamHandler:
@@ -118,3 +119,3 @@ """Helper function to create a stream handler with a specified level and stream.
 
-def _get_effective_level(level: Optional[int]) -> int:
+def _get_effective_level(level: int | None) -> int:
     """Helper function to get the effective level - if the level is None,
@@ -121,0 +122,0 @@ we use the current level of the CLEASEroot logger, otherwise use the provided level

clease/montecarlo/averager.py

@@ -0,1 +1,4 @@
+from typing_extensions import Self
+
+
 class Averager:
@@ -19,3 +22,3 @@ """
 
-    def __iadd__(self, value):
+    def __iadd__(self, value) -> Self:
         """
@@ -58,3 +61,3 @@ += Operator
 
-    def __idiv__(self, number):
+    def __idiv__(self, number) -> Self:
         """Divide by number.
@@ -61,0 +64,0 @@

clease/montecarlo/barrier_models.py

 from abc import ABC, abstractmethod
-from typing import Dict, Union
 
@@ -17,3 +16,3 @@ from ase import Atoms
 
-    def __call__(self, system: Union[Atoms, MCEvaluator], system_changes: SystemChanges) -> float:
+    def __call__(self, system: Atoms | MCEvaluator, system_changes: SystemChanges) -> float:
         """Get the height of the barrier"""
@@ -61,3 +60,3 @@ # Ensure we have a valid evaluator object
         self,
-        dilute_barrier: Dict[str, float],
+        dilute_barrier: dict[str, float],
         alpha: float = 0.5,
@@ -64,0 +63,0 @@ ):

clease/montecarlo/base.py

 from abc import ABC
-from typing import Union
 
@@ -17,3 +16,3 @@ from ase import Atoms
     Args:
-        system (Union[Atoms, MCEvaluator]): Either an ASE Atoms object
+        system (Atoms | MCEvaluator): Either an ASE Atoms object
             with an attached calculator, or a pre-initialized
@@ -25,3 +24,3 @@ :class:`~clease.montecarlo.mc_evaluator.MCEvaluator`
 
-    def __init__(self, system: Union[Atoms, MCEvaluator], temp: float):
+    def __init__(self, system: Atoms | MCEvaluator, temp: float):
         self.evaluator = system
@@ -44,3 +43,3 @@ self.temperature = temp
     @evaluator.setter
-    def evaluator(self, value: Union[Atoms, MCEvaluator]) -> None:
+    def evaluator(self, value: Atoms | MCEvaluator) -> None:
         """Set the evaluator object. If the value is an Atoms object,
@@ -47,0 +46,0 @@ the evaluator is created on basis of the attached calculator object.

clease/montecarlo/bias_potential.py

 # XXX: Some funny imports here. This file needs to be cleaned up some
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -58,3 +58,3 @@ import numpy as np
             dill.dump(self, outfile)
-        print(f"Pseudo binary free energy bias potential written to " f"{fname}")
+        print(f"Pseudo binary free energy bias potential written to {fname}")
 
@@ -179,3 +179,3 @@ @staticmethod
         self._conc_init = pseudo_bin_conc_init
-        super(PseudoBinaryFreeEnergyBias, self).__init__(reac_crd, free_eng)
+        super().__init__(reac_crd, free_eng)
 
@@ -187,5 +187,3 @@ @property
         if not isinstance(self._conc_init, PseudoBinaryConcInitializer):
-            raise TypeError(
-                "pseudo_bin_conc_init has to be of type " "PseudoBinaryConcInitializer!"
-            )
+            raise TypeError("pseudo_bin_conc_init has to be of type PseudoBinaryConcInitializer!")
         return self._conc_init
@@ -198,5 +196,3 @@
         if not isinstance(init, PseudoBinaryConcInitializer):
-            raise TypeError(
-                "pseudo_bin_conc_init has to be of type " "PseudoBinaryConcInitializer!"
-            )
+            raise TypeError("pseudo_bin_conc_init has to be of type PseudoBinaryConcInitializer!")
         self._conc_init = init
@@ -255,3 +251,3 @@
     def __init__(self, cov_range=None, reac_crd=(), free_eng=()):
-        super(CovarianceBiasPotential, self).__init__(reac_crd, free_eng)
+        super().__init__(reac_crd, free_eng)
         self._cov_range = cov_range
@@ -258,0 +254,0 @@

clease/montecarlo/binned_bias_potential.py

@@ -1,2 +0,2 @@
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -3,0 +3,0 @@ import numpy as np

clease/montecarlo/constraints/constrain_swap_by_basis.py

@@ -1,2 +0,2 @@
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -23,3 +23,3 @@ from ase import Atoms
 
-    index_by_basis: List[List[int]]
+    index_by_basis: list[list[int]]
         Indices ordered by basis (same as ``index_by_basis`` parameter in
@@ -26,0 +26,0 @@ the :class:`~clease.settings.settings.ClusterExpansionSettings`

clease/montecarlo/constraints/fixed_indices.py

@@ -1,2 +0,2 @@
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -3,0 +3,0 @@ from clease.datastructures import SystemChanges

clease/montecarlo/gaussian_kernel_bias_potential.py

@@ -1,2 +0,2 @@
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -3,0 +3,0 @@ import numpy as np

clease/montecarlo/kinetic_monte_carlo.py

@@ -0,5 +1,5 @@
+from collections.abc import Sequence
 import logging
 import random
 import time
-from typing import List, Sequence, Tuple, Union
 import warnings
@@ -52,3 +52,3 @@
         self,
-        system: Union[Atoms, MCEvaluator],
+        system: Atoms | MCEvaluator,
         temp: float,
@@ -80,3 +80,3 @@ barrier: BarrierModel,
 
-    def _update_epr(self, current: int, choice: int, swaps: List[int], cum_rates: np.ndarray):
+    def _update_epr(self, current: int, choice: int, swaps: list[int], cum_rates: np.ndarray):
         """
@@ -104,3 +104,3 @@ Update the EPR tracker (if set).
 
-    def _rates(self, vac_idx: int) -> Tuple[List[int], np.ndarray]:
+    def _rates(self, vac_idx: int) -> tuple[list[int], np.ndarray]:
         """
@@ -141,3 +141,3 @@ Return the rates for all possible swaps for a vacancy
 
-    def _mc_step(self, vac_idx: int, step_no: int) -> Tuple[int, MCStep]:
+    def _mc_step(self, vac_idx: int, step_no: int) -> tuple[int, MCStep]:
         """
@@ -206,3 +206,3 @@ Perform an MC step and return the new index of the moving vacancy
             raise ValueError(
-                f"Index {vac_idx} is not a vacancy. " f"Symbol: {self.atoms[vac_idx].symbol}"
+                f"Index {vac_idx} is not a vacancy. Symbol: {self.atoms[vac_idx].symbol}"
             )
@@ -209,0 +209,0 @@

clease/montecarlo/kmc_events.py

 from abc import ABCMeta
-from typing import List
 
@@ -15,3 +14,3 @@ from ase import Atoms
 
-    def get_swaps(self, atoms: Atoms, vac_idx: int) -> List[int]:
+    def get_swaps(self, atoms: Atoms, vac_idx: int) -> list[int]:
         raise NotImplementedError(f"get_swaps has not been implemented for class {self.name}")
@@ -35,3 +34,3 @@
 
-    def get_swaps(self, atoms: Atoms, vac_idx: int) -> List[int]:
+    def get_swaps(self, atoms: Atoms, vac_idx: int) -> list[int]:
         return self.nl[vac_idx]

clease/montecarlo/mc_evaluator.py

 import logging
-from typing import Optional, Union
 
@@ -31,3 +30,3 @@ from ase import Atoms
 
-    def get_energy(self, applied_changes: Optional[SystemChanges] = None) -> float:
+    def get_energy(self, applied_changes: SystemChanges | None = None) -> float:
         """Evaluate the energy of a system.
@@ -81,3 +80,3 @@ If a change is sufficiently local/small, it there, in some situations,
 
-    def keep_system_changes(self, system_changes: Optional[SystemChanges] = None) -> None:
+    def keep_system_changes(self, system_changes: SystemChanges | None = None) -> None:
         """A set of system changes are to be kept. Perform necessary actions to prepare
@@ -124,3 +123,3 @@ for a new evaluation."""
 
-    def get_energy(self, applied_changes: Optional[SystemChanges] = None) -> float:
+    def get_energy(self, applied_changes: SystemChanges | None = None) -> float:
         return self.calc.get_energy()
@@ -156,3 +155,3 @@
 
-    def keep_system_changes(self, system_changes: Optional[SystemChanges] = None) -> None:
+    def keep_system_changes(self, system_changes: SystemChanges | None = None) -> None:
         """A set of system changes are to be kept. Perform necessary actions to prepare
@@ -197,3 +196,3 @@ for a new evaluation."""
 
-def construct_evaluator(system: Union[Atoms, MCEvaluator]) -> MCEvaluator:
+def construct_evaluator(system: Atoms | MCEvaluator) -> MCEvaluator:
     """Helper function for constructing a new evaluator object, either by passing
@@ -203,3 +202,3 @@ in an ASE atoms object or an explicit evaluator object.
     Args:
-        system (Union[Atoms, MCEvaluator]): If the system is an Atoms object,
+        system (Atoms | MCEvaluator): If the system is an Atoms object,
             then a new :class:`~clease.montecarlo.mc_evaluator.MCEvaluator`
@@ -206,0 +205,0 @@ object is created on the basis of the attached calculator.

clease/montecarlo/metadynamics_sampler.py

@@ -145,3 +145,3 @@ from copy import deepcopy
             raise PeakNotAcceptedError(
-                ("Observer does not accept peak as a keyword argument to __call__")
+                "Observer does not accept peak as a keyword argument to __call__"
             )
@@ -175,3 +175,3 @@
             if time.perf_counter() - now > self.log_freq:
-                msg = f"Sweep no. {int(counter/len(self.mc.atoms))} "
+                msg = f"Sweep no. {int(counter / len(self.mc.atoms))} "
                 msg += f"Average visits: {self.progress_info['mean']:.2e}. "
@@ -178,0 +178,0 @@ msg += f"Min/avg: {self.progress_info['minval']:.2e} "

clease/montecarlo/montecarlo.py

 """Monte Carlo method for ase."""
+
 from collections import Counter
+from collections.abc import Iterator
 from datetime import datetime
@@ -8,3 +10,3 @@ import logging
 import time
-from typing import Any, Dict, Iterator, Optional, Union
+from typing import Any
 
@@ -32,3 +34,3 @@ from ase import Atoms
     Args:
-        system (Union[ase.Atoms, MCEvaluator]): Either an ASE Atoms object
+        system (ase.Atoms | MCEvaluator): Either an ASE Atoms object
             with an attached calculator, or a pre-initialized
@@ -47,5 +49,5 @@ :class:`~clease.montecarlo.mc_evaluator.MCEvaluator`
         self,
-        system: Union[Atoms, MCEvaluator],
+        system: Atoms | MCEvaluator,
         temp: float,
-        generator: Optional[TrialMoveGenerator] = None,
+        generator: TrialMoveGenerator | None = None,
     ):
@@ -266,3 +268,3 @@ # We cannot cause an energy calculation trigger in init,
     @property
-    def meta_info(self) -> Dict[str, str]:
+    def meta_info(self) -> dict[str, str]:
         """Return dict with meta info."""
@@ -284,3 +286,3 @@ # Get the timestamp with millisecond precision
 
-    def get_thermodynamic_quantities(self) -> Dict[str, Any]:
+    def get_thermodynamic_quantities(self) -> dict[str, Any]:
         """Compute thermodynamic quantities."""
@@ -310,3 +312,3 @@ quantities = {}
 
-    def _get_obs_averages(self) -> Dict[str, Any]:
+    def _get_obs_averages(self) -> dict[str, Any]:
         """Get average measurements from observers"""
@@ -373,3 +375,3 @@ obs_avgs = {}
 
-    def count_atoms(self) -> Dict[str, int]:
+    def count_atoms(self) -> dict[str, int]:
         """Count the number of each element."""
@@ -376,0 +378,0 @@ return dict(Counter(self.atoms.symbols))

clease/montecarlo/observers/acceptance_rate.py

@@ -1,3 +0,1 @@
-from typing import Dict
-
 from clease.datastructures import SystemChanges
@@ -48,3 +46,3 @@ from clease.montecarlo.observers.mc_observer import MCObserver
 
-    def get_averages(self) -> Dict[str, float]:
+    def get_averages(self) -> dict[str, float]:
         """
@@ -51,0 +49,0 @@ Return dictionary with the rate such that it is added to thermodynaic quantities

clease/montecarlo/observers/buffered_array.py

 from pathlib import Path
-from typing import Optional, Union
 
@@ -19,3 +18,3 @@ import numpy as np
 
-    def __init__(self, size: int = 1000, fname: Optional[Union[str, Path]] = None):
+    def __init__(self, size: int = 1000, fname: str | Path | None = None):
         self._buffer = np.zeros(size)
@@ -22,0 +21,0 @@ self._next = 0

clease/montecarlo/observers/concentration_observer.py

@@ -1,3 +0,1 @@
-from typing import Dict
-
 import ase
@@ -80,3 +78,3 @@
 
-    def get_averages(self) -> Dict[str, float]:
+    def get_averages(self) -> dict[str, float]:
         mean_conc = self.avg_conc.mean
@@ -83,0 +81,0 @@ var_conc = self.avg_conc_sq.mean - mean_conc**2

clease/montecarlo/observers/entropy_prod_rate.py

 from pathlib import Path
-from typing import List, Union
 
@@ -35,3 +34,3 @@ import numpy as np
 
-    def __init__(self, buffer_length: int = 10000, logfile: Union[str, Path] = "epr.txt"):
+    def __init__(self, buffer_length: int = 10000, logfile: str | Path = "epr.txt"):
         self._buffer = BufferedArray(size=buffer_length, fname=logfile)
@@ -59,3 +58,3 @@ self.prev_swap = -1
 
-    def update(self, current: int, choice: int, cumulative_rates: np.ndarray, swaps: List[int]):
+    def update(self, current: int, choice: int, cumulative_rates: np.ndarray, swaps: list[int]):
         """
@@ -62,0 +61,0 @@ Update the buffer

clease/montecarlo/observers/move_obs.py

 from __future__ import annotations
 
-from typing import Iterator, List
+from collections.abc import Iterator
 
@@ -45,3 +45,3 @@ import ase
 
-    def reconstruct(self) -> List[ase.Atoms]:
+    def reconstruct(self) -> list[ase.Atoms]:
         """Rebuild the atoms objects as defined by the observed changes."""
@@ -48,0 +48,0 @@ return list(self.reconstruct_iter())

clease/montecarlo/observers/multi_state_sgc_observer.py

@@ -1,3 +0,1 @@
-from typing import Dict, List, Tuple
-
 from ase.units import kB
@@ -22,3 +20,3 @@ import numpy as np
 
-    def __init__(self, temp: float, chem_pot: Dict[str, float]):
+    def __init__(self, temp: float, chem_pot: dict[str, float]):
         self.temp = temp
@@ -40,3 +38,3 @@ self.chem_pot = chem_pot
         prefix += "_".join(
-            f"{k}{sign_indicator(v)}{int(1000.0*abs(v))}" for k, v in self.chem_pot.items()
+            f"{k}{sign_indicator(v)}{int(1000.0 * abs(v))}" for k, v in self.chem_pot.items()
         )
@@ -96,3 +94,3 @@ return prefix
 
-    def __init__(self, ref_state: SGCState, thermo_states: List[SGCState], calc: Clease):
+    def __init__(self, ref_state: SGCState, thermo_states: list[SGCState], calc: Clease):
         super().__init__()
@@ -110,3 +108,3 @@ self.thermo_states = thermo_states
 
-    def __call__(self, system_changes: List[Tuple[int, str, str]]):
+    def __call__(self, system_changes: list[tuple[int, str, str]]):
         """
@@ -130,3 +128,3 @@ Observers are called after system update, thus the calcualtor already
 
-    def get_averages(self) -> Dict[str, float]:
+    def get_averages(self) -> dict[str, float]:
         """
@@ -133,0 +131,0 @@ Return a dictionary with the calculated averages

clease/montecarlo/observers/snapshot_observer.py

@@ -1,2 +0,2 @@
-from typing import Sequence
+from collections.abc import Sequence
 
@@ -3,0 +3,0 @@ from ase.atoms import Atoms

clease/montecarlo/sgc_montecarlo.py

@@ -1,2 +0,3 @@
-from typing import Any, Dict, Optional, Sequence
+from collections.abc import Sequence
+from typing import Any
 
@@ -37,3 +38,3 @@ from ase import Atoms
         symbols: Sequence[str] = (),
-        generator: Optional[TrialMoveGenerator] = None,
+        generator: TrialMoveGenerator | None = None,
         observe_singlets: bool = False,
@@ -105,3 +106,3 @@ ):
     @chemical_potential.setter
-    def chemical_potential(self, chem_pot: Dict[str, float]):
+    def chemical_potential(self, chem_pot: dict[str, float]):
         eci = self.calc.eci
@@ -122,3 +123,3 @@ if any(key not in eci for key in chem_pot):
 
-    def _include_chemical_potential_in_eci(self, chem_pot: Dict[str, float], eci: Dict[str, float]):
+    def _include_chemical_potential_in_eci(self, chem_pot: dict[str, float], eci: dict[str, float]):
         """
@@ -150,3 +151,3 @@ Including the chemical potentials in the ECIs
 
-    def _reset_eci_to_original(self, eci_with_chem_pot: Dict[str, float]):
+    def _reset_eci_to_original(self, eci_with_chem_pot: dict[str, float]):
         """
@@ -174,3 +175,3 @@ Resets the ECIs to their original value
         call_observers: bool = True,
-        chem_pot: Optional[Dict[str, float]] = None,
+        chem_pot: dict[str, float] | None = None,
     ):
@@ -200,3 +201,3 @@ """
 
-    def singlet2composition(self, avg_singlets: Dict[str, float]):
+    def singlet2composition(self, avg_singlets: dict[str, float]):
         """Convert singlets to composition."""
@@ -229,3 +230,3 @@ bf = self.settings.basis_functions
 
-    def get_thermodynamic_quantities(self, reset_eci: bool = False) -> Dict[str, Any]:
+    def get_thermodynamic_quantities(self, reset_eci: bool = False) -> dict[str, Any]:
         """Compute thermodynamic quantities.
@@ -232,0 +233,0 @@

clease/montecarlo/swap_move_index_tracker.py

@@ -0,3 +1,3 @@
+from collections.abc import Sequence
 import random
-from typing import Dict, List, Optional, Sequence
 
@@ -17,4 +17,4 @@ import ase
 
-    def __init__(self, atoms, indices: Optional[Sequence[int]] = None):
-        self.tracker: Dict[str, List[int]] = {}
+    def __init__(self, atoms, indices: Sequence[int] | None = None):
+        self.tracker: dict[str, list[int]] = {}
         self.index_loc = None
@@ -29,6 +29,6 @@ self._last_move = []
     @staticmethod
-    def _unique_symbols_from_atoms(atoms) -> List[str]:
+    def _unique_symbols_from_atoms(atoms) -> list[str]:
         return sorted(set(atoms.symbols))
 
-    def get_unique_symbols(self) -> List[str]:
+    def get_unique_symbols(self) -> list[str]:
         return list(self.tracker.keys())
@@ -46,3 +46,3 @@
 
-    def _init_tracker(self, atoms: ase.Atoms, indices: Optional[Sequence[int]] = None) -> None:
+    def _init_tracker(self, atoms: ase.Atoms, indices: Sequence[int] | None = None) -> None:
         """Initialize the tracker with the numbers."""
@@ -124,5 +124,5 @@ # Cache the unique symbols for faster access
 
-    def get_two_random_symbols(self) -> List[str]:
+    def get_two_random_symbols(self) -> list[str]:
         """Get two different random unique symbols"""
         # random.sample samples without replacement.
         return random.sample(self.unique_symbols, 2)

clease/montecarlo/trial_move_generator.py

 from abc import ABC, abstractmethod
+from collections.abc import Sequence
 import random
 from random import choice
-from typing import Dict, List, Optional, Sequence, Set, Tuple
 
@@ -133,3 +133,3 @@ from ase import Atoms
 
-    def made_changes(self, changes: Sequence[SystemChange]) -> List[SystemChange]:
+    def made_changes(self, changes: Sequence[SystemChange]) -> list[SystemChange]:
         """
@@ -155,5 +155,3 @@ Extract the subset system changes made by an instance of itself.
 
-    def __init__(
-        self, symbols: Set[str], atoms: Atoms, indices: Optional[List[int]] = None, **kwargs
-    ):
+    def __init__(self, symbols: set[str], atoms: Atoms, indices: list[int] | None = None, **kwargs):
         super().__init__(**kwargs)
@@ -176,3 +174,3 @@ self.symbols = symbols
 
-    def _make_possible_flips(self) -> Dict[str, List[str]]:
+    def _make_possible_flips(self) -> dict[str, list[str]]:
         """Compute a map of possible flips, given a site has a particular symbol."""
@@ -184,3 +182,3 @@ possible = {}
 
-    def get_single_trial_move(self) -> List[SystemChange]:
+    def get_single_trial_move(self) -> list[SystemChange]:
         """Get a random flip of an included site into a different element."""
@@ -208,3 +206,3 @@ pos = choice(self.indices)
 
-    def __init__(self, atoms: Atoms, indices: Optional[List[int]] = None, **kwargs):
+    def __init__(self, atoms: Atoms, indices: list[int] | None = None, **kwargs):
         super().__init__(**kwargs)
@@ -217,7 +215,6 @@ self.indices = indices
             raise TooFewElementsError(
-                "After filtering there are less than two symbol type left. "
-                "Must have at least two."
+                "After filtering there are less than two symbol type left. Must have at least two."
             )
 
-    def get_single_trial_move(self) -> List[SystemChange]:
+    def get_single_trial_move(self) -> list[SystemChange]:
         """
@@ -278,7 +275,7 @@ Create a swap move
     @property
-    def generators(self) -> Tuple[SingleTrialMoveGenerator]:
+    def generators(self) -> tuple[SingleTrialMoveGenerator]:
         return (self.flipper, self.swapper)
 
     @property
-    def weights(self) -> Tuple[float]:
+    def weights(self) -> tuple[float]:
         """The probability weights for each generator"""
@@ -351,3 +348,3 @@ return (self.flip_prob, 1.0 - self.flip_prob)
         atoms: Atoms,
-        indices: Optional[Sequence[Sequence[int]]] = None,
+        indices: Sequence[Sequence[int]] | None = None,
         **kwargs,
@@ -354,0 +351,0 @@ ):

clease/mp_logger.py

 """Logger that can be used together with multiprocessing funcions."""
+
 import logging as lg
@@ -3,0 +4,0 @@ import multiprocessing as mp

clease/plot_post_process.py

@@ -1,3 +0,1 @@
-from typing import List, Optional, Tuple
-
 from matplotlib.figure import Figure
@@ -12,3 +10,3 @@ import matplotlib.pyplot as plt
 def plot_fit(
-    evaluate: Evaluate, plot_args: Optional[dict] = None, interactive: bool = False
+    evaluate: Evaluate, plot_args: dict | None = None, interactive: bool = False
 ) -> Figure:
@@ -67,3 +65,3 @@ """
         0.01,
-        cv_name + f" = {cv:.3f} meV/atom\n" f"RMSE = {rmse:.3f} meV/atom",
+        cv_name + f" = {cv:.3f} meV/atom\nRMSE = {rmse:.3f} meV/atom",
         verticalalignment="bottom",
@@ -96,3 +94,3 @@ horizontalalignment="right",
 def plot_fit_residual(
-    evaluate: Evaluate, plot_args: Optional[dict] = None, interactive: bool = False
+    evaluate: Evaluate, plot_args: dict | None = None, interactive: bool = False
 ) -> Figure:
@@ -159,3 +157,3 @@ """
     evaluate: Evaluate,
-    plot_args: Optional[dict] = None,
+    plot_args: dict | None = None,
     ignore_sizes=(),
@@ -245,3 +243,3 @@ interactive: bool = False,
 
-def plot_cv(evaluate: Evaluate, plot_args: Optional[dict] = None) -> Figure:
+def plot_cv(evaluate: Evaluate, plot_args: dict | None = None) -> Figure:
     """
@@ -365,3 +363,3 @@ Figure object of CV values according to alpha values
 
-def _make_annotations_hull(evaluate: Evaluate) -> Tuple[List[str], List[str]]:
+def _make_annotations_hull(evaluate: Evaluate) -> tuple[list[str], list[str]]:
     """Helper function to make annotations for interactive plots."""
@@ -389,3 +387,3 @@ e_pred = evaluate.get_energy_predict()
 
-def _make_annotations_plot_fit(evaluate: Evaluate) -> Tuple[List[str], List[str]]:
+def _make_annotations_plot_fit(evaluate: Evaluate) -> tuple[list[str], list[str]]:
     """Helper function to make annotations for interactive plots."""
@@ -392,0 +390,0 @@ e_pred = evaluate.get_energy_predict()

clease/regression/bayesian_compressive_sensing.py

@@ -361,3 +361,3 @@ from itertools import product
         bayes.fname = fname
-        with open(fname, "r") as infile:
+        with open(fname) as infile:
             data = json.load(infile)
@@ -442,8 +442,8 @@
                 msg = f"Iter: {iteration} "
-                msg += f"RMSE: {1000.0*self.rmse():.3E} "
-                msg += f"LOOCV (approx.): {1000.0*self.estimate_loocv():.3E} "
+                msg += f"RMSE: {1000.0 * self.rmse():.3E} "
+                msg += f"LOOCV (approx.): {1000.0 * self.estimate_loocv():.3E} "
                 msg += f"Num ECI: {self.num_ecis} "
                 msg += f"Lamb: {self.lamb:.3E} "
                 msg += f"Shape lamb: {self.shape_lamb:.3E} "
-                msg += f"Noise: {np.sqrt(1.0/self.inv_variance):.3E}"
+                msg += f"Noise: {np.sqrt(1.0 / self.inv_variance):.3E}"
                 logger.info(msg)
@@ -450,0 +450,0 @@ now = time.perf_counter()

clease/regression/constrained_ridge.py

@@ -1,3 +0,1 @@
-from typing import Tuple
-
 import numpy as np
@@ -41,3 +39,3 @@
 
-    def kkt_system(self, X: np.ndarray, y: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
+    def kkt_system(self, X: np.ndarray, y: np.ndarray) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -44,0 +42,0 @@ Return the Karush-Kuhn-Tucker (KKT) system of equations. By solving

clease/regression/ga_fit.py

@@ -128,3 +128,3 @@ import logging
         individuals = []
-        with open(self.fname, "r") as infile:
+        with open(self.fname) as infile:
             for line in infile:
@@ -131,0 +131,0 @@ individual = np.zeros(self.num_genes, dtype=np.uint8)

clease/regression/generalized_ridge_regression.py

 import logging
-from typing import Optional
 
@@ -45,3 +44,3 @@ import numpy as np
 
-    def __init__(self, alpha: Optional[np.ndarray] = None) -> None:
+    def __init__(self, alpha: np.ndarray | None = None) -> None:
         super().__init__()
@@ -199,6 +198,4 @@ self.alpha = alpha
             logger.warning(
-                (
-                    "Warning! The effective number of parameters is negative. "
-                    "Try to change the initial guess for alpha."
-                )
+                "Warning! The effective number of parameters is negative. "
+                "Try to change the initial guess for alpha."
             )
@@ -205,0 +202,0 @@ logger.info("Best GCV: %.3f", np.sqrt(res.fun))

clease/regression/physical_ridge.py

@@ -0,4 +1,4 @@
+from collections.abc import Callable, Sequence
 import logging
 import time
-from typing import Callable, Dict, List, Optional, Sequence, Tuple, Union
 
@@ -78,6 +78,6 @@ import numpy as np
         lamb_dia: float = 1e-6,
-        size_decay: Union[str, Callable[[int], float]] = "linear",
-        dia_decay: Union[str, Callable[[int], float]] = "linear",
+        size_decay: str | Callable[[int], float] = "linear",
+        dia_decay: str | Callable[[int], float] = "linear",
         normalize: bool = True,
-        cf_names: Optional[List[str]] = None,
+        cf_names: list[str] | None = None,
     ) -> None:
@@ -109,3 +109,3 @@ super().__init__()
 
-    def fit_data(self, X: np.ndarray, y: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
+    def fit_data(self, X: np.ndarray, y: np.ndarray) -> tuple[np.ndarray, np.ndarray]:
         """
@@ -130,3 +130,3 @@ If normalize is True, a normalized version of the passed data is
     @size_decay.setter
-    def size_decay(self, decay: Union[str, Callable[[int], float]]) -> None:
+    def size_decay(self, decay: str | Callable[[int], float]) -> None:
         self._size_decay = get_size_decay(decay)
@@ -139,6 +139,6 @@
     @dia_decay.setter
-    def dia_decay(self, decay: Union[str, Callable[[int], float]]) -> None:
+    def dia_decay(self, decay: str | Callable[[int], float]) -> None:
         self._dia_decay = get_dia_decay(decay)
 
-    def sizes_from_names(self, names: List[str]) -> None:
+    def sizes_from_names(self, names: list[str]) -> None:
         """
@@ -154,3 +154,3 @@ Extract the sizes from a list of correlation function names
 
-    def diameters_from_names(self, names: List[str]) -> None:
+    def diameters_from_names(self, names: list[str]) -> None:
         """
@@ -238,3 +238,3 @@ Extract the diameters from a list of correltion function names
 
-def get_size_decay(decay: Union[str, Callable[[int], float]]) -> Callable[[int], float]:
+def get_size_decay(decay: str | Callable[[int], float]) -> Callable[[int], float]:
     if isinstance(decay, str):
@@ -255,3 +255,3 @@ if decay == "linear":
 
-def get_dia_decay(decay: Union[str, Callable[[int], float]]) -> Callable[[int], float]:
+def get_dia_decay(decay: str | Callable[[int], float]) -> Callable[[int], float]:
     if isinstance(decay, str):
@@ -274,3 +274,3 @@ if decay == "linear":
     phys_ridge: PhysicalRidge,
-    params: Dict,
+    params: dict,
     X: np.ndarray,
@@ -281,3 +281,3 @@ y: np.ndarray,
     groups: Sequence[int] = (),
-) -> Dict:
+) -> dict:
     """
@@ -355,4 +355,4 @@ Estimate the hyper parameters of the Physical Ridge by random search.
             msg = (
-                f"{i} of {num_trials}. CV: {best_cv*1000.0} meV/atom. "
-                f"MSE: {best_mse*1000.0} meV/atom. Params: {best_param}"
+                f"{i} of {num_trials}. CV: {best_cv * 1000.0} meV/atom. "
+                f"MSE: {best_mse * 1000.0} meV/atom. Params: {best_param}"
             )
@@ -359,0 +359,0 @@ logger.info(msg)

clease/regression/regression.py

 """Collection of classess to perform regression."""
-from typing import List, Optional, Union
 
@@ -57,3 +56,3 @@ import numpy as np
 
-    def get_instance_array(self) -> List[object]:
+    def get_instance_array(self) -> list[object]:
         return [LinearRegression()]
@@ -94,3 +93,3 @@
         self,
-        alpha: Union[float, np.ndarray] = 1e-5,
+        alpha: float | np.ndarray = 1e-5,
         penalize_bias_term: bool = False,
@@ -171,4 +170,4 @@ normalize: bool = True,
         num_alpha: int = 10,
-        scale: Optional[str] = "log",
-    ) -> List[object]:
+        scale: str | None = "log",
+    ) -> list[object]:
         if scale == "log":
@@ -225,3 +224,3 @@ alpha = np.logspace(
         alpha_min: float, alpha_max: float, num_alpha: int = 10, scale: str = "log"
-    ) -> List[object]:
+    ) -> list[object]:
         if scale == "log":
@@ -228,0 +227,0 @@ alpha = np.logspace(

clease/regression/sequential_cluster_ridge.py

@@ -64,5 +64,3 @@ import logging
         for i in range(20):
-            print(
-                f"Num. coeff: {len(coeffs[srt_idx[i]]):9d} " f"CV: {cvs[srt_idx[i]]:9.3f} meV/atom"
-            )
+            print(f"Num. coeff: {len(coeffs[srt_idx[i]]):9d} CV: {cvs[srt_idx[i]]:9.3f} meV/atom")
         print("--------------------------------------------")
@@ -69,0 +67,0 @@

clease/settings/atoms_manager.py

 from __future__ import annotations
 
-from typing import List, Optional, Sequence
+from collections.abc import Sequence
 
@@ -19,3 +19,3 @@ from ase import Atoms
 
-    def __init__(self, atoms: Optional[Atoms] = None) -> None:
+    def __init__(self, atoms: Atoms | None = None) -> None:
         self.atoms = atoms
@@ -28,9 +28,6 @@
     @atoms.setter
-    def atoms(self, other: Optional[Atoms]) -> None:
+    def atoms(self, other: Atoms | None) -> None:
         if other is not None and not isinstance(other, Atoms):
             raise TypeError(
-                (
-                    "Trying to set atoms with a non-atoms object. "
-                    f"Expected atoms, got {type(other)}"
-                )
+                f"Trying to set atoms with a non-atoms object. Expected atoms, got {type(other)}"
             )
@@ -45,3 +42,3 @@ self._atoms = other
 
-    def index_by_tag(self) -> List[List[int]]:
+    def index_by_tag(self) -> list[list[int]]:
         """Return atomic indices that are grouped by their tags.
@@ -60,3 +57,3 @@
 
-    def index_by_symbol(self, symbols: List) -> List[List[int]]:
+    def index_by_symbol(self, symbols: list) -> list[list[int]]:
         """Group atomic indices by its atomic symbols.
@@ -91,3 +88,3 @@
 
-    def unique_elements(self, ignore: Sequence[str] = ()) -> List[str]:
+    def unique_elements(self, ignore: Sequence[str] = ()) -> list[str]:
         """Return a list of symbols of unique elements.
@@ -100,3 +97,3 @@
 
-    def single_element_sites(self, allowed_elements: List[List[str]]) -> List[int]:
+    def single_element_sites(self, allowed_elements: list[list[str]]) -> list[int]:
         """
@@ -103,0 +100,0 @@ Return a list of sites that can only be occupied by a single element

clease/settings/concentration.py

@@ -443,3 +443,3 @@ """Class containing a manager for setting up concentrations of species."""
             if ref_element is None:
-                raise RuntimeError(f"Did not find reference element for symbol" f" {variable}")
+                raise RuntimeError(f"Did not find reference element for symbol {variable}")
 
@@ -535,3 +535,3 @@ # Extract the coefficients of other elements with their
 
-        raise ValueError(f"Did not find any basis containing " f"{variable_symbol}")
+        raise ValueError(f"Did not find any basis containing {variable_symbol}")
 
@@ -538,0 +538,0 @@ def _get_col_of_element_in_basis(self, basis, element):

clease/settings/settings_bulk.py

@@ -6,2 +6,3 @@ """Definitions of Cluster Expansion settings for bulk.
 """
+
 from copy import deepcopy
@@ -38,3 +39,3 @@
 
-        concentration (Union[Concentration, dict]):
+        concentration (Concentration | dict):
             Concentration object or dictionary specifying the basis elements and
@@ -86,3 +87,5 @@ concentration range of constituting species
 
-    settings = ClusterExpansionSettings(prim, concentration, **kwargs)
+    settings = ClusterExpansionSettings(
+        prim, concentration, _allow_direct_instantiation=True, **kwargs
+    )
 
@@ -116,3 +119,3 @@ settings.kwargs.update(
 
-        concentration (Union[Concentration, dict]):
+        concentration (Concentration | dict):
             Concentration object or dictionary specifying the basis elements and
@@ -125,3 +128,3 @@ concentration range of constituting species
 
-        basis (List[float]):
+        basis (list[float]):
             List of scaled coordinates.
@@ -170,3 +173,5 @@ Positions of the unique sites corresponding to symbols given
 
-    settings = ClusterExpansionSettings(prim, concentration, **kwargs)
+    settings = ClusterExpansionSettings(
+        prim, concentration, _allow_direct_instantiation=True, **kwargs
+    )
     settings.kwargs.update(
@@ -173,0 +178,0 @@ {

clease/settings/settings_slab.py

@@ -1,2 +0,2 @@
-from typing import Optional, Sequence, Tuple, Union
+from collections.abc import Sequence
 
@@ -17,6 +17,6 @@ from ase import Atoms
 def CESlab(
-    conventional_cell: Union[Atoms, str],
-    miller: Tuple[int],
+    conventional_cell: Atoms | str,
+    miller: tuple[int],
     concentration: Concentration,
-    size: Optional[Sequence[int]] = (1, 1, 1),
+    size: Sequence[int] | None = (1, 1, 1),
     **kwargs,
@@ -52,3 +52,5 @@ ) -> ClusterExpansionSettings:
     # Slab should always have one cell vector along the z-axis
-    settings = ClusterExpansionSettings(prim, concentration, size=size, **kwargs)
+    settings = ClusterExpansionSettings(
+        prim, concentration, size=size, _allow_direct_instantiation=True, **kwargs
+    )
 
@@ -71,3 +73,3 @@ dict_rep = conventional_cell.todict()
 
-def get_prim_slab_cell(conventional_cell: Union[Atoms, str], miller: Tuple[int]) -> Atoms:
+def get_prim_slab_cell(conventional_cell: Atoms | str, miller: tuple[int]) -> Atoms:
     """
@@ -74,0 +76,0 @@ Returns the primitive cell used for slab CE

clease/settings/settings.py

@@ -6,7 +6,9 @@ """Definition of ClusterExpansionSettings Class.
 """
+
 from __future__ import annotations
 
+from collections.abc import Sequence
 from copy import deepcopy
 import logging
-from typing import Any, Dict, List, Optional, Sequence, Set, Union
+from typing import Any
 
@@ -43,10 +45,17 @@ from ase import Atoms
 
+    .. note::
+        This class should not be instantiated directly. Use one of the factory
+        functions instead:
+        - :func:`clease.settings.CEBulk` for bulk materials with crystal structures
+        - :func:`clease.settings.CECrystal` for bulk materials with space groups
+        - :func:`clease.settings.CESlab` for slab/surface calculations
+
     Args:
         prim (Atoms): The primitive atoms object.
 
-        concentration (Union[Concentration, dict]): Concentration object or
+        concentration (Concentration | dict): Concentration object or
             dictionary specifying the basis elements and
             concentration range of constituting species.
 
-        size (List[int] | None, optional): Size of the supercell
+        size (list[int] | None, optional): Size of the supercell
             (e.g., [2, 2, 2] for 2x2x2 cell).
@@ -97,5 +106,5 @@ ``supercell_factor`` is ignored if both ``size`` and ``supercell_factor``
         prim: Atoms,
-        concentration: Union[Concentration, dict],
-        size: Optional[List[int]] = None,
-        supercell_factor: Optional[int] = 27,
+        concentration: Concentration | dict,
+        size: list[int] | None = None,
+        supercell_factor: int | None = 27,
         db_name: str = "clease.db",
@@ -105,3 +114,13 @@ max_cluster_dia: Sequence[float] = (5.0, 5.0, 5.0),
         basis_func_type="polynomial",
+        _allow_direct_instantiation: bool = False,
     ) -> None:
+        if not _allow_direct_instantiation:
+            raise TypeError(
+                "ClusterExpansionSettings cannot be instantiated directly. "
+                "Please use one of the factory functions instead:\n"
+                "  - CEBulk() for bulk materials with crystal structures\n"
+                "  - CECrystal() for bulk materials with space groups\n"
+                "  - CESlab() for slab/surface calculations\n"
+                "See clease.settings documentation for details."
+            )
         self._include_background_atoms = include_background_atoms
@@ -212,3 +231,3 @@ self._cluster_mng = None
     @property
-    def all_elements(self) -> List[str]:
+    def all_elements(self) -> list[str]:
         return sorted([item for row in self.basis_elements for item in row])
@@ -221,3 +240,3 @@
     @property
-    def unique_elements(self) -> List[str]:
+    def unique_elements(self) -> list[str]:
         return sorted(list(set(deepcopy(self.all_elements))))
@@ -230,3 +249,3 @@
     @property
-    def ref_index_trans_symm(self) -> List[int]:
+    def ref_index_trans_symm(self) -> list[int]:
         return [i[0] for i in self.index_by_sublattice]
@@ -247,3 +266,3 @@
     @property
-    def background_indices(self) -> List[int]:
+    def background_indices(self) -> list[int]:
         """Get indices of the background atoms."""
@@ -259,3 +278,3 @@ # check if any basis consists of only one element type
     @property
-    def non_background_indices(self) -> List[int]:
+    def non_background_indices(self) -> list[int]:
         """Indices of sites which are not background"""
@@ -289,3 +308,3 @@ bkg = set(self.background_indices)
     @property
-    def spin_dict(self) -> Dict[str, float]:
+    def spin_dict(self) -> dict[str, float]:
         return self.basis_func_type.spin_dict
@@ -302,3 +321,3 @@
     @property
-    def multiplicity_factor(self) -> Dict[str, float]:
+    def multiplicity_factor(self) -> dict[str, float]:
         """Return the multiplicity factor of each cluster."""
@@ -309,3 +328,3 @@ num_sites_in_group = [len(x) for x in self.index_by_sublattice]
     @property
-    def all_cf_names(self) -> List[str]:
+    def all_cf_names(self) -> list[str]:
         num_bf = len(self.basis_functions)
@@ -320,3 +339,3 @@ return self.cluster_list.get_all_cf_names(num_bf)
     @property
-    def index_by_basis(self) -> List[List[int]]:
+    def index_by_basis(self) -> list[list[int]]:
         first_symb_in_basis = [x[0] for x in self.basis_elements]
@@ -326,3 +345,3 @@ return self.atoms_mng.index_by_symbol(first_symb_in_basis)
     @property
-    def index_by_sublattice(self) -> List[List[int]]:
+    def index_by_sublattice(self) -> list[list[int]]:
         return self.atoms_mng.index_by_tag()
@@ -339,3 +358,3 @@
     @property
-    def size(self) -> Union[np.ndarray, None]:
+    def size(self) -> np.ndarray | None:
         return self.template_atoms.size
@@ -349,3 +368,3 @@
     @property
-    def supercell_factor(self) -> Union[int, None]:
+    def supercell_factor(self) -> int | None:
         return self.template_atoms.supercell_factor
@@ -370,3 +389,3 @@
 
-    def get_active_sublattices(self) -> List[bool]:
+    def get_active_sublattices(self) -> list[bool]:
         """List of booleans indicating if a (grouped) sublattice is active"""
@@ -378,3 +397,3 @@ unique_no_bkg = self.unique_element_without_background()
     @property
-    def ignored_species_and_conc(self) -> Dict[str, float]:
+    def ignored_species_and_conc(self) -> dict[str, float]:
         """
@@ -396,6 +415,4 @@ Return the ignored species and their concentrations normalised to the total number
                     raise ValueError(
-                        (
-                            "Ignored sublattice contains multiple elements -"
-                            "this does not make any sense"
-                        )
+                        "Ignored sublattice contains multiple elements -"
+                        "this does not make any sense"
                     )
@@ -453,3 +470,3 @@ if elem not in ignored:
 
-    def get_bg_syms(self) -> Set[str]:
+    def get_bg_syms(self) -> set[str]:
         """
@@ -612,3 +629,3 @@ Return the symbols in the basis where there is only one element
 
-    def get_all_figures_as_atoms(self) -> List[Atoms]:
+    def get_all_figures_as_atoms(self) -> list[Atoms]:
         """Get the list of all possible figures, in their
@@ -661,3 +678,3 @@ ASE Atoms representation."""
 
-    def todict(self) -> Dict:
+    def todict(self) -> dict:
         """Return a dictionary representation of the settings class.
@@ -680,3 +697,3 @@
     @classmethod
-    def from_dict(cls, dct: Dict[str, Any]) -> ClusterExpansionSettings:
+    def from_dict(cls, dct: dict[str, Any]) -> ClusterExpansionSettings:
         """Load a new ClusterExpansionSettings class from a dictionary representation.
@@ -711,2 +728,4 @@
 
+        # Allow instantiation from from_dict (legitimate use case)
+        kwargs["_allow_direct_instantiation"] = True
         settings = cls(*args, **kwargs)
@@ -786,3 +805,3 @@
 
-def _get_concentration(concentration: Union[Concentration, dict]) -> Concentration:
+def _get_concentration(concentration: Concentration | dict) -> Concentration:
     """Helper function to format the concentration"""
@@ -789,0 +808,0 @@ if isinstance(concentration, Concentration):

clease/settings/template_atoms.py

 """Class containing a manager for creating template atoms."""
+
+from collections.abc import Iterator
 from contextlib import contextmanager
 from itertools import product
-from typing import Iterator, List, Optional, Union
 
@@ -58,3 +59,3 @@ import ase
     @property
-    def size(self) -> Union[None, np.ndarray]:
+    def size(self) -> None | np.ndarray:
         return self._size
@@ -73,7 +74,7 @@
     @property
-    def supercell_factor(self) -> Union[int, None]:
+    def supercell_factor(self) -> int | None:
         return self._supercell_factor
 
     @supercell_factor.setter
-    def supercell_factor(self, value: Optional[int]) -> None:
+    def supercell_factor(self, value: int | None) -> None:
         if value is not None:
@@ -244,7 +245,7 @@ # Unset the size, since we provided a value for supercell_factor
 
-    def get_all_templates(self) -> List[ase.Atoms]:
+    def get_all_templates(self) -> list[ase.Atoms]:
         """Return a list with all templates."""
         return list(self.iterate_all_templates())
 
-    def iterate_all_templates(self, max_per_size: Optional[int] = None) -> Iterator[ase.Atoms]:
+    def iterate_all_templates(self, max_per_size: int | None = None) -> Iterator[ase.Atoms]:
         """Get all possible templates in an iterator.
@@ -417,3 +418,3 @@
 
-def _to_3x3_matrix(size: Union[List[int], List[List[int]], np.ndarray]) -> np.ndarray:
+def _to_3x3_matrix(size: list[int] | list[list[int]] | np.ndarray) -> np.ndarray:
     """Convert a list of ints (1D) or list of list of ints (2D) into a
@@ -420,0 +421,0 @@ 3x3 transformation matrix, if possible."""

clease/settings/template_filters.py

@@ -160,3 +160,3 @@ from abc import ABC, abstractmethod
         third_vec = set([0, 1, 2]) - set(span)
-        third_vec = list(third_vec)[0]
+        third_vec = next(iter(third_vec))
         d = normal.dot(cell[third_vec, :])
@@ -163,0 +163,0 @@ return abs(d)

clease/sparsifier.py

@@ -1,3 +0,1 @@
-from typing import List, Tuple
-
 from matplotlib import pyplot as plt
@@ -28,3 +26,3 @@ import numpy as np
         self, fitter: LinearRegression, X: np.ndarray, y: np.ndarray
-    ) -> Tuple[List[int], np.ndarray]:
+    ) -> tuple[list[int], np.ndarray]:
         """
@@ -31,0 +29,0 @@ Remove the feature corresponding to the smallest coefficient

clease/structgen/new_struct.py

 """Module for generating new structures for training."""
+
 from copy import deepcopy
@@ -7,3 +8,3 @@ from functools import reduce
 import os
-from typing import Any, Dict, List, Optional, Tuple, Union
+from typing import Any
 
@@ -65,3 +66,3 @@ from ase import Atoms
         settings: ClusterExpansionSettings,
-        generation_number: Optional[int] = None,
+        generation_number: int | None = None,
         struct_per_gen: int = 5,
@@ -99,5 +100,5 @@ check_db: bool = True,
         self,
-        atoms: Optional[Atoms] = None,
-        init_temp: Optional[float] = None,
-        final_temp: Optional[float] = None,
+        atoms: Atoms | None = None,
+        init_temp: float | None = None,
+        final_temp: float | None = None,
         num_temp: int = 5,
@@ -204,3 +205,3 @@ num_steps_per_temp: int = 1000,
         self,
-        eci: Dict[str, float],
+        eci: dict[str, float],
         num_templates: int = 20,
@@ -245,3 +246,3 @@ num_prim_cells: int = 2,
         self,
-        eci: Dict[str, float],
+        eci: dict[str, float],
         num_templates: int = 20,
@@ -315,4 +316,4 @@ num_prim_cells: int = 2,
         self,
-        atoms: Union[Atoms, List[Atoms]],
-        eci: Dict[str, float],
+        atoms: Atoms | list[Atoms],
+        eci: dict[str, float],
         init_temp: float = 2000.0,
@@ -411,3 +412,3 @@ final_temp: float = 1.0,
 
-    def generate_random_structures(self, atoms: Optional[Atoms] = None) -> None:
+    def generate_random_structures(self, atoms: Atoms | None = None) -> None:
         """
@@ -452,3 +453,3 @@ Generate random structures until the number of structures with
 
-    def generate_one_random_structure(self, atoms: Optional[Atoms] = None) -> bool:
+    def generate_one_random_structure(self, atoms: Atoms | None = None) -> bool:
         """
@@ -493,4 +494,4 @@ Generate and insert a random structure to database if a unique
     def _set_initial_structures(
-        self, atoms: Union[Atoms, List[Atoms]], random_composition: bool = False
-    ) -> List[Atoms]:
+        self, atoms: Atoms | list[Atoms], random_composition: bool = False
+    ) -> list[Atoms]:
         structs = []
@@ -549,3 +550,3 @@ if isinstance(atoms, Atoms):
 
-    def generate_initial_pool(self, atoms: Optional[Atoms] = None) -> None:
+    def generate_initial_pool(self, atoms: Atoms | None = None) -> None:
         """
@@ -576,6 +577,4 @@ Generate initial pool of structures.
         logger.info(
-            (
-                "Generating one structure per concentration where the number "
-                "of an element is at max/min"
-            )
+            "Generating one structure per concentration where the number "
+            "of an element is at max/min"
         )
@@ -605,3 +604,3 @@ indx_in_each_basis = []
     def generate_metropolis_trajectory(
-        self, atoms: Optional[Atoms] = None, random_comp: bool = True
+        self, atoms: Atoms | None = None, random_comp: bool = True
     ) -> None:
@@ -664,3 +663,3 @@ """
                 raise NotValidTemplateException(
-                    ("Did not find an atoms with a satisfactory concentration.")
+                    "Did not find an atoms with a satisfactory concentration."
                 )
@@ -671,3 +670,3 @@
     def insert_structures(
-        self, traj_init: str, traj_final: Optional[str] = None, cb=lambda num, tot: None
+        self, traj_init: str, traj_final: str | None = None, cb=lambda num, tot: None
     ) -> None:
@@ -719,9 +718,9 @@ """
         self,
-        init_struct: Union[Atoms, str],
-        final_struct: Optional[Union[Atoms, str]] = None,
-        name: Optional[str] = None,
-        cf: Optional[Dict[str, float]] = None,
-        meta: Optional[Dict[str, Any]] = None,
+        init_struct: Atoms | str,
+        final_struct: Atoms | str | None = None,
+        name: str | None = None,
+        cf: dict[str, float] | None = None,
+        meta: dict[str, Any] | None = None,
         warn_on_skip: bool = True,
-    ) -> Optional[Union[int, Tuple[int, int]]]:
+    ) -> int | tuple[int, int] | None:
         """Insert a structure to the database.
@@ -807,3 +806,3 @@
 
-    def _exists_in_db(self, atoms: Atoms, formula_unit: Optional[str] = None) -> bool:
+    def _exists_in_db(self, atoms: Atoms, formula_unit: str | None = None) -> bool:
         """
@@ -846,3 +845,3 @@ Check to see if the passed atoms already exists in DB.
 
-    def _get_kvp(self, formula_unit: Optional[str] = None) -> Dict:
+    def _get_kvp(self, formula_unit: str | None = None) -> dict:
         """
@@ -849,0 +848,0 @@ Create a dictionary of key-value pairs and return it.

clease/structgen/structure_generator.py

 """Module for generating new structures."""
+
 from abc import ABC
@@ -119,3 +120,3 @@ from copy import deepcopy
                     msg = (
-                        f"T Step: {Ti+1:>{len_max_temp}} of {N}, "
+                        f"T Step: {Ti + 1:>{len_max_temp}} of {N}, "
                         # Print same len as max required
@@ -126,3 +127,3 @@ f"Temp: {temp:>10.3f}, "
                         f"{self.num_steps_per_temp}, "
-                        f"Accept. rate: {acc_rate*100:.2f} %"
+                        f"Accept. rate: {acc_rate * 100:.2f} %"
                     )
@@ -372,3 +373,3 @@ logger.info(msg)
             else:
-                raise IOError(f"'{fname}' not found.")
+                raise OSError(f"'{fname}' not found.")
         else:
@@ -556,4 +557,4 @@ self.o_mv = mean_variance_approx(self.o_cfm)
     # if inverting matrix leads to a singular matrix, reduce the matrix
-    except np.linalg.linalg.LinAlgError:
+    except np.linalg.LinAlgError:
         prec = pinv(cfm.T.dot(cfm))
     return prec

clease/structure_mapper.py

@@ -1,3 +0,1 @@
-from typing import Tuple
-
 from ase.atoms import Atoms, Cell
@@ -78,3 +76,3 @@ from ase.geometry import find_mic
 
-    def snap_to_lattice(self, atoms: Atoms, template: Atoms) -> Tuple[Atoms, TransformInfo]:
+    def snap_to_lattice(self, atoms: Atoms, template: Atoms) -> tuple[Atoms, TransformInfo]:
         """
@@ -81,0 +79,0 @@ Snaps atoms to an ideal lattice. If the number of atoms in the

clease/tools.py

 """A collection of miscellaneous functions used for Cluster Expansion."""
-from collections.abc import Iterable
+
+from collections.abc import Iterable, Sequence
+from collections.abc import Iterable as tIterable
 from itertools import chain, combinations, filterfalse, permutations, product
@@ -7,4 +9,3 @@ import logging
 import re
-from typing import Dict, List, Optional, Sequence, Set, Tuple, Union
-from typing import Iterable as tIterable
+from typing import Protocol
 
@@ -20,3 +21,2 @@ import ase
 from scipy.spatial import cKDTree as KDTree
-from typing_extensions import Protocol
 
@@ -56,3 +56,3 @@ ASE_VERSION = parse(ase.__version__)
 
-def index_by_position(atoms: ase.Atoms) -> List[int]:
+def index_by_position(atoms: ase.Atoms) -> list[int]:
     """Set atomic indices by its position."""
@@ -179,4 +179,4 @@ # add zero to avoid negative zeros
     db_name=None,
-    custom_kvp_init: Optional[dict] = None,
-    custom_kvp_final: Optional[dict] = None,
+    custom_kvp_init: dict | None = None,
+    custom_kvp_final: dict | None = None,
 ):
@@ -289,3 +289,3 @@ """Update the database.
     X: np.ndarray, y: np.ndarray, nsplits: int = 10, groups: Sequence[int] = ()
-) -> List[Dict[str, np.ndarray]]:
+) -> list[dict[str, np.ndarray]]:
     """Split the dataset such that it can be used for k-fold
@@ -350,3 +350,3 @@ cross validation.
 def random_validation_set(
-    num: int = 10, select_cond: Optional[list] = None, db_name: Optional[str] = None
+    num: int = 10, select_cond: list | None = None, db_name: str | None = None
 ):
@@ -372,3 +372,3 @@ """
 
-def exclude_ids(ids: List[int]) -> List[tuple]:
+def exclude_ids(ids: list[int]) -> list[tuple]:
     """
@@ -455,3 +455,3 @@ Construct a select condition based on the ids passed.
         # Opposite facet
-        remaining = list(set([0, 1, 2]) - set(plane))[0]
+        remaining = next(iter(set([0, 1, 2]) - set(plane)))
         vec = cell[remaining, :]
@@ -873,3 +873,3 @@ dist = np.abs(n.dot(x - vec))
 
-def get_extension(fname: Union[str, Path]) -> str:
+def get_extension(fname: str | Path) -> str:
     """
@@ -886,3 +886,3 @@ Return the file extension of a filename
 
-def add_file_extension(fname: Union[str, Path], ext: str) -> str:
+def add_file_extension(fname: str | Path, ext: str) -> str:
     """
@@ -934,3 +934,3 @@ Adds the wanted file extension to a filename. If a file extension
 
-def sort_cf_names(cf_names: tIterable[str]) -> List[str]:
+def sort_cf_names(cf_names: tIterable[str]) -> list[str]:
     """
@@ -952,3 +952,3 @@ Return a sorted list of correlation function names. The names are
 
-def get_ids(select_cond: List[tuple], db_name: str) -> List[int]:
+def get_ids(select_cond: list[tuple], db_name: str) -> list[int]:
     """
@@ -978,3 +978,3 @@ Return ids in the datase that correspond to the passed selection.
 class SQLCursor(Protocol):
-    def execute(self, sql: str, placeholder: Tuple[str]) -> None:
+    def execute(self, sql: str, placeholder: tuple[str]) -> None:
         pass
@@ -986,3 +986,3 @@
 
-def get_attribute(ids: List[int], cur: SQLCursor, key: str, table: str) -> list:
+def get_attribute(ids: list[int], cur: SQLCursor, key: str, table: str) -> list:
     """
@@ -1016,3 +1016,3 @@ Retrieve the value of the given key for the rows with the given ID of
 
-def common_cf_names(ids: Set[int], cur: SQLCursor, table: str) -> Set[str]:
+def common_cf_names(ids: set[int], cur: SQLCursor, table: str) -> set[str]:
     """
@@ -1048,4 +1048,4 @@ Extracts all correlation function names that are present for all
     d: float,
-    A_eq: Optional[np.ndarray] = None,
-    b_eq: Optional[np.ndarray] = None,
+    A_eq: np.ndarray | None = None,
+    b_eq: np.ndarray | None = None,
 ) -> bool:
@@ -1086,5 +1086,5 @@ """
     b_lb: np.ndarray,
-    A_eq: Optional[np.ndarray] = None,
-    b_eq: Optional[np.ndarray] = None,
-) -> Tuple[np.ndarray, np.ndarray]:
+    A_eq: np.ndarray | None = None,
+    b_eq: np.ndarray | None = None,
+) -> tuple[np.ndarray, np.ndarray]:
     """
@@ -1234,3 +1234,3 @@ Remove all redundant constraints from A_lb and b_lb.
     if n_target != len(superatoms):
-        msg = "Number of atoms in supercell: {}, expected: {}".format(n_target, len(superatoms))
+        msg = f"Number of atoms in supercell: {n_target}, expected: {len(superatoms)}"
         raise ase_sc.SupercellError(msg)
@@ -1299,3 +1299,3 @@
     center=(0.5, 0.5, 0.5),
-    cell_T_inv: Optional[np.ndarray] = None,
+    cell_T_inv: np.ndarray | None = None,
     eps: float = 1e-7,
@@ -1302,0 +1302,0 @@ ) -> np.ndarray:

cxx/cython/ce_updater.pxd

@@ -31,3 +31,3 @@ # distutils: language = c++
 
-      object get_singlets()
+      vector[double] get_singlets()
 
@@ -34,0 +34,0 @@ const vector[string]& get_symbols() const

cxx/cython/clease_cxx.pyx

 # distutils: language = c++
 # cython: c_string_type=str, c_string_encoding=ascii
-import numpy as np
-cimport numpy as np  # Initialize the Numpy API
 from libcpp cimport bool
 
-np.import_array()
 include "pyce_updater.pyx"
@@ -14,5 +11,1 @@ include "py_cluster.pyx"
     cpdef bool has_parallel()
-
-# Files that use the Numpy Array API are included here
-cdef extern from "with_numpy/ce_updater.cpp":
-    pass

cxx/include/additional_tools.tpp

 #include <cmath>
 
-template<class key,class value>
-std::ostream& operator <<( std::ostream &out, const std::map<key,value> &map )
-{
-  for ( auto iter=map.begin(); iter != map.end(); ++iter )
-  {
-    out << iter->first << ":" << iter->second << "\n";
-  }
-  return out;
+template <class key, class value>
+std::ostream &operator<<(std::ostream &out, const std::map<key, value> &map) {
+    for (auto iter = map.begin(); iter != map.end(); ++iter) {
+        out << iter->first << ":" << iter->second << "\n";
+    }
+    return out;
 };
 
-template<class T>
-std::ostream& operator << (std::ostream &out, const std::vector<T> &vec )
-{
-  out << "[";
-  for ( unsigned int i=0;i<vec.size(); i++ )
-  {
-    out << vec[i] << " ";
-  }
-  out << "]";
-  return out;
+template <class T>
+std::ostream &operator<<(std::ostream &out, const std::vector<T> &vec) {
+    out << "[";
+    for (unsigned int i = 0; i < vec.size(); i++) {
+        out << vec[i] << " ";
+    }
+    out << "]";
+    return out;
 };
 
-template<class T>
-std::vector<T>& cyclic_permute( std::vector<T> &vec )
-{
-  for ( unsigned int i=0;i<vec.size()-1;i++ )
-  {
-    T prev = vec[0];
-    vec[0] = vec[i+1];
-    vec[i+1] = prev;
-  }
-  return vec;
+template <class T>
+std::vector<T> &cyclic_permute(std::vector<T> &vec) {
+    for (unsigned int i = 0; i < vec.size() - 1; i++) {
+        T prev = vec[0];
+        vec[0] = vec[i + 1];
+        vec[i + 1] = prev;
+    }
+    return vec;
 };
 
-template<class T>
-void keys( std::map<std::string,T> &dict, std::vector<std::string> &names )
-{
-  names.clear();
-  for ( auto iter=dict.begin(); iter != dict.end(); ++iter )
-  {
-    names.push_back(iter->first);
-  }
+template <class T>
+void keys(std::map<std::string, T> &dict, std::vector<std::string> &names) {
+    names.clear();
+    for (auto iter = dict.begin(); iter != dict.end(); ++iter) {
+        names.push_back(iter->first);
+    }
 };
 
-template<class T>
-void set2vector( const std::set<T> &s, std::vector<T> &vec )
-{
-  vec.clear();
-  for ( auto iter = s.begin(); iter != s.end(); ++iter )
-  {
-    vec.push_back(*iter);
-  }
+template <class T>
+void set2vector(const std::set<T> &s, std::vector<T> &vec) {
+    vec.clear();
+    for (auto iter = s.begin(); iter != s.end(); ++iter) {
+        vec.push_back(*iter);
+    }
 };
 
-template<class T, unsigned int N>
-std::ostream& operator <<(std::ostream& out, std::array<T, N> &array)
-{
-  for (unsigned int i=0;i<N;i++ )
-  {
-    out << array[i];
-  }
-  out << "\n";
-  return out;
+template <class T, unsigned int N>
+std::ostream &operator<<(std::ostream &out, std::array<T, N> &array) {
+    for (unsigned int i = 0; i < N; i++) {
+        out << array[i];
+    }
+    out << "\n";
+    return out;
 };
 
-template<class T>
-std::ostream& operator <<(std::ostream& out, std::set<T> &set){
-  for (auto iter=set.begin(); iter != set.end(); ++iter){
-    out << *iter << " ";
-  }
-  out << "\n";
-  return out;
+template <class T>
+std::ostream &operator<<(std::ostream &out, std::set<T> &set) {
+    for (auto iter = set.begin(); iter != set.end(); ++iter) {
+        out << *iter << " ";
+    }
+    out << "\n";
+    return out;
 }
 
-template<class T>
-bool is_in_vector(const T& value, const std::vector<T> &vector){
-  for (const T& elem : vector){
-    if (value == elem) return true;
-  }
-  return false;
+template <class T>
+bool is_in_vector(const T &value, const std::vector<T> &vector) {
+    for (const T &elem : vector) {
+        if (value == elem) return true;
+    }
+    return false;
 }
 
-template<class T>
-void insert_in_set(const std::vector<T> &vec, std::set<T> &unique){
-  for (const T& elem : vec){
-    unique.insert(elem);
-  }
+template <class T>
+void insert_in_set(const std::vector<T> &vec, std::set<T> &unique) {
+    for (const T &elem : vec) {
+        unique.insert(elem);
+    }
 }

cxx/include/ce_updater.hpp

@@ -81,4 +81,3 @@ #ifndef CE_UPDATER_H
     /** Returns the value of the singlets */
-    void get_singlets(PyObject *npy_array) const;
-    PyObject *get_singlets() const;
+    std::vector<double> get_singlets() const;
 
@@ -85,0 +84,0 @@ /** Updates the CF */

cxx/include/matrix.tpp

 
-template<class T>
-Matrix<T>::Matrix( unsigned int n_rows, unsigned int n_cols ):n_rows(n_rows), n_cols(n_cols)
-{
-  data = new T[n_rows*n_cols];
+template <class T>
+Matrix<T>::Matrix(unsigned int n_rows, unsigned int n_cols) : n_rows(n_rows), n_cols(n_cols) {
+    data = new T[n_rows * n_cols];
 };
 
-template<class T>
-Matrix<T>::~Matrix()
-{
-  delete [] data;
+template <class T>
+Matrix<T>::~Matrix() {
+    delete[] data;
 };
 
 template <class T>
-T& Matrix<T>::operator()( unsigned int i, unsigned int j )
-{
-  return data[j*n_rows+i];
+T &Matrix<T>::operator()(unsigned int i, unsigned int j) {
+    return data[j * n_rows + i];
 };
 
-template<class T>
-const T& Matrix<T>::operator()( unsigned int i, unsigned int j ) const
-{
-  return data[j*n_rows+i];
+template <class T>
+const T &Matrix<T>::operator()(unsigned int i, unsigned int j) const {
+    return data[j * n_rows + i];
 };
 
-template<class T>
-void Matrix<T>::set_size( unsigned int nr, unsigned int nc )
-{
-  n_rows = nr;
-  n_cols = nc;
-  data = new T[n_rows*n_cols];
+template <class T>
+void Matrix<T>::set_size(unsigned int nr, unsigned int nc) {
+    n_rows = nr;
+    n_cols = nc;
+    data = new T[n_rows * n_cols];
 };
 
-template<class T>
-void swap( Matrix<T> &first, const Matrix<T> &second )
-{
-  first.n_rows = second.n_rows;
-  first.n_cols = second.n_cols;
-  first.data = new T[second.n_rows*second.n_cols];
-  for ( unsigned int i=0;i<first.n_rows*first.n_cols;i++ )
-  {
-    first.data[i] = second.data[i];
-  }
+template <class T>
+void swap(Matrix<T> &first, const Matrix<T> &second) {
+    first.n_rows = second.n_rows;
+    first.n_cols = second.n_cols;
+    first.data = new T[second.n_rows * second.n_cols];
+    for (unsigned int i = 0; i < first.n_rows * first.n_cols; i++) {
+        first.data[i] = second.data[i];
+    }
 }
 
-template<class T>
-Matrix<T>::Matrix( const Matrix<T> &other )
-{
-  swap(*this,other);
+template <class T>
+Matrix<T>::Matrix(const Matrix<T> &other) {
+    swap(*this, other);
 }
 
-template<class T>
-Matrix<T>& Matrix<T>::operator=(const Matrix<T> &other )
-{
-  swap(*this,other);
-  return *this;
+template <class T>
+Matrix<T> &Matrix<T>::operator=(const Matrix<T> &other) {
+    swap(*this, other);
+    return *this;
 }

doc/source/index.rst

@@ -10,2 +10,11 @@ .. CLEASE documentation master file, created by
 
+.. note::
+   **This documentation has moved!**
+
+   The updated version of the CLEASE documentation is now available at:
+
+   https://www.clease.org
+
+   Please visit the new site for the latest documentation, tutorials, and resources.
+
 Cluster expansion (CE) is a widely used method for studying thermondynamic
@@ -12,0 +21,0 @@ properties of disordered materials. CLEASE is a cluster expansion code which strives

doc/source/releasenotes.rst

@@ -7,2 +7,11 @@ .. _releasenotes:
 
+1.2.0
+=====
+* **Python 3.10 is now the minimum supported version**
+* Removed requirement `numpy < 2`
+* No longer requires Numpy to build the C++ extension
+* `ClusterExpansionSettings` can no longer be instantiated directly,
+  unless overridden by the `_allow_direct_instantiation` keyword.
+
+
 1.1.0
@@ -9,0 +18,0 @@ ======

PKG-INFO

@@ -1,4 +0,4 @@
-Metadata-Version: 2.1
+Metadata-Version: 2.4
 Name: clease
-Version: 1.1.0
+Version: 1.2.0
 Summary: CLuster Expansion in Atomistic Simulation Environment
@@ -14,13 +14,10 @@ Home-page: https://gitlab.com/computationalmaterials/clease/
 Classifier: License :: OSI Approved :: Mozilla Public License 2.0 (MPL 2.0)
-Classifier: Programming Language :: Python :: 3.7
-Classifier: Programming Language :: Python :: 3.8
-Classifier: Programming Language :: Python :: 3.9
-Classifier: Programming Language :: Python :: 3.10
+Classifier: Programming Language :: Python :: 3
 Classifier: Topic :: Scientific/Engineering :: Physics
 Classifier: Topic :: Scientific/Engineering :: Chemistry
-Requires-Python: >=3.7
+Requires-Python: >=3.10
 Description-Content-Type: text/markdown
 License-File: LICENSE.md
 Requires-Dist: ase>=3.22
-Requires-Dist: numpy<2
+Requires-Dist: numpy
 Requires-Dist: cython
@@ -52,5 +49,4 @@ Requires-Dist: matplotlib
 Requires-Dist: twine; extra == "dev"
-Requires-Dist: black>=22.1.0; extra == "dev"
-Requires-Dist: clang-format>=14.0.3; extra == "dev"
-Requires-Dist: ruff; extra == "dev"
+Requires-Dist: clang-format==21.1.5; extra == "dev"
+Requires-Dist: ruff==0.14.5; extra == "dev"
 Requires-Dist: pyclean>=2.0.0; extra == "dev"
@@ -62,21 +58,22 @@ Requires-Dist: pytest-cov; extra == "dev"
 Provides-Extra: all
-Requires-Dist: black>=22.1.0; extra == "all"
+Requires-Dist: pyclean>=2.0.0; extra == "all"
+Requires-Dist: pre-commit; extra == "all"
+Requires-Dist: ipython; extra == "all"
+Requires-Dist: sphinx; extra == "all"
+Requires-Dist: pytest; extra == "all"
+Requires-Dist: twine; extra == "all"
+Requires-Dist: cython; extra == "all"
+Requires-Dist: ruff==0.14.5; extra == "all"
+Requires-Dist: pytest-benchmark[histogram]>=3.4.1; extra == "all"
 Requires-Dist: mock; extra == "all"
-Requires-Dist: twine; extra == "all"
-Requires-Dist: ruff; extra == "all"
+Requires-Dist: clang-format==21.1.5; extra == "all"
 Requires-Dist: sphinx_rtd_theme; extra == "all"
-Requires-Dist: pre-commit; extra == "all"
-Requires-Dist: pyclean>=2.0.0; extra == "all"
 Requires-Dist: pip; extra == "all"
-Requires-Dist: pytest-cov; extra == "all"
 Requires-Dist: clease-gui; extra == "all"
+Requires-Dist: tox>=4; extra == "all"
+Requires-Dist: pytest-mock; extra == "all"
 Requires-Dist: build; extra == "all"
-Requires-Dist: sphinx; extra == "all"
-Requires-Dist: pytest-mock; extra == "all"
-Requires-Dist: pytest; extra == "all"
-Requires-Dist: cython; extra == "all"
-Requires-Dist: ipython; extra == "all"
-Requires-Dist: clang-format>=14.0.3; extra == "all"
-Requires-Dist: tox>=4; extra == "all"
-Requires-Dist: pytest-benchmark[histogram]>=3.4.1; extra == "all"
+Requires-Dist: pytest-cov; extra == "all"
+Dynamic: license-file
+Dynamic: provides-extra
 
@@ -94,2 +91,10 @@ # CLEASE
 
+# Documentation
+
+   The updated version of the CLEASE documentation is now available at:
+
+   https://www.clease.org
+
+   Please visit the new site for the latest documentation, tutorials, and resources.
+
 # Installation
@@ -96,0 +101,0 @@

pyproject.toml

 [build-system]
 # setuptools 61.0 broke, see !488
-requires = ["setuptools >=58, !=61.0", "wheel", "Cython", "numpy<2"]
+requires = ["setuptools >=62", "wheel", "Cython"]
 build-backend = "setuptools.build_meta"
 
-[tool.black]
-line-length = 100
-target-version = ["py37"]
-
 [tool.ruff]
-select = ["E", "F", "NPY", "RUF", "I"]
-exclude = ["tests", "__init__.py"]
-ignore = ["NPY002", "RUF012"]
-# Same as Black.
 line-length = 100
-target-version = "py37"
+target-version = "py310"
+extend-exclude = ["tests", "__init__.py", "clease/gui/*", "__pycache__"]
 
-[tool.ruff.isort]
+[tool.ruff.lint]
+select = ["E", "F", "NPY", "RUF", "I", "UP"]
+ignore = [
+    "E402",  # Allow module level import not at top of file
+    "N802",  # Allow class names to start with an underscore
+    "N803",  # Allow function names to start with an underscore
+    "N806",   # Allow function names to start with an underscore
+    "NPY002", # Checks for the use of legacy np.random function calls.
+    "NPY201", # Checks for deprecated function calls removed in NumPy 2.0.
+    "RUF012", # Mutable class attributes should be annotated with typing.ClassVar
+    "RUF015", # Prefer next({iterable}) over single element slice
+    "RUF022", #  `__all__` is not sorted
+]
+
+[tool.ruff.lint.per-file-ignores]
+"__init__.py" = ["F401", "F403", "F405"]
+
+[tool.ruff.lint.isort]
 force-sort-within-sections = true
+
+[tool.cibuildwheel]
+# Documentation: https://cibuildwheel.pypa.io/en/stable/options/#options
+# cibuildwheel automatically reads the supported Python versions from setup.cfg
+
+environment = "CLEASE_OMP='-fopenmp'"
+skip = ["cp3??t-*"]  # Disable free-threading builds
+build = "cp*"  # Must be CPython
+archs = ["auto"]

README.md

@@ -12,2 +12,10 @@ # CLEASE
 
+# Documentation
+
+   The updated version of the CLEASE documentation is now available at:
+
+   https://www.clease.org
+
+   Please visit the new site for the latest documentation, tutorials, and resources.
+
 # Installation
@@ -14,0 +22,0 @@

setup.cfg

@@ -22,6 +22,3 @@ [metadata]
 	License :: OSI Approved :: Mozilla Public License 2.0 (MPL 2.0)
-	Programming Language :: Python :: 3.7
-	Programming Language :: Python :: 3.8
-	Programming Language :: Python :: 3.9
-	Programming Language :: Python :: 3.10
+	Programming Language :: Python :: 3
 	Topic :: Scientific/Engineering :: Physics
@@ -32,7 +29,7 @@ Topic :: Scientific/Engineering :: Chemistry
 packages = find:
-python_requires = >=3.7
+python_requires = >=3.10
 include_package_data = true
 install_requires = 
 	ase>=3.22
-	numpy<2
+	numpy
 	cython
@@ -46,3 +43,3 @@ matplotlib
 	attrs>=21.4.0
-	scipy>=1.7.0  # Last minor version which allows python 3.7
+	scipy>=1.7.0
 	packaging  # Replaces deprecated distutils
@@ -49,0 +46,0 @@ threadpoolctl  # For controlling parallelization of NumPy

setup.py

@@ -32,8 +32,2 @@ import os
 
-def get_npy_include_folder():
-    import numpy as np
-
-    return np.get_include()
-
-
 def build_ext(ext_module):
@@ -56,4 +50,3 @@ from Cython.Build import cythonize
 
-    # C++ sources, excluding the ones which rely on NumPy
-    # those are included directly from the Cython source
+    # C++ sources
     sources = [str(f) for f in cpp_sources]
@@ -95,3 +88,2 @@
         cxx_inc_folder,
-        get_npy_include_folder(),
         cxx_src_folder,
@@ -129,5 +121,4 @@ py_include,
         "twine",
-        "black>=22.1.0",  # Style formatting
-        "clang-format>=14.0.3",  # c++ style formatting
-        "ruff",
+        "clang-format==21.1.5",  # c++ style formatting
+        "ruff==0.14.5",  # Python style and formatting
         "pyclean>=2.0.0",  # For removing __pycache__ and .pyc files
@@ -134,0 +125,0 @@ "pytest-cov",

cxx/src/with_numpy/ce_updater.cpp

-#include "ce_updater.hpp"
-
-#include <algorithm>
-#include <cassert>
-#include <iostream>
-#include <iterator>
-#include <sstream>
-#include <stdexcept>
-
-#include "additional_tools.hpp"
-#include "atomic_numbers.hpp"
-#include "basis_function.hpp"
-#include "cluster_name.hpp"
-#include "config.hpp"
-
-using namespace std;
-
-CEUpdater::CEUpdater(){};
-CEUpdater::~CEUpdater(){};
-
-void CEUpdater::init(PyObject *py_atoms, PyObject *settings, PyObject *corrFunc, PyObject *pyeci,
-                     PyObject *cluster_list) {
-    atoms = py_atoms;
-    if (settings == nullptr) {
-        throw invalid_argument("Settings object is nullptr!");
-    }
-
-#ifdef PRINT_DEBUG
-    cout << "Getting symbols from settings object\n";
-#endif
-    PyObject *py_ignore_bck = get_attr(settings, "ignore_background_atoms");
-    ignore_background_indices = PyObject_IsTrue(py_ignore_bck);
-    Py_DECREF(py_ignore_bck);
-
-    unsigned int n_atoms = PyObject_Length(atoms);
-    if (n_atoms < 0) {
-        throw invalid_argument("Could not retrieve the length of the atoms object!");
-    }
-
-    // Read the atomic symbols
-    std::vector<std::string> symbols = get_symbols_from_atoms(py_atoms);
-
-    trans_symm_group.resize(n_atoms);
-    set<string> unique_symbols;
-
-    // Extract unique symbols from settings
-    PyObject *py_unique_symb = get_attr(settings, "unique_elements");
-    for (unsigned int i = 0; i < list_size(py_unique_symb); i++) {
-        unique_symbols.insert(py2string(PyList_GetItem(py_unique_symb, i)));
-    }
-    Py_DECREF(py_unique_symb);
-
-    insert_in_set(symbols, unique_symbols);
-    symbols_with_id = std::make_unique<Symbols>(symbols, unique_symbols);
-
-    // Build read the translational sites
-    PyObject *py_trans_symm_group = get_attr(settings, "index_by_sublattice");
-
-    if (py_trans_symm_group == nullptr) {
-        throw invalid_argument("index_by_sublattice is nullptr!");
-    }
-
-#ifdef PRINT_DEBUG
-    cout << "Reading background indices\n";
-#endif
-    // Read the backgound indices from settings
-    PyObject *bkg_indx = get_attr(settings, "background_indices");
-    read_background_indices(bkg_indx);
-    Py_DECREF(bkg_indx);
-
-    count_non_bkg_sites();
-    build_trans_symm_group(py_trans_symm_group);
-    Py_DECREF(py_trans_symm_group);
-
-#ifdef PRINT_DEBUG
-    cout << "Getting cluster names from atoms object\n";
-#endif
-
-    // Read cluster names
-    create_cname_with_dec(corrFunc);
-
-#ifdef PRINT_DEBUG
-    cout << "Cluster names with decoration number created...\n";
-#endif
-
-    PyObject *py_num_elements = get_attr(settings, "num_unique_elements");
-
-    if (py_num_elements == nullptr) {
-        throw invalid_argument("num_unique_elements is nullptr!");
-    }
-    int num_bfs = py2int(py_num_elements) - 1;
-    Py_DECREF(py_num_elements);
-
-    if (cluster_list == nullptr) {
-        throw invalid_argument("cluster_list is nullptr!");
-    }
-    // unsigned int num_trans_symm = list_size(cluster_info);
-    unsigned int num_clusters = PySequence_Size(cluster_list);
-
-#ifdef PRINT_DEBUG
-    cout << "Parsing cluster list...\n";
-#endif
-
-    PyObject *py_no_si = get_attr(cluster_list, "assume_no_self_interactions");
-    assume_no_self_interactions = PyObject_IsTrue(py_no_si);
-    Py_DECREF(py_no_si);
-
-#ifdef PRINT_DEBUG
-    std::cout << "Assuming no self-interaction?: " << assume_no_self_interactions << std::endl;
-#endif
-
-    for (unsigned int i = 0; i < num_clusters; i++) {
-        PyObject *py_cluster = PySequence_GetItem(cluster_list, i);
-
-        Cluster new_clst(py_cluster);
-        PyObject *py_cluster_name = get_attr(py_cluster, "name");
-        string cluster_name = py2string(py_cluster_name);
-        Py_DECREF(py_cluster_name);
-        Py_DECREF(py_cluster);
-
-        new_clst.construct_equivalent_deco(num_bfs);
-        clusters.append(new_clst);
-
-        int norm_fact = new_clst.get().size() * trans_symm_group_count[new_clst.symm_group];
-        if (normalisation_factor.find(cluster_name) == normalisation_factor.end()) {
-            normalisation_factor[cluster_name] = norm_fact;
-        } else {
-            normalisation_factor[cluster_name] += norm_fact;
-        }
-    }
-#ifdef PRINT_DEBUG
-    cout << "Finished reading cluster_info\n";
-#endif
-
-#ifdef PRINT_DEBUG
-    cout << "Reading basis functions from settings object\n";
-#endif
-
-    PyObject *bfs = get_attr(settings, "basis_functions");
-    if (bfs == NULL) {
-        status = Status_t::INIT_FAILED;
-        return;
-    }
-
-    // Reading basis functions from python object
-    PyObject *key;
-    PyObject *value;
-    unsigned int n_bfs = list_size(bfs);
-    bf_list basis_func_raw;
-    for (unsigned int i = 0; i < n_bfs; i++) {
-        Py_ssize_t pos = 0;
-        map<string, double> new_entry;
-        PyObject *bf_dict = PyList_GetItem(bfs, i);
-        while (PyDict_Next(bf_dict, &pos, &key, &value)) {
-            new_entry[py2string(key)] = PyFloat_AS_DOUBLE(value);
-        }
-        basis_func_raw.push_back(new_entry);
-    }
-
-    this->basis_functions = std::make_unique<BasisFunction>(basis_func_raw, *symbols_with_id);
-
-#ifdef PRINT_DEBUG
-    cout << "Reading translation matrix from settings\n";
-#endif
-    // Retrieve the TransMatrix object
-    PyObject *trans_mat_obj = get_attr(settings, "trans_matrix");
-    if (trans_mat_obj == NULL) {
-        status = Status_t::INIT_FAILED;
-        return;
-    }
-
-    // Get the internal trans_matrix object from within the TransMatrix object
-    PyObject *trans_mat_orig = get_attr(trans_mat_obj, "trans_matrix");
-
-    read_trans_matrix(trans_mat_orig);
-    Py_DECREF(trans_mat_obj);
-    Py_DECREF(trans_mat_orig);
-
-    // Read the ECIs, and parse the names.
-    this->set_eci(pyeci);
-#ifdef PRINT_DEBUG
-    cout << "Parsing correlation function\n";
-#endif
-
-    vector<string> flattened_cnames;
-    flattened_cf_names(flattened_cnames);
-
-    history = std::make_unique<CFHistoryTracker>(eci.get_names());
-    history->insert(corrFunc, nullptr);
-
-    // Store the singlets names
-    for (unsigned int i = 0; i < flattened_cnames.size(); i++) {
-        std::string name = flattened_cnames[i];
-        // Fetch the pre-parsed version of the name.
-        const ParsedName parsed = this->m_parsed_names[i];
-        if (parsed.size == 1) {
-            singlets.push_back(name);
-        }
-    }
-
-    status = Status_t::READY;
-    clear_history();
-#ifdef PRINT_DEBUG
-    cout << "CEUpdater initialized sucessfully!\n";
-#endif
-
-    // Verify that the ECIs given corresponds to a correlation function
-    if (!all_eci_corresponds_to_cf()) {
-        throw invalid_argument("All ECIs does not correspond to a correlation function!");
-    }
-}
-
-void CEUpdater::parse_eci_names() {
-    std::size_t num = eci.size();
-    this->m_parsed_names.clear();
-    this->m_parsed_names.reserve(num);
-
-    for (unsigned int i = 0; i < num; i++) {
-        std::string name = eci.name(i);
-        ClusterName c_name = ClusterName(name);
-        ParsedName parsed = c_name.get_parsed();
-        this->m_parsed_names.emplace_back(parsed);
-    }
-}
-
-double CEUpdater::get_energy() {
-    double energy = 0.0;
-    cf &corr_func = history->get_current();
-    energy = eci.dot(corr_func);
-    return energy * symbols_with_id->size();
-}
-
-double CEUpdater::spin_product_one_atom(int ref_indx, const Cluster &cluster,
-                                        const vector<int> &dec, int ref_id) const {
-    double sp = 0.0;
-
-    // Note: No duplication factor is included here, since this method is used for calculating the
-    // CF from scratch (which will account for self-interactions with no issue), and not updating
-    // (which must account for self-interactions via the duplication factor).
-    // i.e. spin_product_one_atom_delta must accounts for the self-interaction.
-
-    const vector<vector<int>> &indx_list = cluster.get();
-    unsigned int num_indx = indx_list.size();
-    unsigned int n_memb = indx_list[0].size();
-
-    // Cache the relevant row from the trans matrix.
-    int *trans_matrix_row = trans_matrix.get_row(ref_indx);
-
-    for (unsigned int i = 0; i < num_indx; i++) {
-        double sp_temp = 1.0;
-
-        // Use pointer arithmetics in the inner most loop
-        const int *indices = &indx_list[i][0];
-
-        for (unsigned int j = 0; j < n_memb; j++) {
-            int trans_index = trans_matrix.lookup_in_row(trans_matrix_row, indices[j]);
-            int id = (trans_index == ref_indx) ? ref_id : symbols_with_id->id(trans_index);
-            sp_temp *= basis_functions->get(dec[j], id);
-        }
-        sp += sp_temp;
-    }
-    return sp;
-}
-
-int CEUpdater::get_original_index(int ref_indx) const {
-    int *trans_matrix_row = trans_matrix.get_row(ref_indx);
-    int *allowed_lu = trans_matrix.get_allowed_lookup_values();
-    for (unsigned int j = 0; j < trans_matrix.get_num_non_zero(); j++) {
-        int col = allowed_lu[j];
-        int indx = trans_matrix.lookup_in_row(trans_matrix_row, col);
-        if (indx == ref_indx) {
-            return col;
-        }
-    }
-    std::stringstream err;
-    err << "Did not locate original index for ref index: " << ref_indx;
-    throw std::runtime_error(err.str());
-}
-
-double CEUpdater::spin_product_one_atom_delta_no_si(const SpinProductCache &sp_cache,
-                                                    const Cluster &cluster,
-                                                    const deco_t &deco) const {
-    /* Note: This function assumes no self-interaction within a cluster. */
-
-    // Figure out how many times we need to iterate
-    unsigned int num_indx = cluster.get_num_figures();  // Outer loop count
-    // Assume 1 ref site in a figure, so we iterate 1 less
-    unsigned int n_non_ref = cluster.get_size() - 1;  // Inner loop count
-
-    int *tm_row = sp_cache.trans_matrix_row;
-
-    /* There are n_non_ref sites per each ref site, so the non_ref_site_ptr
-    iterates faster than the ref_site_ptr. Figures are placed contiguously
-    in a 1D array.*/
-    const ClusterSite *non_ref_site_ptr = &cluster.get_non_ref_sites()[0];
-    const int *ref_site_ptr = &cluster.get_ref_cluster_sites()[0];
-
-    // Keep track of the change in spin-product
-    double sp_delta = 0.0;
-
-    // Iterate each figure in the cluster. 1 reference site is assumed per cluster
-    for (unsigned int i = 0; i < num_indx; i++, ++ref_site_ptr) {
-        /* Calculate the spin product for both new and the old (ref)
-        The constant term to the spin product from the sites which didn't change.*/
-        const int dec_ref = deco[*ref_site_ptr];
-        double new_bf = basis_functions->get(dec_ref, sp_cache.new_symb_id);
-        double old_bf = basis_functions->get(dec_ref, sp_cache.old_symb_id);
-        double sp_change = new_bf - old_bf;
-
-        /* Iterate the remaining non-reference sites, as we already took care of the reference
-        site (assuming no self-interaction)*/
-        for (unsigned int j = 0; j < n_non_ref; j++, ++non_ref_site_ptr) {
-            const ClusterSite &site = *non_ref_site_ptr;
-            const int dec_j = deco[site.cluster_index];
-
-            const int trans_index = trans_matrix.lookup_in_row(tm_row, site.lattice_index);
-            sp_change *= basis_functions->get(dec_j, symbols_with_id->id(trans_index));
-        }
-
-        sp_delta += sp_change;
-    }
-    return sp_delta;
-}
-
-double CEUpdater::spin_product_one_atom_delta(const SpinProductCache &sp_cache,
-                                              const Cluster &cluster, const deco_t &deco) const {
-    // Keep track of the change in spin-product
-    double sp_delta = 0.0;
-
-    // List of figures in the cluster
-    const vector<vector<int>> &indx_list = cluster.get();
-    // Account for the self-interaction, in case we updated 2 sites with 1 change
-    const std::vector<double> &dup_factors = cluster.get_duplication_factors();
-    unsigned int num_indx = indx_list.size();
-    unsigned int n_memb = indx_list[0].size();
-
-    int *tm_row = sp_cache.trans_matrix_row;
-
-    // Iterate each site in the cluster
-    for (unsigned int i = 0; i < num_indx; i++) {
-        // Use pointer arithmetics in the inner most loop
-        const int *indices = &indx_list[i][0];
-
-        // Calculate the spin product for both new and the old (ref)
-        double sp_temp_new = 1.0, sp_temp_ref = 1.0;
-        // The constant term to the spin product from the sites which didn't change.
-        double sp_const = dup_factors[i];
-
-        for (unsigned int j = 0; j < n_memb; j++) {
-            const int site_index = indices[j];
-            const int dec_j = deco[j];
-            if (site_index == sp_cache.original_index) {
-                // This site is the reference index.
-                // Look up the BF values directly for the new and old symbols
-                sp_temp_new *= basis_functions->get(dec_j, sp_cache.new_symb_id);
-                sp_temp_ref *= basis_functions->get(dec_j, sp_cache.old_symb_id);
-            } else {
-                // Look up the symbol of the non-reference site, which hasn't changed.
-                const int trans_index = trans_matrix.lookup_in_row(tm_row, site_index);
-                sp_const *= basis_functions->get(dec_j, symbols_with_id->id(trans_index));
-            }
-        }
-        // The change in spin-product is the difference in SP between the site(s) which
-        // changed, multiplied by the constant SP from the other un-changed sites (since
-        // these contribute equally before and after the change).
-        sp_delta += (sp_temp_new - sp_temp_ref) * sp_const;
-    }
-    return sp_delta;
-}
-
-void CEUpdater::update_cf(PyObject *single_change) {
-    SymbolChange symb_change = SymbolChange(single_change);
-    update_cf(symb_change);
-}
-
-void CEUpdater::py_changes2_symb_changes(PyObject *all_changes,
-                                         vector<SymbolChange> &symb_changes) {
-    unsigned int size = list_size(all_changes);
-    for (unsigned int i = 0; i < size; i++) {
-        SymbolChange symb_change = SymbolChange(PyList_GetItem(all_changes, i));
-        symb_changes.push_back(symb_change);
-    }
-}
-
-SpinProductCache CEUpdater::build_sp_cache(const SymbolChange &symb_change) const {
-    int ref_indx = symb_change.indx;
-    // Look up the untranslated index of the reference index.
-    int orig_indx = this->get_original_index(ref_indx);
-    // Cache the relevant row from the trans matrix.
-    int *trans_matrix_row = this->trans_matrix.get_row(ref_indx);
-
-    unsigned int old_symb_id = symbols_with_id->get_symbol_id(symb_change.old_symb);
-    unsigned int new_symb_id = symbols_with_id->get_symbol_id(symb_change.new_symb);
-
-    SpinProductCache sp_cache = {ref_indx, orig_indx, trans_matrix_row, new_symb_id, old_symb_id};
-    return sp_cache;
-}
-
-cf &CEUpdater::get_next_cf(SymbolChange &symb_change) {
-    SymbolChange *symb_change_track;
-    cf *next_cf_ptr = nullptr;
-    history->get_next(&next_cf_ptr, &symb_change_track);
-    cf &next_cf = *next_cf_ptr;
-
-    symb_change_track->indx = symb_change.indx;
-    symb_change_track->old_symb = symb_change.old_symb;
-    symb_change_track->new_symb = symb_change.new_symb;
-    symb_change_track->track_indx = symb_change.track_indx;
-    return next_cf;
-}
-
-void CEUpdater::update_cf(SymbolChange &symb_change) {
-    if (symb_change.old_symb == symb_change.new_symb) {
-        return;
-    }
-
-    if (is_background_index[symb_change.indx]) {
-        throw runtime_error("Attempting to move a background atom!");
-    }
-
-    cf &current_cf = history->get_current();
-    cf &next_cf = get_next_cf(symb_change);
-
-    symbols_with_id->set_symbol(symb_change.indx, symb_change.new_symb);
-
-    /* The following prepares a range of properties which will be
-     useful for all of the clusters, so we don't compute more
-     than we have to inside the main ECI loop */
-    SpinProductCache sp_cache = this->build_sp_cache(symb_change);
-
-    if (atoms != nullptr) {
-        set_symbol_in_atoms(atoms, symb_change.indx, symb_change.new_symb);
-    }
-
-    int symm = trans_symm_group[symb_change.indx];
-    const std::vector<ClusterCache> &clusters_cache = m_cluster_by_symm[symm];
-
-    // Loop over all ECIs
-    // As work load for different clusters are different due to a different
-    // multiplicity factor, we need to apply a dynamic schedule
-#ifdef HAS_OMP
-    bool is_par = this->cf_update_num_threads > 1;
-#pragma omp parallel for if (is_par) num_threads(this->cf_update_num_threads) schedule(dynamic)
-#endif
-    for (unsigned int i = 0; i < eci.size(); i++) {
-        // The pre-parsed version of the cluster name.
-        const ParsedName &parsed = this->m_parsed_names[i];
-
-        // 0-body
-        if (parsed.size == 0) {
-            next_cf[i] = current_cf[i];
-            continue;
-        }
-
-        // Singlet
-        if (parsed.size == 1) {
-            unsigned int dec = parsed.dec_num;
-            double new_bf = basis_functions->get(dec, sp_cache.new_symb_id);
-            double old_bf = basis_functions->get(dec, sp_cache.old_symb_id);
-            next_cf[i] = current_cf[i] + (new_bf - old_bf) / num_non_bkg_sites;
-            continue;
-        }
-
-        // n-body
-        const ClusterCache &cluster_cache = clusters_cache[i];
-        const Cluster *cluster_ptr = cluster_cache.cluster_ptr;
-
-        if (cluster_ptr == nullptr) {
-            // This cluster was not present in the symmetry group.
-            next_cf[i] = current_cf[i];
-            continue;
-        }
-        // The cluster is in the symmetry group, so calculate the spin product
-        // change for this cluster.
-        const Cluster &cluster = *cluster_ptr;
-        const equiv_deco_t &equiv_deco = *cluster_cache.equiv_deco_ptr;
-
-        double delta_sp = 0.0;
-        // Calculate the change (delta) in spin product
-        for (const deco_t &deco : equiv_deco) {
-            if (this->assume_no_self_interactions) {
-                // Faster version for large cells with no self interaction
-                delta_sp += spin_product_one_atom_delta_no_si(sp_cache, cluster, deco);
-            } else {
-                // Safe fall-back version
-                delta_sp += spin_product_one_atom_delta(sp_cache, cluster, deco);
-            }
-        }
-
-        delta_sp *= cluster_cache.normalization;
-        next_cf[i] = current_cf[i] + delta_sp;
-    }
-}
-
-void CEUpdater::undo_changes() {
-    unsigned int buf_size = history->history_size();
-    undo_changes(buf_size - 1);
-}
-
-void CEUpdater::undo_changes(int num_steps) {
-    int buf_size = history->history_size();
-
-    if (num_steps > buf_size - 1) {
-        throw invalid_argument("Can't reset history beyond the buffer size!");
-    }
-
-    SymbolChange *last_changes;
-    for (int i = 0; i < num_steps; i++) {
-        history->pop(&last_changes);
-        symbols_with_id->set_symbol(last_changes->indx, last_changes->old_symb);
-
-        if (atoms != nullptr) {
-            set_symbol_in_atoms(atoms, last_changes->indx, last_changes->old_symb);
-        }
-    }
-}
-
-double CEUpdater::calculate(PyObject *system_changes) {
-    unsigned int size = list_size(system_changes);
-    if (size == 0) {
-        return get_energy();
-    } else if (size == 1) {
-        for (unsigned int i = 0; i < size; i++) {
-            update_cf(PyList_GetItem(system_changes, i));
-        }
-        return get_energy();
-    }
-
-    if (size % 2 == 0) {
-        bool sequence_arbitrary_moves = false;
-        vector<swap_move> sequence;
-        for (unsigned int i = 0; i < size / 2; i++) {
-            swap_move changes;
-            changes[0] = SymbolChange(PyList_GetItem(system_changes, 2 * i));
-            changes[1] = SymbolChange(PyList_GetItem(system_changes, 2 * i + 1));
-
-            if (!is_swap_move(changes)) {
-                sequence_arbitrary_moves = true;
-                break;
-            }
-            sequence.push_back(changes);
-        }
-
-        if (!sequence_arbitrary_moves) {
-            return calculate(sequence);
-        }
-    }
-
-    // Last option is that this is a sequence of arbitrary moves
-    vector<SymbolChange> changes(size);
-    for (unsigned int i = 0; i < size; i++) {
-        changes[i] = SymbolChange(PyList_GetItem(system_changes, i));
-    }
-    return calculate(changes);
-}
-
-double CEUpdater::calculate(swap_move &system_changes) {
-    if (symbols_with_id->id(system_changes[0].indx) ==
-        symbols_with_id->id(system_changes[1].indx)) {
-        cout << system_changes[0] << endl;
-        cout << system_changes[1] << endl;
-        throw runtime_error(
-            "This version of the calculate function assumes that the provided update is swapping "
-            "two atoms\n");
-    }
-
-    if (symbols_with_id->get_symbol(system_changes[0].indx) != system_changes[0].old_symb) {
-        throw runtime_error("The atom position tracker does not match the current state\n");
-    } else if (symbols_with_id->get_symbol(system_changes[1].indx) != system_changes[1].old_symb) {
-        throw runtime_error("The atom position tracker does not match the current state\n");
-    }
-
-    // Update correlation function
-    update_cf(system_changes[0]);
-    update_cf(system_changes[1]);
-
-    return get_energy();
-}
-
-void CEUpdater::clear_history() {
-    history->clear();
-}
-
-void CEUpdater::flattened_cf_names(vector<string> &flattened) {
-    flattened = eci.get_names();
-
-    // Sort the cluster names for consistency
-    sort(flattened.begin(), flattened.end());
-}
-
-PyObject *CEUpdater::get_cf() {
-    PyObject *cf_dict = PyDict_New();
-    cf &corrfunc = history->get_current();
-
-    for (unsigned int i = 0; i < corrfunc.size(); i++) {
-        PyObject *pyvalue = PyFloat_FromDouble(corrfunc[i]);
-        PyDict_SetItemString(cf_dict, corrfunc.name(i).c_str(), pyvalue);
-        Py_DECREF(pyvalue);
-    }
-    return cf_dict;
-}
-
-void CEUpdater::set_symbols(const vector<string> &new_symbs) {
-    if (new_symbs.size() != symbols_with_id->size()) {
-        throw runtime_error(
-            "The number of atoms in the updater cannot be changed via the set_symbols function\n");
-    }
-    symbols_with_id->set_symbols(new_symbs);
-}
-
-void CEUpdater::cluster_by_symm_group() {
-    m_cluster_by_symm.clear();
-
-    // Find unique symmetry groups
-    std::set<int> unique;
-    insert_in_set(this->trans_symm_group, unique);
-
-    for (const int symm : unique) {
-        if (symm == -1) {
-            // Background symmetry group
-            continue;
-        }
-        // 1 ClusterCache per ECI value
-        std::vector<ClusterCache> cluster_cache;
-        cluster_cache.reserve(this->m_parsed_names.size());
-        for (const ParsedName &parsed : this->m_parsed_names) {
-            ClusterCache cache;
-            if (parsed.size == 0 || parsed.size == 1 ||
-                !clusters.is_in_symm_group(parsed.prefix, symm)) {
-                /* Either 0- or 1-body cluster, or cluster is not in this
-                symmetry group. */
-                cluster_cache.push_back(cache);
-                continue;
-            }
-            Cluster *cluster_ptr = clusters.get_ptr(parsed.prefix, symm);
-            equiv_deco_t *equiv_ptr = cluster_ptr->get_equiv_deco_ptr(parsed.dec_str);
-            // Calculate the normalization of the resulting cluster functions
-            double normalization = static_cast<double>(cluster_ptr->size);
-            normalization /= equiv_ptr->size();
-            normalization /= normalisation_factor.at(parsed.prefix);
-
-            // Populate the new cache object
-            cache.cluster_ptr = cluster_ptr;
-            cache.equiv_deco_ptr = equiv_ptr;
-            cache.normalization = normalization;
-
-            cluster_cache.push_back(cache);
-        }
-        m_cluster_by_symm.insert({symm, cluster_cache});
-    }
-}
-void CEUpdater::set_eci(PyObject *pyeci) {
-    PyObject *key;
-    PyObject *value;
-
-    // Read the ECIs
-    Py_ssize_t pos = 0;
-    std::map<std::string, double> temp_eci;
-
-    while (PyDict_Next(pyeci, &pos, &key, &value)) {
-        temp_eci[py2string(key)] = PyFloat_AS_DOUBLE(value);
-    }
-    this->eci.init(temp_eci);
-    // Pre-parse the names of the clusters.
-    this->parse_eci_names();
-
-    // If status is not READY, then we're still initializing, and CF's are missing.
-    if (this->status == Status_t::READY && !all_eci_corresponds_to_cf()) {
-        throw invalid_argument("All ECIs has to correspond to a correlation function!");
-    }
-    // Update the cluster pointers to match the order with the ECI's.
-    cluster_by_symm_group();
-}
-
-bool CEUpdater::all_decoration_nums_equal(const vector<int> &dec_nums) const {
-    for (unsigned int i = 1; i < dec_nums.size(); i++) {
-        if (dec_nums[i] != dec_nums[0]) {
-            return false;
-        }
-    }
-    return true;
-}
-
-void CEUpdater::get_singlets(PyObject *npy_obj) const {
-    PyObject *npy_array = PyArray_FROM_OTF(npy_obj, NPY_DOUBLE, NPY_OUT_ARRAY);
-    unsigned int npy_array_size = PyArray_SIZE(npy_array);
-    if (npy_array_size < singlets.size()) {
-        string msg("The passed Numpy array is too small to hold all the singlets terms!\n");
-        stringstream ss;
-        ss << "Minimum size: " << singlets.size() << ". Given size: " << npy_array_size;
-        msg += ss.str();
-        Py_DECREF(npy_array);
-        throw runtime_error(msg);
-    }
-    cf &cfs = history->get_current();
-    for (unsigned int i = 0; i < singlets.size(); i++) {
-        double *ptr = static_cast<double *>(PyArray_GETPTR1(npy_array, i));
-        *ptr = cfs[singlets[i]];
-    }
-    Py_DECREF(npy_array);
-}
-
-PyObject *CEUpdater::get_singlets() const {
-    npy_intp dims[1] = {static_cast<npy_intp>(singlets.size())};
-    PyObject *npy_array = PyArray_SimpleNew(1, dims, NPY_DOUBLE);
-    get_singlets(npy_array);
-    return npy_array;
-}
-
-void CEUpdater::create_cname_with_dec(PyObject *cf) {
-    if (!PyDict_Check(cf)) {
-        throw invalid_argument("Correlation functons has to be dictionary!");
-    }
-    Py_ssize_t pos = 0;
-    PyObject *key;
-    PyObject *value;
-    while (PyDict_Next(cf, &pos, &key, &value)) {
-        string new_key = py2string(key);
-        ClusterName c_name = ClusterName(new_key);
-        unsigned int size = c_name.get_size();
-#ifdef PRINT_DEBUG
-        cout << "Read CF: " << new_key << endl;
-#endif
-        if ((size == 0 || size == 1)) {
-            cname_with_dec[new_key] = new_key;
-        } else {
-            std::string prefix = c_name.get_prefix();
-            cname_with_dec[prefix] = new_key;
-        }
-    }
-}
-
-void CEUpdater::build_trans_symm_group(PyObject *py_trans_symm_group) {
-    // Fill the symmetry group array with -1 indicating an invalid value
-    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
-        trans_symm_group[i] = -1;
-    }
-
-    unsigned int py_list_size = list_size(py_trans_symm_group);
-    for (unsigned int i = 0; i < py_list_size; i++) {
-        PyObject *sublist = PyList_GetItem(py_trans_symm_group, i);
-        unsigned int n_sites = list_size(sublist);
-        for (unsigned int j = 0; j < n_sites; j++) {
-            int indx = py2int(PyList_GetItem(sublist, j));
-            if (trans_symm_group[indx] != -1) {
-                throw runtime_error(
-                    "One site appears to be present in more than one translation symmetry "
-                    "group!");
-            }
-            trans_symm_group[indx] = i;
-        }
-    }
-
-    // Check that all sites belongs to one translational symmetry group
-    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
-        if ((trans_symm_group[i] == -1) && !is_background_index[i]) {
-            stringstream msg;
-            msg << "Site " << i << " has not been assigned to any translational symmetry group!";
-            throw runtime_error(msg.str());
-        }
-    }
-
-    // Count the number of atoms in each symmetry group
-    trans_symm_group_count.resize(py_list_size);
-    fill(trans_symm_group_count.begin(), trans_symm_group_count.end(), 0);
-    for (unsigned int i = 0; i < trans_symm_group.size(); i++) {
-        if (trans_symm_group[i] >= 0) {
-            trans_symm_group_count[trans_symm_group[i]] += 1;
-        }
-    }
-}
-
-bool CEUpdater::all_eci_corresponds_to_cf() {
-    cf &corrfunc = history->get_current();
-    return eci.names_are_equal(corrfunc);
-}
-
-double CEUpdater::calculate(vector<swap_move> &sequence) {
-    if (sequence.size() >= history->max_history / 2) {
-        throw invalid_argument(
-            "The length of sequence of swap move exceeds the buffer size for the history "
-            "tracker");
-    }
-
-    for (unsigned int i = 0; i < sequence.size(); i++) {
-        calculate(sequence[i]);
-    }
-    return get_energy();
-}
-
-double CEUpdater::calculate(vector<SymbolChange> &sequence) {
-    for (auto &change : sequence) {
-        update_cf(change);
-    }
-    return get_energy();
-}
-
-void CEUpdater::read_trans_matrix(PyObject *py_trans_mat) {
-    bool is_list = PyList_Check(py_trans_mat);
-
-#ifdef PRINT_DEBUG
-    cout << "read_trans_matrix: Extracting unique indices" << endl;
-#endif
-
-    set<int> unique_indx;
-    clusters.unique_indices(unique_indx);
-    vector<int> unique_indx_vec;
-    set2vector(unique_indx, unique_indx_vec);
-
-    // Compute the max index that is ever going to be checked
-    unsigned int max_indx = clusters.max_index();
-    if (is_list) {
-        unsigned int size = list_size(py_trans_mat);
-
-#ifdef PRINT_DEBUG
-        cout << "read_trans_matrix: Updating size of trans_matrix" << endl;
-#endif
-
-        trans_matrix.set_size(size, unique_indx_vec.size(), max_indx);
-
-#ifdef PRINT_DEBUG
-        cout << "read_trans_matrix: Setting lookup values in trans_matrix" << endl;
-#endif
-        trans_matrix.set_lookup_values(unique_indx_vec);
-
-#ifdef PRINT_DEBUG
-        cout << "read_trans_matrix: Reading translation matrix from list of dictionaries" << endl;
-#endif
-        unsigned int n_elements_insterted = 0;
-        for (unsigned int i = 0; i < size; i++) {
-            // Background atoms are ignored (and should never be accessed)
-            if (is_background_index[i] && ignore_background_indices) {
-                continue;
-            }
-
-            PyObject *dict = PyList_GetItem(py_trans_mat, i);
-            for (unsigned int j = 0; j < unique_indx_vec.size(); j++) {
-                int col = unique_indx_vec[j];
-                PyObject *value = PyDict_GetItem(dict, int2py(col));
-
-                if (value == NULL) {
-                    stringstream ss;
-                    ss << "Requested value " << col << " is not a key in the dictionary!";
-                    throw invalid_argument(ss.str());
-                }
-                trans_matrix(i, col) = py2int(value);
-                n_elements_insterted++;
-            }
-        }
-#ifdef PRINT_DEBUG
-        cout << "Inserted " << n_elements_insterted << " into the translation matrix\n";
-#endif
-    } else {
-        PyObject *trans_mat = PyArray_FROM_OTF(py_trans_mat, NPY_INT32, NPY_ARRAY_IN_ARRAY);
-
-        npy_intp *size = PyArray_DIMS(trans_mat);
-        trans_matrix.set_size(size[0], unique_indx_vec.size(), max_indx);
-        trans_matrix.set_lookup_values(unique_indx_vec);
-
-#ifdef PRINT_DEBUG
-        cout << "Dimension of translation matrix stored: " << size[0] << " "
-             << unique_indx_vec.size() << endl;
-#endif
-
-        if (max_indx + 1 > size[1]) {
-            stringstream ss;
-            ss << "Something is wrong with the translation matrix passed.\n";
-            ss << "Shape of translation matrix (" << size[0] << "," << size[1] << ")\n";
-            ss << "Maximum index encountered in the cluster lists: " << max_indx << endl;
-            throw invalid_argument(ss.str());
-        }
-        for (unsigned int i = 0; i < size[0]; i++)
-            for (unsigned int j = 0; j < unique_indx_vec.size(); j++) {
-                int col = unique_indx_vec[j];
-                trans_matrix(i, col) = *static_cast<int *>(PyArray_GETPTR2(trans_mat, i, col));
-            }
-        Py_DECREF(trans_mat);
-    }
-}
-
-void CEUpdater::sort_indices(int indices[], const vector<int> &order, unsigned int n_indices) {
-    // This function is called many times
-    // profiling (with YEP) revealed that
-    // [] operator of the vector used quite a bit of time
-    // Therefore we here use raw C-arrays or pointer arithmetics
-    int sorted[4];
-    const int *ptr = &order[0];
-    for (unsigned int i = 0; i < n_indices; i++) {
-        sorted[i] = indices[*(ptr + i)];
-    }
-    memcpy(indices, sorted, n_indices * sizeof(int));
-}
-
-bool CEUpdater::is_swap_move(const swap_move &move) const {
-    return (move[0].old_symb == move[1].new_symb) && (move[1].new_symb == move[1].old_symb);
-}
-
-void CEUpdater::read_background_indices(PyObject *bkg_indices) {
-    // Fill array with false
-    is_background_index.resize(symbols_with_id->size());
-    fill(is_background_index.begin(), is_background_index.end(), false);
-
-    // Set to true if index is in bkg_indices
-    int size = list_size(bkg_indices);
-    for (int i = 0; i < size; i++) {
-        PyObject *py_indx = PyList_GetItem(bkg_indices, i);
-        int indx = py2int(py_indx);
-        is_background_index[indx] = true;
-    }
-}
-
-void CEUpdater::count_non_bkg_sites() {
-    // Count and store the number of non-background sites
-    num_non_bkg_sites = 0;
-    for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
-        if (!is_background_index[atom_no] || !ignore_background_indices) {
-            num_non_bkg_sites += 1;
-        }
-    }
-}
-
-void CEUpdater::get_changes(const std::vector<std::string> &new_symbols,
-                            std::vector<unsigned int> &changed_sites) const {
-    if (new_symbols.size() != symbols_with_id->size()) {
-        throw invalid_argument("Size of passed atoms does not match!");
-    }
-
-    for (unsigned int i = 0; i < new_symbols.size(); i++) {
-        unsigned int symb_id = symbols_with_id->id(i);
-        if (symbols_with_id->get_symbol_id(new_symbols[i]) != symb_id) {
-            changed_sites.push_back(i);
-        }
-    }
-}
-
-void CEUpdater::calculate_cf_from_scratch(const vector<string> &cf_names, map<string, double> &cf) {
-    cf.clear();
-
-    // Initialise all cluster names
-    for (const string &name : cf_names) {
-        cf[name] = 0.0;
-    }
-
-    // Loop over all clusters
-    for (const string &name : cf_names) {
-        ClusterName c_name = ClusterName(name);
-        unsigned int cluster_size = c_name.get_size();
-
-        // Handle empty cluster
-        if (cluster_size == 0) {
-            cf[name] = 1.0;
-            continue;
-        }
-
-        // Handle singlet cluster
-        if (cluster_size == 1) {
-            unsigned int dec = c_name.get_dec_num();
-            double new_value = 0.0;
-            // Normalise with respect to the actual number of atoms included
-            for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
-                if (!is_background_index[atom_no] || !ignore_background_indices) {
-                    new_value += basis_functions->get(dec, symbols_with_id->id(atom_no));
-                }
-            }
-            cf[name] = new_value / num_non_bkg_sites;
-            continue;
-        }
-
-        // Handle the rest of the clusters
-        std::string prefix, dec_str;
-        c_name.get_prefix_and_dec_str(prefix, dec_str);
-
-        double sp = 0.0;
-        double count = 0;
-        for (unsigned int atom_no = 0; atom_no < symbols_with_id->size(); atom_no++) {
-            int symm = trans_symm_group[atom_no];
-
-            if ((!clusters.is_in_symm_group(prefix, symm)) ||
-                (is_background_index[atom_no] && ignore_background_indices)) {
-                continue;
-            }
-
-            const Cluster &cluster = clusters.get(prefix, symm);
-
-            const equiv_deco_t &equiv_deco = cluster.get_equiv_deco(dec_str);
-            unsigned int ref_id = symbols_with_id->id(atom_no);
-            double sp_temp = 0.0;
-            for (const vector<int> &deco : equiv_deco) {
-                sp_temp += spin_product_one_atom(atom_no, cluster, deco, ref_id);
-            }
-            sp += sp_temp / equiv_deco.size();
-            count += cluster.get().size();
-        }
-
-        if (count == 0) {
-            cf[name] = 0.0;
-        } else {
-            cf[name] = sp / count;
-        }
-    }
-
-    history->get_current().init(cf);
-}
-
-void CEUpdater::set_atoms(PyObject *py_atoms) {
-    unsigned int num_atoms = PySequence_Length(py_atoms);
-
-    if (num_atoms != symbols_with_id->size()) {
-        throw invalid_argument("Length of passed atoms object is different from current");
-    }
-
-    std::vector<std::string> symbols = get_symbols_from_atoms(py_atoms);
-
-    this->atoms = py_atoms;
-    symbols_with_id->set_symbols(symbols);
-}

Alert delta unavailable

Currently unable to show alert delta for PyPI packages.

Worsened metrics

Total package byte prevSize

1932976

-3.75%

Lines of code

18988

-0.11%