cell_manager.hh

Go to the documentation of this file.

#ifndef UTOPIA_CORE_CELL_MANAGER_HH
#define UTOPIA_CORE_CELL_MANAGER_HH
 
#include <algorithm>
#include <type_traits>
#include <string_view>
#include <unordered_set>
 
#include <armadillo>
 
#include "logging.hh"
#include "types.hh"
#include "exceptions.hh"
#include "cell.hh"
#include "grids.hh"
#include "apply.hh"
#include "select.hh"
 
 
namespace Utopia {
 
template<class CellTraits, class Model>
class CellManager {
public:
    using Self = CellManager<CellTraits, Model>;
 
    using Cell = Utopia::Cell<CellTraits>;
 
    using Entity = Cell;
 
    using CellState = typename CellTraits::State;
 
    using Space = typename Model::Space;
 
    static constexpr DimType dim = Model::Space::dim;
 
    using GridType = Grid<Space>;
 
    using SpaceVec = SpaceVecType<dim>;
 
    using MultiIndex = typename GridType::MultiIndex;
 
    using RNG = typename Model::RNG;
 
    using Config = typename Model::Config;
 
    using NBFuncCell = std::function<CellContainer<Cell>(const Cell&)>;
 
 
    using RuleFunc =
        typename std::function<CellState(const std::shared_ptr<Cell>&)>;
 
 
    using VoidRuleFunc =
        typename std::function<void(const std::shared_ptr<Cell>&)>;
 
 
private:
    // -- Members -------------------------------------------------------------
    const std::shared_ptr<spdlog::logger> _log;
 
    const Config _cfg;
 
    const std::shared_ptr<RNG> _rng;
 
    const std::shared_ptr<Space> _space;
 
    const std::shared_ptr<GridType> _grid;
 
    CellContainer<Cell> _cells;
 
    std::vector<CellContainer<Cell>> _cell_neighbors;
 
    NBFuncCell _nb_func;
 
public:
 
    bool _empty_nb_warning_emitted = false;
 
 
    // -- Constructors --------------------------------------------------------
 
    CellManager (const Model& model,
                 const Config& custom_cfg = {})
    :
        _log(model.get_logger()),
        _cfg(setup_cfg(model, custom_cfg)),
        _rng(model.get_rng()),
        _space(model.get_space()),
        _grid(setup_grid()),
        _cells(setup_cells()),
        _cell_neighbors()
    {
        // Set default value for _nb_func
        _nb_func = _nb_compute_each_time_empty;
 
        // Use setup function to set up neighborhood from configuration
        setup_nb_funcs();
 
        _log->info("CellManager is all set up.");
    }
 
 
 
    CellManager (const Model& model,
                 const CellState initial_state,
                 const Config& custom_cfg = {})
    :
        _log(model.get_logger()),
        _cfg(setup_cfg(model, custom_cfg)),
        _rng(model.get_rng()),
        _space(model.get_space()),
        _grid(setup_grid()),
        _cells(setup_cells(initial_state)),
        _cell_neighbors()
    {
        // Set default value for _nb_func
        _nb_func = _nb_compute_each_time_empty;
 
        // Use setup function to set up neighborhood from configuration
        setup_nb_funcs();
 
        _log->info("CellManager is all set up.");
    }
 
 
    const auto& log () const {
        return _log;
    }
 
    const Config& cfg () const {
        return _cfg;
    }
 
    const std::shared_ptr<RNG>& rng () const {
        return _rng;
    }
 
    const std::shared_ptr<Space>& space () const {
        return _space;
    }
 
    const std::shared_ptr<GridType>& grid () const {
        return _grid;
    }
 
    const CellContainer<Cell>& cells () const {
        return _cells;
    }
 
    const CellContainer<Cell>& entities () const {
        return cells();
    }
 
 
    // -- Public interface ----------------------------------------------------
    // .. Position-related ....................................................
 
 
    MultiIndex midx_of(const Cell& cell) const {
        return _grid->midx_of(cell.id());
    }
 
 
    MultiIndex midx_of(const std::shared_ptr<Cell>& cell) const {
        return _grid->midx_of(cell->id());
    }
 
    SpaceVec barycenter_of(const Cell& cell) const {
        return _grid->barycenter_of(cell.id());
    }
 
    SpaceVec barycenter_of(const std::shared_ptr<Cell>& cell) const {
        return _grid->barycenter_of(cell->id());
    }
 
    SpaceVec extent_of(const Cell& cell) const {
        return _grid->extent_of(cell.id());
    }
 
