text:
/*
 * Copyright (c) Facebook, Inc. and its affiliates.
 *
 * This source code is licensed under the MIT license found in the
 * LICENSE file in the root directory of this source tree.
 */

#include "Reachability.h"

#include <boost/bimap/bimap.hpp>
#include <boost/bimap/unordered_set_of.hpp>
#include <boost/range/adaptor/map.hpp>

#include "BinarySerialization.h"
#include "DexUtil.h"
#include "Pass.h"
#include "ReachableClasses.h"
#include "Resolver.h"
#include "Timer.h"
#include "Walkers.h"

using namespace reachability;

namespace bs = binary_serialization;
namespace mog = method_override_graph;

namespace {

static ReachableObject SEED_SINGLETON{};

DexMethod* resolve(const DexMethodRef* method, const DexClass* cls) {
  if (!cls) return nullptr;
  for (auto const& m : cls->get_vmethods()) {
    if (method::signatures_match(method, m)) {
      return m;
    }
  }
  for (auto const& m : cls->get_dmethods()) {
    if (method::signatures_match(method, m)) {
      return m;
    }
  }
  {
    auto const& superclass = type_class(cls->get_super_class());
    auto const resolved = resolve(method, superclass);
    if (resolved) {
      return resolved;
    }
  }
  for (auto const& interface : cls->get_interfaces()->get_type_list()) {
    auto const resolved = resolve(method, type_class(interface));
    if (resolved) {
      return resolved;
    }
  }
  return nullptr;
}

} // namespace

