quantifier.cpp

Go to the documentation of this file.

#include "faulttree/quantifier.hpp"
 
#include "bncore/graph/graph.hpp"
#include "bncore/inference/compiler.hpp"
#include "bncore/inference/engine.hpp"
 
#include <sstream>
#include <vector>
 
namespace faulttree {
 
namespace {
 
// State indices for every variable in the compiled network. Index 0 is the
// failure state, matching the crew BBN's convention.
constexpr std::size_t kFailed = 0;
constexpr std::size_t kWorking = 1;
 
bool gateFails(const Node &node, std::size_t failed_inputs, std::size_t n) {
  switch (node.type) {
    case GateType::And: return failed_inputs == n;
    case GateType::Or: return failed_inputs > 0;
    case GateType::KofN:
      return failed_inputs >= static_cast<std::size_t>(node.k);
    case GateType::Not: return failed_inputs == 0;
    case GateType::Basic: break;
  }
  throw FaultTreeError("gate '" + node.name + "' has no deterministic table");
}
 
} // namespace
 
Quantifier::Quantifier(const FaultTree &tree) : tree_(tree) {
  tree_.validate();
  build();
}
 
Quantifier::~Quantifier() = default;
 
void Quantifier::build() {
  graph_ = std::make_unique<bncore::Graph>();
  var_of_node_.clear();
 
  // One binary variable per fault tree node.
  for (const auto &node : tree_.nodes()) {
    const auto id = graph_->add_variable(node.name, {"Failed", "Working"});
    var_of_node_.emplace(node.name, static_cast<std::size_t>(id));
  }
 
  // Each input is a parent of the gate it feeds. Insertion order here fixes
  // the parent order BNCore uses to index conditional tables, so the tables
  // below must enumerate inputs in this same order.
  for (const auto &node : tree_.nodes()) {
    for (const auto &input : node.children) {
      graph_->add_edge(static_cast<bncore::NodeId>(var_of_node_.at(input)),
                       static_cast<bncore::NodeId>(var_of_node_.at(node.name)));
    }
  }
 
  for (const auto &node : tree_.nodes()) {
    const auto id = static_cast<bncore::NodeId>(var_of_node_.at(node.name));
 
    if (node.type == GateType::Basic) {
      graph_->set_cpt(id, {node.probability, 1.0 - node.probability});
      continue;
    }
 
    const std::size_t n = node.children.size();
    if (n > kMaxGateInputs) {
      std::ostringstream out;
      out << "gate '" << node.name << "' has " << n << " inputs, above the "
          << kMaxGateInputs << " this quantifier will tabulate";
      throw FaultTreeError(out.str());
    }
 
    // BNCore lays a conditional table out row-major over
    // [parent 0, ..., parent n-1, self], so `self` varies fastest and parent 0
    // is the most significant bit of the parent configuration.
    const std::size_t combinations = std::size_t{1} << n;
    std::vector<double> cpt(combinations * 2);
    for (std::size_t combo = 0; combo < combinations; ++combo) {
      std::size_t failed_inputs = 0;
      for (std::size_t i = 0; i < n; ++i) {
        const std::size_t state = (combo >> (n - 1 - i)) & 1u;
        if (state == kFailed) { ++failed_inputs; }
      }
      const bool fails = gateFails(node, failed_inputs, n);
      cpt[combo * 2 + kFailed] = fails ? 1.0 : 0.0;
      cpt[combo * 2 + kWorking] = fails ? 0.0 : 1.0;
    }
    graph_->set_cpt(id, cpt);
  }
 
  junction_tree_ = bncore::JunctionTreeCompiler::compile(*graph_);
  engine_ = std::make_unique<bncore::BatchExecutionEngine>(*junction_tree_, 1, 1);
}
 
double Quantifier::probabilityOf(const std::string &name) const {
  auto it = var_of_node_.find(name);
  if (it == var_of_node_.end()) {
    throw FaultTreeError("no fault tree node named '" + name + "'");
  }
 
  const std::size_t num_vars = tree_.size();
  std::vector<int> evidence(num_vars, -1); // no hard evidence
  double marginal[2] = {0.0, 0.0};
  engine_->evaluate(evidence.data(), 1, num_vars, marginal, it->second);
  return marginal[kFailed];
}
 
double Quantifier::topEventProbability() const {
  return probabilityOf(tree_.top());
}
 
void Quantifier::setBasicEventProbability(const std::string &name,
                                          double probability) {
  tree_.setBasicEventProbability(name, probability);
  // Only a number changed, so the junction tree structure still holds: update
  // the prior in place and drop the cached workspace so the next query rebuilds
  // its potentials from the new table.
  graph_->set_cpt(static_cast<bncore::NodeId>(var_of_node_.at(name)),
                  {probability, 1.0 - probability});
  engine_->invalidate_workspace_cache();
}
 
std::vector<std::string> Quantifier::basicEventNames() const {
  return tree_.basicEventNames();
}
 
bool Quantifier::isDefined(const std::string &name) const {
  return tree_.hasNode(name);
}
 
} // namespace faulttree