    SpaceVec extent_of(const std::shared_ptr<Cell>& cell) const {
        return _grid->extent_of(cell->id());
    }
 
 
    std::vector<SpaceVec> vertices_of(const Cell& cell) const {
        return _grid->vertices_of(cell.id());
    }
 
 
    std::vector<SpaceVec> vertices_of(const std::shared_ptr<Cell>& cell) const
    {
        return _grid->vertices_of(cell->id());
    }
 
 
    const std::shared_ptr<Cell>& cell_at(const SpaceVec& pos) const {
        return _cells[_grid->cell_at(pos)];
    }
 
    // .. Cell Selection ......................................................
 
    CellContainer<Cell> boundary_cells(const std::string& select="all") const {
        if (_space->periodic) {
            _log->warn("Selecting boundary cells (mode '{}') of a periodic "
                       "space will always return an empty container!", select);
        }
 
        return entity_pointers_from_ids(_grid->boundary_cells(select));
    }
 
 
    template<SelectionMode mode, class... Args>
    CellContainer<Cell> select_cells(Args&&... args) const {
        return select_entities<mode>(*this, std::forward<Args>(args)...);
    }
 
 
    CellContainer<Cell> select_cells(const Config& sel_cfg) const {
        return select_entities(*this, sel_cfg);
    }
 
 
    // .. Setting the cell states .............................................
 
 
    template<class ElementT=double, class SetterFunc>
    void set_cell_states(const std::string& hdf5_file,
                         const std::string& dset_path,
                         const SetterFunc& setter_func)
    {
        static_assert(dim == 2, "Loading cell state is only supported in 2D!");
        static_assert(CellTraits::mode == Update::manual,
                      "Setting cell states for cell update modes other than "
                      "Update::manual is currently not supported.");
        // TODO Once passing of additional arguments to rule functions becomes
        //      possible, remove this restriction and use apply_rule below!
 
        _log->debug("Setting cell states using HDF5 data ...");
        _log->debug("  File:          {}", hdf5_file);
        _log->debug("  Dataset path:  {}", dset_path);
 
        // Load the array data
        arma::Mat<ElementT> data;
        data.load(arma::hdf5_name(hdf5_file, dset_path,
                                  arma::hdf5_opts::trans));
        // NOTE Need to load transposed; otherwise the index-access below would
        //      be reversed. This is to account for the different ordering
        //      between armadillo and CellManager coordinates.
        //      Thus: n_rows corresponds to N_x, n_cols to N_y!
 
        // Check that data was loaded
        if (not data.size()) {
            throw std::runtime_error("Failed loading HDF5 data! Is the file "
                                     "used by another program?");
        }
 
        // Check against grid shape
        const auto grid_shape = _grid->shape();
        if (data.n_rows != grid_shape[0] or data.n_cols != grid_shape[1]) {
            throw std::invalid_argument("Shape mismatch between loaded data ("
                + std::to_string(data.n_rows) + ", "
                + std::to_string(data.n_cols) + ") and grid ("
                + std::to_string(grid_shape[0]) + ", "
                + std::to_string(grid_shape[1]) + ")!");
        }
 
        // Go over all cells, extract multi index, and let the setter function
        // set the cell state.
        for (const auto& cell : _cells) {
            const auto midx = this->midx_of(cell);
            setter_func(cell, data(midx[0], midx[1]));
        }
 
        _log->debug("Cell states set successfully.");
    }
 