namespace reachability {

std::ostream& operator<<(std::ostream& os, const ReachableObject& obj) {
  switch (obj.type) {
  case ReachableObjectType::ANNO:
    return os << show_deobfuscated(obj.anno->type());
  case ReachableObjectType::CLASS:
    return os << show_deobfuscated(obj.cls);
  case ReachableObjectType::FIELD:
    return os << show_deobfuscated(obj.field);
  case ReachableObjectType::METHOD:
    return os << show_deobfuscated(obj.method);
  case ReachableObjectType::SEED: {
    if (obj.keep_reason) {
      return os << *obj.keep_reason;
    }
    return os << "<SEED>";
  }
  }
}

bool RootSetMarker::is_canary(const DexClass* cls) {
  return strstr(cls->get_name()->c_str(), "Canary");
}

bool RootSetMarker::should_mark_cls(const DexClass* cls) {
  return root(cls) || is_canary(cls);
}

void RootSetMarker::mark_all_as_seed(const Scope& scope) {
  walk::parallel::classes(scope, [&](const DexClass* cls) {
    TRACE(REACH, 3, "Visiting seed: %s", SHOW(cls));
    push_seed(cls);

    for (auto const& f : cls->get_ifields()) {
      TRACE(REACH, 3, "Visiting seed: %s", SHOW(f));
      push_seed(f);
    }
    for (auto const& f : cls->get_sfields()) {
      TRACE(REACH, 3, "Visiting seed: %s", SHOW(f));
      push_seed(f);
    }
    for (auto const& m : cls->get_dmethods()) {
      TRACE(REACH, 3, "Visiting seed: %s", SHOW(m));
      push_seed(m);
    }
    for (auto const& m : cls->get_vmethods()) {
      TRACE(REACH, 3, "Visiting seed: %s (root)", SHOW(m));
      push_seed(m);
    }
  });
}

/*
 * Initializes the root set by marking and pushing nodes onto the work queue.
 * Also conditionally marks class member seeds.
 */
void RootSetMarker::mark(const Scope& scope) {
  walk::parallel::classes(scope, [&](const DexClass* cls) {
    if (should_mark_cls(cls)) {
      TRACE(REACH, 3, "Visiting seed: %s", SHOW(cls));
      push_seed(cls);
    }
    for (auto const& f : cls->get_ifields()) {
      if (root(f) || is_volatile(f)) {
        TRACE(REACH, 3, "Visiting seed: %s", SHOW(f));
        push_seed(f);
      }
    }
    for (auto const& f : cls->get_sfields()) {
      if (root(f)) {
        TRACE(REACH, 3, "Visiting seed: %s", SHOW(f));
        push_seed(f);
      }
    }
    for (auto const& m : cls->get_dmethods()) {
      if (root(m)) {
        TRACE(REACH, 3, "Visiting seed: %s", SHOW(m));
        push_seed(m);
      }
    }
    for (auto const& m : cls->get_vmethods()) {
      if (root(m)) {
        TRACE(REACH, 3, "Visiting seed: %s (root)", SHOW(m));
        push_seed(m);
      }
    }
  });

  mark_external_method_overriders();
}

void RootSetMarker::push_seed(const DexClass* cls) {
  if (!cls) return;
  record_is_seed(cls);
  m_reachable_objects->mark(cls);
  m_root_set->emplace(cls);
}

void RootSetMarker::push_seed(const DexField* field) {
  if (!field) return;
  m_cond_marked->fields.insert(field);
}

void RootSetMarker::push_seed(const DexMethod* method) {
  if (!method) return;
  m_cond_marked->methods.insert(method);
}

template <class Seed>
void RootSetMarker::record_is_seed(Seed* seed) {
  if (m_record_reachability) {
    redex_assert(seed != nullptr);
    m_reachable_objects->record_is_seed(seed);
  }
}

/*
 * Mark as seeds all methods that override or implement an external method.
 */
void RootSetMarker::mark_external_method_overriders() {
  std::unordered_set<const DexMethod*> visited;
  for (auto& pair : m_method_override_graph.nodes()) {
    auto method = pair.first;
    if (!method->is_external() || visited.count(method)) {
      continue;
    }
    const auto& overriding_methods =
        mog::get_overriding_methods(m_method_override_graph, method);
    for (auto* overriding : overriding_methods) {
      // Avoid re-visiting methods found in overriding sets since we would
      // already have conditionally marked all their children.
      visited.emplace(overriding);
      if (!overriding->is_external()) {
        TRACE(REACH, 3, "Visiting seed: %s (implements %s)", SHOW(overriding),
              SHOW(method));
        push_seed(overriding);
      }
    }
  }
}

/*
 * Marks :obj and pushes its immediately reachable neighbors onto the local
 * task queue of the current worker.
 */
void TransitiveClosureMarker::visit(const ReachableObject& obj) {
  switch (obj.type) {
  case ReachableObjectType::CLASS:
    visit_cls(obj.cls);
    break;
  case ReachableObjectType::FIELD:
    visit(obj.field);
    break;
  case ReachableObjectType::METHOD:
    visit(obj.method);
    break;
  case ReachableObjectType::ANNO:
  case ReachableObjectType::SEED:
    always_assert_log(false, "Unexpected ReachableObject type");
    break;
  }
}

template <class Parent, class InputIt>
void TransitiveClosureMarker::push(const Parent* parent,
                                   InputIt begin,
                                   InputIt end) {
  for (auto it = begin; it != end; ++it) {
    push(parent, *it);
  }
}

template <class Parent>
void TransitiveClosureMarker::push(const Parent* parent, const DexType* type) {
  type = type::get_element_type_if_array(type);
  push(parent, type_class(type));
}

template <class Parent>
void TransitiveClosureMarker::push(const Parent* parent, const DexClass* cls) {
  if (!cls) {
    return;
  }
  record_reachability(parent, cls);
  if (m_reachable_objects->marked(cls)) {
    return;
  }
  m_reachable_objects->mark(cls);
  m_worker_state->push_task(ReachableObject(cls));
}

template <class Parent>
void TransitiveClosureMarker::push(const Parent* parent,
                                   const DexFieldRef* field) {
  if (!field) {
    return;
  }
  record_reachability(parent, field);
  if (m_reachable_objects->marked(field)) {
    return;
  }
  auto f = field->as_def();
  if (f) {
    gather_and_push(f);
  }
  m_reachable_objects->mark(field);
  m_worker_state->push_task(ReachableObject(field));
}

template <class Parent>
void TransitiveClosureMarker::push(const Parent* parent,
                                   const DexMethodRef* method) {
  if (!method) {
    return;
  }
  record_reachability(parent, method);
  if (m_reachable_objects->marked(method)) {
    return;
  }
  m_reachable_objects->mark(method);
  m_worker_state->push_task(ReachableObject(method));
}

void TransitiveClosureMarker::push_cond(const DexMethod* method) {
  if (!method || m_reachable_objects->marked(method)) return;
  TRACE(REACH, 4, "Conditionally marking method: %s", SHOW(method));
  auto clazz = type_class(method->get_class());
  m_cond_marked->methods.insert(method);
  // If :clazz has been marked, we cannot count on visit(DexClass*) to move
  // the conditionally-marked methods into the actually-marked ones -- we have
  // to do it ourselves. Note that we must do this check after adding :method
  // to m_cond_marked to avoid a race condition where we add to m_cond_marked
  // after visit(DexClass*) has finished moving its contents over to
  // m_reachable_objects.
  if (m_reachable_objects->marked(clazz)) {
    push(clazz, method);
  }
}

template <class T>
static References generic_gather(T t) {
  References refs;
  t->gather_strings(refs.strings);
  t->gather_types(refs.types);
  t->gather_fields(refs.fields);
  t->gather_methods(refs.methods);
  return refs;
}

References TransitiveClosureMarker::gather(const DexAnnotation* anno) const {
  return generic_gather(anno);
}

References TransitiveClosureMarker::gather(const DexMethod* method) const {
  return generic_gather(method);
}

References TransitiveClosureMarker::gather(const DexField* field) const {
  return generic_gather(field);
}

bool TransitiveClosureMarker::has_class_forname(DexMethod* meth) {
  auto code = meth->get_code();
  if (!code || !m_class_forname) {
    return false;
  }
  for (auto& mie : InstructionIterable(code)) {
    auto insn = mie.insn;
    if (insn->has_method() && insn->get_method() == m_class_forname) {
      return true;
    }
  }
  return false;
}

void TransitiveClosureMarker::gather_and_push(DexMethod* meth) {
  auto* type = meth->get_class();
  auto* cls = type_class(type);
  auto refs = gather(meth);
  bool check_strings = m_ignore_sets.keep_class_in_string;
  if (!check_strings && refs.strings.size() > 0 && has_class_forname(meth)) {
    check_strings = true;
  }
  if (m_ignore_sets.string_literals.count(type)) {
    ++m_stats->num_ignore_check_strings;
    check_strings = false;
  }
  if (cls && check_strings) {
    for (const auto& ignore_anno_type : m_ignore_sets.string_literal_annos) {
      if (has_anno(cls, ignore_anno_type)) {
        ++m_stats->num_ignore_check_strings;
        check_strings = false;
        break;
      }
    }
  }
  if (check_strings) {
    push_typelike_strings(meth, refs.strings);
  }
  push(meth, refs.types.begin(), refs.types.end());
  push(meth, refs.fields.begin(), refs.fields.end());
  push(meth, refs.methods.begin(), refs.methods.end());
}

template <typename T>
void TransitiveClosureMarker::gather_and_push(T t) {
  auto refs = gather(t);
  push_typelike_strings(t, refs.strings);
  push(t, refs.types.begin(), refs.types.end());
  push(t, refs.fields.begin(), refs.fields.end());
  push(t, refs.methods.begin(), refs.methods.end());
}

template <class Parent>
void TransitiveClosureMarker::push_typelike_strings(
    const Parent* parent, const std::vector<DexString*>& strings) {
  for (auto const& str : strings) {
    auto internal = java_names::external_to_internal(str->c_str());
    auto type = DexType::get_type(internal.c_str());
    if (!type) {
      continue;
    }
    push(parent, type);
  }
}

void TransitiveClosureMarker::visit_cls(const DexClass* cls) {
  TRACE(REACH, 4, "Visiting class: %s", SHOW(cls));
  for (auto& m : cls->get_dmethods()) {
    if (method::is_clinit(m)) {
      push(cls, m);
    } else if (method::is_init(m)) {
      // Push the parameterless constructor, in case it's constructed via
      // .class or Class.forName()
      if (m->get_proto()->get_args()->get_type_list().size() == 0) {
        push(cls, m);
      }
    }
  }
  push(cls, type_class(cls->get_super_class()));
  for (auto const& t : cls->get_interfaces()->get_type_list()) {
    push(cls, t);
  }
  const DexAnnotationSet* annoset = cls->get_anno_set();
  if (annoset) {
    for (auto const& anno : annoset->get_annotations()) {
      if (m_ignore_sets.system_annos.count(anno->type())) {
        TRACE(REACH,
              5,
              "Stop marking from %s by system anno: %s",
              SHOW(cls),
              SHOW(anno->type()));
        continue;
      }
      record_reachability(cls, anno);
      gather_and_push(anno);
    }
  }
  for (auto const& m : cls->get_ifields()) {
    if (m_cond_marked->fields.count(m)) {
      push(cls, m);
    }
  }
  for (auto const& m : cls->get_sfields()) {
    if (m_cond_marked->fields.count(m)) {
      push(cls, m);
    }
  }
  for (auto const& m : cls->get_dmethods()) {
    if (m_cond_marked->methods.count(m)) {
      push(cls, m);
    }
  }
  for (auto const& m : cls->get_vmethods()) {
    if (m_cond_marked->methods.count(m)) {
      push(cls, m);
    }
  }
}

void TransitiveClosureMarker::visit(const DexFieldRef* field) {
  TRACE(REACH, 4, "Visiting field: %s", SHOW(field));
  if (!field->is_concrete()) {
    auto const& realfield =
        resolve_field(field->get_class(), field->get_name(), field->get_type());
    push(field, realfield);
  }
  push(field, field->get_class());
  push(field, field->get_type());
}

void TransitiveClosureMarker::visit(const DexMethodRef* method) {
  TRACE(REACH, 4, "Visiting method: %s", SHOW(method));
  auto resolved_method = resolve(method, type_class(method->get_class()));
  if (resolved_method != nullptr) {
    TRACE(REACH, 5, "    Resolved to: %s", SHOW(resolved_method));
    push(method, resolved_method);
    gather_and_push(resolved_method);
  }
  push(method, method->get_class());
  push(method, method->get_proto()->get_rtype());
  for (auto const& t : method->get_proto()->get_args()->get_type_list()) {
    push(method, t);
  }
  auto m = method->as_def();
  if (!m) {
    return;
  }
  // If we're keeping an interface or virtual method, we have to keep its
  // implementations and overriding methods respectively.
  if (m->is_virtual() || !m->is_concrete()) {
    const auto& overriding_methods =
        mog::get_overriding_methods(m_method_override_graph, m);
    for (auto* overriding : overriding_methods) {
      push_cond(overriding);
    }
  }
}

template <class Parent, class Object>
void TransitiveClosureMarker::record_reachability(Parent* parent,
                                                  Object* object) {
  if (m_record_reachability) {
    redex_assert(parent != nullptr && object != nullptr);
    m_reachable_objects->record_reachability(parent, object);
  }
}

std::unique_ptr<ReachableObjects> compute_reachable_objects(
    const DexStoresVector& stores,
    const IgnoreSets& ignore_sets,
    int* num_ignore_check_strings,
    bool record_reachability,
    bool should_mark_all_as_seed,
    std::unique_ptr<const mog::Graph>* out_method_override_graph) {
  Timer t("Marking");
  auto scope = build_class_scope(stores);
  auto reachable_objects = std::make_unique<ReachableObjects>();
  ConditionallyMarked cond_marked;
  auto method_override_graph = mog::build_graph(scope);

  ConcurrentSet<ReachableObject, ReachableObjectHash> root_set;
  RootSetMarker root_set_marker(*method_override_graph,
                                record_reachability,
                                &cond_marked,
                                reachable_objects.get(),
                                &root_set);
  if (should_mark_all_as_seed) {
    root_set_marker.mark_all_as_seed(scope);
  } else {
    root_set_marker.mark(scope);
  }

  size_t num_threads = redex_parallel::default_num_threads();
  auto stats_arr = std::make_unique<Stats[]>(num_threads);
  MarkWorkQueue work_queue(
      [&](MarkWorkerState* worker_state, const ReachableObject& obj) {
        TransitiveClosureMarker transitive_closure_marker(
            ignore_sets, *method_override_graph, record_reachability,
            &cond_marked, reachable_objects.get(), worker_state,
            &stats_arr[worker_state->worker_id()]);
        transitive_closure_marker.visit(obj);
        return nullptr;
      },
      num_threads);
  for (const auto& obj : root_set) {
    work_queue.add_item(obj);
  }
  work_queue.run_all();

  if (num_ignore_check_strings != nullptr) {
    for (size_t i = 0; i < num_threads; ++i) {
      *num_ignore_check_strings += stats_arr[i].num_ignore_check_strings;
    }
  }

  if (out_method_override_graph) {
    *out_method_override_graph = std::move(method_override_graph);
  }

  return reachable_objects;
}

void ReachableObjects::record_reachability(const DexMethodRef* member,
                                           const DexClass* cls) {
  // Each class member trivially retains its containing class; let's filter out
  // this uninteresting information from our diagnostics.
  if (member->get_class() == cls->get_type()) {
    return;
  }
  m_retainers_of.update(ReachableObject(cls),
                        [&](const ReachableObject&, ReachableObjectSet& set,
                            bool /* exists */) { set.emplace(member); });
}

void ReachableObjects::record_reachability(const DexFieldRef* member,
                                           const DexClass* cls) {
  if (member->get_class() == cls->get_type()) {
    return;
  }
  m_retainers_of.update(ReachableObject(cls),
                        [&](const ReachableObject&, ReachableObjectSet& set,
                            bool /* exists */) { set.emplace(member); });
}

template <class Object>
void ReachableObjects::record_reachability(Object* parent, Object* object) {
  if (parent == object) {
    return;
  }
  m_retainers_of.update(ReachableObject(object),
                        [&](const ReachableObject&, ReachableObjectSet& set,
                            bool /* exists */) { set.emplace(parent); });
}

template <class Parent, class Object>
void ReachableObjects::record_reachability(Parent* parent, Object* object) {
  m_retainers_of.update(ReachableObject(object),
                        [&](const ReachableObject&, ReachableObjectSet& set,
                            bool /* exists */) { set.emplace(parent); });
}

template <class Seed>
void ReachableObjects::record_is_seed(Seed* seed) {
  redex_assert(seed != nullptr);
  const auto& keep_reasons = seed->rstate.keep_reasons();
  m_retainers_of.update(
      ReachableObject(seed),
      [&](const ReachableObject&, ReachableObjectSet& set, bool /* exists */) {
        for (const auto& reason : keep_reasons) {
          // -keepnames rules are irrelevant when analyzing reachability
          if (reason->type == keep_reason::KEEP_RULE &&
              reason->keep_rule->allowshrinking) {
            continue;
          }
          set.emplace(reason);
        }
      });
}

/*
 * Remove unmarked classes / methods / fields. and add all swept objects to
 * :removed_symbols.
 */
template <class Container, class FnPtr>
static void sweep_if_unmarked(const ReachableObjects& reachables,
                              FnPtr erase_hook,
                              Container* c,
                              ConcurrentSet<std::string>* removed_symbols) {
  auto p = [&](const auto& m) {
    if (reachables.marked_unsafe(m) == 0) {
      TRACE(RMU, 2, "Removing %s", SHOW(m));
      if (erase_hook) {
        erase_hook(m);
      }
      return false;
    }
    return true;
  };
  const auto it = std::partition(c->begin(), c->end(), p);
  if (removed_symbols) {
    for (auto i = it; i != c->end(); i++) {
      removed_symbols->insert(show_deobfuscated(*i));
    }
  }
  c->erase(it, c->end());
}

void sweep(DexStoresVector& stores,
           const ReachableObjects& reachables,
           ConcurrentSet<std::string>* removed_symbols) {
  Timer t("Sweep");
  for (auto& dex : DexStoreClassesIterator(stores)) {
    sweep_if_unmarked(reachables, (void (*)(DexClass*))(nullptr), &dex,
                      removed_symbols);
    walk::parallel::classes(dex, [&](DexClass* cls) {
      sweep_if_unmarked(reachables, DexField::erase_field, &cls->get_ifields(),
                        removed_symbols);
      sweep_if_unmarked(reachables, DexField::erase_field, &cls->get_sfields(),
                        removed_symbols);
      sweep_if_unmarked(reachables, DexMethod::erase_method,
                        &cls->get_dmethods(), removed_symbols);
      sweep_if_unmarked(reachables, DexMethod::erase_method,
                        &cls->get_vmethods(), removed_symbols);
    });
  }
}

ObjectCounts count_objects(const DexStoresVector& stores) {
  ObjectCounts counts;
  for (auto const& dex : DexStoreClassesIterator(stores)) {
    counts.num_classes += dex.size();
    for (auto const& cls : dex) {
      counts.num_fields += cls->get_ifields().size();
      counts.num_fields += cls->get_sfields().size();
      counts.num_methods += cls->get_dmethods().size();
      counts.num_methods += cls->get_vmethods().size();
    }
  }
  return counts;
}

// Graph serialization helpers
namespace {

void write_reachable_object(std::ostream& os, const ReachableObject& obj) {
  bs::write<uint8_t>(os, static_cast<uint8_t>(obj.type));
  std::ostringstream ss;
  ss << obj;
  const auto& s = ss.str();
  bs::write<uint32_t>(os, s.size());
  os << s;
}

} // namespace

void dump_graph(std::ostream& os, const ReachableObjectGraph& retainers_of) {
  auto compare = [](const ReachableObject& lhs, const ReachableObject& rhs) {
    if (lhs.type != rhs.type) {
      return lhs.type < rhs.type;
    }
    switch (lhs.type) {
    case ReachableObjectType::CLASS:
      return compare_dexclasses(lhs.cls, rhs.cls);
    case ReachableObjectType::FIELD:
      return compare_dexfields(lhs.field, rhs.field);
    case ReachableObjectType::METHOD:
      return compare_dexmethods(lhs.method, rhs.method);

    case ReachableObjectType::ANNO:
    case ReachableObjectType::SEED: {
      // Pretty terrible, optimize.
      std::ostringstream oss1;
      oss1 << lhs;
      std::ostringstream oss2;
      oss2 << rhs;
      return oss1.str() < oss2.str();
    }
    }
    __builtin_unreachable();
  };

  bs::write_header(os, /* version */ 1);
  bs::GraphWriter<ReachableObject, ReachableObjectHash> gw(
      write_reachable_object,
      [&](const ReachableObject& obj) -> std::vector<ReachableObject> {
        if (!retainers_of.count(obj)) {
          return {};
        }
        const auto& preds = retainers_of.at(obj);
        std::vector<ReachableObject> preds_vec(preds.begin(), preds.end());
        // Gotta sort the reachables or the output is nondeterministic.
        std::sort(preds_vec.begin(), preds_vec.end(), compare);
        return preds_vec;
      });

  // Gotta sort the keys or the output is nondeterministic.
  auto key_adaptor = boost::adaptors::keys(retainers_of);
  std::vector<ReachableObject> keys(key_adaptor.begin(), key_adaptor.end());
  std::sort(keys.begin(), keys.end(), compare);
  gw.write(os, keys);
}

template void TransitiveClosureMarker::push<DexClass>(const DexClass* parent,
                                                      const DexType* type);
} // namespace reachability

