text:
/* ----------------------------------------------------------------------------

 * GTSAM Copyright 2010, Georgia Tech Research Corporation,
 * Atlanta, Georgia 30332-0415
 * All Rights Reserved
 * Authors: Frank Dellaert, et al. (see THANKS for the full author list)

 * See LICENSE for the license information

 * -------------------------------------------------------------------------- */

/**
 * @file    Rot3.cpp
 * @brief   Rotation, common code between Rotation matrix and Quaternion
 * @author  Alireza Fathi
 * @author  Christian Potthast
 * @author  Frank Dellaert
 * @author  Richard Roberts
 */

#include <gtsam/geometry/Rot3.h>
#include <gtsam/geometry/SO3.h>
#include <boost/math/constants/constants.hpp>

#include <cmath>
#include <random>

using namespace std;

namespace gtsam {

/* ************************************************************************* */
void Rot3::print(const std::string& s) const {
  gtsam::print((Matrix)matrix(), s);
}

/* ************************************************************************* */
Rot3 Rot3::Random(std::mt19937& rng) {
  // TODO allow any engine without including all of boost :-(
  Unit3 axis = Unit3::Random(rng);
  uniform_real_distribution<double> randomAngle(-M_PI, M_PI);
  double angle = randomAngle(rng);
  return AxisAngle(axis, angle);
}

/* ************************************************************************* */
Rot3 Rot3::AlignPair(const Unit3& axis, const Unit3& a_p, const Unit3& b_p) {
  // if a_p is already aligned with b_p, return the identity rotation
  if (std::abs(a_p.dot(b_p)) > 0.999999999) {
    return Rot3();
  }

  // Check axis was not degenerate cross product
  const Vector3 z = axis.unitVector();
  if (z.hasNaN())
    throw std::runtime_error("AlignSinglePair: axis has Nans");

  // Now, calculate rotation that takes b_p to a_p
  const Matrix3 P = I_3x3 - z * z.transpose();  // orthogonal projector
  const Vector3 a_po = P * a_p.unitVector();    // point in a orthogonal to axis
  const Vector3 b_po = P * b_p.unitVector();    // point in b orthogonal to axis
  const Vector3 x = a_po.normalized();          // x-axis in axis-orthogonal plane, along a_p vector
  const Vector3 y = z.cross(x);                 // y-axis in axis-orthogonal plane
  const double u = x.dot(b_po);                 // x-coordinate for b_po
  const double v = y.dot(b_po);                 // y-coordinate for b_po
  double angle = std::atan2(v, u);
  return Rot3::AxisAngle(z, -angle);
}

/* ************************************************************************* */
Rot3 Rot3::AlignTwoPairs(const Unit3& a_p, const Unit3& b_p,  //
                         const Unit3& a_q, const Unit3& b_q) {
  // there are three frames in play:
  // a: the first frame in which p and q are measured
  // b: the second frame in which p and q are measured
  // i: intermediate, after aligning first pair

  // First, find rotation around that aligns a_p and b_p
  Rot3 i_R_b = AlignPair(a_p.cross(b_p), a_p, b_p);

  // Rotate points in frame b to the intermediate frame,
  // in which we expect the point p to be aligned now
  Unit3 i_q = i_R_b * b_q;
  assert(assert_equal(a_p, i_R_b * b_p, 1e-6));

  // Now align second pair: we need to align i_q to a_q
  Rot3 a_R_i = AlignPair(a_p, a_q, i_q);
  assert(assert_equal(a_p, a_R_i * a_p, 1e-6));
  assert(assert_equal(a_q, a_R_i * i_q, 1e-6));

  // The desired rotation is the product of both
  Rot3 a_R_b = a_R_i * i_R_b;
  return a_R_b;
}

/* ************************************************************************* */
bool Rot3::equals(const Rot3 & R, double tol) const {
  return equal_with_abs_tol(matrix(), R.matrix(), tol);
}

/* ************************************************************************* */
Point3 Rot3::operator*(const Point3& p) const {
  return rotate(p);
}

/* ************************************************************************* */
Unit3 Rot3::rotate(const Unit3& p,
    OptionalJacobian<2,3> HR, OptionalJacobian<2,2> Hp) const {
  Matrix32 Dp;
  Unit3 q = Unit3(rotate(p.point3(Hp ? &Dp : 0)));
  if (Hp) *Hp = q.basis().transpose() * matrix() * Dp;
  if (HR) *HR = -q.basis().transpose() * matrix() * p.skew();
  return q;
}

/* ************************************************************************* */
Unit3 Rot3::unrotate(const Unit3& p,
    OptionalJacobian<2,3> HR, OptionalJacobian<2,2> Hp) const {
  Matrix32 Dp;
  Unit3 q = Unit3(unrotate(p.point3(Dp)));
  if (Hp) *Hp = q.basis().transpose() * matrix().transpose () * Dp;
  if (HR) *HR = q.basis().transpose() * q.skew();
  return q;
}

/* ************************************************************************* */
Unit3 Rot3::operator*(const Unit3& p) const {
  return rotate(p);
}

/* ************************************************************************* */
// see doc/math.lyx, SO(3) section
Point3 Rot3::unrotate(const Point3& p, OptionalJacobian<3,3> H1,
    OptionalJacobian<3,3> H2) const {
  const Matrix3& Rt = transpose();
  Point3 q(Rt * p); // q = Rt*p
  const double wx = q.x(), wy = q.y(), wz = q.z();
  if (H1)
    *H1 << 0.0, -wz, +wy, +wz, 0.0, -wx, -wy, +wx, 0.0;
  if (H2)
    *H2 = Rt;
  return q;
}

/* ************************************************************************* */
Point3 Rot3::column(int index) const{
  if(index == 3)
    return r3();
  else if(index == 2)
    return r2();
  else if(index == 1)
    return r1(); // default returns r1
  else
    throw invalid_argument("Argument to Rot3::column must be 1, 2, or 3");
}

/* ************************************************************************* */
Vector3 Rot3::xyz() const {
  Matrix3 I;Vector3 q;
  boost::tie(I,q)=RQ(matrix());
  return q;
}

/* ************************************************************************* */
Vector3 Rot3::ypr() const {
  Vector3 q = xyz();
  return Vector3(q(2),q(1),q(0));
}

/* ************************************************************************* */
Vector3 Rot3::rpy() const {
  return xyz();
}

/* ************************************************************************* */
Vector Rot3::quaternion() const {
  gtsam::Quaternion q = toQuaternion();
  Vector v(4);
  v(0) = q.w();
  v(1) = q.x();
  v(2) = q.y();
  v(3) = q.z();
  return v;
}

/* ************************************************************************* */
Matrix3 Rot3::ExpmapDerivative(const Vector3& x) {
  return SO3::ExpmapDerivative(x);
}

/* ************************************************************************* */
Matrix3 Rot3::LogmapDerivative(const Vector3& x)    {
  return SO3::LogmapDerivative(x);
}

/* ************************************************************************* */
pair<Matrix3, Vector3> RQ(const Matrix3& A) {

  double x = -atan2(-A(2, 1), A(2, 2));
  Rot3 Qx = Rot3::Rx(-x);
  Matrix3 B = A * Qx.matrix();

  double y = -atan2(B(2, 0), B(2, 2));
  Rot3 Qy = Rot3::Ry(-y);
  Matrix3 C = B * Qy.matrix();

  double z = -atan2(-C(1, 0), C(1, 1));
  Rot3 Qz = Rot3::Rz(-z);
  Matrix3 R = C * Qz.matrix();

  Vector xyz = Vector3(x, y, z);
  return make_pair(R, xyz);
}

/* ************************************************************************* */
ostream &operator<<(ostream &os, const Rot3& R) {
  os << "\n";
  os << '|' << R.r1().x() << ", " << R.r2().x() << ", " << R.r3().x() << "|\n";
  os << '|' << R.r1().y() << ", " << R.r2().y() << ", " << R.r3().y() << "|\n";
  os << '|' << R.r1().z() << ", " << R.r2().z() << ", " << R.r3().z() << "|\n";
  return os;
}

/* ************************************************************************* */
Rot3 Rot3::slerp(double t, const Rot3& other) const {
  return interpolate(*this, other, t);
}

/* ************************************************************************* */

} // namespace gtsam

max_stars_repo_name: ori-drs/gtsam

max_issues_repo_name: ori-drs/gtsam

max_forks_repo_name: ori-drs/gtsam


text:
#ifndef STAN_MATH_PRIM_SCAL_FUN_LB_FREE_HPP
#define STAN_MATH_PRIM_SCAL_FUN_LB_FREE_HPP

#include <stan/math/prim/scal/fun/identity_free.hpp>
#include <stan/math/prim/scal/err/check_greater_or_equal.hpp>
#include <boost/math/tools/promotion.hpp>
#include <cmath>
#include <limits>

namespace stan {
namespace math {

/**
 * Return the unconstrained value that produces the specified
 * lower-bound constrained value.
 *
 * If the lower bound is negative infinity, it is ignored and
 * the function reduces to identity_free(y).
 *
 * @tparam T type of scalar
 * @tparam L type of lower bound
 * @param[in] y input scalar
 * @param[in] lb lower bound
 * @return unconstrained value that produces the input when
 * constrained
 * @throw std::domain_error if y is lower than the lower bound
 */
template <typename T, typename L>
inline typename boost::math::tools::promote_args<T, L>::type lb_free(
    const T& y, const L& lb) {
  using std::log;
  if (lb == -std::numeric_limits<double>::infinity())
    return identity_free(y);
  check_greater_or_equal("lb_free", "Lower bounded variable", y, lb);
  return log(y - lb);
}

}  // namespace math
}  // namespace stan
#endif

ext: hpp

max_stars_repo_name: jrmie/math

max_issues_repo_name: jrmie/math

max_forks_repo_name: jrmie/math


text:
// Copyright 2020 Matthias Heinz
#include "negele/matrix.h"

#include <lapack.h>

#include <cmath>
#include <iostream>
#include <vector>

#include <boost/numeric/odeint.hpp>

#include "negele/constants.h"

namespace boost {
namespace numeric {
namespace odeint {

template <>
struct is_resizeable<negele::matrix::MomentumSpaceMatrix> {  // declare
                                                             // resizeability
  typedef boost::true_type type;
  const static bool value = type::value;
};

template <>
struct same_size_impl<negele::matrix::MomentumSpaceMatrix,
                      negele::matrix::MomentumSpaceMatrix> {  // define how to
                                                              // check size
  static bool same_size(const negele::matrix::MomentumSpaceMatrix& v1,
                        const negele::matrix::MomentumSpaceMatrix& v2) {
    return v1.same_size(v2);
  }
};

template <>
struct resize_impl<negele::matrix::MomentumSpaceMatrix,
                   negele::matrix::MomentumSpaceMatrix> {  // define how to
                                                           // resize
  static void resize(negele::matrix::MomentumSpaceMatrix& v1,
                     const negele::matrix::MomentumSpaceMatrix& v2) {
    v1.resize(v2);
  }
};

}  // namespace odeint
}  // namespace numeric
}  // namespace boost

negele::matrix::SquareMatrix::SquareMatrix(int dim, double value)
    : m_dim(dim), m_matrix(dim * dim, value) {}

negele::matrix::SquareMatrix::SquareMatrix(int dim,
                                           const std::vector<double>& mels)
    : m_dim(dim), m_matrix(mels) {
  if (m_matrix.size() != m_dim * m_dim) {
    std::cout << "Matrix initialized with incompatible dimensions: \n"
              << "\tsingle-axis dimension given: " << m_dim << "\n"
              << "\tactual size of matrix: " << m_matrix.size() << "\n";
    exit(EXIT_FAILURE);
  }
}

negele::matrix::SquareMatrix negele::matrix::SquareMatrix::from_diag(
    const std::vector<double>& diag_mels) {
  int dim = diag_mels.size();
  std::vector<double> mels(dim * dim, 0.0);
  for (int i = 0; i < dim; i++) {
    mels[i * dim + i] = diag_mels[i];
  }

  return negele::matrix::SquareMatrix(dim, mels);
}

negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::operator+=(
    const negele::matrix::SquareMatrix& other) {
  if (m_dim != other.m_dim) {
    std::cout << "Matrices have incompatible dimensions: " << m_dim << ", and "
              << other.m_dim << "\n";
    exit(EXIT_FAILURE);
  }
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_matrix[i] += other.m_matrix[i];
  }
  return *this;
}

negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::operator-=(
    const negele::matrix::SquareMatrix& other) {
  if (m_dim != other.m_dim) {
    std::cout << "Matrices have incompatible dimensions: " << m_dim << ", and "
              << other.m_dim << "\n";
    exit(EXIT_FAILURE);
  }
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_matrix[i] -= other.m_matrix[i];
  }
  return *this;
}
negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::operator*=(
    double factor) {
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_matrix[i] *= factor;
  }
  return *this;
}
negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::operator/=(
    double factor) {
  if (factor == 0.0) {
    std::cout << "Attempted division by 0.\n";
    exit(EXIT_FAILURE);
  }
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_matrix[i] /= factor;
  }
  return *this;
}

negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::operator/=(
    const negele::matrix::SquareMatrix& other) {
  if (m_dim != other.m_dim) {
    std::cout << "Matrices have incompatible dimensions: " << m_dim << ", and "
              << other.m_dim << "\n";
    exit(EXIT_FAILURE);
  }
  for (int i = 0; i < m_dim * m_dim; i++) {
    if (other.m_matrix[i] != 0.0) {
      m_matrix[i] /= other.m_matrix[i];
    }
  }
  return *this;
}

negele::matrix::SquareMatrix& negele::matrix::SquareMatrix::abs() {
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_matrix[i] = fabs(m_matrix[i]);
  }
  return *this;
}

double negele::matrix::SquareMatrix::max() const {
  double val = 0.0;
  if (m_dim > 0) {
    val = m_matrix[0];
  }
  for (int i = 0; i < m_dim * m_dim; i++) {
    if (m_matrix[i] > val) {
      val = m_matrix[i];
    }
  }
  return val;
}