    // .. Neighborhood-related ................................................
 
 
    const NBMode& nb_mode () const {
        return _grid->nb_mode();
    }
 
 
    auto nb_size () const {
        return _grid->nb_size();
    }
 
 
    CellContainer<Cell> neighbors_of(const Cell& cell) const {
        return _nb_func(cell);
    }
 
 
    CellContainer<Cell> neighbors_of(const std::shared_ptr<Cell>& cell) const {
        return _nb_func(*cell);
    }
 
 
    void select_neighborhood(const Config& nb_cfg) {
        // Extract the desired values
        if (not nb_cfg["mode"]) {
            throw KeyError("mode", nb_cfg, "Could not select neighborhood!");
        }
        const auto nb_mode = get_as<std::string>("mode", nb_cfg);
        const bool compute = get_as<bool>("compute_and_store", nb_cfg, true);
 
        // Call the string-based selection function, passing through the whole
        // config node. The fact that the above two keys are also present in
        // the node is not a problem.
        select_neighborhood(nb_mode, compute, nb_cfg);
    }
 
 
    void select_neighborhood(const std::string& nb_mode,
                             const bool compute_and_store = true,
                             const Config& nb_params = {})
    {
        // Check if the string is valid
        if (not nb_mode_map.count(nb_mode)) {
            throw std::invalid_argument("Got unexpected neighborhood mode '"
                + nb_mode + "'! Available modes: empty, vonNeumann, Moore, "
                "hexagonal.");
        }
 
        // Translate string; pass all other arguments through
        select_neighborhood(nb_mode_map.at(nb_mode),
                            compute_and_store, nb_params);
    }
 
 
    void select_neighborhood(const NBMode& nb_mode,
                             const bool compute_and_store = true,
                             const Config& nb_params = {})
    {
        // Only change the neighborhood if it is different to the one already
        // set in the grid or if it is set to be empty
        if ((nb_mode != _grid->nb_mode()) or (nb_mode == NBMode::empty)) {
            _log->info("Selecting '{}' neighborhood ...",
                       nb_mode_to_string(nb_mode));
 
            // Tell the grid which mode to use
            _grid->select_neighborhood(nb_mode, nb_params);
 
            // Adjust function object that the public interface calls
            if (nb_mode == NBMode::empty) {
                // Issue a warning alongside the neighborhood calculation
                _nb_func = _nb_compute_each_time_empty;
            }
            else {
                // Compute the cell neighbors each time
                _nb_func = _nb_compute_each_time;
            }
 
            // Clear the no-longer valid neighborhood relationships
            if (_cell_neighbors.size() > 0) {
                _cell_neighbors.clear();
                _log->debug("Cleared cell neighborhood cache.");
            }
 
            _log->debug("Successfully selected '{}' neighborhood (size: {}).",
                        nb_mode_to_string(_grid->nb_mode()), nb_size());
        }
        else {
            _log->debug("Neighborhood was already set to '{}'; not changing.",
                        nb_mode_to_string(_grid->nb_mode()));
        }
 
        // Still allow to compute the neighbors, regardless of all the above
        if (compute_and_store) {
            compute_cell_neighbors();
        }
    }
 
 
    void compute_cell_neighbors() {
        _log->info("Computing and storing '{}' neighbors of all {} cells ...",
                   nb_mode_to_string(_grid->nb_mode()), _cells.size());
 
        // Clear cell neighbors container and pre-allocate space
        _cell_neighbors.clear();
        _cell_neighbors.reserve(_cells.size());
 
        // Compute cell neighbors and store them
        for (const auto& cell : _cells) {
            _cell_neighbors.push_back(neighbors_of(cell));
        }
 
        // Change access function to access the storage directly. Done.
        _nb_func = _nb_from_cache;
        _log->info("Computed and stored cell neighbors.");
    }
 
 
    // -- Public Helpers ------------------------------------------------------
 
 
    template<class IndexContainer>
    CellContainer<Cell> entity_pointers_from_ids(IndexContainer&& ids) const {
        CellContainer<Cell> ret;
        ret.reserve(ids.size());
 
        for (const auto& id : ids) {
            ret.emplace_back(std::shared_ptr<Cell>(_cells[id]));
        }
 
        return ret;
    }
 
private:
    // -- Private Helpers -----------------------------------------------------
    // ...
 
 
    // -- std::functions to call from neighbors_of ----------------------------
 
    NBFuncCell _nb_from_cache = [this](const Cell& cell) {
        return this->_cell_neighbors[cell.id()];
    };
 
    NBFuncCell _nb_compute_each_time = [this](const Cell& cell) {
        return
            this->entity_pointers_from_ids(
                this->_grid->neighbors_of(cell.id()));
    };
 
    NBFuncCell _nb_compute_each_time_empty = [this](const Cell& cell) {
        if (not this->_empty_nb_warning_emitted) {
            this->_log->warn("No neighborhood selected! Calls to the "
                "CellManager::neighbors_of method will always return an empty "
                "container. There will be no further warning.");
            this->_empty_nb_warning_emitted = true;
        }
        return
            this->entity_pointers_from_ids(
                this->_grid->neighbors_of(cell.id()));
    };
 
 
    // -- Setup functions -----------------------------------------------------
 
 
    Config setup_cfg(const Model& model, const Config& custom_cfg) const {
        Config cfg;
 
        if (custom_cfg.size() > 0) {
            _log->debug("Using custom config for cell manager setup ...");
            cfg = custom_cfg;
        }
        else {
            _log->debug("Using '{}' model's configuration for cell manager "
                        "setup ... ", model.get_name());
 
            if (not model.get_cfg()["cell_manager"]) {
                throw std::invalid_argument("Missing config entry "
                    "'cell_manager' in model configuration! Either specify "
                    "that key or pass a custom configuration node to the "
                    "CellManager constructor.");
            }
            cfg = model.get_cfg()["cell_manager"];
        }
        return cfg;
    }
 