max_stars_repo_name: rybalkinsd/redex

max_issues_repo_name: rybalkinsd/redex

max_forks_repo_name: rybalkinsd/redex


text:
/**
 * @author      : adamlamson (adamlamson@LamsonMacbookPro)
 * @file        : two_step_prob
 * @created     : Wednesday Jan 29, 2020 14:05:35 MST
 */

#ifndef TWO_STEP_PROB_HPP

#define TWO_STEP_PROB_HPP
#include <boost/math/tools/roots.hpp>
#include <cassert>
#include <cmath>

// Solve the for the time that maximizes the probability of two events occuring
template <class T>
class MaxTimeOfTwoStepProbFunctor {
  public:
    // Constructor just stores values find root of function.
    MaxTimeOfTwoStepProbFunctor(T const &rate1, T const &rate2,
                                T const &total_time)
        : k_ab(rate1), k_bc(rate2), tot_time(total_time) {}

    T operator()(T const &t) {
        T fx = k_ab * (exp(k_bc * (tot_time - t)) - 1.) -
               k_bc * (exp(k_ab * t) - 1.);
        return fx;
    }

  private:
    T k_ab;
    T k_bc;
    T tot_time;
};

template <class T>
T max_time_of_two_step_prob(T const &rate1, T const &rate2,
                            T const &total_time) {

    T guess = .5 * total_time;
    T factor = 2; // How big steps to take when searching.

    const boost::uintmax_t maxit = 20; // Limit to maximum iterations.
    boost::uintmax_t it =
        maxit; // Initally our chosen max iterations, but updated with actual.
    bool is_rising =
        false; // So if guess is too low, then try increasing guess.
    int digits = std::numeric_limits<T>::digits; // Maximum possible binary
                                                 // digits accuracy for type T.
    // Some fraction of digits is used to control how accurate o try to make the
    // result.
    int get_digits =
        digits - 3; // We have to have a non-zero interval at each step, so
                    // maximum accuracy is digits - 1.  But we also have to
                    // allow for inaccuracy in f(x), otherwise the last few
                    // iterations just thrash around.
    boost::math::tools::eps_tolerance<T> rel_tol(
        get_digits); // Set the tolerance.
    std::pair<T, T> t = boost::math::tools::bracket_and_solve_root(
        MaxTimeOfTwoStepProbFunctor<T>(rate1, rate2, total_time), guess, factor,
        is_rising, rel_tol, it);
    return t.first + (t.second - t.first) * .5;
};

template <class T>
T two_step_prob(T const &rate1, T const &rate2, T const &total_time,
                T const &t) {
    assert(total_time >= t);
    T prob = (1. - exp(-1. * rate1 * t)) *
             (1. - exp(-1. * rate2 * (total_time - t)));
    return prob;
};

template <class T>
T two_step_max_prob(T const &rate1, T const &rate2, T const &total_time) {
    T t_max = max_time_of_two_step_prob(rate1, rate2, total_time);
    T prob = two_step_prob(rate1, rate2, total_time, t_max);
    return prob;
};

#endif /* end of include guard TWO_STEP_PROB_HPP */

ext: hpp


text:
#ifndef LBV_EIGEN_ROS_CONVERSIONS_HPP
#define LBV_EIGEN_ROS_CONVERSIONS_HPP

#include <geometry_msgs/Pose.h>
#include <geometry_msgs/Transform.h>

#include <Eigen/Dense>
#include <Eigen/StdVector>

////////////////////////////////////////////////////////////////////////////////

template <typename Output>
Output ConvertTo(geometry_msgs::Point const& geo);

template <typename Output>
Output ConvertTo(geometry_msgs::Vector3 const& eig);

template <typename Output>
Output ConvertTo(Eigen::Vector3d const& eig);

template <typename Output>
Output ConvertTo(Eigen::Vector3f const& eig);

template <>
inline Eigen::Vector3f ConvertTo<Eigen::Vector3f>(geometry_msgs::Point const& geo) {
  return Eigen::Vector3f(geo.x, geo.y, geo.z);
}

template <>
inline Eigen::Vector3f ConvertTo<Eigen::Vector3f>(geometry_msgs::Vector3 const& eig) {
  return Eigen::Vector3f(eig.x, eig.y, eig.z);
}

template <>
inline geometry_msgs::Point ConvertTo<geometry_msgs::Point>(Eigen::Vector3f const& eig) {
  geometry_msgs::Point geo;
  geo.x = eig.x();
  geo.y = eig.y();
  geo.z = eig.z();
  return geo;
}

template <>
inline geometry_msgs::Vector3 ConvertTo<geometry_msgs::Vector3>(Eigen::Vector3f const& eig) {
  geometry_msgs::Vector3 geo;
  geo.x = eig.x();
  geo.y = eig.y();
  geo.z = eig.z();
  return geo;
}

template <>
inline Eigen::Vector3d ConvertTo<Eigen::Vector3d>(geometry_msgs::Point const& geo) {
  return Eigen::Vector3d(geo.x, geo.y, geo.z);
}

template <>
inline Eigen::Vector3d ConvertTo<Eigen::Vector3d>(geometry_msgs::Vector3 const& eig) {
  return Eigen::Vector3d(eig.x, eig.y, eig.z);
}

template <>
inline geometry_msgs::Point ConvertTo<geometry_msgs::Point>(Eigen::Vector3d const& eig) {
  geometry_msgs::Point geo;
  geo.x = eig.x();
  geo.y = eig.y();
  geo.z = eig.z();
  return geo;
}

template <>
inline geometry_msgs::Vector3 ConvertTo<geometry_msgs::Vector3>(Eigen::Vector3d const& eig) {
  geometry_msgs::Vector3 geo;
  geo.x = eig.x();
  geo.y = eig.y();
  geo.z = eig.z();
  return geo;
}

////////////////////////////////////////////////////////////////////////////////

template <typename Output>
Output ConvertTo(geometry_msgs::Quaternion const& quat);

template <typename Output>
Output ConvertTo(Eigen::Quaterniond const& quat);

template <typename Output>
Output ConvertTo(Eigen::Matrix3d const& rot);

template <typename Output>
Output ConvertTo(Eigen::Quaternionf const& quat);

template <typename Output>
Output ConvertTo(Eigen::Matrix3f const& rot);

template <>
inline Eigen::Quaternionf ConvertTo<Eigen::Quaternionf>(geometry_msgs::Quaternion const& quat) {
  return Eigen::Quaternionf(quat.w, quat.x, quat.y, quat.z);
}

template <>
inline geometry_msgs::Quaternion ConvertTo<geometry_msgs::Quaternion>(Eigen::Quaternionf const& quat) {
  geometry_msgs::Quaternion geo_quat;
  geo_quat.x = quat.x();
  geo_quat.y = quat.y();
  geo_quat.z = quat.z();
  geo_quat.w = quat.w();
  return geo_quat;
}

template <>
inline geometry_msgs::Quaternion ConvertTo<geometry_msgs::Quaternion>(Eigen::Matrix3f const& rot) {
  return ConvertTo<geometry_msgs::Quaternion>(Eigen::Quaternionf(rot));
}

template <>
inline Eigen::Quaterniond ConvertTo<Eigen::Quaterniond>(geometry_msgs::Quaternion const& quat) {
  return Eigen::Quaterniond(quat.w, quat.x, quat.y, quat.z);
}

template <>
inline geometry_msgs::Quaternion ConvertTo<geometry_msgs::Quaternion>(Eigen::Quaterniond const& quat) {
  geometry_msgs::Quaternion geo_quat;
  geo_quat.x = quat.x();
  geo_quat.y = quat.y();
  geo_quat.z = quat.z();
  geo_quat.w = quat.w();
  return geo_quat;
}

template <>
inline geometry_msgs::Quaternion ConvertTo<geometry_msgs::Quaternion>(Eigen::Matrix3d const& rot) {
  return ConvertTo<geometry_msgs::Quaternion>(Eigen::Quaterniond(rot));
}

////////////////////////////////////////////////////////////////////////////////

template <typename Output>
Output ConvertTo(geometry_msgs::Pose const& pose);

template <typename Output>
Output ConvertTo(Eigen::Affine3d const& transform);

template <typename Output>
Output ConvertTo(Eigen::Isometry3d const& transform);

template <typename Output>
Output ConvertTo(Eigen::Translation3d const& translation);

template <>
inline Eigen::Isometry3d ConvertTo<Eigen::Isometry3d>(geometry_msgs::Pose const& pose) {
  Eigen::Translation3d const trans(pose.position.x, pose.position.y, pose.position.z);
  Eigen::Quaterniond const quat(pose.orientation.w, pose.orientation.x, pose.orientation.y, pose.orientation.z);
  return trans * quat;
}

template <>
inline Eigen::Affine3d ConvertTo<Eigen::Affine3d>(geometry_msgs::Pose const& pose) {
  Eigen::Translation3d const trans(pose.position.x, pose.position.y, pose.position.z);
  Eigen::Quaterniond const quat(pose.orientation.w, pose.orientation.x, pose.orientation.y, pose.orientation.z);
  return trans * quat;
}

template <>
inline geometry_msgs::Transform ConvertTo<geometry_msgs::Transform>(Eigen::Affine3d const& transform) {
  geometry_msgs::Transform geo_transform;
  geo_transform.translation.x = transform.translation().x();
  geo_transform.translation.y = transform.translation().y();
  geo_transform.translation.z = transform.translation().z();
  const Eigen::Quaterniond quat(transform.rotation());
  geo_transform.rotation.x = quat.x();
  geo_transform.rotation.y = quat.y();
  geo_transform.rotation.z = quat.z();
  geo_transform.rotation.w = quat.w();
  return geo_transform;
}

template <>
inline geometry_msgs::Transform ConvertTo<geometry_msgs::Transform>(Eigen::Isometry3d const& transform) {
  geometry_msgs::Transform geo_transform;
  geo_transform.translation.x = transform.translation().x();
  geo_transform.translation.y = transform.translation().y();
  geo_transform.translation.z = transform.translation().z();
  const Eigen::Quaterniond quat(transform.rotation());
  geo_transform.rotation.x = quat.x();
  geo_transform.rotation.y = quat.y();
  geo_transform.rotation.z = quat.z();
  geo_transform.rotation.w = quat.w();
  return geo_transform;
}

