engine_adaptive.h

Go to the documentation of this file.

#pragma once
#include "gridfire/engine/engine_abstract.h"
#include "gridfire/engine/views/engine_view_abstract.h"
#include "gridfire/screening/screening_abstract.h"
#include "gridfire/screening/screening_types.h"
#include "gridfire/network.h"
 
#include "fourdst/composition/atomicSpecies.h"
#include "fourdst/config/config.h"
#include "fourdst/logging/logging.h"
 
#include "gridfire/engine/procedures/priming.h"
#include "gridfire/engine/procedures/construction.h"
 
#include "quill/Logger.h"
 
namespace gridfire {
    class AdaptiveEngineView final : public DynamicEngine, public EngineView<DynamicEngine> {
    public:
        explicit AdaptiveEngineView(DynamicEngine& baseEngine);

        fourdst::composition::Composition update(const NetIn &netIn) override;
 
        bool isStale(const NetIn& netIn) override;

        [[nodiscard]] const std::vector<fourdst::atomic::Species>& getNetworkSpecies() const override;

        [[nodiscard]] std::expected<StepDerivatives<double>, expectations::StaleEngineError> calculateRHSAndEnergy(
            const std::vector<double> &Y_culled,
            const double T9,
            const double rho
        ) const override;

        void generateJacobianMatrix(
            const std::vector<double> &Y_dynamic,
            const double T9,
            const double rho
        ) const override;

        [[nodiscard]] double getJacobianMatrixEntry(
            const int i_culled,
            const int j_culled
        ) const override;

        void generateStoichiometryMatrix() override;

        [[nodiscard]] int getStoichiometryMatrixEntry(
            const int speciesIndex_culled,
            const int reactionIndex_culled
        ) const override;

        [[nodiscard]] double calculateMolarReactionFlow(
            const reaction::Reaction &reaction,
            const std::vector<double> &Y_culled,
            double T9,
            double rho
        ) const override;

        [[nodiscard]] const reaction::LogicalReactionSet& getNetworkReactions() const override;
 
        void setNetworkReactions(const reaction::LogicalReactionSet& reactions) override;

        [[nodiscard]] std::expected<std::unordered_map<fourdst::atomic::Species, double>, expectations::StaleEngineError> getSpeciesTimescales(
            const std::vector<double> &Y_culled,
            double T9,
            double rho
        ) const override;
 
        [[nodiscard]] std::expected<std::unordered_map<fourdst::atomic::Species, double>, expectations::StaleEngineError> getSpeciesDestructionTimescales(
            const std::vector<double> &Y,
            double T9,
            double rho
        ) const override;

        [[nodiscard]] const DynamicEngine& getBaseEngine() const override { return m_baseEngine; }

        void setScreeningModel(screening::ScreeningType model) override;

        [[nodiscard]] screening::ScreeningType getScreeningModel() const override;
 
        [[nodiscard]] int getSpeciesIndex(const fourdst::atomic::Species &species) const override;
 
        [[nodiscard]] std::vector<double> mapNetInToMolarAbundanceVector(const NetIn &netIn) const override;
 
        [[nodiscard]] PrimingReport primeEngine(const NetIn &netIn) override;
    private:
        using Config = fourdst::config::Config;
        using LogManager = fourdst::logging::LogManager;
        Config& m_config = Config::getInstance();
        quill::Logger* m_logger = LogManager::getInstance().getLogger("log");

        DynamicEngine& m_baseEngine;

        std::vector<fourdst::atomic::Species> m_activeSpecies;
        reaction::LogicalReactionSet m_activeReactions;

        std::vector<size_t> m_speciesIndexMap;
        std::vector<size_t> m_reactionIndexMap;

        bool m_isStale = true;
 
    private:
        struct ReactionFlow {
            const reaction::LogicalReaction* reactionPtr;
            double flowRate;
        };
    private:
        [[nodiscard]] std::vector<size_t> constructSpeciesIndexMap() const;

        [[nodiscard]] std::vector<size_t> constructReactionIndexMap() const;

        [[nodiscard]] std::vector<double> mapCulledToFull(const std::vector<double>& culled) const;

        [[nodiscard]] std::vector<double> mapFullToCulled(const std::vector<double>& full) const;

        [[nodiscard]] size_t mapCulledToFullSpeciesIndex(size_t culledSpeciesIndex) const;

        [[nodiscard]] size_t mapCulledToFullReactionIndex(size_t culledReactionIndex) const;

        void validateState() const;

        std::vector<ReactionFlow> calculateAllReactionFlows(
            const NetIn& netIn,
            std::vector<double>& out_Y_Full
        ) const;
        [[nodiscard]] std::unordered_set<fourdst::atomic::Species> findReachableSpecies(
            const NetIn& netIn
        ) const;
        [[nodiscard]] std::vector<const reaction::LogicalReaction*> cullReactionsByFlow(
            const std::vector<ReactionFlow>& allFlows,
            const std::unordered_set<fourdst::atomic::Species>& reachableSpecies,
            const std::vector<double>& Y_full,
            double maxFlow
        ) const;
 
        typedef std::pair<std::unordered_set<const reaction::LogicalReaction*>, std::unordered_set<fourdst::atomic::Species>> RescueSet;
        [[nodiscard]] RescueSet rescueEdgeSpeciesDestructionChannel(
            const std::vector<double>& Y_full,
            const double T9,
            const double rho,
            const std::vector<fourdst::atomic::Species>& activeSpecies,
            const reaction::LogicalReactionSet& activeReactions
        ) const;
        void finalizeActiveSet(
            const std::vector<const reaction::LogicalReaction*>& finalReactions
        );
    };
}