    std::shared_ptr<GridType> setup_grid() const {
        // Check if the required parameter nodes are available
        // TODO Throw KeyErrors here instead!
        if (not _cfg["grid"]) {
            throw std::invalid_argument("Missing entry 'grid' in the "
                "configuration node supplied to the CellManager! Check that "
                "the model configuration includes such an entry.");
        }
        else if (not _cfg["grid"]["structure"]) {
            throw std::invalid_argument("Missing required grid configuration "
                                        "entry 'structure'!");
        }
 
        // Get the structure parameter
        auto structure = get_as<std::string>("structure", _cfg["grid"]);
 
        _log->info("Setting up grid discretization with '{}' cells ...",
                   structure);
 
        // Create the respective grids, distinguishing by structure
        if (structure == "triangular") {
            using GridSpec = TriangularGrid<Space>;
            return std::make_shared<GridSpec>(_space, _cfg["grid"]);
        }
        else if (structure == "square") {
            using GridSpec = SquareGrid<Space>;
            return std::make_shared<GridSpec>(_space, _cfg["grid"]);
        }
        else if (structure == "hexagonal") {
            using GridSpec = HexagonalGrid<Space>;
            return std::make_shared<GridSpec>(_space, _cfg["grid"]);
        }
        else {
            throw std::invalid_argument("Invalid value for grid "
                "'structure' argument: '" + structure + "'! Allowed "
                "values: 'square', 'hexagonal', 'triangular'");
        }
    }
 
    CellContainer<Cell> setup_cells(const CellState& initial_state) const {
        CellContainer<Cell> cont;
 
        // Construct all the cells
        for (IndexType i=0; i<_grid->num_cells(); i++) {
            cont.emplace_back(std::make_shared<Cell>(i, initial_state));
        }
 
        // Done. Shrink it.
        cont.shrink_to_fit();
        _log->info("Populated cell container with {:d} cells.", cont.size());
 
        return cont;
    }
 
 
    CellContainer<Cell> setup_cells() const {
        // Distinguish depending on constructor.
        // Is the default constructor to be used?
        if constexpr (CellTraits::use_default_state_constructor) {
            static_assert(
                std::is_default_constructible<CellState>(),
                "CellTraits were configured to use the default constructor to "
                "create cell states, but the CellState is not "
                "default-constructible! Either implement such a constructor, "
                "unset the flag in the CellTraits, or pass an explicit "
                "initial cell state to the CellManager."
            );
 
            _log->info("Setting up cells using default constructor ...");
 
            // Create the initial state (same for all cells)
            return setup_cells(CellState());
        }
 
        // Is there a constructor available that allows passing the RNG?
        else if constexpr (std::is_constructible<CellState,
                                                 const Config&,
                                                 const std::shared_ptr<RNG>&
                                                 >())
        {
            _log->info("Setting up cells using config constructor (with RNG) "
                       "...");
 
            // Extract the configuration parameter
            if (not _cfg["cell_params"]) {
                throw std::invalid_argument("CellManager is missing the "
                    "configuration entry 'cell_params' to set up the cells' "
                    "initial states!");
            }
            const auto cell_params = _cfg["cell_params"];
 
            // The cell container to be populated
            CellContainer<Cell> cont;
 
            // Populate the container, creating the cell state anew each time
            for (IndexType i=0; i<_grid->num_cells(); i++) {
                cont.emplace_back(
                    std::make_shared<Cell>(i, CellState(cell_params, _rng))
                );
            }
            // Done. Shrink it.
            cont.shrink_to_fit();
            _log->info("Populated cell container with {:d} cells.",
                       cont.size());
            return cont;
        }
 
        // As default, require a Config constructor
        else {
            static_assert(
                std::is_constructible<CellState, const Config&>(),
                "CellManager::CellState needs to be constructible using "
                "const Config& as only argument. Either implement "
                "such a constructor, pass an explicit initial cell state to "
                "the CellManager, or set the CellTraits such that a default "
                "constructor is to be used."
            );
 
            _log->info("Setting up cells using config constructor ...");
 
            // Extract the configuration parameter
            if (not _cfg["cell_params"]) {
                throw std::invalid_argument("CellManager is missing the "
                    "configuration entry 'cell_params' to set up the cells' "
                    "initial states!");
            }
 
            // Create the initial state (same for all cells)
            return setup_cells(CellState(_cfg["cell_params"]));
        }
        // This point is never reached.
    }
 
    void setup_nb_funcs() {
        // If there is a neighborhood key, use it to set up the neighborhood
        if (_cfg["neighborhood"]) {
            _log->debug("Setting up neighborhood from config entry ...");
            select_neighborhood(_cfg["neighborhood"]);
            return;
        }
        // else: Use empty.
 
        _log->debug("No neighborhood configuration given; using empty.");
        select_neighborhood(NBMode::empty);
    }
};
 
 
// end group CellManager
} // namespace Utopia
 
#endif // UTOPIA_CORE_CELL_MANAGER_HH