template <>
inline geometry_msgs::Transform ConvertTo<geometry_msgs::Transform>(Eigen::Translation3d const& translation) {
  geometry_msgs::Transform geo_transform;
  geo_transform.translation.x = translation.translation().x();
  geo_transform.translation.y = translation.translation().y();
  geo_transform.translation.z = translation.translation().z();
  geo_transform.rotation.x = 0;
  geo_transform.rotation.y = 0;
  geo_transform.rotation.z = 0;
  geo_transform.rotation.w = 1;
  return geo_transform;
}

////////////////////////////////////////////////////////////////////////////////

template <typename Output>
Output ConvertTo(geometry_msgs::Transform const& transform);

template <>
inline Eigen::Isometry3d ConvertTo<Eigen::Isometry3d>(geometry_msgs::Transform const& transform) {
  Eigen::Translation3d const trans(transform.translation.x, transform.translation.y, transform.translation.z);
  Eigen::Quaterniond const quat(transform.rotation.w, transform.rotation.x, transform.rotation.y, transform.rotation.z);
  return trans * quat;
}

template <>
inline Eigen::Affine3d ConvertTo<Eigen::Affine3d>(geometry_msgs::Transform const& transform) {
  Eigen::Translation3d const trans(transform.translation.x, transform.translation.y, transform.translation.z);
  Eigen::Quaterniond const quat(transform.rotation.w, transform.rotation.x, transform.rotation.y, transform.rotation.z);
  return trans * quat;
}

////////////////////////////////////////////////////////////////////////////////

template <typename Input, typename Output>
std::vector<Output, Eigen::aligned_allocator<Output>> ConvertTo(std::vector<Input> const& geom) {
  std::vector<Output, Eigen::aligned_allocator<Output>> eig(geom.size());
  for (size_t idx = 0; idx < geom.size(); ++idx) {
    eig[idx] = ConvertTo<Output>(geom[idx]);
  }
  return eig;
}

#endif  // LBV_EIGEN_ROS_CONVERSIONS_HPP


text:
/**
 * @file
 * @author University of Warwick
 * @version 1.0
 *
 * @section LICENSE
 *
 * @section DESCRIPTION
 *
 * Tests for the concrete methods of the MeshConfig class.
 */

#define BOOST_TEST_MODULE MeshConfig

// Boost setup
#include <boost/test/unit_test.hpp>
#include <boost/test/output_test_stream.hpp>
#include <stdexcept>

#include "Error.h"
#include "MeshConfig.h"
#include "PartitionerConfig.h"
#include "PartitionerNaiveConfig.h"
#include "MeshSourceStructGenConfig.h"
#include "UnstructuredMeshInterface.h"
#include "MeshSourceFileConfig.h"
#include "CupCfdAoSMesh.h"

// Setup
BOOST_AUTO_TEST_CASE(setup)
{
    int argc = boost::unit_test::framework::master_test_suite().argc;
    char ** argv = boost::unit_test::framework::master_test_suite().argv;

    MPI_Init(&argc, &argv);
}

// Namespaces
namespace utf = boost::unit_test;

using namespace cupcfd::geometry::mesh;

// === Constructor ===

// Test 1: Test Constructor with provided PartitionerConfigNaive is built correctly.
BOOST_AUTO_TEST_CASE(constructor_test1, * utf::tolerance(0.00001))
{
	// Setup PartitionerConfig

	// Constructor

	// Delete PartitionerConfig

	// Check Constructor Copy still exists
}

// Test 2: Test Constructor with provided PartitionerConfigParmetis is built correctly.
BOOST_AUTO_TEST_CASE(constructor_test2, * utf::tolerance(0.00001))
{
	// Setup PartitionerConfig

	// Constructor

	// Delete PartitionerConfig

	// Check Constructor Copy still exists
}

// Test 3: Test Constructor with Referenced Config copies values correctly
BOOST_AUTO_TEST_CASE(constructor_test3, * utf::tolerance(0.00001))
{

}

// === set/getPartitionerConfig
BOOST_AUTO_TEST_CASE(getPartitionerConfig_test1, * utf::tolerance(0.00001))
{

}

// === operator= ===
// Test 1: Test values are correctly copied from one config to another
BOOST_AUTO_TEST_CASE(operator_equals_test1, * utf::tolerance(0.00001))
{

}

// === buildUnstructuredMesh ===
// Test 1: Test correct build of a mesh from a HDF5 file format
BOOST_AUTO_TEST_CASE(buildUnstructuredMesh_test1, * utf::tolerance(0.00001))
{
	cupcfd::error::eCodes status;
    cupcfd::comm::Communicator comm(MPI_COMM_WORLD);

	// Setup the partitioner config
	// Use naive so the outcome is predictable
    cupcfd::partitioner::PartitionerNaiveConfig<int,int> naiveConfig;
    cupcfd::partitioner::PartitionerConfig<int,int> partConfig(naiveConfig);

	// Setup the source config
	MeshSourceFileConfig<int, double> meshSourceConfig(MESH_FILE_FORMAT_HDF5, "../tests/geometry/mesh/data/MeshHDF5.hdf5");

	// Setup the config to use for building
	MeshConfig<int, double, int> meshConfig(partConfig, meshSourceConfig);

	CupCfdAoSMesh<int,double,int> * mesh;
	status = meshConfig.buildUnstructuredMesh(&mesh, comm);
	BOOST_CHECK_EQUAL(status, cupcfd::error::E_SUCCESS);

	// Naive means we should have 1 cell per rank
	if(comm.rank == 0)
	{
		// The graphs will use labels, not internal indexes
		// The local indexes they map to should still be zero-indexed however
		int localCellCmp[1] = {1};
		int localNodes[1];

		int ghostCellCmp[1] = {2};
		int ghostNodes[1];

		// === Check Local and Ghost Cells ===
		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLONodes, 1);
		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLGhNodes, 1);

		mesh->cellConnGraph->getLocalNodes(localNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(localNodes, localNodes + 1, localCellCmp, localCellCmp + 1);

		mesh->cellConnGraph->getGhostNodes(ghostNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(ghostNodes, ghostNodes + 1, ghostCellCmp, ghostCellCmp + 1);

		// === Check Internal IDs ===
		BOOST_CHECK_EQUAL(mesh->getCellID(1), 0);
		BOOST_CHECK_EQUAL(mesh->getCellID(2), 1);

		// === Check Cell Data ===
		// Local Cell Data
		int nLocalGVerticesCmp[1] = {8};		// Global Vertex Counts per cell
		int nLocalSVerticesCmp[1] = {8};		// Stored Vertex Counts per cell
		int nLocalGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nLocalSFaceCmp[1] = {6};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(localCellCmp[i])),  nLocalGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(localCellCmp[i])),  nLocalSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(localCellCmp[i])),  nLocalGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(localCellCmp[i])),  nLocalSFaceCmp[i]);
		}

		// Ghost Cell Data
		int nGhostGVerticesCmp[1] = {8};		// Global Vertex Counts
		int nGhostSVerticesCmp[1] = {4};		// Stored Vertex Counts
		int nGhostGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nGhostSFaceCmp[1] = {1};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostSFaceCmp[i]);
		}

		// Check Face Data

		// Check Vertex Data

		// Check Boundary Data

		// Check Region Data
	}
	else if(comm.rank == 1)
	{
		int localCellCmp[1] = {2};
		int localNodes[1];

		int ghostCellCmp[3] = {1, 3, 4};
		int ghostNodes[3];

		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLONodes, 1);
		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLGhNodes, 3);

		mesh->cellConnGraph->getLocalNodes(localNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(localNodes, localNodes + 1, localCellCmp, localCellCmp + 1);

		mesh->cellConnGraph->getGhostNodes(ghostNodes, 3);
		BOOST_CHECK_EQUAL_COLLECTIONS(ghostNodes, ghostNodes + 3, ghostCellCmp, ghostCellCmp + 3);

		// === Check Internal IDs ===
		BOOST_CHECK_EQUAL(mesh->getCellID(2), 0);
		BOOST_CHECK_EQUAL(mesh->getCellID(1), 1);
		BOOST_CHECK_EQUAL(mesh->getCellID(3), 2);
		BOOST_CHECK_EQUAL(mesh->getCellID(4), 3);

		// === Check Cell Data ===
		// Local Cell Data
		int nLocalGVerticesCmp[1] = {8};		// Global Vertex Counts
		int nLocalSVerticesCmp[1] = {8};		// Stored Vertex Counts
		int nLocalGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nLocalSFaceCmp[1] = {6};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(localCellCmp[i])),  nLocalGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(localCellCmp[i])),  nLocalSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(localCellCmp[i])),  nLocalGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(localCellCmp[i])),  nLocalSFaceCmp[i]);
		}

		// Ghost Cell Data
		int nGhostGVerticesCmp[3] = {8, 6, 8};		// Global Vertex Counts
		int nGhostSVerticesCmp[3] = {4, 4, 4};		// Stored Vertex Counts
		int nGhostGFaceCmp[3] = {6, 5, 6};			// Global Face Counts per cell
		int nGhostSFaceCmp[3] = {1, 1, 1};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostSFaceCmp[i]);
		}

		// Check Face Data

		// Check Vertex Data

		// Check Boundary Data

		// Check Region Data
	}
	else if(comm.rank == 2)
	{
		int localCellCmp[1] = {3};
		int localNodes[1];

		int ghostCellCmp[1] = {2};
		int ghostNodes[1];

		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLONodes, 1);
		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLGhNodes, 1);

		mesh->cellConnGraph->getLocalNodes(localNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(localNodes, localNodes + 1, localCellCmp, localCellCmp + 1);

		mesh->cellConnGraph->getGhostNodes(ghostNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(ghostNodes, ghostNodes + 1, ghostCellCmp, ghostCellCmp + 1);

		// === Check Internal IDs ===
		BOOST_CHECK_EQUAL(mesh->getCellID(3), 0);
		BOOST_CHECK_EQUAL(mesh->getCellID(2), 1);

		// === Check Cell Data ===
		// Local Cell Data
		int nLocalGVerticesCmp[1] = {6};		// Global Vertex Counts
		int nLocalSVerticesCmp[1] = {6};		// Stored Vertex Counts
		int nLocalGFaceCmp[1] = {5};			// Global Face Counts per cell
		int nLocalSFaceCmp[1] = {5};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(localCellCmp[i])),  nLocalGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(localCellCmp[i])),  nLocalSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(localCellCmp[i])),  nLocalGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(localCellCmp[i])),  nLocalSFaceCmp[i]);
		}

		// Ghost Cell Data
		int nGhostGVerticesCmp[1] = {8};		// Global Vertex Counts
		int nGhostSVerticesCmp[1] = {4};		// Stored Vertex Counts
		int nGhostGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nGhostSFaceCmp[1] = {1};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostSFaceCmp[i]);
		}

		// Check Face Data

		// Check Vertex Data

		// Check Boundary Data

		// Check Region Data
	}
	else if(comm.rank == 3)
	{
		int localCellCmp[1] = {4};
		int localNodes[1];

		int ghostCellCmp[1] = {2};
		int ghostNodes[1];

		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLONodes, 1);
		BOOST_CHECK_EQUAL(mesh->cellConnGraph->nLGhNodes, 1);

		mesh->cellConnGraph->getLocalNodes(localNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(localNodes, localNodes + 1, localCellCmp, localCellCmp + 1);

		mesh->cellConnGraph->getGhostNodes(ghostNodes, 1);
		BOOST_CHECK_EQUAL_COLLECTIONS(ghostNodes, ghostNodes + 1, ghostCellCmp, ghostCellCmp + 1);

		// === Check Internal IDs ===
		BOOST_CHECK_EQUAL(mesh->getCellID(4), 0);
		BOOST_CHECK_EQUAL(mesh->getCellID(2), 1);

		// === Check Cell Data ===
		// Local Cell Data
		int nLocalGVerticesCmp[1] = {8};		// Global Vertex Counts
		int nLocalSVerticesCmp[1] = {8};		// Stored Vertex Counts
		int nLocalGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nLocalSFaceCmp[1] = {6};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(localCellCmp[i])),  nLocalGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(localCellCmp[i])),  nLocalSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(localCellCmp[i])),  nLocalGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(localCellCmp[i])),  nLocalSFaceCmp[i]);
		}

		// Ghost Cell Data
		int nGhostGVerticesCmp[1] = {8};		// Global Vertex Counts
		int nGhostSVerticesCmp[1] = {4};		// Stored Vertex Counts
		int nGhostGFaceCmp[1] = {6};			// Global Face Counts per cell
		int nGhostSFaceCmp[1] = {1};			// Stored Face Counts per cell

		for(int i = 0; i < 1; i++)
		{
			BOOST_CHECK_EQUAL(mesh->getCellNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostGVerticesCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNVertices(mesh->getCellID(ghostCellCmp[i])),  nGhostSVerticesCmp[i]);

			BOOST_CHECK_EQUAL(mesh->getCellNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostGFaceCmp[i]);
			BOOST_CHECK_EQUAL(mesh->getCellStoredNFaces(mesh->getCellID(ghostCellCmp[i])),  nGhostSFaceCmp[i]);
		}

		// Check Face Data

		// Check Vertex Data

		// Check Boundary Data

		// Check Region Data
	}
}

