setter.hpp

Go to the documentation of this file.

#pragma once
#include "sigma/uncertain/uncertain.hpp"
#include <vector>

 
namespace sigma::detail_ {

template<typename UncertainType>
class Setter {
public:
    using my_t = Setter<UncertainType>;

    using uncertain_t = UncertainType;

    using value_t = typename uncertain_t::value_t;

    using dep_sd_t = typename uncertain_t::dep_sd_t;

    using dep_sd_ptr = typename uncertain_t::dep_sd_ptr;

    using deps_map_t = typename uncertain_t::deps_map_t;

    using deps_vector_t = typename std::vector<dep_sd_ptr>;

    Setter(uncertain_t& u) : m_x_(u) {}

    void update_mean(value_t mean) {
        if(std::abs(m_x_.m_mean_ - mean) >= m_x_.threshold())
            m_x_.m_mean_ = mean;
    }

    void recalculate_sd() {
        m_x_.m_sd_ = 0.0;
        for(const auto& [dep, deriv] : m_x_.m_deps_) {
            m_x_.m_sd_ += std::pow(*dep.get() * deriv, 2.0);
        }
        m_x_.m_sd_ = std::sqrt(m_x_.m_sd_);
    }

    void update_derivatives(value_t dxda) {
        // If the derivative is effectively zero, we can just clear the
        // dependencies and set the standard deviation to zero.
        if(std::abs(dxda) < m_x_.threshold() || dxda == 0.0) {
            m_x_.m_deps_.clear();
            m_x_.m_sd_ = 0.0;
            return;
        }
        // If the derivative is one, we don't need to do anything since the
        // derivatives are unchanged.
        if(dxda == 1.0) return;
        // Update the derivatives in place, tracking any that are reduced below
        // our threshold and will need to be removed from the dependencies map.
        m_removed_deps_.clear();
        for(const auto& [dep, deriv] : m_x_.m_deps_) {
            m_x_.m_deps_[dep] *= dxda;
            if(std::abs(m_x_.m_deps_[dep]) < m_x_.threshold()) {
                m_removed_deps_.emplace_back(dep);
            }
        }
        // If no dependencies need to be removed, we can just update the
        // standard deviation by multiplying by the absolute value of the
        // partial derivative. If the partial derivate is negative one, we don't
        // need to update the standard deviation since it is unchanged. If any
        // dependencies were reduced below our threshold, we need to remove them
        // from the dependencies map and recalculate the standard deviation.
        if(m_removed_deps_.size() == 0) {
            if(dxda != -1.0) m_x_.m_sd_ *= std::abs(dxda);
        } else {
            for(const auto& dep : m_removed_deps_) m_x_.m_deps_.erase(dep);
            recalculate_sd();
        }
    }

    void update_derivatives(const deps_map_t& deps, value_t dxda) {
        // If the derivative is zero, we can skip the update since it won't
        // change anything.
        if(std::abs(dxda) < m_x_.threshold() || dxda == 0.0) { return; }
        // Update the map of contributions and their derivatives, keeping track
        // of any that are reduced to zero and will need to be removed. If the
        // uncertainty or contribution of a new dependency is below the
        // threshold, we can skip adding it to the map.
        m_removed_deps_.clear();
        size_t n_updated = 0;
        for(const auto& [dep, deriv] : deps) {
            if(m_x_.m_deps_.count(dep)) {
                m_x_.m_deps_[dep] += dxda * deriv;
                if(std::abs(m_x_.m_deps_[dep]) < m_x_.threshold() ||
                   m_x_.m_deps_[dep] == 0.0) {
                    m_removed_deps_.emplace_back(dep);
                }
                n_updated++;
            } else {
                if(std::abs(dxda * deriv) < m_x_.threshold()) continue;
                if(std::abs(*dep.get()) < m_x_.threshold()) continue;
                m_x_.m_deps_.emplace(dep, dxda * deriv);
            }
        }
        // If more than half of the existing dependencies are being updated,
        // it's more efficient to update the derivatives and then recalculate
        // the standard deviation from scratch. Otherwise, we can update the
        // standard deviation in place by removing the contribution of the
        // zeroed and updated dependencies and then adding the new contribution
        // of the updated dependencies.
        if(n_updated > (m_x_.m_deps_.size() / 2)) {
            for(const auto& dep : m_removed_deps_) m_x_.m_deps_.erase(dep);
            recalculate_sd();
        } else {
            // Return to the variance so we can update the sum.
            m_x_.m_sd_ = std::pow(m_x_.m_sd_, 2.0);
            // Step through the updated dependencies and update the variance.
            // accordingly.
            for(const auto& [dep, deriv] : deps) {
                // Skip any dependencies that aren't in the map at this point.
                if(!m_x_.m_deps_.count(dep)) continue;
                // Calculate the previous value of the derivative for this
                // dependency.
                auto old_deriv = m_x_.m_deps_[dep] - dxda * deriv;
                // As long as the dependency hasn't been reduce to zero, we need
                // to add its contribution to the variance.
                if(m_x_.m_deps_[dep] != 0.0) {
                    m_x_.m_sd_ += std::pow(*dep.get() * m_x_.m_deps_[dep], 2.0);
                }
                // As long as the dependency isn't new, we need to remove its
                // previous contribution to the variance.
                if(old_deriv != 0.0) {
                    m_x_.m_sd_ -= std::pow(*dep.get() * old_deriv, 2.0);
                }
            }
            m_x_.m_sd_ = std::sqrt(m_x_.m_sd_);
            for(const auto& dep : m_removed_deps_) m_x_.m_deps_.erase(dep);
        }
    }
 
private:
    uncertain_t& m_x_;

    deps_vector_t m_removed_deps_{};
};
 
} // namespace sigma::detail_