negele::matrix::SquareMatrix negele::matrix::operator+(
    const SquareMatrix& a, const negele::matrix::SquareMatrix& b) {
  negele::matrix::SquareMatrix c(a);
  c += b;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::operator-(
    const negele::matrix::SquareMatrix& a,
    const negele::matrix::SquareMatrix& b) {
  negele::matrix::SquareMatrix c(a);
  c -= b;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::operator*(
    double factor, const negele::matrix::SquareMatrix& mat) {
  negele::matrix::SquareMatrix c(mat);
  c *= factor;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::operator*(
    const negele::matrix::SquareMatrix& mat, double factor) {
  negele::matrix::SquareMatrix c(mat);
  c *= factor;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::operator/(
    const negele::matrix::SquareMatrix& mat, double factor) {
  negele::matrix::SquareMatrix c(mat);
  c /= factor;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::operator/(
    const negele::matrix::SquareMatrix& a,
    const negele::matrix::SquareMatrix& b) {
  negele::matrix::SquareMatrix c(a);
  c /= b;
  return c;
}
negele::matrix::SquareMatrix negele::matrix::abs(
    const negele::matrix::SquareMatrix& a) {
  negele::matrix::SquareMatrix c(a);
  c.abs();
  return c;
}

negele::matrix::SquareMatrix negele::matrix::operator*(
    const negele::matrix::SquareMatrix& a,
    const negele::matrix::SquareMatrix& b) {
  if (a.m_dim != b.m_dim) {
    std::cout << "Matrices have incompatible dimensions: " << a.m_dim
              << ", and " << b.m_dim << "\n";
    exit(EXIT_FAILURE);
  }
  auto new_mels = negele::matrix::mat_mul(a.m_matrix, b.m_matrix, a.m_dim);
  return negele::matrix::SquareMatrix(a.m_dim, new_mels);
}

std::vector<double> negele::matrix::SquareMatrix::eigenvalues() const {
  std::vector<double> evs(m_dim, 0.0);
  std::vector<double> matrix_copy(m_matrix);
  std::vector<double> work(20 * m_dim, 0.0);
  int info = 0;
  char uplo = 'U';
  char jobz = 'N';
  int lwork = work.size();
  LAPACK_dsyev(&jobz, &uplo, &m_dim, matrix_copy.data(), &m_dim, evs.data(),
               work.data(), &lwork, &info);

  return evs;
}

negele::matrix::MomentumSpaceMatrix::MomentumSpaceMatrix(
    const std::vector<double>& pts, const std::vector<double>& weights,
    const std::vector<double>& mels)
    : m_dim(pts.size()),
      m_pts(pts),
      m_weights(weights),
      m_mels_with_weights(pts.size(), 0.0),
      m_kin_e(pts.size(), 0.0) {
  std::vector<double> weights_sqrt;
  for (const auto& w : m_weights) {
    weights_sqrt.push_back(sqrt(w));
  }

  auto weights_mat = negele::matrix::SquareMatrix::from_diag(weights_sqrt);
  negele::matrix::SquareMatrix mels_wo_weights(m_dim, mels);

  m_mels_with_weights = weights_mat * mels_wo_weights * weights_mat;

  std::vector<double> kin_diag;
  for (const auto& p : pts) {
    kin_diag.push_back(p * p * HBARC * HBARC / (2 * MASS));
  }
  m_kin_e = negele::matrix::SquareMatrix::from_diag(kin_diag);
}

void negele::matrix::MomentumSpaceMatrix::clear_mels() {
  m_mels_with_weights *= 0.0;
}
double negele::matrix::MomentumSpaceMatrix::max() const {
  return negele::matrix::abs(m_mels_with_weights).max();
}

void negele::matrix::MomentumSpaceMatrix::resize(
    const negele::matrix::MomentumSpaceMatrix& other) {
  m_dim = other.m_dim;
  m_pts = std::vector<double>(other.m_pts);
  m_weights = std::vector<double>(other.m_weights);

  m_mels_with_weights = negele::matrix::SquareMatrix(m_dim, 0.0);
  m_kin_e = negele::matrix::SquareMatrix(m_dim, 0.0);
}

bool negele::matrix::MomentumSpaceMatrix::same_size(
    const negele::matrix::MomentumSpaceMatrix& other) const {
  if (m_dim != other.m_dim) {
    return false;
  }
  for (int i = 0; i < m_dim; i++) {
    if (m_pts[i] != other.m_pts[i]) {
      return false;
    }
    if (m_weights[i] != other.m_weights[i]) {
      return false;
    }
  }
  if (m_mels_with_weights.size() != other.m_mels_with_weights.size()) {
    return false;
  }
  if (m_kin_e.size() != other.m_kin_e.size()) {
    return false;
  }
  return true;
}

negele::matrix::SquareMatrix
negele::matrix::MomentumSpaceMatrix::matrix_elements() const {
  std::vector<double> weights_inv_sqrt;
  for (const auto& w : m_weights) {
    weights_inv_sqrt.push_back(1.0 / sqrt(w));
  }
  auto weights_inv_mat =
      negele::matrix::SquareMatrix::from_diag(weights_inv_sqrt);
  return weights_inv_mat * m_mels_with_weights * weights_inv_mat;
}

std::vector<double> negele::matrix::mat_mul(const std::vector<double>& a,
                                            const std::vector<double>& b,
                                            int dim) {
  std::vector<double> new_matrix(dim * dim, 0.0);

  for (int i = 0; i < dim; i++) {
    for (int j = 0; j < dim; j++) {
      for (int k = 0; k < dim; k++) {
        new_matrix[i * dim + j] += a[i * dim + k] * b[k * dim + j];
      }
    }
  }

  return new_matrix;
}

negele::matrix::MomentumSpaceMatrix& negele::matrix::MomentumSpaceMatrix::
operator+=(const negele::matrix::MomentumSpaceMatrix& other) {
  m_mels_with_weights += other.m_mels_with_weights;
  return *this;
}

negele::matrix::MomentumSpaceMatrix& negele::matrix::MomentumSpaceMatrix::
operator+=(const double other) {
  for (int i = 0; i < m_dim * m_dim; i++) {
    m_mels_with_weights[i] += copysign(other, m_mels_with_weights[i]);
  }
  return *this;
}

negele::matrix::MomentumSpaceMatrix& negele::matrix::MomentumSpaceMatrix::
operator*=(const double factor) {
  m_mels_with_weights *= factor;
  return *this;
}

negele::matrix::MomentumSpaceMatrix negele::matrix::operator/(
    const negele::matrix::MomentumSpaceMatrix& a,
    const negele::matrix::MomentumSpaceMatrix& b) {
  negele::matrix::MomentumSpaceMatrix c(a);
  c.m_mels_with_weights /= b.m_mels_with_weights;
  return c;
}

negele::matrix::MomentumSpaceMatrix negele::matrix::abs(
    const negele::matrix::MomentumSpaceMatrix& a) {
  negele::matrix::MomentumSpaceMatrix c(a);
  c.m_mels_with_weights.abs();
  return c;
}

ext: cpp

max_stars_repo_name: cheshyre/negele-srg-solver-cpp

max_issues_repo_name: cheshyre/negele-srg-solver-cpp

max_forks_repo_name: cheshyre/negele-srg-solver-cpp


text:
/*==============================================================================
    Copyright (c) 2001-2010 Joel de Guzman
    Copyright (c) 2010 Thomas Heller

    Distributed under the Boost Software License, Version 1.0. (See accompanying
    file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
==============================================================================*/
#ifndef BOOST_PHOENIX_STATEMENT_WHILE_HPP
#define BOOST_PHOENIX_STATEMENT_WHILE_HPP

#include <boost/phoenix/core/limits.hpp>
#include <boost/phoenix/core/call.hpp>
#include <boost/phoenix/core/expression.hpp>
#include <boost/phoenix/core/meta_grammar.hpp>

BOOST_PHOENIX_DEFINE_EXPRESSION(
    (boost)(phoenix)(while_)
  , (meta_grammar) // Cond
    (meta_grammar) // Do
)

namespace boost { namespace phoenix
{
    struct while_eval
    {
        typedef void result_type;

        template <typename Cond, typename Do, typename Context>
        result_type
        operator()(Cond const& cond, Do const& do_, Context & ctx) const
        {
            while(boost::phoenix::eval(cond, ctx))
            {
                boost::phoenix::eval(do_, ctx);
            }
        }
    };

    template <typename Dummy>
    struct default_actions::when<rule::while_, Dummy>
        : call<while_eval, Dummy>
    {};

    template <typename Cond>
    struct while_gen
    {
        while_gen(Cond const& cond) : cond(cond) {}

        template <typename Do>
        typename expression::while_<Cond, Do>::type const
        operator[](Do const& do_) const
        {
            return expression::while_<Cond, Do>::make(cond, do_);
        }

        Cond const& cond;
    };

    template <typename Cond>
    inline
    while_gen<Cond> const
    while_(Cond const& cond)
    {
        return while_gen<Cond>(cond);
    }

}}

#endif

max_stars_repo_name: brinkqiang/dmspirit


text:
#include <crs_perception/model_to_point_cloud.hpp>

#include <pcl/filters/voxel_grid.h>
#include <vtkCellArray.h>
#include <vtkPolyDataMapper.h>
#include <vtkTriangle.h>
#include <vtkTriangleFilter.h>
#include <vtkPLYReader.h>
#include <vtkSTLReader.h>
#include <vtkOBJReader.h>

#include <map>
#include <boost/format.hpp>
#include <boost/filesystem.hpp>

enum class Extensions : int
{
  STL = 0,
  PLY,
  OBJ
};
const static std::map<std::string, Extensions> EXTENSIONS = { { ".stl", Extensions::STL },
                                                              { ".ply", Extensions::PLY },
                                                              { ".obj", Extensions::OBJ } };

namespace crs_perception
{
bool ModelToPointCloud::convertToPCL(const std::string& file_path,
                                     pcl::PointCloud<pcl::PointXYZ>::Ptr point_cloud,
                                     std::string& err_msg)
{
  namespace fs = boost::filesystem;

  fs::path file_path_obj(file_path);
  std::string ext = fs::path(file_path).extension().string();
  if (!std::any_of(EXTENSIONS.begin(), EXTENSIONS.end(), [ext](const decltype(EXTENSIONS)::value_type& kv) {
        return kv.first == ext;
      }))
  {
    err_msg = boost::str(boost::format("The extension %s is not supported") % ext);
    return false;
  }

  // todo(ayoungs): handle multiple file formats?
  // todo(ayoungs): handle error for bad file
  vtkSmartPointer<vtkPolyData> polydata1 = vtkSmartPointer<vtkPolyData>::New();
  switch (EXTENSIONS.at(ext))
  {
    case Extensions::STL:
    {
      vtkSmartPointer<vtkSTLReader> readerQuery = vtkSmartPointer<vtkSTLReader>::New();
      readerQuery->SetFileName(file_path.c_str());
      readerQuery->Update();
      polydata1 = readerQuery->GetOutput();
    }
    break;

    case Extensions::PLY:
    {
      vtkSmartPointer<vtkPLYReader> reader = vtkSmartPointer<vtkPLYReader>::New();
      reader->SetFileName(file_path.c_str());
      reader->Update();
      polydata1 = reader->GetOutput();
    }
    break;

    case Extensions::OBJ:
    {
      vtkSmartPointer<vtkOBJReader> reader = vtkSmartPointer<vtkOBJReader>::New();
      reader->SetFileName(file_path.c_str());
      reader->Update();
      polydata1 = reader->GetOutput();
    }
    break;
  }

  // make sure that the polygons are triangles!
  vtkSmartPointer<vtkTriangleFilter> triangleFilter = vtkSmartPointer<vtkTriangleFilter>::New();
  triangleFilter->SetInputData(polydata1);
  triangleFilter->Update();

  vtkSmartPointer<vtkPolyDataMapper> triangleMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
  triangleMapper->SetInputConnection(triangleFilter->GetOutputPort());
  triangleMapper->Update();
  polydata1 = triangleMapper->GetInput();

  pcl::PointCloud<pcl::PointXYZ>::Ptr cloud = boost::make_shared<pcl::PointCloud<pcl::PointXYZ> >();
  uniformSampling(polydata1, num_samples_, *cloud);

  // Voxelgrid
  pcl::VoxelGrid<pcl::PointXYZ> grid_;
  grid_.setInputCloud(cloud);
  grid_.setLeafSize(leaf_size_, leaf_size_, leaf_size_);

  grid_.filter(*point_cloud);

  return true;
}

void ModelToPointCloud::uniformSampling(vtkSmartPointer<vtkPolyData> polydata,
                                        std::size_t n_samples,
                                        pcl::PointCloud<pcl::PointXYZ>& cloud_out)
{
  polydata->BuildCells();
  vtkSmartPointer<vtkCellArray> cells = polydata->GetPolys();

  double p1[3], p2[3], p3[3], totalArea = 0;
  std::vector<double> cumulativeAreas(cells->GetNumberOfCells(), 0);
  vtkIdType npts = 0, *ptIds = nullptr;
  std::size_t cellId = 0;
  for (cells->InitTraversal(); cells->GetNextCell(npts, ptIds); cellId++)
  {
    polydata->GetPoint(ptIds[0], p1);
    polydata->GetPoint(ptIds[1], p2);
    polydata->GetPoint(ptIds[2], p3);
    totalArea += vtkTriangle::TriangleArea(p1, p2, p3);
    cumulativeAreas[cellId] = totalArea;
  }

  cloud_out.points.resize(n_samples);
  cloud_out.width = static_cast<std::uint32_t>(n_samples);
  cloud_out.height = 1;

  for (std::size_t i = 0; i < n_samples; i++)
  {
    Eigen::Vector3f p;
    Eigen::Vector3f n(0, 0, 0);
    Eigen::Vector3f c(0, 0, 0);
    randPSurface(polydata, &cumulativeAreas, totalArea, p);
    cloud_out.points[i].x = p[0];
    cloud_out.points[i].y = p[1];
    cloud_out.points[i].z = p[2];
  }
}

inline void ModelToPointCloud::randPSurface(vtkSmartPointer<vtkPolyData> polydata,
                                            std::vector<double>* cumulativeAreas,
                                            double totalArea,
                                            Eigen::Vector3f& p)
{
  float r = static_cast<float>(uniformDeviate(rand()) * totalArea);

  std::vector<double>::iterator low = std::lower_bound(cumulativeAreas->begin(), cumulativeAreas->end(), r);
  vtkIdType el = vtkIdType(low - cumulativeAreas->begin());

  double A[3], B[3], C[3];
  vtkIdType npts = 0;
  vtkIdType* ptIds = nullptr;
  polydata->GetCellPoints(el, npts, ptIds);
  polydata->GetPoint(ptIds[0], A);
  polydata->GetPoint(ptIds[1], B);
  polydata->GetPoint(ptIds[2], C);
  float r1 = static_cast<float>(uniformDeviate(rand()));
  float r2 = static_cast<float>(uniformDeviate(rand()));
  randomPointTriangle(float(A[0]),
                      float(A[1]),
                      float(A[2]),
                      float(B[0]),
                      float(B[1]),
                      float(B[2]),
                      float(C[0]),
                      float(C[1]),
                      float(C[2]),
                      r1,
                      r2,
                      p);
}

inline void ModelToPointCloud::randomPointTriangle(float a1,
                                                   float a2,
                                                   float a3,
                                                   float b1,
                                                   float b2,
                                                   float b3,
                                                   float c1,
                                                   float c2,
                                                   float c3,
                                                   float r1,
                                                   float r2,
                                                   Eigen::Vector3f& p)
{
  float r1sqr = std::sqrt(r1);
  float OneMinR1Sqr = (1 - r1sqr);
  float OneMinR2 = (1 - r2);
  a1 *= OneMinR1Sqr;
  a2 *= OneMinR1Sqr;
  a3 *= OneMinR1Sqr;
  b1 *= OneMinR2;
  b2 *= OneMinR2;
  b3 *= OneMinR2;
  c1 = r1sqr * (r2 * c1 + b1) + a1;
  c2 = r1sqr * (r2 * c2 + b2) + a2;
  c3 = r1sqr * (r2 * c3 + b3) + a3;
  p[0] = c1;
  p[1] = c2;
  p[2] = c3;
}

inline double ModelToPointCloud::uniformDeviate(int seed)
{
  double ran = seed * (1.0 / (RAND_MAX + 1.0));
  return ran;
}

}  // namespace crs_perception


text:
/// \file   test_quantity.cpp
///
/// \brief
///
/// \authors    Maarten P. Scholl
/// \date       2019-08-15
/// \copyright  Copyright 2017-2019 The Institute for New Economic Thinking,
///             Oxford Martin School, University of Oxford
///
///             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.
///
///             You may obtain instructions to fulfill the attribution
///             requirements in CITATION.cff
///
#define BOOST_TEST_DYN_LINK
#define BOOST_TEST_MODULE quantity

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

#include <esl/quantity.hpp>

BOOST_AUTO_TEST_SUITE(ESL)

BOOST_AUTO_TEST_CASE(quantity_constructor)
{
    {
        esl::quantity q_default_;

        BOOST_CHECK_EQUAL(q_default_.amount, 0);
        BOOST_CHECK_EQUAL(q_default_.basis, 1);
    }

    BOOST_CHECK_THROW(esl::quantity(13, 0), std::logic_error);

    {
        esl::quantity one_(1);

        BOOST_CHECK_EQUAL(one_.amount, 1);
        BOOST_CHECK_EQUAL(one_.basis, 1);
    }

    {
        esl::quantity threefifty_(350, 100);

        BOOST_CHECK_EQUAL(threefifty_.amount, 350);
        BOOST_CHECK_EQUAL(threefifty_.basis, 100);

        esl::quantity again_(threefifty_);

        BOOST_CHECK_EQUAL(again_.amount, 350);
        BOOST_CHECK_EQUAL(again_.basis, 100);
    }
}

BOOST_AUTO_TEST_CASE(quantity_arithmetic_guards)
{
    // test exception on adding incompatible quantities
    esl::quantity threefifty_(350, 100);
#if EXCEPTION_ON_BASIS_CHANGE
    BOOST_CHECK_THROW(threefifty_ + esl::quantity(1, 42), std::logic_error);
#else
    BOOST_CHECK_EQUAL(threefifty_ + esl::quantity(1, 42),
                      esl::quantity(148, 42));
#endif

    // test exception on subtracting with negative result
    BOOST_CHECK_THROW(threefifty_ - esl::quantity(351, 100), std::logic_error);
    BOOST_CHECK_THROW(threefifty_ - esl::quantity(176, 50), std::logic_error);
    BOOST_CHECK_THROW(threefifty_ - esl::quantity(701, 200), std::logic_error);

    // test conversion on subtraction
#if EXCEPTION_ON_BASIS_CHANGE
    BOOST_CHECK_THROW(threefifty_ - esl::quantity(1, 42), std::logic_error);
#else
    BOOST_CHECK_EQUAL(threefifty_ - esl::quantity(1, 42),
                      esl::quantity(146, 42));

    BOOST_CHECK_EQUAL(threefifty_ - esl::quantity(3, 1),
                      esl::quantity(50, 100));
#endif

#if EXCEPTION_ON_BASIS_CHANGE
    {
        esl::quantity again_(threefifty_);
        BOOST_CHECK_THROW(again_ += esl::quantity(1, 42), std::logic_error);
    }
    {
        esl::quantity again_(threefifty_);
        BOOST_CHECK_THROW(again_ -= esl::quantity(1, 42), std::logic_error);
    }
#endif
}

BOOST_AUTO_TEST_CASE(quantity_arithmetic)
{
    esl::quantity threefifty_(350, 100);

    // .divide integer
    {
        auto constexpr_d_ = threefifty_.divide<4>();
        BOOST_CHECK_EQUAL(constexpr_d_[0], esl::quantity(88, 100));
        BOOST_CHECK_EQUAL(constexpr_d_[1], esl::quantity(88, 100));
        BOOST_CHECK_EQUAL(constexpr_d_[2], esl::quantity(87, 100));
        BOOST_CHECK_EQUAL(constexpr_d_[3], esl::quantity(87, 100));
    }

    // / integer
    {
        auto dynamic_d_ = threefifty_ / 3;
        auto q_117_     = esl::quantity(117, 100);
        auto q_116_     = esl::quantity(116, 100);
        BOOST_CHECK_EQUAL(dynamic_d_[0], q_117_);
        BOOST_CHECK_EQUAL(dynamic_d_[1], q_117_);
        BOOST_CHECK_EQUAL(dynamic_d_[2], q_116_);
    }

    // +
    BOOST_CHECK_EQUAL(threefifty_ + esl::quantity(123, 100),
                      esl::quantity(473, 100));
    BOOST_CHECK_EQUAL(threefifty_ + esl::quantity(0, 100),
                      esl::quantity(350, 100));
    // -
    BOOST_CHECK_EQUAL(threefifty_ - esl::quantity(123, 100),
                      esl::quantity(227, 100));
    // * integer
    BOOST_CHECK_EQUAL(threefifty_ * 3, esl::quantity(1050, 100));
    // * quantity
    BOOST_CHECK_EQUAL(threefifty_ * esl::quantity(50, 100),
                      esl::quantity(175, 100));
    BOOST_CHECK_EQUAL(threefifty_ * esl::quantity(200, 100),
                      esl::quantity(700, 100));
    BOOST_CHECK_EQUAL(threefifty_ * esl::quantity(0, 100),
                      esl::quantity(0, 100));
    BOOST_CHECK_EQUAL(threefifty_ * esl::quantity(100, 100),
                      esl::quantity(350, 100));
    BOOST_CHECK_EQUAL(threefifty_ * esl::quantity(10, 100),
                      esl::quantity(35, 100));
}

BOOST_AUTO_TEST_CASE(quantity_logic)
{
    esl::quantity threefifty_(350, 100);

    BOOST_CHECK_EQUAL(threefifty_, esl::quantity(350, 100));
    BOOST_CHECK_NE(threefifty_, esl::quantity(123, 100));

    BOOST_CHECK_LE(threefifty_, esl::quantity(456, 100));
    BOOST_CHECK_GE(threefifty_, esl::quantity(123, 100));

    BOOST_CHECK_EQUAL(threefifty_, 3.5);

    BOOST_CHECK_EQUAL(esl::quantity(142'857'142'857, 100'000'000'000),
                      10.0 / 7);
}

BOOST_AUTO_TEST_SUITE_END()  // ESL


text:
#include <string>

#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>

#include <boost/foreach.hpp>

typedef CGAL::Simple_cartesian<double> K;
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
typedef Mesh::Vertex_index vertex_descriptor;
typedef Mesh::Face_index face_descriptor;

int main()
{

  Mesh m;
  vertex_descriptor v0 = m.add_vertex(K::Point_3(0,2,0));
  vertex_descriptor v1 = m.add_vertex(K::Point_3(2,2,0));
  vertex_descriptor v2 = m.add_vertex(K::Point_3(0,0,0));
  vertex_descriptor v3 = m.add_vertex(K::Point_3(2,0,0));
  vertex_descriptor v4 = m.add_vertex(K::Point_3(1,1,0));
  m.add_face(v3, v1, v4);
  m.add_face(v0, v4, v1);
  m.add_face(v0, v2, v4);
  m.add_face(v2, v3, v4);

  // give each vertex a name, the default is empty
  Mesh::Property_map<vertex_descriptor,std::string> name;
  bool created;
  boost::tie(name, created) = m.add_property_map<vertex_descriptor,std::string>("v:name","");
  assert(created);
  // add some names to the vertices
  name[v0] = "hello";
  name[v2] = "world";

  {
    // You get an existing property, and created will be false
    Mesh::Property_map<vertex_descriptor,std::string> name;
    bool created;
    boost::tie(name, created) = m.add_property_map<vertex_descriptor,std::string>("v:name", "");
    assert(! created);
  }

  //  You can't get a property that does not exist
  Mesh::Property_map<face_descriptor,std::string> gnus;
  bool found;
  boost::tie(gnus, found) = m.property_map<face_descriptor,std::string>("v:gnus");
  assert(! found);

  // retrieve the point property for which exists a convenience function
  Mesh::Property_map<vertex_descriptor, K::Point_3> location = m.points();
  BOOST_FOREACH( vertex_descriptor vd, m.vertices()) {
    std::cout << name[vd] << " @ " << location[vd] << std::endl;
  }

  // delete the string property again
  m.remove_property_map(name);

  return 0;
}

max_stars_repo_name: hlzz/dotfiles

max_issues_repo_name: hlzz/dotfiles

max_forks_repo_name: hlzz/dotfiles


text:
/**@file	tests/UT_NOK_missing_mandatory-4.cpp
 * @author	Luc Hermitte <EMAIL:luc{dot}hermitte{at}gmail{dot}com>
 *
 * Distributed under the Boost Software License, Version 1.0.
 * See accompanying file LICENSE_1_0.txt or copy at
 * http://www.boost.org/LICENSE_1_0.txt
 */
#define BOOST_TEST_MODULE "Named parameters related tests"
#define BOOST_TEST_DYN_LINK

#include "named-parameters.hpp"
#include <boost/test/unit_test.hpp>
#include <vector>

using namespace na::literals;

// ===[ copy_no_default1 ]==================================== {{{1
template <typename ...Args>
void f_copy_with_default4(Args&& ...args)
{
    auto a = na::get("a"_na = 1    , std::forward<Args>(args)...);
    auto b = na::get("b"_na        , std::forward<Args>(args)...);
    auto c = na::get("c"_na = "bar", std::forward<Args>(args)...);
    auto d = na::get("d"_na = 24.5 , std::forward<Args>(args)...);
}

BOOST_AUTO_TEST_CASE(copy_no_default5)
{
    auto       v = std::vector<int>{ 1, 2, 3, 4, 5};
    auto const x = std::vector<int>{ 11, 12, 13, 14, 15};
    f_copy_with_default4(
            "a"_na=42,
            // "b"_na=v,
            "c"_na="foo",
            "d"_na=12.5
            );
    BOOST_CHECK_EQUAL(v.size(), 5); // check v was copied
}

// =========================================================== }}}1
// vim:set fdm=marker:

max_stars_repo_name: LucHermitte/NamedParameter

max_issues_repo_name: LucHermitte/NamedParameter

max_forks_repo_name: LucHermitte/NamedParameter


text:
/**
 *  @file
 *  @author     2012-2013 Stefan Radomski (stefan.radomski@cs.tu-darmstadt.de)
 *  @copyright  Simplified BSD
 *
 *  @cond
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the FreeBSD license as published by the FreeBSD
 *  project.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 *
 *  You should have received a copy of the FreeBSD license along with this
 *  program. If not, see <http://www.opensource.org/licenses/bsd-license>.
 *  @endcond
 */

#include <boost/algorithm/string.hpp>

#include "uscxml/Common.h"
#include "NullDataModel.h"

#ifdef BUILD_AS_PLUGINS
#include <Pluma/Connector.hpp>
#endif

namespace uscxml {

#ifdef BUILD_AS_PLUGINS
PLUMA_CONNECTOR
bool pluginConnect(pluma::Host& host) {
	host.add( new NullDataModelProvider() );
	return true;
}
#endif

NullDataModel::NullDataModel() {
}

std::shared_ptr<DataModelImpl> NullDataModel::create(DataModelCallbacks* callbacks) {
	std::shared_ptr<NullDataModel> dm(new NullDataModel());
	dm->_callbacks = callbacks;
	dm->setup();
	return dm;
}

void NullDataModel::setup() {
}

NullDataModel::~NullDataModel() {
}

Data NullDataModel::getAsData(const std::string& content) {
	Data data = Data::fromJSON(content);
	if (data.empty()) {
		data = Data(content, Data::VERBATIM);
	}
	return data;
}

/**
 * The boolean expression language consists of the In predicate only. It has the
 * form 'In(id)', where id is the id of a state in the enclosing state machine.
 * The predicate must return 'true' if and only if that state is in the current
 * state configuration.
 */
bool NullDataModel::evalAsBool(const XERCESC_NS::DOMElement* scriptNode, const std::string& expr) {
	std::string trimmedExpr = expr;
	boost::trim(trimmedExpr);
	if (!boost::istarts_with(trimmedExpr, "in"))
		return false;

	// find string in between brackets
	size_t start = trimmedExpr.find_first_of("(");
	size_t end = trimmedExpr.find_last_of(")");
	if (start == std::string::npos || end == std::string::npos || start >= end)
		return false;
	start++;

	// split at comma
	std::stringstream ss(trimmedExpr.substr(start, end - start));
	std::list<std::string> stateExprs;
	std::string item;
	while(std::getline(ss, item, ',')) {
		stateExprs.push_back(item);
	}

	for (std::list<std::string>::const_iterator stateIter = stateExprs.begin(); stateIter != stateExprs.end(); stateIter++) {
		// remove ticks
		size_t start = stateIter->find_first_of("'");
		size_t end = stateIter->find_last_of("'");

		std::string stateName;
		if (start != std::string::npos && end != std::string::npos && start < end) {
			start++;
			stateName = stateIter->substr(start, end - start);
		} else {
			stateName = *stateIter;
		}

		if (_callbacks->isInState(stateName)) {
			continue;
		}
		return false;
	}
	return true;
}

}

max_stars_repo_name: stil/uscxml

max_issues_repo_name: stil/uscxml

max_forks_repo_name: stil/uscxml


text:
#include "program_options_parser.hpp"

#include <fstream>
#include <iostream>
#include <limits>
#include <map>
#include <vector>
#include <iomanip>

#include "logger/logger.hpp"
#include "program_options.hpp"
#include "time_stepping/time_integrator.hpp"
#include "utils/registry.hpp"

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

namespace po = boost::program_options;

template <typename T>
std::ostream& operator<<(std::ostream& os, const std::vector<T>& vec)
{
  for (size_t i = 0; i < vec.size(); ++i) {
    os << vec[i] << " ";
  }
  return os;
}

namespace std
{
std::ostream& operator<<(std::ostream& os, const std::vector<double>& vec)
{
  for (auto i : vec) {
    os << i << " ";
  }
  return os;
}
}

namespace std
{
std::ostream& operator<<(std::ostream& os, const std::vector<std::string>& vec)
{
  for (auto i : vec) {
    os << i << " ";
  }
  return os;
}
}

std::istream& operator>>(std::istream& in, LogLevel& level)
{
  std::string token;
  in >> token;
  if (boost::iequals(token, "trace")) {
    level = LogLevel::trace;
  } else if (boost::iequals(token, "debug")) {
    level = LogLevel::debug;
  } else if (boost::iequals(token, "info")) {
    level = LogLevel::info;
  } else if (boost::iequals(token, "warning")) {
    level = LogLevel::warning;
  } else if (boost::iequals(token, "error")) {
    level = LogLevel::error;
  } else if (boost::iequals(token, "fatal")) {
    level = LogLevel::fatal;
  } else {
    in.setstate(std::ios_base::failbit);
  }
  return in;
}

ProgramOptionsParser::ProgramOptionsParser()
{
  base.add_options()("help", "Produces this help message.");

  file.add_options()(
      "config_file",
      po::value<std::string>(&Options::get().config_file_),
      "Input configuration options file. Any command line options will override config file "
      "settings.")(
      "solution_output_file",
      po::value<std::string>(&Options::get().solution_output_file_),
      "Solution output file name, without file extension.")(
      "initial_condition_file",
      po::value<std::string>(&Options::get().initial_condition_file_)->default_value("none"),
      "Filename containing and initial solution. Overrides input values if not 'none'.")(
      "log_file",
      po::value<std::string>(&Options::get().log_file_)->default_value("maDGiCart"),
      "Solution output file name, without file extension.")(
      "log_frequency",
      po::value<int>(&Options::get().log_frequency_)->default_value(20),
      "Logger update frequency.")(
      "save_every",
      po::value<int>(&Options::get().save_every_)->default_value(0),
      "If not zero, save a solution every N time step.");

  physics.add_options()(
      "equation_type",
      po::value<std::string>(&Options::get().equation_type_)->default_value("cahn-hilliard"),
      "Governing equations to solve.")(
      "m",
      po::value<double>(&Options::get().ch_m_)->default_value(0.0),
      "Macroscopic average volume fraction of concentration.")(
      "eps2",
      po::value<double>(&Options::get().ch_eps2_)->default_value(2.419753086419753e-05),
      "epsilon^2, the biharmonic coefficient.")(
      "sigma",
      po::value<double>(&Options::get().ch_sigma_)->default_value(291.4488109293305),
      "sigma, the coefficient of the linearly stabilizing term about m.");

  discretization.add_options()(
      "domain_x_begin",
      po::value<double>(&Options::get().domain_x_begin_)->default_value(0),
      "Minimum x value for Cartesian domain")(
      "domain_x_end",
      po::value<double>(&Options::get().domain_x_end_)->default_value(1),
      "Maximum x value for Cartesian domain")(
      "domain_resolution",
      po::value<int>(&Options::get().domain_resolution_)->default_value(128),
      "Cartesian mesh with [domain_resolution]^N nodes")(
      "dimension",
      po::value<int>(&Options::get().dimension_)->default_value(2),
      "Spatial dimension, 2 or 3.")(
      "bc_x",
      po::value<std::string>(&Options::get().bc_x_)->default_value("periodic"),
      "Boundary condition in x-direction, periodic or neumann. Only for 3D.")(
      "bc_y",
      po::value<std::string>(&Options::get().bc_y_)->default_value("periodic"),
      "Boundary condition in y-direction, periodic or neumann. Only for 3D.")(
      "bc_z",
      po::value<std::string>(&Options::get().bc_z_)->default_value("periodic"),
      "Boundary condition in z-direction, periodic or neumann. Only for 3D.");

  time_stepping.add_options()(
      "time_integrator",
      po::value<std::string>(&Options::get().time_integrator_)->default_value("rk3ssp"),
      "Time integration method.")(
      "time_step_size",
      po::value<double>(&Options::get().time_step_size_)->default_value(1.e-6),
      "Unsteady time step size. For adaptive methods this is the initial step size.")(
      "initial_time",
      po::value<double>(&Options::get().initial_time_)->default_value(0),
      "Initial time for unsteady.")(
      "final_time",
      po::value<double>(&Options::get().final_time_)->default_value(1),
      "Stopping time for unsteady.")(
      "max_time_steps",
      po::value<unsigned>(&Options::get().max_time_steps_)->default_value(10000),
      "Maximum number of unsteady time steps.")(
      "use_cfl_time_step",
      po::value<bool>(&Options::get().use_cfl_time_step_)->default_value(false),
      "Use CFL-based time stepping instead of dimensional time step.")(
      "cfl",
      po::value<double>(&Options::get().cfl_)->default_value(1.0),
      "CFL number for CFL-based time stepping.")(
      "use_adaptive_time_step",
      po::value<bool>(&Options::get().use_adaptive_time_step_)->default_value(false),
      "Use error-adaptive time stepping.")(
      "min_time_step_size",
      po::value<double>(&Options::get().min_time_step_size_)
          ->default_value(std::numeric_limits<double>::epsilon()*100),
      "For adaptive time stepping, the minimum time step size.")(
      "max_time_step_size",
      po::value<double>(&Options::get().max_time_step_size_)
          ->default_value(std::numeric_limits<double>::max()/100),
      "For adaptive time stepping, the maximum time step size.")(
      "time_rel_err_tol",
      po::value<double>(&Options::get().time_rel_err_tol_)->default_value(1e-2),
      "Relative error tolerance for adaptive time stepping.")(
      "converged_rel_tol",
      po::value<double>(&Options::get().converged_rel_tol_)->default_value(0),
      "Relative residual to consider solution converged. 0 indicates use max time or max timesteps.")(
      "converged_abs_tol",
      po::value<double>(&Options::get().converged_abs_tol_)->default_value(0),
      "Absolute residual to consider solution converged. 0 indicates use max time or max timesteps.");

}

void ProgramOptionsParser::parseInputOptions(const std::vector<std::string>& cmd_line)
{
  vm.clear();

  base.add(file);
  base.add(discretization);
  base.add(physics);
  base.add(time_stepping);
  visible.add(base);

  po::positional_options_description p;
  p.add("config_file", -1);
  po::store(
      po::command_line_parser(cmd_line)
          .options(base)
          .style(po::command_line_style::unix_style ^ po::command_line_style::allow_short)
          .positional(p)
          .run(),
      vm);

  po::notify(vm);
  std::ifstream ifs(Options::get().config_file_.c_str());
  po::store(po::parse_config_file(ifs, base), vm);
  po::notify(vm);

  if (vm.count("help")) {
    std::cout << visible << "\n";
    printRegistryContents(std::cout);
    std::stringstream ss;
    std::cout << getConfigurationFileTemplate();
    exit(0);
    return;
  }
}

std::ostream& operator<<(std::ostream& os, const boost::any& a)
{
  // must determine the type of the boost::any value stored by program options
  // in order to use the << operator.
  using boost::any_cast;
  try {
    std::string str = any_cast<std::string>(a);
    os << str;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    double dbl = any_cast<double>(a);
    os << dbl;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    int in = any_cast<int>(a);
    os << in;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    bool b = any_cast<bool>(a);
    os << b;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    unsigned uin = any_cast<unsigned>(a);
    os << uin;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    std::vector<double> vec = any_cast<std::vector<double> >(a);
    os << vec;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  try {
    std::vector<std::string> vec = any_cast<std::vector<std::string> >(a);
    os << vec;
    return os;
  } catch (boost::bad_any_cast&) {
  }
  return os;
}

void ProgramOptionsParser::printRegistryContents(std::ostream& out) const
{
  out << FactoryRegistry<TimeIntegrator>::get() << std::endl;
}

std::string ProgramOptionsParser::getConfigurationFileTemplate() const
{
  std::stringstream out;

  out.setf(std::ios_base::boolalpha);

  out << "\n# ****** File Options ****** #\n\n";
  for (const auto& i : file.options()) {
    std::string varname = i->long_name();
    out << std::setw(30) << varname << " = " << vm[varname].value() << std::endl;
  }

  out << "\n# ****** Physics Options ****** #\n\n";
  for (const auto& i : physics.options()) {
    std::string varname = i->long_name();
    out << std::setw(30) << varname << " = " << vm[varname].value() << std::endl;
  }

  out << "\n# ****** Discretization Options ****** #\n\n";
  for (const auto& i : discretization.options()) {
    std::string varname = i->long_name();
    out << std::setw(30) << varname << " = " << vm[varname].value() << std::endl;
  }

  out << "\n# ****** Time Stepping Options ****** #\n\n";
  for (const auto& i : time_stepping.options()) {
    std::string varname = i->long_name();
    out << std::setw(30) << varname << " = " << vm[varname].value() << std::endl;
  }

  out << "\n";
  return out.str();
}

max_stars_repo_name: mlohry/maDGiCart-CH

max_issues_repo_name: mlohry/maDGiCart-CH

max_forks_repo_name: mlohry/maDGiCart-CH


text:
/*
   For more information, please see: http://software.sci.utah.edu

   The MIT License

   Copyright (c) 2015 Scientific Computing and Imaging Institute,
   University of Utah.

   License for the specific language governing rights and limitations under
   Permission is hereby granted, free of charge, to any person obtaining a
   copy of this software and associated documentation files (the "Software"),
   to deal in the Software without restriction, including without limitation
   the rights to use, copy, modify, merge, publish, distribute, sublicense,
   and/or sell copies of the Software, and to permit persons to whom the
   Software is furnished to do so, subject to the following conditions:

   The above copyright notice and this permission notice shall be included
   in all copies or substantial portions of the Software.

   THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
   OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
   FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
   THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
   LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
   FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
   DEALINGS IN THE SOFTWARE.
*/

#include <Testing/ModuleTestBase/ModuleTestBase.h>
#include <Modules/Factory/HardCodedModuleFactory.h>
#include <Dataflow/Engine/Controller/NetworkEditorController.h>
#include <Modules/Legacy/Fields/CreateLatVol.h>
#include <Dataflow/State/SimpleMapModuleState.h>
#include <Dataflow/Engine/Scheduler/DesktopExecutionStrategyFactory.h>
#include <Core/Python/PythonInterpreter.h>
#include <boost/filesystem.hpp>

using namespace SCIRun;
using namespace Core;
using namespace Testing;
using namespace Modules::Factory;
using namespace Modules::Fields;
using namespace Dataflow::Networks;
using namespace Dataflow::Engine;
using namespace ReplacementImpl;
using namespace Dataflow::Engine;
using namespace Dataflow::State;
using namespace Algorithms;

class PythonControllerFunctionalTests : public ModuleTest
{
public:
  PythonControllerFunctionalTests()
  {
    PythonInterpreter::Instance().initialize(false, "Engine_Python_Tests", boost::filesystem::current_path().string());
    PythonInterpreter::Instance().importSCIRunLibrary();
  }
};

TEST_F(PythonControllerFunctionalTests, CanAddModule)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);
  //TODO: expose API directly on NEC?
  //controller.runPython("addModule(\"CreateLatVol\")");

  ASSERT_EQ(1, controller.getNetwork()->nmodules());
}

TEST_F(PythonControllerFunctionalTests, CanAddMultipleModule)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);
  PythonInterpreter::Instance().run_string(command);

  ASSERT_EQ(2, controller.getNetwork()->nmodules());
}

TEST_F(PythonControllerFunctionalTests, CanChangeModuleState)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  ModuleStateFactoryHandle sf(new SimpleMapModuleStateFactory);
  NetworkEditorController controller(mf, sf, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "m = scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);

  ASSERT_EQ(1, controller.getNetwork()->nmodules());
  auto mod = controller.getNetwork()->module(0);
  ASSERT_TRUE(mod != nullptr);
  EXPECT_EQ(16, mod->get_state()->getValue(CreateLatVol::XSize).toInt());
  command = "scirun_set_module_state(m, \"XSize\", 14)";
  PythonInterpreter::Instance().run_string(command);
  EXPECT_EQ(14, mod->get_state()->getValue(CreateLatVol::XSize).toInt());
}

TEST_F(PythonControllerFunctionalTests, CanConnectModules)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  PythonInterpreter::Instance().run_string("m1 = scirun_add_module(\"CreateLatVol\")");
  PythonInterpreter::Instance().run_string("m2 = scirun_add_module(\"CreateLatVol\")");

  ASSERT_EQ(2, controller.getNetwork()->nmodules());

  ASSERT_EQ(0, controller.getNetwork()->nconnections());

  PythonInterpreter::Instance().run_string("scirun_connect_modules(m1, 0, m2, 0)");
  ASSERT_EQ(1, controller.getNetwork()->nconnections());
}

//TODO: this test is unstable
TEST_F(PythonControllerFunctionalTests, DISABLED_CanExecuteNetwork)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  ModuleStateFactoryHandle sf(new SimpleMapModuleStateFactory);
  ExecutionStrategyFactoryHandle exe(new DesktopExecutionStrategyFactory(boost::none));
  NetworkEditorController controller(mf, sf, exe, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  PythonInterpreter::Instance().run_string("m1 = scirun_add_module(\"CreateLatVol\")");
  ASSERT_TRUE(controller.getNetwork()->module(0)->executionState().currentState() == ModuleExecutionState::NotExecuted);
  PythonInterpreter::Instance().run_string("m2 = scirun_add_module(\"CreateLatVol\")");
  PythonInterpreter::Instance().run_string("scirun_connect_modules(m1, 0, m2, 0)");
  PythonInterpreter::Instance().run_string("scirun_execute_all()");
 // boost::this_thread::sleep(boost::posix_time::milliseconds(500));
  ASSERT_TRUE(controller.getNetwork()->module(0)->executionState().currentState() == ModuleExecutionState::Completed);
  //TODO: how do i assert on
}

TEST_F(PythonControllerFunctionalTests, CanAddModuleWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);
  //TODO: expose API directly on NEC?
  //controller.runPython("addModule(\"CreateLatVol\")");

  ASSERT_EQ(1, controller.getNetwork()->nmodules());
}