// Test 2: Test correct build of a mesh from a structured gen format
BOOST_AUTO_TEST_CASE(buildUnstructuredMesh_test2, * utf::tolerance(0.00001))
{
	cupcfd::error::eCodes status;
    cupcfd::comm::Communicator comm(MPI_COMM_WORLD);

	// Setup the partitioner config
	// Use naive so the outcome is predictable
    cupcfd::partitioner::PartitionerNaiveConfig<int,int> naiveConfig;
	cupcfd::partitioner::PartitionerConfig<int,int> partConfig(naiveConfig);

	// Setup the source config
	MeshSourceStructGenConfig<int, double> meshSourceConfig(20, 20, 20, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0);

	// Setup the config to use for building
	MeshConfig<int, double,int> meshConfig(partConfig, meshSourceConfig);

	CupCfdAoSMesh<int,double,int> * mesh;
	status = meshConfig.buildUnstructuredMesh(&mesh, comm);

	// ToDo: Check Mesh Values

	BOOST_CHECK_EQUAL(status, cupcfd::error::E_SUCCESS);
}

BOOST_AUTO_TEST_CASE(cleanup)
{
    MPI_Finalize();
}

ext: cpp

max_stars_repo_name: thorbenlouw/CUP-CFD

max_issues_repo_name: thorbenlouw/CUP-CFD

max_forks_repo_name: thorbenlouw/CUP-CFD


text:
/***************************************************************************
 * Copyright 1998-2018 by authors (see AUTHORS.txt)                        *
 *                                                                         *
 *   This file is part of LuxCoreRender.                                   *
 *                                                                         *
 * Licensed under the Apache License, Version 2.0 (the "License");         *
 * you may not use this file except in compliance with the License.        *
 * You may obtain a copy of the License at                                 *
 *                                                                         *
 *     http://www.apache.org/licenses/LICENSE-2.0                          *
 *                                                                         *
 * Unless required by applicable law or agreed to in writing, software     *
 * distributed under the License is distributed on an "AS IS" BASIS,       *
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.*
 * See the License for the specific language governing permissions and     *
 * limitations under the License.                                          *
 ***************************************************************************/

#if !defined(LUXRAYS_DISABLE_OPENCL)

#include <iosfwd>
#include <limits>

#include <boost/lexical_cast.hpp>
#include <boost/algorithm/string.hpp>

#include "slg/engines/pathoclbase/compiledscene.h"
#include "slg/kernels/kernels.h"

#include "slg/cameras/orthographic.h"
#include "slg/cameras/perspective.h"
#include "slg/cameras/stereo.h"
#include "slg/cameras/environment.h"

using namespace std;
using namespace luxrays;
using namespace slg;

void CompiledScene::CompileCamera() {
	wasCameraCompiled = true;

	//SLG_LOG("Compile Camera");

	//--------------------------------------------------------------------------
	// Camera definition
	//--------------------------------------------------------------------------

	const Camera *sceneCamera = scene->camera;

	// Initialize CameraBase

	camera.base.yon = sceneCamera->clipYon;
	camera.base.hither = sceneCamera->clipHither;
	camera.base.shutterOpen = sceneCamera->shutterOpen;
	camera.base.shutterClose = sceneCamera->shutterClose;
	camera.base.volumeIndex = sceneCamera->volume ?
		scene->matDefs.GetMaterialIndex(sceneCamera->volume) : NULL_INDEX;

	if (sceneCamera->motionSystem) {
		if (sceneCamera->motionSystem->interpolatedTransforms.size() > CAMERA_MAX_INTERPOLATED_TRANSFORM)
			throw runtime_error("Too many interpolated transformations in camera motion system: " +
					ToString(sceneCamera->motionSystem->interpolatedTransforms.size()));

		for (u_int i = 0; i < sceneCamera->motionSystem->interpolatedTransforms.size(); ++i) {
			const InterpolatedTransform &it = sceneCamera->motionSystem->interpolatedTransforms[i];

			// Here, I assume that luxrays::ocl::InterpolatedTransform and
			// luxrays::InterpolatedTransform are the same
			camera.base.interpolatedTransforms[i] = *((const luxrays::ocl::InterpolatedTransform *)&it);
		}

		camera.base.motionSystem.interpolatedTransformFirstIndex = 0;
		camera.base.motionSystem.interpolatedTransformLastIndex = sceneCamera->motionSystem->interpolatedTransforms.size() - 1;
	} else {
		camera.base.motionSystem.interpolatedTransformFirstIndex = NULL_INDEX;
		camera.base.motionSystem.interpolatedTransformLastIndex = NULL_INDEX;
	}

	// Initialize Camera specific data

	switch (sceneCamera->GetType()) {
		case Camera::ORTHOGRAPHIC: {
			const OrthographicCamera *orthoCamera = (OrthographicCamera *)sceneCamera;
			cameraType = slg::ocl::ORTHOGRAPHIC;

			memcpy(camera.base.rasterToCamera.m.m, orthoCamera->GetRasterToCamera().m.m, 4 * 4 * sizeof(float));
			memcpy(camera.base.cameraToWorld.m.m, orthoCamera->GetCameraToWorld().m.m, 4 * 4 * sizeof(float));

			camera.ortho.projCamera.lensRadius = orthoCamera->lensRadius;
			camera.ortho.projCamera.focalDistance = orthoCamera->focalDistance;

			enableCameraOculusRiftBarrel = false;
			if (orthoCamera->enableClippingPlane) {
				enableCameraClippingPlane = true;
				ASSIGN_VECTOR(camera.ortho.projCamera.clippingPlaneCenter, orthoCamera->clippingPlaneCenter);
				ASSIGN_VECTOR(camera.ortho.projCamera.clippingPlaneNormal, orthoCamera->clippingPlaneNormal);
			} else
				enableCameraClippingPlane = false;
			break;
		}
		case Camera::PERSPECTIVE: {
			const PerspectiveCamera *perspCamera = (PerspectiveCamera *)sceneCamera;
			cameraType = slg::ocl::PERSPECTIVE;

			memcpy(camera.base.rasterToCamera.m.m, perspCamera->GetRasterToCamera().m.m, 4 * 4 * sizeof(float));
			memcpy(camera.base.cameraToWorld.m.m, perspCamera->GetCameraToWorld().m.m, 4 * 4 * sizeof(float));

			camera.persp.projCamera.lensRadius = perspCamera->lensRadius;
			camera.persp.projCamera.focalDistance = perspCamera->focalDistance;

			camera.persp.screenOffsetX = perspCamera->screenOffsetX;
			camera.persp.screenOffsetY = perspCamera->screenOffsetY;

			enableCameraOculusRiftBarrel = perspCamera->enableOculusRiftBarrel;
			if (perspCamera->enableClippingPlane) {
				enableCameraClippingPlane = true;
				ASSIGN_VECTOR(camera.persp.projCamera.clippingPlaneCenter, perspCamera->clippingPlaneCenter);
				ASSIGN_VECTOR(camera.persp.projCamera.clippingPlaneNormal, perspCamera->clippingPlaneNormal);
			} else
				enableCameraClippingPlane = false;
			break;
		}
		case Camera::STEREO: {
			const StereoCamera *stereoCamera = (StereoCamera *)sceneCamera;
			cameraType = slg::ocl::STEREO;

			camera.stereo.perspCamera.projCamera.lensRadius = stereoCamera->lensRadius;
			camera.stereo.perspCamera.projCamera.focalDistance = stereoCamera->focalDistance;

			memcpy(camera.stereo.leftEyeRasterToCamera.m.m, stereoCamera->GetRasterToCamera(0).m.m, 4 * 4 * sizeof(float));
			memcpy(camera.stereo.leftEyeCameraToWorld.m.m, stereoCamera->GetCameraToWorld(0).m.m, 4 * 4 * sizeof(float));
			memcpy(camera.stereo.rightEyeRasterToCamera.m.m, stereoCamera->GetRasterToCamera(1).m.m, 4 * 4 * sizeof(float));
			memcpy(camera.stereo.rightEyeCameraToWorld.m.m, stereoCamera->GetCameraToWorld(1).m.m, 4 * 4 * sizeof(float));

			enableCameraOculusRiftBarrel = stereoCamera->enableOculusRiftBarrel;
			if (stereoCamera->enableClippingPlane) {
				enableCameraClippingPlane = true;
				ASSIGN_VECTOR(camera.stereo.perspCamera.projCamera.clippingPlaneCenter, stereoCamera->clippingPlaneCenter);
				ASSIGN_VECTOR(camera.stereo.perspCamera.projCamera.clippingPlaneNormal, stereoCamera->clippingPlaneNormal);
			} else
				enableCameraClippingPlane = false;
			break;
		}
		case Camera::ENVIRONMENT: {
			const EnvironmentCamera *envCamera = (EnvironmentCamera *)sceneCamera;
			cameraType = slg::ocl::ENVIRONMENT;

			memcpy(camera.base.rasterToCamera.m.m, envCamera->GetRasterToCamera().m.m, 4 * 4 * sizeof(float));
			memcpy(camera.base.cameraToWorld.m.m, envCamera->GetCameraToWorld().m.m, 4 * 4 * sizeof(float));

			camera.env.projCamera.lensRadius = envCamera->lensRadius;
			camera.env.projCamera.focalDistance = envCamera->focalDistance;

			enableCameraClippingPlane = false;
			break;
		}
		default:
			throw runtime_error("Unknown camera type in CompiledScene::CompileCamera(): " + boost::lexical_cast<string>(sceneCamera->GetType()));
	}
}

#endif

max_stars_repo_name: Alpistinho/LuxCore

max_issues_repo_name: Alpistinho/LuxCore

max_forks_repo_name: Alpistinho/LuxCore


text:
/**
 *  @file
 *  @copyright defined in eos/LICENSE.txt
 */
#include <appbase/application.hpp>

#include <eosio/http_plugin/http_plugin.hpp>
#include <eosio/wallet_plugin/wallet_plugin.hpp>
#include <eosio/wallet_api_plugin/wallet_api_plugin.hpp>

#include <fc/log/logger_config.hpp>
#include <fc/exception/exception.hpp>

#include <boost/exception/diagnostic_information.hpp>

#ifndef WIN32
#include <pwd.h>
#else
#pragma comment(lib, "crypt32.lib")
#endif

using namespace appbase;
using namespace eosio;

bfs::path determine_home_directory()
{
	bfs::path home;
#ifndef WIN32
	struct passwd* pwd = getpwuid(getuid());
	if (pwd) {
		home = pwd->pw_dir;
	}
	else {
		home = getenv("PATH");
	}
#else
	char tmp_dir[100] = { 0 };
	GetEnvironmentVariable("HOME", tmp_dir, sizeof(tmp_dir));
	home = tmp_dir;
#endif
   if(home.empty())
      home = "./";
   return home;
}

int main(int argc, char** argv)
{
   try {
      bfs::path home = determine_home_directory();
      app().set_default_data_dir(home / "eosio-wallet");
      app().set_default_config_dir(home / "eosio-wallet");
      app().register_plugin<wallet_api_plugin>();
      if(!app().initialize<wallet_plugin, wallet_api_plugin, http_plugin>(argc, argv))
         return -1;
      auto& http = app().get_plugin<http_plugin>();
      http.add_handler("/v1/keosd/stop", [](string, string, url_response_callback cb) { cb(200, "{}"); /*std::raise(SIGTERM);*/ } );
      app().startup();
      app().exec();
   } catch (const fc::exception& e) {
      elog("${e}", ("e",e.to_detail_string()));
   } catch (const boost::exception& e) {
      elog("${e}", ("e",boost::diagnostic_information(e)));
   } catch (const std::exception& e) {
      elog("${e}", ("e",e.what()));
   } catch (...) {
      elog("unknown exception");
   }
   return 0;
}


text:
#include "tf_misc.h"

#include <algorithm>
#include <boost/algorithm/string.hpp>
#include <boost/format.hpp>