TEST_F(PythonControllerFunctionalTests, CanAddMultipleModulesWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);
  PythonInterpreter::Instance().run_string(command);

  ASSERT_EQ(2, controller.getNetwork()->nmodules());
}

TEST_F(PythonControllerFunctionalTests, CanGetModuleStateWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  ModuleStateFactoryHandle sf(new SimpleMapModuleStateFactory);
  NetworkEditorController controller(mf, sf, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "m = scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);

  ASSERT_EQ(1, controller.getNetwork()->nmodules());
  auto mod = controller.getNetwork()->module(0);
  ASSERT_TRUE(mod != nullptr);
  EXPECT_EQ(16, mod->get_state()->getValue(CreateLatVol::XSize).toInt());
  command = "scirun_get_module_state(m, \"XSize\")";
  PythonInterpreter::Instance().run_string(command);

  ////???? need to get value back!!! how??
}

TEST_F(PythonControllerFunctionalTests, CanChangeModuleStateWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  ModuleStateFactoryHandle sf(new SimpleMapModuleStateFactory);
  NetworkEditorController controller(mf, sf, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  std::string command = "m = scirun_add_module(\"CreateLatVol\")";
  PythonInterpreter::Instance().run_string(command);

  ASSERT_EQ(1, controller.getNetwork()->nmodules());
  auto mod = controller.getNetwork()->module(0);
  ASSERT_TRUE(mod != nullptr);
  EXPECT_EQ(16, mod->get_state()->getValue(CreateLatVol::XSize).toInt());
  command = "scirun_set_module_state(m, \"XSize\", 14)";
  PythonInterpreter::Instance().run_string(command);
  EXPECT_EQ(14, mod->get_state()->getValue(CreateLatVol::XSize).toInt());
 // FAIL() << "todo";
}

TEST_F(PythonControllerFunctionalTests, CanConnectModulesWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  PythonInterpreter::Instance().run_string("m1 = scirun_add_module(\"CreateLatVol\")");
  PythonInterpreter::Instance().run_string("m2 = scirun_add_module(\"CreateLatVol\")");

  ASSERT_EQ(2, controller.getNetwork()->nmodules());

  ASSERT_EQ(0, controller.getNetwork()->nconnections());

  PythonInterpreter::Instance().run_string("scirun_connect_modules(m1, 0, m2, 0)");
  ASSERT_EQ(1, controller.getNetwork()->nconnections());
}

TEST_F(PythonControllerFunctionalTests, CanDisconnectModulesWithStaticFunction)
{
  ModuleFactoryHandle mf(new HardCodedModuleFactory);
  NetworkEditorController controller(mf, nullptr, nullptr, nullptr, nullptr, nullptr, nullptr);
  initModuleParameters(false);

  ASSERT_EQ(0, controller.getNetwork()->nmodules());

  PythonInterpreter::Instance().run_string("m1 = scirun_add_module(\"CreateLatVol\")");
  PythonInterpreter::Instance().run_string("m2 = scirun_add_module(\"CreateLatVol\")");

  ASSERT_EQ(2, controller.getNetwork()->nmodules());

  ASSERT_EQ(0, controller.getNetwork()->nconnections());

  PythonInterpreter::Instance().run_string("c = scirun_connect_modules(m1, 0, m2, 0)");
  ASSERT_EQ(1, controller.getNetwork()->nconnections());

  PythonInterpreter::Instance().run_string("scirun_disconnect_modules(m1, 0, m2, 0)");

  ASSERT_EQ(0, controller.getNetwork()->nconnections());
}

ext: cc

max_stars_repo_name: Nahusa/SCIRun


text:
#ifndef BOOST_THREAD_PERMIT_PTHREAD_HPP
#define BOOST_THREAD_PERMIT_PTHREAD_HPP
// Distributed under the Boost Software License, Version 1.0. (See
// accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt)
// (C) Copyright 2007-10 Anthony Williams
// (C) Copyright 2011-2012 Vicente J. Botet Escriba
// (C) Copyright 2014 Niall Douglas

#include <boost/assert.hpp>
#include <boost/throw_exception.hpp>
#include <pthread.h>
#include <boost/thread/cv_status.hpp>
#include <boost/thread/mutex.hpp>
#include <boost/thread/lock_types.hpp>
#include <boost/thread/thread_time.hpp>
#include <boost/thread/pthread/timespec.hpp>
#include <boost/thread/detail/delete.hpp>
#include <boost/date_time/posix_time/posix_time_duration.hpp>
#include <boost/thread/pthread/pthread_mutex_scoped_lock.hpp>
#if defined BOOST_THREAD_USES_DATETIME
#include <boost/thread/xtime.hpp>
#endif
#ifdef BOOST_THREAD_USES_CHRONO
#include <boost/chrono/system_clocks.hpp>
#include <boost/chrono/ceil.hpp>
#endif
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
#include <boost/thread/pthread/thread_data.hpp>
#endif

// Set up the C11 pthreads permit implementation to use an internal C++ namespace
// to prevent it colliding with anything else.
#ifdef PTHREAD_PERMIT_H
#error The C11 pthreads permit implementation is already known to this compilation unit. I cannot proceed.
#endif
#define PTHREAD_PERMIT_USE_BOOST
#define PTHREAD_PERMIT_CXX_NAMESPACE_BEGIN namespace boost { namespace c_permit {
#define PTHREAD_PERMIT_CXX_NAMESPACE_END } }
#define PTHREAD_PERMIT_MANGLEAPI(api) pthread_##api
#define PTHREAD_PERMIT_HEADER_ONLY 1
#include "../permit/pthread_permit.h"

#include <boost/config/abi_prefix.hpp>

namespace boost
{
    namespace detail
    {
        template<bool consuming> struct permit_impl_selector
        {
            typedef boost::c_permit::pthread_permit1_t pthread_permit_t;
            static  int pthread_permit_init     (pthread_permit_t *permit, bool initial)                                    { return boost::c_permit::pthread_permit1_init     (permit, initial); }
            static void pthread_permit_destroy  (pthread_permit_t *permit)                                                  {        boost::c_permit::pthread_permit1_destroy  (permit); }
            static  int pthread_permit_grant    (pthread_permit_t *permit)                                                  { return boost::c_permit::pthread_permit1_grant    (permit); }
            static void pthread_permit_revoke   (pthread_permit_t *permit)                                                  {        boost::c_permit::pthread_permit1_revoke   (permit); }
            static  int pthread_permit_wait     (pthread_permit_t *permit, pthread_mutex_t *mtx)                            { return boost::c_permit::pthread_permit1_wait_locked_grant(permit, mtx); }
            static  int pthread_permit_timedwait(pthread_permit_t *permit, pthread_mutex_t *mtx, const struct timespec *ts) { return boost::c_permit::pthread_permit1_timedwait_locked_grant(permit, mtx, ts); }
        };
        template<> struct permit_impl_selector<false>
        {
            typedef boost::c_permit::pthread_permitnc_t pthread_permit_t;
            static  int pthread_permit_init     (pthread_permit_t *permit, bool initial)                                    { return boost::c_permit::pthread_permitnc_init     (permit, initial); }
            static void pthread_permit_destroy  (pthread_permit_t *permit)                                                  {        boost::c_permit::pthread_permitnc_destroy  (permit); }
            static  int pthread_permit_grant    (pthread_permit_t *permit)                                                  { return boost::c_permit::pthread_permitnc_grant    (permit); }
            static void pthread_permit_revoke   (pthread_permit_t *permit)                                                  {        boost::c_permit::pthread_permitnc_revoke   (permit); }
            static  int pthread_permit_wait     (pthread_permit_t *permit, pthread_mutex_t *mtx)                            { return boost::c_permit::pthread_permitnc_wait_locked_grant(permit, mtx); }
            static  int pthread_permit_timedwait(pthread_permit_t *permit, pthread_mutex_t *mtx, const struct timespec *ts) { return boost::c_permit::pthread_permitnc_timedwait_locked_grant(permit, mtx, ts); }
        };
    }

    template<bool consuming=true> class permit : private detail::permit_impl_selector<consuming>
    {
      typedef typename detail::permit_impl_selector<consuming>::pthread_permit_t pthread_permit_t;
    private:
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
        pthread_mutex_t internal_mutex;
#endif
        pthread_permit_t perm;

    public:
    //private: // used by boost::thread::try_join_until

        inline bool do_wait_until(
            unique_lock<mutex>& lock,
            struct timespec const &timeout);

        bool do_wait_for(
            unique_lock<mutex>& lock,
            struct timespec const &timeout)
        {
          return do_wait_until(lock, boost::detail::timespec_plus(timeout, boost::detail::timespec_now()));
        }

    public:
      BOOST_THREAD_NO_COPYABLE(permit)
        permit(bool initial_state=false)
        {
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
            int const res=pthread_mutex_init(&internal_mutex,NULL);
            if(res)
            {
                boost::throw_exception(thread_resource_error(res, "boost::permit::permit() constructor failed in pthread_mutex_init"));
            }
#endif
            int const res2=pthread_permit_init(&perm,initial_state);
            if(res2)
            {
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
                BOOST_VERIFY(!pthread_mutex_destroy(&internal_mutex));
#endif
                boost::throw_exception(thread_resource_error(res2, "boost::permit::permit() constructor failed in pthread_permit_init"));
            }
        }
        ~permit()
        {
            int ret;
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
            do {
              ret = pthread_mutex_destroy(&internal_mutex);
            } while (ret == EINTR);
            BOOST_ASSERT(!ret);
#endif
            pthread_permit_destroy(&perm);
        }

        void wait(unique_lock<mutex>& m);

        template<typename predicate_type>
        void wait(unique_lock<mutex>& m,predicate_type pred)
        {
            while(!pred()) wait(m);
        }

#if defined BOOST_THREAD_USES_DATETIME
        inline bool timed_wait(
            unique_lock<mutex>& m,
            boost::system_time const& abs_time)
        {
#if defined BOOST_THREAD_WAIT_BUG
            struct timespec const timeout=detail::to_timespec(abs_time + BOOST_THREAD_WAIT_BUG);
            return do_wait_until(m, timeout);
#else
            struct timespec const timeout=detail::to_timespec(abs_time);
            return do_wait_until(m, timeout);
#endif
        }
        bool timed_wait(
            unique_lock<mutex>& m,
            xtime const& abs_time)
        {
            return timed_wait(m,system_time(abs_time));
        }

        template<typename duration_type>
        bool timed_wait(
            unique_lock<mutex>& m,
            duration_type const& wait_duration)
        {
            return timed_wait(m,get_system_time()+wait_duration);
        }

        template<typename predicate_type>
        bool timed_wait(
            unique_lock<mutex>& m,
            boost::system_time const& abs_time,predicate_type pred)
        {
            while (!pred())
            {
                if(!timed_wait(m, abs_time))
                    return pred();
            }
            return true;
        }

        template<typename predicate_type>
        bool timed_wait(
            unique_lock<mutex>& m,
            xtime const& abs_time,predicate_type pred)
        {
            return timed_wait(m,system_time(abs_time),pred);
        }

        template<typename duration_type,typename predicate_type>
        bool timed_wait(
            unique_lock<mutex>& m,
            duration_type const& wait_duration,predicate_type pred)
        {
            return timed_wait(m,get_system_time()+wait_duration,pred);
        }
#endif

#ifdef BOOST_THREAD_USES_CHRONO

        template <class Duration>
        cv_status
        wait_until(
                unique_lock<mutex>& lock,
                const chrono::time_point<chrono::system_clock, Duration>& t)
        {
          using namespace chrono;
          typedef time_point<system_clock, nanoseconds> nano_sys_tmpt;
          wait_until(lock,
                        nano_sys_tmpt(ceil<nanoseconds>(t.time_since_epoch())));
          return system_clock::now() < t ? cv_status::no_timeout :
                                             cv_status::timeout;
        }

        template <class Clock, class Duration>
        cv_status
        wait_until(
                unique_lock<mutex>& lock,
                const chrono::time_point<Clock, Duration>& t)
        {
          using namespace chrono;
          system_clock::time_point     s_now = system_clock::now();
          typename Clock::time_point  c_now = Clock::now();
          wait_until(lock, s_now + ceil<nanoseconds>(t - c_now));
          return Clock::now() < t ? cv_status::no_timeout : cv_status::timeout;
        }

        template <class Clock, class Duration, class Predicate>
        bool
        wait_until(
                unique_lock<mutex>& lock,
                const chrono::time_point<Clock, Duration>& t,
                Predicate pred)
        {
            while (!pred())
            {
                if (wait_until(lock, t) == cv_status::timeout)
                    return pred();
            }
            return true;
        }

        template <class Rep, class Period>
        cv_status
        wait_for(
                unique_lock<mutex>& lock,
                const chrono::duration<Rep, Period>& d)
        {
          using namespace chrono;
          system_clock::time_point s_now = system_clock::now();
          steady_clock::time_point c_now = steady_clock::now();
          wait_until(lock, s_now + ceil<nanoseconds>(d));
          return steady_clock::now() - c_now < d ? cv_status::no_timeout :
                                                   cv_status::timeout;

        }

        template <class Rep, class Period, class Predicate>
        bool
        wait_for(
                unique_lock<mutex>& lock,
                const chrono::duration<Rep, Period>& d,
                Predicate pred)
        {
          return wait_until(lock, chrono::steady_clock::now() + d, boost::move(pred));

//          while (!pred())
//          {
//              if (wait_for(lock, d) == cv_status::timeout)
//                  return pred();
//          }
//          return true;
        }
#endif

#define BOOST_THREAD_DEFINES_PERMIT_NATIVE_HANDLE
        typedef pthread_permit_t* native_handle_type;
        native_handle_type native_handle()
        {
            return &perm;
        }

        void grant() BOOST_NOEXCEPT;
        void revoke() BOOST_NOEXCEPT;

        void notify_one() BOOST_NOEXCEPT { BOOST_STATIC_ASSERT_MSG(consuming, "Use permit<true> if you wish to notify_one()"); grant(); }
        void notify_all() BOOST_NOEXCEPT { BOOST_STATIC_ASSERT_MSG(!consuming, "Use permit<true> if you wish to notify_all()"); grant(); revoke(); }

#ifdef BOOST_THREAD_USES_CHRONO
        inline cv_status wait_until(
            unique_lock<mutex>& lk,
            chrono::time_point<chrono::system_clock, chrono::nanoseconds> tp)
        {
            using namespace chrono;
            nanoseconds d = tp.time_since_epoch();
            timespec ts = boost::detail::to_timespec(d);
            if (do_wait_until(lk, ts)) return cv_status::no_timeout;
            else return cv_status::timeout;
        }
#endif
    };

#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
    namespace this_thread
    {
        void BOOST_THREAD_DECL interruption_point();
    }
#endif

    namespace thread_permit_detail
    {
        template<typename MutexType>
        struct lock_on_exit
        {
            MutexType* m;

            lock_on_exit():
                m(0)
            {}

            void activate(MutexType& m_)
            {
                m_.unlock();
                m=&m_;
            }
            ~lock_on_exit()
            {
                if(m)
                {
                    m->lock();
                }
           }
        };
    }

    template<bool consuming> inline void permit<consuming>::wait(unique_lock<mutex>& m)
    {
#if defined BOOST_THREAD_THROW_IF_PRECONDITION_NOT_SATISFIED
        if(! m.owns_lock())
        {
            boost::throw_exception(condition_error(-1, "boost::permit::wait() failed precondition mutex not owned"));
        }
#endif
        int res=0;
        {
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
            thread_permit_detail::lock_on_exit<unique_lock<mutex> > guard;
            detail::interruption_checker check_for_interruption(&internal_mutex,pthread_permit_get_internal_cond(&perm));
            guard.activate(m);
            //do {
              res = pthread_permit_wait(&perm,&internal_mutex);
            //} while (res == EINTR);
#else
            //boost::pthread::pthread_mutex_scoped_lock check_for_interruption(&internal_mutex);
            pthread_mutex_t* the_mutex = m.mutex()->native_handle();
            //do {
              res = pthread_permit_wait(&perm,the_mutex);
            //} while (res == EINTR);
#endif
        }
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
        this_thread::interruption_point();
#endif
        if(res)
        {
            boost::throw_exception(condition_error(res, "boost::permit::wait failed in pthread_permit_wait"));
        }
    }

    template<bool consuming> inline bool permit<consuming>::do_wait_until(
                unique_lock<mutex>& m,
                struct timespec const &timeout)
    {
#if defined BOOST_THREAD_THROW_IF_PRECONDITION_NOT_SATISFIED
        if (!m.owns_lock())
        {
            boost::throw_exception(condition_error(EPERM, "boost::permit::do_wait_until() failed precondition mutex not owned"));
        }
#endif
        thread_permit_detail::lock_on_exit<unique_lock<mutex> > guard;
        int cond_res;
        {
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
            detail::interruption_checker check_for_interruption(&internal_mutex,pthread_permit_get_internal_cond(&perm));
            guard.activate(m);
            cond_res=pthread_permit_timedwait(&perm,&internal_mutex,&timeout);
#else
            //boost::pthread::pthread_mutex_scoped_lock check_for_interruption(&internal_mutex);
            pthread_mutex_t* the_mutex = m.mutex()->native_handle();
            cond_res=pthread_permit_timedwait(&perm,the_mutex,&timeout);
#endif
        }
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
        this_thread::interruption_point();
#endif
        if(cond_res==ETIMEDOUT)
        {
            return false;
        }
        if(cond_res)
        {
            boost::throw_exception(condition_error(cond_res, "boost::permit::do_wait_until failed in pthread_permit_timedwait"));
        }
        return true;
    }

    template<bool consuming> inline void permit<consuming>::grant() BOOST_NOEXCEPT
    {
        BOOST_VERIFY(!pthread_permit_grant(&perm));
    }

    template<bool consuming> inline void permit<consuming>::revoke() BOOST_NOEXCEPT
    {
        pthread_permit_revoke(&perm);
    }

    template<bool consuming=true> class permit_any : private detail::permit_impl_selector<consuming>
    {
        typedef typename detail::permit_impl_selector<consuming>::pthread_permit_t pthread_permit_t;
        pthread_mutex_t internal_mutex;
        pthread_permit_t perm;

    public:
        BOOST_THREAD_NO_COPYABLE(permit_any)
        permit_any(bool initial_state=false)
        {
            int const res=pthread_mutex_init(&internal_mutex,NULL);
            if(res)
            {
                boost::throw_exception(thread_resource_error(res, "boost::permit_any::permit_any() failed in pthread_mutex_init"));
            }
            int const res2=pthread_permit_init(&perm, initial_state);
            if(res2)
            {
                BOOST_VERIFY(!pthread_mutex_destroy(&internal_mutex));
                boost::throw_exception(thread_resource_error(res, "boost::permit_any::permit_any() failed in pthread_permit_init"));
            }
        }
        ~permit_any()
        {
            BOOST_VERIFY(!pthread_mutex_destroy(&internal_mutex));
            pthread_permit_destroy(&perm);
        }

        template<typename lock_type>
        void wait(lock_type& m)
        {
            int res=0;
            {
                thread_permit_detail::lock_on_exit<lock_type> guard;
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
                detail::interruption_checker check_for_interruption(&internal_mutex,pthread_permit_get_internal_cond(&perm));
#else
            boost::pthread::pthread_mutex_scoped_lock check_for_interruption(&internal_mutex);
#endif
                guard.activate(m);
                res=pthread_permit_wait(&perm,&internal_mutex);
            }
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
            this_thread::interruption_point();
#endif
            if(res)
            {
                boost::throw_exception(condition_error(res, "boost::permit_any::wait() failed in pthread_permit_wait"));
            }
        }

        template<typename lock_type,typename predicate_type>
        void wait(lock_type& m,predicate_type pred)
        {
            while(!pred()) wait(m);
        }

#if defined BOOST_THREAD_USES_DATETIME
        template<typename lock_type>
        bool timed_wait(lock_type& m,boost::system_time const& abs_time)
        {
            struct timespec const timeout=detail::to_timespec(abs_time);
            return do_wait_until(m, timeout);
        }
        template<typename lock_type>
        bool timed_wait(lock_type& m,xtime const& abs_time)
        {
            return timed_wait(m,system_time(abs_time));
        }

        template<typename lock_type,typename duration_type>
        bool timed_wait(lock_type& m,duration_type const& wait_duration)
        {
            return timed_wait(m,get_system_time()+wait_duration);
        }

        template<typename lock_type,typename predicate_type>
        bool timed_wait(lock_type& m,boost::system_time const& abs_time, predicate_type pred)
        {
            while (!pred())
            {
                if(!timed_wait(m, abs_time))
                    return pred();
            }
            return true;
        }

        template<typename lock_type,typename predicate_type>
        bool timed_wait(lock_type& m,xtime const& abs_time, predicate_type pred)
        {
            return timed_wait(m,system_time(abs_time),pred);
        }

        template<typename lock_type,typename duration_type,typename predicate_type>
        bool timed_wait(lock_type& m,duration_type const& wait_duration,predicate_type pred)
        {
            return timed_wait(m,get_system_time()+wait_duration,pred);
        }
#endif
#ifdef BOOST_THREAD_USES_CHRONO
        template <class lock_type,class Duration>
        cv_status
        wait_until(
                lock_type& lock,
                const chrono::time_point<chrono::system_clock, Duration>& t)
        {
          using namespace chrono;
          typedef time_point<system_clock, nanoseconds> nano_sys_tmpt;
          wait_until(lock,
                        nano_sys_tmpt(ceil<nanoseconds>(t.time_since_epoch())));
          return system_clock::now() < t ? cv_status::no_timeout :
                                             cv_status::timeout;
        }

        template <class lock_type, class Clock, class Duration>
        cv_status
        wait_until(
                lock_type& lock,
                const chrono::time_point<Clock, Duration>& t)
        {
          using namespace chrono;
          system_clock::time_point     s_now = system_clock::now();
          typename Clock::time_point  c_now = Clock::now();
          wait_until(lock, s_now + ceil<nanoseconds>(t - c_now));
          return Clock::now() < t ? cv_status::no_timeout : cv_status::timeout;
        }

        template <class lock_type, class Clock, class Duration, class Predicate>
        bool
        wait_until(
                lock_type& lock,
                const chrono::time_point<Clock, Duration>& t,
                Predicate pred)
        {
            while (!pred())
            {
                if (wait_until(lock, t) == cv_status::timeout)
                    return pred();
            }
            return true;
        }

        template <class lock_type, class Rep, class Period>
        cv_status
        wait_for(
                lock_type& lock,
                const chrono::duration<Rep, Period>& d)
        {
          using namespace chrono;
          system_clock::time_point s_now = system_clock::now();
          steady_clock::time_point c_now = steady_clock::now();
          wait_until(lock, s_now + ceil<nanoseconds>(d));
          return steady_clock::now() - c_now < d ? cv_status::no_timeout :
                                                   cv_status::timeout;

        }

        template <class lock_type, class Rep, class Period, class Predicate>
        bool
        wait_for(
                lock_type& lock,
                const chrono::duration<Rep, Period>& d,
                Predicate pred)
        {
          return wait_until(lock, chrono::steady_clock::now() + d, boost::move(pred));

//          while (!pred())
//          {
//              if (wait_for(lock, d) == cv_status::timeout)
//                  return pred();
//          }
//          return true;
        }

        template <class lock_type>
        cv_status wait_until(
            lock_type& lk,
            chrono::time_point<chrono::system_clock, chrono::nanoseconds> tp)
        {
            using namespace chrono;
            nanoseconds d = tp.time_since_epoch();
            timespec ts = boost::detail::to_timespec(d);
            if (do_wait_until(lk, ts)) return cv_status::no_timeout;
            else return cv_status::timeout;
        }
#endif

        void grant() BOOST_NOEXCEPT
        {
            BOOST_VERIFY(!pthread_permit_grant(&perm));
        }

        void revoke() BOOST_NOEXCEPT
        {
            BOOST_VERIFY(!pthread_permit_revoke(&perm));
        }
        void notify_one() BOOST_NOEXCEPT { BOOST_STATIC_ASSERT_MSG(consuming, "Use permit<true> if you wish to notify_one()"); grant(); }
        void notify_all() BOOST_NOEXCEPT { BOOST_STATIC_ASSERT_MSG(!consuming, "Use permit<true> if you wish to notify_all()"); grant(); revoke(); }
    private: // used by boost::thread::try_join_until

        template <class lock_type>
        inline bool do_wait_until(
          lock_type& m,
          struct timespec const &timeout)
        {
          int res=0;
          {
              thread_permit_detail::lock_on_exit<lock_type> guard;
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
              detail::interruption_checker check_for_interruption(&internal_mutex,pthread_permit_get_internal_cond(&perm));
#else
            boost::pthread::pthread_mutex_scoped_lock check_for_interruption(&internal_mutex);
#endif
              guard.activate(m);
              res=pthread_permit_timedwait(&perm,&internal_mutex,&timeout);
          }
#if defined BOOST_THREAD_PROVIDES_INTERRUPTIONS
          this_thread::interruption_point();
#endif
          if(res==ETIMEDOUT)
          {
              return false;
          }
          if(res)
          {
              boost::throw_exception(condition_error(res, "boost::permit_any::do_wait_until() failed in pthread_permit_timedwait"));
          }
          return true;
        }

    };

    typedef permit<true> permit_c;
    typedef permit_any<true> permit_c_any;
    typedef permit<false> permit_nc;
    typedef permit_any<false> permit_nc_any;

}

#include <boost/config/abi_suffix.hpp>

#endif


text:
// vim:ts=2:sw=2:expandtab:autoindent:filetype=cpp:
/*
The MIT License

Copyright (c) 2008, 2009 Flusspferd contributors (see "CONTRIBUTORS" or
                                       http://flusspferd.org/contributors.txt)

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
THE SOFTWARE.
*/

#include "flusspferd/load_core.hpp"
#include "flusspferd/modules.hpp"
#include "flusspferd/evaluate.hpp"
#include "flusspferd/properties_functions.hpp"
#include "flusspferd/binary.hpp"
#include "flusspferd/encodings.hpp"
#include "flusspferd/system.hpp"
#include "flusspferd/getopt.hpp"
#include "flusspferd/io/io.hpp"
#include "flusspferd/io/filesystem-base.hpp"
#include "flusspferd/create.hpp"

#include <boost/filesystem.hpp>
#include <boost/spirit/home/phoenix/core.hpp>
#include <boost/spirit/home/phoenix/bind.hpp>
#include <boost/spirit/home/phoenix/operator.hpp>
namespace phoenix = boost::phoenix;
namespace args = phoenix::arg_names;

using namespace flusspferd;

namespace flusspferd { namespace detail {
  extern char const * const json2;
}}

void flusspferd::load_core(object const &scope_, std::string const &argv0) {
  root_object scope(scope_);

  // Initalize boost's copy of cwd as early as possible
  boost::filesystem::initial_path<boost::filesystem::path>();

  flusspferd::gc();//FIXME

  flusspferd::load_require_function(scope);

  flusspferd::gc();//FIXME

  flusspferd::load_properties_functions(scope);

  flusspferd::object require_fn = scope.get_property_object("require");

  // Create the top level |module| and |exports| properties.
  scope.define_property("module", require_fn.get_property("main"), dont_enumerate);

  root_object exports(create_object());
  scope.define_property("exports", exports, dont_enumerate);

  root_object preload(require_fn.get_property_object("preload"));

  flusspferd::create_native_method(
    preload, "binary",
    &flusspferd::load_binary_module);

  flusspferd::create_native_method(
    preload, "encodings",
    &flusspferd::load_encodings_module);

  flusspferd::create_native_method(
    preload, "io",
    &flusspferd::io::load_io_module);

  flusspferd::create_native_method(
    preload, "system",
    &flusspferd::load_system_module);

  flusspferd::create_native_method(
    preload, "getopt",
    &flusspferd::load_getopt_module);

  flusspferd::create_native_method(
    preload, "filesystem-base",
    &flusspferd::load_filesystem_base_module);
  // alias fs-base -> filesystem-base
  require_fn.get_property_object("alias").set_property("fs-base", "filesystem-base");

  // Curry argv[0] into the preload function
  flusspferd::create_native_method(
    preload, "flusspferd",
    boost::function<void (object)>((
      phoenix::bind(&flusspferd::load_flusspferd_module, args::arg1, argv0)
    )));

  if (!scope_.has_own_property("JSON")) {
    flusspferd::evaluate_in_scope(
      detail::json2,
      std::strlen(detail::json2),
      0x0, // Filename
      0,   // Linenumber
      scope_);
  }
}

max_stars_repo_name: ngoyal/flusspferd

max_issues_repo_name: ngoyal/flusspferd

max_forks_repo_name: ngoyal/flusspferd


text:
#include <htslib/sam.h>
#include <iostream>
#include <fstream>
#include <cstring>
#include <cstdint>
#include <algorithm>
#include <unordered_map>
#include <stdexcept>
#include <boost/graph/adjacency_list.hpp>
#include <boost/graph/connected_components.hpp>

using namespace std;

typedef boost::adjacency_list<boost::vecS, boost::vecS, boost::undirectedS> Graph;

class BamRecord {
  using deleter = void(*)(bam1_t*);
  private:
    unique_ptr<bam1_t, deleter> unq;
    static unique_ptr<bam1_t, deleter> alloc() {
      if (auto p = bam_init1()) {
        return unique_ptr<bam1_t, deleter>(p, bam_destroy1);
      } else {
        throw bad_alloc();
      }
    }
  public:
    bool eof{false};
    BamRecord() : unq(alloc()) {};
    bam1_t* ptr() {
      return unq.get();
    }
    string qname() {
      return string(bam_get_qname(unq.get()));
    }
    int32_t ref_id() const {
      return unq.get()->core.tid;
    }
    int32_t mate_ref_id() const {
      return unq.get()->core.mtid;
    }
};

class BamHeader {
  using deleter = void(*)(bam_hdr_t*);
  private:
    shared_ptr<bam_hdr_t> sptr;
    static shared_ptr<bam_hdr_t> alloc() {
      if (auto p = bam_hdr_init()) {
        return shared_ptr<bam_hdr_t>(p, bam_hdr_destroy);
      } else {
        throw bad_alloc();
      }
    }
  public:
    BamHeader() {};
    BamHeader(bam_hdr_t *header) : sptr(header, bam_hdr_destroy) {};
    bam_hdr_t* ptr() const {
      return sptr.get();
    }
    bool const operator!(){
      return !sptr;
    }
    int n_refs() const {
      return ptr()->n_targets;
    }
};

enum class BamMode { Read, Write };

class BamFile {
  using deleter = int(*)(htsFile*);
  private:
    unique_ptr<samFile, deleter> unq{nullptr, hts_close};
    BamHeader header;
    void read_header() {
      if (!header) {
        bam_hdr_t *hdr = sam_hdr_read(unq.get());
        if (!hdr) {
          throw runtime_error("Could not read header");
        }
        header = BamHeader(hdr);
      }
    }
  public:
    BamFile(const char *path, BamMode mode) {
      const char *mode_str;
      switch (mode) {
        case BamMode::Read: mode_str = "r"; break;
        case BamMode::Write: mode_str = "wb"; break;
      }
      samFile *file = sam_open(path, mode_str);
      if (!file) {
        stringstream msg;
        msg << "Could not open input file: " << strerror(errno);
        throw runtime_error(msg.str());
      }
      unq = unique_ptr<samFile, deleter>(file, hts_close);
    };
    samFile* ptr() {
      return unq.get();
    }
    BamHeader get_header() {
      read_header();
      return header;
    }
    BamRecord next_record() {
      BamRecord rec;

      // ensure that the header is read
      read_header();

      int r = sam_read1(unq.get(), header.ptr(), rec.ptr());
      if (r == -1)
        rec.eof = true;
      if (r < -1)
        throw runtime_error("The file appears corrupted");
      return rec;
    }
};

class ConnectedComponents {
  public:
    // component_of[i] is the component number of the element i
    vector<int> component_of;

    ConnectedComponents(const Graph &g) : component_of(num_vertices(g)) {
      connected_components(g, &component_of[0]);
    }
    void write(const char *path, const BamHeader &bh) {
      ofstream out(path);
      out << "seqid\tseqname\tcomponent\n";
      bam_hdr_t *hdr = bh.ptr();
      for (unsigned i = 0; i < component_of.size(); i++) {
        out << i << "\t" << hdr->target_name[i] << "\t" << component_of[i] << "\n";
      }
    }
};

Graph construct_graph(BamFile &file) {
  unordered_map<string,vector<int>> map;
  BamRecord rec;

  while (rec = file.next_record(), !rec.eof) {
    int ref_id = rec.ref_id(),
        mate_ref_id = rec.mate_ref_id();
    if (ref_id >= 0) { // ref_id == -1 => unmapped
      map[rec.qname()].push_back(ref_id);
    }
    if (mate_ref_id >= 0) { // mate_ref_id == -1 => no mate
      map[rec.qname()].push_back(mate_ref_id);
    }
  }

  Graph g(file.get_header().n_refs());
  for (auto pair : map) {
    auto vec = pair.second;
    for (size_t i = 1; i < vec.size(); i++) {
      add_edge(vec[0], vec[i], g);
    }
  }

  return g;
}

int main(int argc, char **argv) {
  if (argc != 3) {
    cerr << "USAGE: bamcc input.bam output.tsv" << endl;
    exit(1);
  }

  try {
    BamFile bam(argv[1], BamMode::Read);
    Graph g = construct_graph(bam);
    ConnectedComponents comps(g);
    comps.write(argv[2], bam.get_header());
  } catch (const std::exception &exc)
  {
    std::cerr << exc.what() << endl;
    exit(1);
  }
}

max_stars_repo_name: feuerbach/bamcc

max_issues_repo_name: feuerbach/bamcc

max_forks_repo_name: feuerbach/bamcc


text:
#include <radioacek/remote_player.h>
#include <radioacek/remote_player_collection.h>
#include <vector>
#include <boost/regex.hpp>

using namespace std;
namespace radioacek {
bool RemotePlayerCollection::parse_start(
        const string to_parse,
        string& host,
        string& arguments,
        int& port) {
    boost::smatch result;
    if (boost::regex_match(to_parse, result, boost::regex(start_regex))) {
        host = result[1];
        arguments = result[2];
        string filename = result[3];
        if(filename == string("-"))
            throw ParseException("Can't use cout remotely");
        port = std::stoi(result[4]);
        return true;
    }
    else {
        return false;
    }
}
bool RemotePlayerCollection::parse_ppqt(
        const string to_parse,
        string& command_name,
        int& id) {
    boost::smatch result;
    if (boost::regex_match(to_parse, result, boost::regex(play_pause_quit_title_regex))) {
        command_name = result[1];
        id = std::stoi(result[2]);
        return true;
    }
    else {
        return false;
    }
}

bool RemotePlayerCollection::parse_at(
        const string to_parse,
        int& hour,
        int& minutes,
        int& time_to_die,
        string& host,
        string& arguments,
        int& port) {
    boost::smatch result;
    if (boost::regex_match(to_parse, result, boost::regex(at_regex))) {
        hour = std::stoi(result[1]);
        minutes = std::stoi(result[2]);
        time_to_die = std::stoi(result[3]);
        host = result[4];
        arguments = result[5];
        port = std::stoi(result[7]);
        return true;
    }
    else {
        return false;
    }
}

/* Cleans inactive players */
int RemotePlayerCollection::insert(
        string server,
        int port,
        string command,
        std::shared_ptr<TCPConnection> telnet,
        int hour,
        int minutes,
        int time_to_die) {
    while ((unsigned int)it != rplayers.size()) {

        if (!rplayers[it]->is_alive()) {
            rplayers[it] = std::make_shared<RemotePlayer>(
                    server,
                    port,
                    command,
                    it,
                    telnet,
                    hour,
                    minutes,
                    time_to_die);
            break;
        }
        else ++it;
    }
    int ret = it;
    if ((unsigned int)it == rplayers.size()) {
        rplayers.push_back((std::make_shared<RemotePlayer>(
                        server,
                        port,
                        command,
                        it,
                        telnet,
                        hour,
                        minutes,
                        time_to_die)));
        it = 0;
    }
    return ret;

}

bool RemotePlayerCollection::active(unsigned int i) {
    if(i < rplayers.size())
        return rplayers[i]->is_alive();
    else {
        return false;
    }
}

std::shared_ptr<RemotePlayer> RemotePlayerCollection::get(unsigned int i) {
    if (i > rplayers.size())
        throw RemotePlayerException(i, "Wrong index");
    return (rplayers[i]);
}

void RemotePlayerCollection::cleanup() {
    for(auto player : rplayers) {
        if(!player->opened && player->to_start.passed()) {
            player->open_player.unlock();
        }
        if(player->is_suicidal() && player->to_die.passed() && !player->is_canceled())
            player->quit();
    }
}
}


text:
#include "testsuite.h"
#include <blitz/array.h>

BZ_USING_NAMESPACE(blitz)

int main()
{
    Array<int,2> A(4,4), B(3,4), C(4,4), D(3,4);
    Array<int,3> E(4,4,4);

    A = 0;
    B = 0;
    C = 0;
    D = 0;
    E = 0;

    beginCheckAssert();
    A = D;
    endCheckAssert();

    A = C;

    beginCheckAssert();
    A = B + D;
    endCheckAssert();

    return 0;
}

max_stars_repo_name: fraguela/depspawn

max_issues_repo_name: MSV-Project/IBAMR

max_forks_repo_name: MSV-Project/IBAMR