static std::string web_encode(str_ref v)
{
	std::string d;
	d.reserve(v.size() << 1);
	while (v)
	{
		switch (v.front())
		{
		case '"':
			d += """;
			break;
		case '&':
			d += "&";
			break;
		case '<':
			d += "<";
			break;
		default:
			d += v.front();
		}
		v.pop_front();
	}
	return d;
}

static std::string web_link(str_ref title, str_ref link)
{
	return "" + web_encode(title.empty() ? link : title) + "";
}

std::string encode_field(str_ref v, bool add_br)
{
	std::string r;
	r.reserve(v.size() << 1);
	while (v)
	{
		if (boost::istarts_with(v, "ftp://")
			|| boost::istarts_with(v, "http://")
			|| boost::istarts_with(v, "https://")
			|| boost::istarts_with(v, "mailto:"))
		{
			size_t p = 0;
			while (p < v.size() && !isspace(v[p] & 0xff) && v[p] != '\"' && v[p] != '<' && v[p] != '>' && v[p] != '[' && v[p] != ']')
				p++;
			if (v[p - 1] == '!' || v[p - 1] == ',' || v[p - 1] == '.' || v[p - 1] == '?')
				p--;
			if (v[p - 1] == ')')
				p--;
			str_ref url = v.substr(0, p);
			if (boost::istarts_with(v, "ftp."))
				r += web_link(url, "ftp://" + url.s());
			else if (boost::istarts_with(v, "www."))
				r += web_link(url, "http://" + url.s());
			else
				r += web_link(boost::istarts_with(v, "mailto:") ? url.substr(7) : url, url);
			while (p--)
				v.pop_front();
		}
		else
		{
			switch (v.front())
			{
			case '\n':
				r += add_br ? "
" : " ";
				break;
			case '\r':
				break;
			case '&':
				r += "&";
				break;
			case '<':
				r += "<";
				break;
			default:
				r += v.front();
			}
			v.pop_front();
		}
	}
	return r;
}

std::string encode_text(str_ref v, bool add_quote_class)
{
	std::string r;
	r.reserve(v.size() << 1);
	while (v)
	{
		str_ref line = read_until(v, '\n');
		r += add_quote_class && boost::istarts_with(line, "> ") ? "" + encode_field(line) + "" : encode_field(line);
		r += "
";
	}
	return r;
}

std::string trim_field(const std::string& v)
{
	return boost::find_format_all_copy(boost::trim_copy(v), boost::token_finder(boost::is_space(), boost::token_compress_on), boost::const_formatter(" "));
}

std::string trim_text(const std::string& v)
{
	std::string r;
	bool copy_white = false;
	for (size_t i = 0; i < v.size(); )
	{
		size_t p = v.find('\n', i);
		if (p == std::string::npos)
			p = v.size();
		std::string line = trim_field(v.substr(i, p - i));
		if (line.empty())
			copy_white = true;
		else
		{
			if (copy_white)
			{
				if (!r.empty())
					r += '\n';
				copy_white = false;
			}
			r += line + '\n';
		}
		i = p + 1;
	}
	return r;
}

enum bb_t
{
	bb_literal,
	bb_none,
	bb_bold,
	bb_bold_close,
	bb_center,
	bb_center_close,
	bb_color,
	bb_color_close,
	bb_quote,
	bb_quote_close,
	bb_strike,
	bb_strike_close,
	bb_underline,
	bb_underline_close,
	bb_url,
	bb_unknown,
	bb_video,
	bb_end,
};

static bool operator==(str_ref a, const char* b)
{
	return a.size() == strlen(b) && !memcmp(a.data(), b, a.size());
}

bb_t get_next(str_ref& s, str_ref& a0)
{
	if (!s)
		return bb_end;
	if (s.front() != '[')
	{
		auto a = std::find(s.begin(), s.end(), '[');
		if (a == s.end())
		{
			a0 = s;
			s.clear();
		}
		else
		{
			a0 = str_ref(s.begin(), a);
			s.set_begin(a);
		}
		return bb_literal;
	}
	auto a = std::find(s.begin(), s.end(), ']');
	if (a == s.end())
	{
		a0 = s;
		s.clear();
		return bb_literal;
	}
	str_ref tag = { &s[1], a };
	s.set_begin(a + 1);
	a0.clear();
	if (tag == "b")
		return bb_bold;
	if (tag == "/b")
		return bb_bold_close;
	if (tag == "center")
		return bb_center;
	if (tag == "/center")
		return bb_center_close;
	if (boost::starts_with(tag, "color="))
	{
		a0 = tag.substr(6);
		return bb_color;
	}
	if (tag == "/color")
		return bb_color_close;
	if (boost::starts_with(tag, "font=") || tag == "/font")
		return bb_none;
	if (tag == "i" || tag == "/i")
		return bb_none;
	if (tag == "img" || tag == "IMG" || tag == "/img" || tag == "/IMG")
		return bb_none;
	if (boost::starts_with(tag, "img="))
	{
		a0 = tag.substr(4);
		return bb_literal;
	}
	if (tag == "q" || tag == "quote")
		return bb_quote;
	if (boost::starts_with(tag, "quote="))
	{
		a0 = tag.substr(6);
		return bb_quote;
	}
	if (tag == "/q" || tag == "/quote")
		return bb_quote_close;
	if (tag == "s")
		return bb_strike;
	if (tag == "/s")
		return bb_strike_close;
	if (boost::starts_with(tag, "size=") || tag == "/size")
		return bb_none;
	if (tag == "u")
		return bb_underline;
	if (tag == "/u")
		return bb_underline_close;
	if (boost::starts_with(tag, "url="))
	{
		a0 = tag.substr(4);
		return bb_url;
	}
	if (tag == "/url")
		return bb_none;
	if (boost::starts_with(tag, "video="))
	{
		a0 = tag.substr(6);
		return bb_video;
	}
	a0 = tag;
	return bb_unknown;
}

std::string bbformat(str_ref s)
{
	std::string d;
	str_ref a0;
	while (1)
	{
		switch (get_next(s, a0))
		{
		case bb_literal:
			d += encode_field(a0, true);
			break;
		case bb_none:
			break;
		case bb_bold:
			d += "";
			break;
		case bb_bold_close:
			d += "";
			break;
		case bb_center:
			d += "";
			break;
		case bb_center_close:
			d += "";
			break;
		case bb_color:
			d += ""; // escape a0
			break;
		case bb_color_close:
			d += "";
			break;
		case bb_quote:
			d += "
";
			if (a0)
				d += "" + encode_field(a0) + " wrote:";
			break;
		case bb_quote_close:
			d += "
";
			break;
		case bb_strike:
			d += "<s>";
			break;
		case bb_strike_close:
			d += "</s>";
			break;
		case bb_underline:
			d += "";
			break;
		case bb_underline_close:
			d += "";
			break;
		case bb_url:
			d += encode_field(a0, true) + " ";
			break;
		case bb_video:
			d += encode_field(a0, true) + " ";
			break;
		case bb_unknown:
			d += "[" + encode_field(a0) + "]";
			break;
		case bb_end:
			return d;
		default:
			assert(false);
		}
	}
}


text:
/**
 * @file advectionfv2d_main.cc
 * @brief NPDE homework AdvectionFV2D code
 * @author Philipp Egg
 * @date 21.06.2020
 * @copyright Developed at ETH Zurich
 */

#include <lf/assemble/assemble.h>
#include <lf/base/base.h>
#include <lf/geometry/geometry.h>
#include <lf/io/io.h>
#include <lf/mesh/mesh.h>
#include <lf/mesh/test_utils/test_meshes.h>
#include <lf/mesh/utils/utils.h>
#include <lf/refinement/refinement.h>

#include <Eigen/Core>
#include <array>
#include <cmath>
#include <cstdlib>
#include <fstream>
#include <iostream>
#include <memory>
#include <string>
#include <vector>

#include "advectionfv2d.h"

// Use this function to plot your solution
void write_vtk(const lf::assemble::DofHandler &dofh,
               const Eigen::VectorXd &solution, std::string name) {
  std::shared_ptr<const lf::mesh::Mesh> mesh_p = dofh.Mesh();
  lf::io::VtkWriter vtk_writer(mesh_p, name + ".vtk");
  auto cell_data_ref =
      lf::mesh::utils::make_CodimMeshDataSet<double>(mesh_p, 0);
  for (const lf::mesh::Entity *cell : mesh_p->Entities(0)) {
    int row = dofh.GlobalDofIndices(*cell)[0];
    cell_data_ref->operator()(*cell) = solution[row];
  }
  vtk_writer.WriteCellData(name, *cell_data_ref);
}

int main() {
  /* SAM_LISTING_BEGIN_1 */
  // Define velocity field beta
  // Note that the problem description requires ||B|| <= 1
  auto beta = [](Eigen::Vector2d x) -> Eigen::Vector2d {
    return Eigen::Vector2d(-x[1], x[0]) / std::sqrt(2.0);
  };

  // Functor for initial bump
  auto u0 = [](Eigen::Vector2d x) -> double {
    Eigen::Vector2d x0(0.8, 0.2);
    double d = 0.2;
    double dist = (x - x0).norm();
    if (dist < d) {
      double cos = std::cos(M_PI / (2.0 * d) * dist);
      return cos * cos;
    } else {
      return 0.0;
    }
  };

  // Task 8-8.o
  // Generate a mesh hierarchy
  double T = 1.0;

  //////////////////////////////////////////////////////////////////////////////
  // TODO inconsistancy: g vs. G
  //////////////////////////////////////////////////////////////////////////////
  auto mesh_p = lf::mesh::test_utils::GenerateHybrid2DTestMesh(0, 1.0 / 3.0);

  auto mesh_seq_p{
      lf::refinement::GenerateMeshHierarchyByUniformRefinemnt(mesh_p, 6)};

  std::vector<int> vector_num_cells;
  std::vector<double> vector_l2error;

  // Iterate over mesh levels starting from thrid refinement
  int num_meshes = mesh_seq_p->NumLevels();
  for (int level = 3; level < num_meshes; ++level) {
    std::cout << "Computing L2Error for level: " << level << std::endl;

#if SOLUTION
    // Get the current mesh
    auto cur_mesh = mesh_seq_p->getMesh(level);

    // Create a DOF Hander for the current mesh
    const lf::assemble::UniformFEDofHandler cur_dofh(
        cur_mesh, {{lf::base::RefEl::kPoint(), 0},
                   {lf::base::RefEl::kSegment(), 0},
                   {lf::base::RefEl::kTria(), 1},
                   {lf::base::RefEl::kQuad(), 1}});
    int N = cur_dofh.NumDofs();

    // Compute cell normals
    std::shared_ptr<lf::mesh::utils::CodimMeshDataSet<
        Eigen::Matrix<double, 2, Eigen::Dynamic>>>
        normal_vectors = AdvectionFV2D::computeCellNormals(cur_dofh.Mesh());

    // Compute adjecent cells
    std::shared_ptr<lf::mesh::utils::CodimMeshDataSet<
        std::array<const lf::mesh::Entity *, 4>>>
        adjacentCells = AdvectionFV2D::getAdjacentCellPointers(cur_dofh.Mesh());

    // Get approximate solution from simulation
    Eigen::VectorXd mu_approx = AdvectionFV2D::simulateAdvection(
        cur_dofh, beta, u0, adjacentCells, normal_vectors, T);
    write_vtk(cur_dofh, mu_approx, "approx" + std::to_string(level));

    // Get exact solution at barycenters of cells
    Eigen::VectorXd mu_exact = AdvectionFV2D::refSolution(cur_dofh, u0, T);
    write_vtk(cur_dofh, mu_exact, "exact" + std::to_string(level));

    // Compute L2 error in barycenter
    double l2_error = 0;
    for (const lf::mesh::Entity *cell : cur_mesh->Entities(0)) {
      const lf::geometry::Geometry *geo_p = cell->Geometry();
      double area = lf::geometry::Volume(*geo_p);
      int idx = cur_dofh.GlobalDofIndices(*cell)[0];
      l2_error += std::pow((mu_approx[idx] - mu_exact[idx]), 2) * area;
    }
    l2_error = std::sqrt(l2_error);
#else
    //====================
    // Your code goes here
    // Compute the number N of DOFs and the L2-error,
    // and replace the lines below:
    int N = 1;
    double l2_error = 1.0;
    // If you want, you can use the function write_vtk(...),
    // defined in this file, to plot your solution.
    //====================
#endif

    vector_num_cells.push_back(N);
    vector_l2error.push_back(l2_error);
  }

  // Write Output file of Task 8-8.o
  std::ofstream csv_file;
  csv_file.open("advectionfv2d.csv");
  for (int i = 0; i < vector_num_cells.size(); ++i) {
    std::cout << "Cells: " << vector_num_cells.at(i)
              << " | L2Error: " << vector_l2error.at(i) << std::endl;
    csv_file << vector_num_cells.at(i) << "," << vector_l2error.at(i) << "\n";
  }
  csv_file.close();

  std::system("python3 " CURRENT_SOURCE_DIR
              "/advectionfv2d.py " CURRENT_BINARY_DIR
              "/advectionfv2d.csv " CURRENT_BINARY_DIR "/solution.eps");
  /* SAM_LISTING_END_1 */

  /* SAM_LISTING_BEGIN_2 */
  // Task 8-8.q
  // Compute threshold for fourth refinement level
  int level = 4;
  auto mesh = mesh_seq_p->getMesh(level);

#if SOLUTION
  // Create a DOF Hander for the current mesh
  const lf::assemble::UniformFEDofHandler dofh(
      mesh, {{lf::base::RefEl::kPoint(), 0},
             {lf::base::RefEl::kSegment(), 0},
             {lf::base::RefEl::kTria(), 1},
             {lf::base::RefEl::kQuad(), 1}});

  int threshold = AdvectionFV2D::findCFLthreshold(dofh, beta, T);
  int cfl_thres = int((T / AdvectionFV2D::computeHmin(mesh) + 1));
#else
  //====================
  // Your code goes here
  // Replace the two variables below:
  int threshold = 0.0;  // Threshold computed by findCFLthreshold(...)
  int cfl_thres = 0.0;  // Threshold obtained from CLF using computeHmin(mesh)
                        //====================
#endif

  std::cout << "Threshold for level " << level << " is " << threshold
            << " | Threshold from CFL is " << cfl_thres << std::endl;
  /* SAM_LISTING_END_2 */

  return 0;
}

ext: cc


text:
// Copyright DEWETRON GmbH 2019

#include "odkapi_update_channels_xml.h"
#include "odkuni_xpugixml.h"

#include <boost/test/unit_test.hpp>
#include <boost/version.hpp>

#include <cmath>

BOOST_AUTO_TEST_SUITE(update_plugin_channels_test)

BOOST_AUTO_TEST_CASE(ParseXML)
{
    const char* const xml_content =
        R"xxx(<?xml version='1.0' encoding='UTF-8'?>
<UpdatePluginChannels protocol_version="1.0">

    <Channel local_id="4294967294"/>

    <Channel local_id="2" default_name="Ch1" domain="UnitTests" deletable="true" valid="true">
        <DataFormat sample_occurrence="SYNC" sample_format="double" sample_dimension="11" />
        <SimpleTimebase frequency="10000"/>
        <Parent local_id="4294967294"/>
        <Property name = "Frequency">
            <ScalarValue>
                <Value>10000</Value>
                <Unit>Hz</Unit>
            </ScalarValue>
        </Property>
    </Channel>

    <ListTopology>
        <Group name="Test">
            <Channel local_id="4294967294">
                <Channel local_id="2"/>
            </Channel>
        </Group>
    </ListTopology>

</UpdatePluginChannels>
            )xxx"
            ;

    using namespace odk;
    UpdateChannelsTelegram telegram;
    BOOST_CHECK(telegram.parse(xml_content));

    BOOST_REQUIRE_EQUAL(telegram.m_channels.size(), 2);

    BOOST_CHECK_EQUAL(telegram.m_channels[0].m_local_id, 4294967294);
    BOOST_CHECK(telegram.m_channels[0].m_dataformat_info.m_sample_occurrence == ChannelDataformat::SampleOccurrence::INVALID);
    BOOST_CHECK(telegram.m_channels[0].m_dataformat_info.m_sample_format == ChannelDataformat::SampleFormat::INVALID);
    BOOST_CHECK_EQUAL(telegram.m_channels[0].m_dataformat_info.m_sample_dimension, -1);
    BOOST_CHECK(telegram.m_channels[0].m_timebase.m_type == Timebase::TimebaseType::NONE);
    BOOST_CHECK(std::isnan(telegram.m_channels[0].m_timebase.m_frequency));
    BOOST_CHECK(!telegram.m_channels[0].m_deletable);
    BOOST_CHECK(!telegram.m_channels[0].m_valid);
    auto default_channel = UpdateChannelsTelegram::PluginChannelInfo(4294967294);
    BOOST_CHECK(telegram.m_channels[0] == default_channel);

    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_local_id, 2);
    BOOST_CHECK(telegram.m_channels[1].m_dataformat_info.m_sample_occurrence == ChannelDataformat::SampleOccurrence::SYNC);
    BOOST_CHECK(telegram.m_channels[1].m_dataformat_info.m_sample_format == ChannelDataformat::SampleFormat::DOUBLE);
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_dataformat_info.m_sample_dimension, 11);
    BOOST_CHECK(telegram.m_channels[1].m_timebase.m_type == Timebase::TimebaseType::SIMPLE);
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_timebase.m_frequency, 10000.0);
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_default_name, "Ch1");
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_domain, "UnitTests");
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_local_parent_id, 4294967294);
    BOOST_CHECK(telegram.m_channels[1].m_deletable);
    BOOST_CHECK(telegram.m_channels[1].m_valid);

    BOOST_REQUIRE_EQUAL(telegram.m_channels[1].m_channel_config.m_properties.size(), 1);
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_channel_config.m_properties[0].getName(), "Frequency");
    BOOST_CHECK_EQUAL(telegram.m_channels[1].m_channel_config.m_properties[0].getType(), odk::Property::SCALAR);
    BOOST_CHECK(telegram.m_channels[1].m_channel_config.m_properties[0].getScalarValue() == odk::Scalar(10000.0, "Hz"));

    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children.size(), 1);
    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children[0].m_group_name, "Test");
    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children[0].m_children.size(), 1);
    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children[0].m_children[0].m_channel_id, 4294967294);

    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children[0].m_children[0].m_children.size(), 1);
    BOOST_REQUIRE_EQUAL(telegram.m_list_topology.m_children[0].m_children[0].m_children[0].m_channel_id, 2);
}

BOOST_AUTO_TEST_CASE(CreateRequestXML)
{
    using namespace odk;
    UpdateChannelsTelegram input;
    input.addChannel(2)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::SYNC,
            ChannelDataformat::SampleFormat::SINT32)
        .m_channel_config.addProperty(odk::Property("IP", "127.0.0.1"));

    input.addChannel(5)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::ASYNC,
            ChannelDataformat::SampleFormat::DOUBLE, 5)
        .setSimpleTimebase(12345)
        .setDeletable(true)
        .setLocalParent(2)
        .setDefaultName("DoubleVector")
        .setDomain("Tests");

    input.addChannel(8)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::SINGLE_VALUE,
            ChannelDataformat::SampleFormat::COMPLEX_FLOAT, 10)
        .setLocalParent(5);

    input.addChannel(11)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::NEVER,
            ChannelDataformat::SampleFormat::NONE, 0)
        .setLocalParent(2);

    UpdateChannelsTelegram expected_output;
    expected_output.addChannel(2)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::SYNC,
            ChannelDataformat::SampleFormat::SINT32)
        .m_channel_config.addProperty(odk::Property("IP", "127.0.0.1"));

    expected_output.addChannel(5)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::ASYNC,
            ChannelDataformat::SampleFormat::DOUBLE, 5)
        .setSimpleTimebase(12345)
        .setDeletable(true)
        .setLocalParent(2)
        .setDefaultName("DoubleVector")
        .setDomain("Tests");

    expected_output.addChannel(11)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::NEVER,
            ChannelDataformat::SampleFormat::NONE, 0)
        .setLocalParent(2);

    expected_output.addChannel(8)
        .setSampleFormat(
            ChannelDataformat::SampleOccurrence::SINGLE_VALUE,
            ChannelDataformat::SampleFormat::COMPLEX_FLOAT, 10)
        .setLocalParent(5);

    auto xml = input.generate();

    odk::UpdateChannelsTelegram output;

    BOOST_CHECK(!(input == output));
    BOOST_CHECK(output.parse(xml.c_str()));
    BOOST_CHECK(expected_output == output);

    BOOST_CHECK(getRootChannels(output) == std::vector<std::uint32_t>(std::initializer_list<std::uint32_t>({ 2u })));
    BOOST_CHECK(getChildrenOfChannel(output, 2, false) == std::vector<std::uint32_t>({5, 11}));
    BOOST_CHECK(getChildrenOfChannel(output, 5, false) == std::vector<std::uint32_t>(std::initializer_list<std::uint32_t>({ 8u })));
    BOOST_CHECK(getChildrenOfChannel(output, 8, false) == std::vector<std::uint32_t>());
    BOOST_CHECK(getChildrenOfChannel(output, 2, true) == std::vector<std::uint32_t>({5, 8, 11}));
}

BOOST_AUTO_TEST_CASE(CreateRequestXML2)
{
    using namespace odk;
    UpdateChannelsTelegram input;
    auto& g = input.m_list_topology.appendGroup("Plugin Group");
    auto& gc = g.appendChannel(1);
    gc.appendChannel(2);
    g.appendChannel(11);
    BOOST_CHECK_EQUAL(input.m_list_topology.m_children[0].m_group_name, "Plugin Group");
    BOOST_CHECK_EQUAL(input.m_list_topology.m_children[0].m_children[0].m_channel_id, 1);

    auto xml = input.generate();

    odk::UpdateChannelsTelegram output;

    BOOST_CHECK(input != output);
    BOOST_CHECK(output.parse(xml.c_str()));
    BOOST_CHECK(input == output);
}

BOOST_AUTO_TEST_SUITE_END()

max_stars_repo_name: HpLightcorner/OXYGEN-SDK

max_issues_repo_name: HpLightcorner/OXYGEN-SDK

max_forks_repo_name: HpLightcorner/OXYGEN-SDK


text:
/*ckwg +29
 * Copyright 2011-2014 by Kitware, Inc.
 * All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are met:
 *
 *  * Redistributions of source code must retain the above copyright notice,
 *    this list of conditions and the following disclaimer.
 *
 *  * Redistributions in binary form must reproduce the above copyright notice,
 *    this list of conditions and the following disclaimer in the documentation
 *    and/or other materials provided with the distribution.
 *
 *  * Neither name of Kitware, Inc. nor the names of any contributors may be used
 *    to endorse or promote products derived from this software without specific
 *    prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS ``AS IS''
 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
 * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS OR CONTRIBUTORS BE LIABLE FOR
 * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
 * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

/**
 * \file
 *
 * \brief Implementation of \link maptk::config_block configuration \endlink object
 */

#include "config_block.h"

#include <boost/algorithm/string/case_conv.hpp>
#include <boost/algorithm/string/predicate.hpp>
#include <boost/bind.hpp>
#include <boost/foreach.hpp>
#include <boost/none.hpp>

#include <algorithm>
#include <iterator>
#include <sstream>

namespace maptk
{

config_block_key_t const config_block::block_sep = config_block_key_t(":");
config_block_key_t const config_block::global_value = config_block_key_t("_global");

/// private helper method for determining key path prefixes
static bool does_not_begin_with(config_block_key_t const& key,
                                config_block_key_t const& name);
/// private helper method to strip a block name from a key path
static config_block_key_t strip_block_name(config_block_key_t const& subblock,
                                           config_block_key_t const& key);

/// Create an empty configuration.
config_block_sptr
config_block
::empty_config(config_block_key_t const& name)
{
  // remember, config_block_sptr is a boost shared pointer
  return config_block_sptr(new config_block(name, config_block_sptr()));
}

/// Destructor
config_block
::~config_block()
{
}

/// Get the name of this \c config_block instance.
config_block_key_t
config_block
::get_name()
{
  return this->m_name;
}

/// Get a subblock from the configuration.
config_block_sptr
config_block
::subblock(config_block_key_t const& key) const
{
  config_block_sptr conf = empty_config(key);

  BOOST_FOREACH (config_block_key_t const& key_name, available_values())
  {
    if (does_not_begin_with(key_name, key))
    {
      continue;
    }

    config_block_key_t const stripped_key_name = strip_block_name(key, key_name);

    conf->set_value(stripped_key_name,
                    m_get_value(key_name),
                    get_description(key_name));
  }

  return conf;
}

/// Get a subblock view into the configuration.
config_block_sptr
config_block
::subblock_view(config_block_key_t const& key)
{
  return config_block_sptr(new config_block(key, shared_from_this()));
}

config_block_description_t
config_block
::get_description(config_block_key_t const& key) const
{
  if (m_parent)
  {
    return m_parent->get_description(m_name + block_sep + key);
  }

  store_t::const_iterator i = m_descr_store.find(key);
  if (i == m_descr_store.end())
  {
    throw no_such_configuration_value_exception(key);
  }

  return i->second;
}

/// Remove a value from the configuration.
void
config_block
::unset_value(config_block_key_t const& key)
{
  if (m_parent)
  {
    m_parent->unset_value(m_name + block_sep + key);
  }
  else
  {
    if (is_read_only(key))
    {
      config_block_value_t const current_value = get_value<config_block_value_t>(key, config_block_value_t());

      throw unset_on_read_only_value_exception(key, current_value);
    }

    store_t::iterator const i = m_store.find(key);
    store_t::iterator const j = m_descr_store.find(key);

    // value and descr stores managed in parallel, so if key doesn't exist in
    // value store, there will be no parallel value in the descr store.
    if (i == m_store.end())
    {
      throw no_such_configuration_value_exception(key);
    }

    m_store.erase(i);
    m_descr_store.erase(j);
  }
}

/// Query if a value is read-only.
bool
config_block
::is_read_only(config_block_key_t const& key) const
{
  return (0 != m_ro_list.count(key));
}

/// Set the value within the configuration as read-only.
void
config_block
::mark_read_only(config_block_key_t const& key)
{
  m_ro_list.insert(key);
}

/// Merge the values in \p config_block into the current config.
void
config_block
::merge_config(config_block_sptr const& conf)
{
  config_block_keys_t const keys = conf->available_values();

  BOOST_FOREACH (config_block_key_t const& key, keys)
  {
    config_block_value_t const& val = conf->get_value<config_block_value_t>(key);
    config_block_description_t const& descr = conf->get_description(key);

    m_set_value(key, val, descr);
  }
}

///Return the values available in the configuration.
config_block_keys_t
config_block
::available_values() const
{
  config_block_keys_t keys;

  if (m_parent)
  {
    config_block_keys_t parent_keys = m_parent->available_values();

    config_block_keys_t::iterator const i = std::remove_if(parent_keys.begin(), parent_keys.end(), boost::bind(does_not_begin_with, _1, m_name));

    parent_keys.erase(i, parent_keys.end());

    std::transform(parent_keys.begin(), parent_keys.end(), std::back_inserter(keys), boost::bind(strip_block_name, m_name, _1));
  }
  else
  {
    BOOST_FOREACH (store_t::value_type const& value, m_store)
    {
      config_block_key_t const& key = value.first;

      keys.push_back(key);
    }
  }

  return keys;
}

/// Check if a value exists for \p key.
bool
config_block
::has_value(config_block_key_t const& key) const
{
  if (m_parent)
  {
    return m_parent->has_value(m_name + block_sep + key);
  }

  return (0 != m_store.count(key));
}

/// Internal constructor
config_block
::config_block(config_block_key_t const& name, config_block_sptr parent)
  : m_parent(parent)
  , m_name(name)
  , m_store()
  , m_descr_store()
  , m_ro_list()
{
}

/// private helper method to extract a value for a key
boost::optional<config_block_value_t>
config_block
::find_value(config_block_key_t const& key) const
{
  if (!has_value(key))
  {
    return boost::none;
  }

  return m_get_value(key);
}

/// private value getter function
config_block_value_t
config_block
::m_get_value(config_block_key_t const& key) const
{
  if (m_parent)
  {
    return m_parent->m_get_value(m_name + block_sep + key);
  }

  store_t::const_iterator i = m_store.find(key);

  if (i == m_store.end())
  {
    return config_block_value_t();
  }

  return i->second;
}

/// private key/value setter
void
config_block
::m_set_value(config_block_key_t const& key,
              config_block_value_t const& value,
              config_block_description_t const& descr)
{
  if (m_parent)
  {
    m_parent->set_value(m_name + block_sep + key, value, descr);
  }
  else
  {
    if (is_read_only(key))
    {
      config_block_value_t const current_value = get_value<config_block_value_t>(key, config_block_value_t());

      throw set_on_read_only_value_exception(key, current_value, value);
    }

    m_store[key] = value;

    // Only assign the description given if there is no stored description
    // for this key, or the given description is non-zero.
    if (m_descr_store.count(key) == 0 || descr.size() > 0)
    {
      m_descr_store[key] = descr;
    }
  }
}

/// Type-specific casting handling, bool specialization
template <>
bool
config_block_cast(config_block_value_t const& value)
{
  static config_block_value_t const true_string = config_block_value_t("true");
  static config_block_value_t const false_string = config_block_value_t("false");
  static config_block_value_t const yes_string = config_block_value_t("yes");
  static config_block_value_t const no_string = config_block_value_t("no");

  config_block_value_t const value_lower = boost::to_lower_copy(value);

  if (value_lower == true_string || value_lower == yes_string)
  {
    return true;
  }
  else if (value_lower == false_string || value_lower == no_string)
  {
    return false;
  }

  return config_block_cast_default<bool, config_block_value_t>(value);
}

/// private helper method for determining key path prefixes
/**
 * \param key   The key string to check.
 * \param name  The prefix string to check for. Should not include a trailing
 *              block separator.
 * \returns True if the given key does not begin with the given name and is
 *          not a global variable.
 */
bool
does_not_begin_with(config_block_key_t const& key, config_block_key_t const& name)
{
  static config_block_key_t const global_start = config_block::global_value + config_block::block_sep;

  return (!boost::starts_with(key, name + config_block::block_sep) &&
          !boost::starts_with(key, global_start));
}

/// private helper method to strip a block name from a key path
/**
 * Conditionally strip the given subblock name from the given key path. If the
 * given key doesn't start with the given subblock, the given key is returned
 * as is.
 *
 * \param subblock  The subblock string to strip if present.
 * \param key       The key to conditionally strip from.
 * \returns The stripped key name.
 */
config_block_key_t
strip_block_name(config_block_key_t const& subblock, config_block_key_t const& key)
{
  if (!boost::starts_with(key, subblock + config_block::block_sep))
  {
    return key;
  }

  return key.substr(subblock.size() + config_block::block_sep.size());
}

}

ext: cxx

max_stars_repo_name: efernandez/maptk

max_issues_repo_name: efernandez/maptk

max_forks_repo_name: efernandez/maptk