text:
#ifndef STAN_MATH_TORSTEN_EVENTS_MANAGER_HPP
#define STAN_MATH_TORSTEN_EVENTS_MANAGER_HPP

#include <stan/math/torsten/dsolve/pk_vars.hpp>
#include <stan/math/torsten/ev_history.hpp>
#include <stan/math/torsten/ev_record.hpp>
#include <stan/math/torsten/event.hpp>
#include <boost/math/tools/promotion.hpp>
#include <Eigen/Dense>
#include <string>
#include <vector>

namespace torsten {
  template <typename T_event_record, typename T_params>
  struct EventsManager;

  template <typename T0, typename T1, typename T2, typename T3, typename T4,
            template<class...> class theta_container, typename... tuple_pars_t, typename... Ts>
  struct EventsManager<NONMENEventsRecord<T0, T1, T2, T3>,
                       NonEventParameters<T0, T4, theta_container, std::tuple<tuple_pars_t...>, Ts...> > {
    using param_t = NonEventParameters<T0, T4, theta_container, std::tuple<tuple_pars_t...>, Ts...>;
    using ER = NONMENEventsRecord<T0, T1, T2, T3>;
    using T_scalar = typename stan::return_type_t<typename ER::T_scalar, T4, tuple_pars_t..., Ts...>;
    using T_time   = typename stan::return_type_t<typename ER::T_time, typename param_t::lag_t>;
    using T_rate   = typename stan::return_type_t<typename ER::T_rate, typename param_t::biovar_t>;
    using T_amt    = typename stan::return_type_t<typename ER::T_amt, typename param_t::biovar_t>;
    using T_par    = T4;
    using T_par_rate = T2;
    using T_par_ii   = T3;

    param_t params;
    EventHistory<T0, T1, T2, T3, typename param_t::lag_t> event_his;

    int nKeep;
    int ncmt;

    static int nCmt(const ER& rec) {
      return rec.ncmt;
    }

    /*
     * the index in the result/input where subject @c id begins.
     */
    static int begin(int id, const ER& rec) {
      return rec.begin_.at(id);
    }

    /*
     * For population models, we need generate events using
     * ragged arrays.
     */
    template<typename... Tss>
    EventsManager(const ER& rec,
                  const std::vector<theta_container<T4>>& theta,
                  const std::vector<std::vector<Tss>>&... non_event_params) :
      EventsManager(0, rec, theta, non_event_params...)
    {}

    /*
     * For population models, we need generate events using
     * ragged arrays.
     */
    template<typename... Tss>
    EventsManager(int id, const ER& rec,
                  const std::vector<theta_container<T4>>& theta,
                  const std::vector<std::vector<Tss>>&... non_event_params) :
      params(id, rec, theta, non_event_params...),
      event_his(rec.ncmt, rec.begin_[id], rec.len_[id], rec.time_, rec.amt_, rec.rate_, rec.ii_, rec.evid_, rec.cmt_, rec.addl_, rec.ss_)
    {
      ncmt = rec.ncmt;

      attach_event_parameters();
      insert_lag_dose();
      event_his.generate_rates(ncmt);
      attach_event_parameters();

      nKeep = event_his.num_event_times;
    }

    template<typename... Tss>
    EventsManager(int id, const ER& rec,
                  int ibegin_theta, int isize_theta,
                  int ibegin_biovar, int isize_biovar,
                  int ibegin_tlag, int isize_tlag,
                  const std::vector<theta_container<T4>>& theta,
                  const std::vector<std::vector<Tss>>&... non_event_params) :
      params(id, rec, ibegin_theta, isize_theta, ibegin_biovar, isize_biovar, ibegin_tlag, isize_tlag, theta, non_event_params...),
      event_his(rec.ncmt, rec.begin_[id], rec.len_[id], rec.time_, rec.amt_, rec.rate_, rec.ii_, rec.evid_, rec.cmt_, rec.addl_, rec.ss_)
    {
      ncmt = rec.ncmt;

      attach_event_parameters();
      insert_lag_dose();
      event_his.generate_rates(ncmt);
      attach_event_parameters();

      nKeep = event_his.num_event_times;
    }

    /**
     * Implement absorption lag times by modifying the times of the dosing events.
     * Two cases: parameters are either constant or vary with each event.
     * Function sorts events at the end of the procedure.
     * The old event is set with a special EVID = 9 and it introduces no action.
     */
    void insert_lag_dose() {
      // reverse loop so we don't process same lagged events twice
      int nEvent = event_his.size();
      int iEvent = nEvent - 1;
      while (iEvent >= 0) {
        if (event_his.is_dosing(iEvent)) {
          if (params.lag_time(iEvent, event_his.cmt(iEvent) - 1) > ER::lag_time_min) {
            event_his.insert_event(iEvent);
            event_his.gen_time.back() += params.lag_time(iEvent, event_his.cmt(iEvent) - 1);
            event_his.idx[iEvent][2] = 9;
          }
        }
        iEvent--;
      }
      event_his.sort_state_time();
    }

    const EventHistory<T0, T1, T2, T3, typename param_t::lag_t>& events() const {
      return event_his;
    }

    void attach_event_parameters() {
      int nEvent = event_his.size();
      assert(nEvent > 0);
      int len_Parameters = params.size();  // numbers of events for which parameters are determined
      assert(len_Parameters > 0);

      if (!params.is_ordered()) params.sort();
      params.pars.resize(nEvent);

      int iEvent = 0;
      for (int i = 0; i < len_Parameters - 1; ++i) {
        while (event_his.isnew(iEvent)) iEvent++;  // skip new events
        assert(params.get_par_time(i) == event_his.time(iEvent));  // compare time of "old' events to time of parameters.
        iEvent++;
      }

      if (len_Parameters == 1)  {
        for (int i = 0; i < nEvent; ++i) {
          params.set_par_time(i, stan::math::value_of(event_his.time(i)));
          params.set_par_array(i, params.get_par_array(0));
          event_his.idx[i][3] = 0;
        }
      } else {  // parameters are event dependent.
        std::vector<double> times(nEvent, 0);
        for (int i = 0; i < nEvent; ++i) times[i] = params.pars[i].first;
        iEvent = 0;

        using par_t = typename param_t::par_t;
        par_t newParameter;
        int j = 0;
        typename std::vector<par_t>::const_iterator lower = params.pars.begin();
        typename std::vector<par_t>::const_iterator it_param_end = params.pars.begin() + len_Parameters;
        for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
          if (event_his.isnew(iEvent)) {
            // Find the index corresponding to the time of the new event in the
            // times vector.
            const double t = stan::math::value_of(event_his.time(iEvent));
            lower = std::lower_bound(lower, it_param_end, t,
                                     [](const par_t& t1, const double& t2) {return t1.first < t2;});
            newParameter = lower == (it_param_end) ? params.pars[len_Parameters-1] : *lower;
            newParameter.first = t;
            params.pars[len_Parameters + j] = newParameter;
            event_his.idx[iEvent][3] = 0; /**< item is "new" no more, set "isnew" false */
            j++;
          }
        }
      }
      params.sort();
    }

    /**
     * Get model param for certain subject
     *
     * @param i subject id
     *
     * @return model param
     */
    inline const theta_container<T4>& theta(int i) const {
      return params.theta(i);
    }

    template<size_t Is>
    inline auto& get_model_array_1d_param(int i) const {
      return params.template get_model_array_1d_param<Is>(i);
    }

    /**
     * Get dosing rate adjusted with bioavailability for a subject
     *
     * @param i subject id
     *
     * @return dosing rate
     */
    inline std::vector<T_rate> fractioned_rates(int i) const {
      const int n = event_his.rates[0].second.size();
      const std::vector<T2>& r = event_his.rates[event_his.rate_index[i]].second;
      std::vector<T_rate> res(r.size());
      for (size_t j = 0; j < r.size(); ++j) {
        res[j] = r[j] * params.bioavailability(i, j);
      }
      return res;
    }

    /**
     * Get dosing amount adjusted with bioavailability for a subject
     *
     * @param i subject id
     *
     * @return dosing amount
     */
    inline T_amt fractioned_amt(int i) const {
      // return bioavailability(i, cmt(i) - 1) * amt(i);
      return params.bioavailability(i, event_his.cmt(i) - 1) * event_his.amt(i);
    }

    Event<T_time, T3, T_amt, T_rate, T2> event(int i) const {
      int id;
      switch (event_his.evid(i)) {
      case 2:                   // "other" type given "-cmt" indicates turn-off/reset
        if (event_his.cmt(i) < 0) {
            id = 5;
        } else {
          id = 0;
        }
        break;
      case 3:                   // reset
        id = 1;
        break;
      case 4:                   // reset + dosing
        if (event_his.is_ss_dosing(i)) {
          // since it's reset, SS reset is irrelevant
          id = 4;
        } else {
          id = 2;
        }
        break;
      case 8:                   // mrgsolve: "evid=9" overwrite cmt
        if (event_his.cmt(i) > 0) {
          id = 6;
        }
        break;
      default:
        if (event_his.is_ss_dosing(i)) {
          if (event_his.ss(i) == 2) {
            id = 3;
          } else {
            id = 4;
          }
        } else {
          id = 0;
        }
      }

      T_time t0, t1;
      if (event_his.is_ss_dosing(i)) {
        // t0 = event_his.time(i);
        // t1 = event_his.ii(i);
        t0 = i == 0 ? event_his.time(0) : event_his.time(i-1);
        t1 = event_his.time(i);
      } else {
        t0 = i == 0 ? event_his.time(0) : event_his.time(i-1);
        t1 = event_his.time(i);
      }
      PKRec<T_amt> amt = PKRec<T_amt>::Zero(ncmt);
      if (event_his.is_bolus_dosing(i) || event_his.is_ss_dosing(i)) {
        amt(event_his.cmt(i) - 1) = fractioned_amt(i);
      }
      std::vector<T_rate> rate(fractioned_rates(i));
      return {id, t0, t1, event_his.ii(i),
              amt, rate, event_his.rate(i), event_his.cmt(i)};
    }

    /*
     * number of events for a sinlge individual
     */
    template<typename... Tss>
    static int num_events(const ER& rec,
                          const std::vector<theta_container<T4>>& theta,
                          const std::vector<std::vector<Tss>>&... non_event_params) {
      return num_events(0, rec, theta, non_event_params...);
    }

    static int num_events(int id, const ER& rec,
                          const std::vector<theta_container<T4>>& theta) {
      int res;
      int n = rec.len_[id];
      for (int i = rec.begin_[id]; i < rec.begin_[id] + rec.len_[id]; ++i) {
        if (rec.evid_[i] == 1 || rec.evid_[i] == 4) {      // is dosing event
          if (rec.addl_[i] > 0 && rec.ii_[i] > 0) {        // has addl doses
            if (rec.rate_[i] > 0 && rec.amt_[i] > 0) {
              n++;                               // end event for original IV dose
              n += 2 * rec.addl_[i];                  // end event for addl IV dose
            } else {
              n += rec.addl_[i];
            }
          } else if (rec.rate_[i] > 0 && rec.amt_[i] > 0) {
            n++;                                 // end event for IV dose
          }
        }
      }
      res = n;
      return res;
    }

    template<typename T5>
    static int num_events(int id, const ER& rec,
                          const std::vector<theta_container<T4>>& theta,
                          const std::vector<std::vector<T5> >& biovar) {
      return num_events(id, rec, theta);
    }

    template<typename T5, typename T6, typename... Tss>
    static int num_events(int id, const ER& rec,
                          const std::vector<theta_container<T4>>& theta,
                          const std::vector<std::vector<T5> >& biovar,
                          const std::vector<std::vector<T6> >& tlag,
                          const std::vector<std::vector<Tss>>&... non_event_params) {
      int res;
      bool has_lag = rec.has_positive_param(id, tlag);

      if (!has_lag) {
        res = num_events(id, rec, theta);
      } else if (rec.len_param(id, tlag) == 1) {
        int n = rec.len_[id];
        std::vector<std::tuple<double, int>> dose;
        dose.reserve(rec.len_[id]);
        for (int i = rec.begin_[id]; i < rec.begin_[id] + rec.len_[id]; ++i) {
          if (rec.evid_[i] == 1 || rec.evid_[i] == 4) {      // is dosing event
            if (tlag[rec.begin_param(id, tlag)][rec.cmt_[i] - 1] > 0.0) {       // tlag dose
              n++;
            }
            if (rec.addl_[i] > 0 && rec.ii_[i] > 0) {        // has addl doses
              if (rec.rate_[i] > 0 && rec.amt_[i] > 0) {
                n++;                               // end ev for IV dose
                n += 2 * rec.addl_[i];                  // end ev for addl IV dose
              } else {
                n += rec.addl_[i];
              }
              if (tlag[rec.begin_param(id, tlag)][rec.cmt_[i] - 1] > 0.0) {     // tlag dose
                n += rec.addl_[i];
              }
            } else if (rec.rate_[i] > 0 && rec.amt_[i] > 0) {
              n++;                                 // end event for IV dose
            }
          }
        }
        res = n;
      } else {
        // FIXME not to use brute force
        res = EventsManager(id, rec, theta, biovar, tlag, non_event_params...).events().size();
      }

      return res;
    }
  };

}

#endif

ext: hpp


text:
/**
 * \copyright
 * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
 *            Distributed under a Modified BSD License.
 *              See accompanying file LICENSE.txt or
 *              http://www.opengeosys.org/LICENSE.txt
 */

#include <limits>

#include <gtest/gtest.h>

#include <Eigen/Eigen>

#include "BaseLib/ConfigTree.h"

#include "MaterialLib/FractureModels/LinearElasticIsotropic.h"
#include "MaterialLib/FractureModels/MohrCoulomb.h"

#include "ParameterLib/ConstantParameter.h"

using namespace MaterialLib::Fracture;

static const double eps_sigma = 1e6*1e-5;
static const double eps_C = 1e10*1e-5;

TEST(MaterialLib_Fracture, LinearElasticIsotropic)
{
    ParameterLib::ConstantParameter<double> const kn("", 1e11);
    ParameterLib::ConstantParameter<double> const ks("", 1e9);
    LinearElasticIsotropic<2>::MaterialProperties const mp{kn, ks};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    LinearElasticIsotropic<2> fractureModel{penalty_aperture_cutoff,
                                            tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<2>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector2d const w_prev = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma0 = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma_prev = sigma0;
    Eigen::Vector2d w(-1e-5, -1e-5);

    // Result vectors, not initialized.
    Eigen::Vector2d sigma;
    Eigen::Matrix2d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-1e4, sigma[0], eps_sigma);
    ASSERT_NEAR(-1e6, sigma[1], eps_sigma);
    ASSERT_NEAR(1e9, C(0,0), eps_C);
    ASSERT_NEAR(0, C(0,1), eps_C);
    ASSERT_NEAR(0, C(1,0), eps_C);
    ASSERT_NEAR(1e11, C(1,1), eps_C);

    w << -1e-5, 1e-5;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(0, sigma[0], eps_sigma);
    ASSERT_NEAR(0, sigma[1], eps_sigma);
    ASSERT_NEAR(0, C(0,0), eps_C);
    ASSERT_NEAR(0, C(0,1), eps_C);
    ASSERT_NEAR(0, C(1,0), eps_C);
    ASSERT_NEAR(0, C(1,1), eps_C);
}

TEST(MaterialLib_Fracture, MohrCoulomb2D_elastic)
{
    ParameterLib::ConstantParameter<double> const kn("", 1e11);
    ParameterLib::ConstantParameter<double> const ks("", 1e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e16);
    MohrCoulomb<2>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<2> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<2>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector2d const w_prev = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma0 = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma_prev = sigma0;
    Eigen::Vector2d w(-1e-5, -1e-5);

    // Result vectors, not initialized.
    Eigen::Vector2d sigma;
    Eigen::Matrix2d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-1e4, sigma[0], eps_sigma);
    ASSERT_NEAR(-1e6, sigma[1], eps_sigma);
    ASSERT_NEAR(1e9, C(0, 0), eps_C);
    ASSERT_NEAR(0, C(0, 1), eps_C);
    ASSERT_NEAR(0, C(1, 0), eps_C);
    ASSERT_NEAR(1e11, C(1, 1), eps_C);

    w << -1e-5, 1e-5;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(0, sigma[0], eps_sigma);
    ASSERT_NEAR(0, sigma[1], eps_sigma);
    ASSERT_NEAR(0, C(0, 0), eps_C);
    ASSERT_NEAR(0, C(0, 1), eps_C);
    ASSERT_NEAR(0, C(1, 0), eps_C);
    ASSERT_NEAR(0, C(1, 1), eps_C);
}

TEST(MaterialLib_Fracture, MohrCoulomb2D_negative_t)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<2>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<2> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<2>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector2d const w_prev = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma0(-3.46e6, -2e6);
    Eigen::Vector2d const sigma_prev = sigma0;
    Eigen::Vector2d const w(-1.08e-5, -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector2d sigma;
    Eigen::Matrix2d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-3.50360e6, sigma[0], eps_sigma);
    ASSERT_NEAR(-2.16271e6, sigma[1], eps_sigma);
    ASSERT_NEAR(1.10723e+09, C(0,0), eps_C);
    ASSERT_NEAR(1.26558e+10, C(0,1), eps_C);
    ASSERT_NEAR(4.13226e+09, C(1,0), eps_C);
    ASSERT_NEAR(4.72319e+10, C(1,1), eps_C);
    ASSERT_NEAR(1.06608E+5, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb2D_positive_t)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<2>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<2> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<2>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector2d const w_prev = Eigen::Vector2d::Zero();
    Eigen::Vector2d const sigma0(0, -2e6);
    Eigen::Vector2d const sigma_prev = sigma0;
    Eigen::Vector2d const w(20.08e-5, -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector2d sigma;
    Eigen::Matrix2d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(253972.0, sigma[0], eps_sigma);
    ASSERT_NEAR(-2.94784e+06, sigma[1], eps_sigma);
    ASSERT_NEAR(1.10723e+09, C(0,0), eps_C);
    ASSERT_NEAR(-1.26558e+10, C(0,1), eps_C);
    ASSERT_NEAR(-4.13226e+09, C(1,0), eps_C);
    ASSERT_NEAR(4.72319e+10, C(1,1), eps_C);
    ASSERT_NEAR(446608.0, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb3D_negative_t1)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<3>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<3> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<3>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector3d const w_prev = Eigen::Vector3d::Zero();
    Eigen::Vector3d const sigma0(-3.46e6, 0, -2e6);
    Eigen::Vector3d const sigma_prev = sigma0;
    Eigen::Vector3d const w(-1.08e-5, 0, -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector3d sigma;
    Eigen::Matrix3d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-3.50360e6, sigma[0], eps_sigma);
    ASSERT_NEAR(0.0, sigma[1], eps_sigma);
    ASSERT_NEAR(-2.16271e6, sigma[2], eps_sigma);
    ASSERT_NEAR(1.10723e+09, C(0,0), eps_C);
    ASSERT_NEAR(0.0, C(0,1), eps_C);
    ASSERT_NEAR(1.26558e+10, C(0,2), eps_C);
    ASSERT_NEAR(0.0, C(1,0), eps_C);
    ASSERT_NEAR(20.e9, C(1,1), eps_C);
    ASSERT_NEAR(0.0, C(1,2), eps_C);
    ASSERT_NEAR(4.13226e+09, C(2,0), eps_C);
    ASSERT_NEAR(0.0, C(2,1), eps_C);
    ASSERT_NEAR(4.72319e+10, C(2,2), eps_C);
    ASSERT_NEAR(1.06608E+5, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb3D_positive_t1)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<3>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<3> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<3>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector3d const w_prev = Eigen::Vector3d::Zero();
    Eigen::Vector3d const sigma0(0, 0, -2e6);
    Eigen::Vector3d const sigma_prev = sigma0;
    Eigen::Vector3d const w(20.08e-5, 0, -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector3d sigma;
    Eigen::Matrix3d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(253972.0, sigma[0], eps_sigma);
    ASSERT_NEAR(0.0, sigma[1], eps_sigma);
    ASSERT_NEAR(-2.94784e+06, sigma[2], eps_sigma);
    ASSERT_NEAR(1.10723e+09, C(0,0), eps_C);
    ASSERT_NEAR(0.0, C(0,1), eps_C);
    ASSERT_NEAR(-1.26558e+10, C(0,2), eps_C);
    ASSERT_NEAR(0.0, C(1,0), eps_C);
    ASSERT_NEAR(20.e9, C(1,1), eps_C);
    ASSERT_NEAR(0.0, C(1,2), eps_C);
    ASSERT_NEAR(-4.13226e+09, C(2,0), eps_C);
    ASSERT_NEAR(0.0, C(2,1), eps_C);
    ASSERT_NEAR(4.72319e+10, C(2,2), eps_C);
    ASSERT_NEAR(446608.0, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb3D_positive_t2)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<3>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<3> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<3>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector3d const w_prev = Eigen::Vector3d::Zero();
    Eigen::Vector3d const sigma0(0, 0, -2e6);
    Eigen::Vector3d const sigma_prev = sigma0;
    Eigen::Vector3d const w(0, 20.08e-5, -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector3d sigma;
    Eigen::Matrix3d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(0.0, sigma[0], eps_sigma);
    ASSERT_NEAR(253972.0, sigma[1], eps_sigma);
    ASSERT_NEAR(-2.94784e+06, sigma[2], eps_sigma);
    ASSERT_NEAR(20.e9, C(0,0), eps_C);
    ASSERT_NEAR(0.0, C(0,1), eps_C);
    ASSERT_NEAR(0.0, C(0,2), eps_C);
    ASSERT_NEAR(0.0, C(1,0), eps_C);
    ASSERT_NEAR(1.10723e+09, C(1,1), eps_C);
    ASSERT_NEAR(-1.26558e+10, C(1,2), eps_C);
    ASSERT_NEAR(0.0, C(2,0), eps_C);
    ASSERT_NEAR(-4.13226e+09, C(2,1), eps_C);
    ASSERT_NEAR(4.72319e+10, C(2,2), eps_C);
    ASSERT_NEAR(446608.0, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb3D_negative_t1t2)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<3>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<3> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<3>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector3d const w_prev = Eigen::Vector3d::Zero();
    Eigen::Vector3d const sigma0(-3.46e6 / std::sqrt(2), -3.46e6 / std::sqrt(2),
                                 -2e6);
    Eigen::Vector3d const sigma_prev = sigma0;
    Eigen::Vector3d const w(-1.08e-5/std::sqrt(2), -1.08e-5/std::sqrt(2), -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector3d sigma;
    Eigen::Matrix3d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-3.50360e6/std::sqrt(2), sigma[0], eps_sigma);
    ASSERT_NEAR(-3.50360e6/std::sqrt(2), sigma[1], eps_sigma);
    ASSERT_NEAR(-2.16271e6, sigma[2], eps_sigma);
    ASSERT_NEAR(1.05536e+10, C(0,0), eps_C);
    ASSERT_NEAR(-9.44638e+09, C(0,1), eps_C);
    ASSERT_NEAR(8.94897e+09, C(0,2), eps_C);
    ASSERT_NEAR(-9.44638e+09, C(1,0), eps_C);
    ASSERT_NEAR(1.05536e+10, C(1,1), eps_C);
    ASSERT_NEAR(8.94897e+09, C(1,2), eps_C);
    ASSERT_NEAR(2.92195e+09, C(2,0), eps_C);
    ASSERT_NEAR(2.92195e+09, C(2,1), eps_C);
    ASSERT_NEAR(4.72319e+10, C(2,2), eps_C);
    ASSERT_NEAR(1.06608E+5, state->getShearYieldFunctionValue(), eps_sigma);
}

TEST(MaterialLib_Fracture, MohrCoulomb3D_negative_t1_positive_t2)
{
    ParameterLib::ConstantParameter<double> const kn("", 50e9);
    ParameterLib::ConstantParameter<double> const ks("", 20e9);
    ParameterLib::ConstantParameter<double> const phi("", 15);
    ParameterLib::ConstantParameter<double> const psi("", 5);
    ParameterLib::ConstantParameter<double> const c("", 3e6);
    MohrCoulomb<3>::MaterialProperties const mp{kn, ks, phi, psi, c};

    double const aperture0 = 1e-5;
    double const penalty_aperture_cutoff = aperture0;
    bool const tension_cutoff = true;
    MohrCoulomb<3> fractureModel{penalty_aperture_cutoff, tension_cutoff, mp};
    std::unique_ptr<FractureModelBase<3>::MaterialStateVariables> state(
        fractureModel.createMaterialStateVariables());

    Eigen::Vector3d const w_prev = Eigen::Vector3d::Zero();
    Eigen::Vector3d const sigma0(-3.46e6 / std::sqrt(2), 3.46e6 / std::sqrt(2),
                                 -2e6);
    Eigen::Vector3d const sigma_prev = sigma0;
    Eigen::Vector3d const w(-1.08e-5/std::sqrt(2), 1.08e-5/std::sqrt(2), -0.25e-5);

    // Result vectors, not initialized.
    Eigen::Vector3d sigma;
    Eigen::Matrix3d C;

    ParameterLib::SpatialPosition x;
    fractureModel.computeConstitutiveRelation(0, x, aperture0, sigma0, w_prev,
                                              w, sigma_prev, sigma, C, *state);

    ASSERT_NEAR(-3.50360e6/std::sqrt(2), sigma[0], eps_sigma);
    ASSERT_NEAR(3.50360e6/std::sqrt(2), sigma[1], eps_sigma);
    ASSERT_NEAR(-2.16271e6, sigma[2], eps_sigma);
    ASSERT_NEAR(1.05536e+10, C(0,0), eps_C);
    ASSERT_NEAR(9.44638e+09, C(0,1), eps_C);
    ASSERT_NEAR(8.94897e+09, C(0,2), eps_C);
    ASSERT_NEAR(9.44638e+09, C(1,0), eps_C);
    ASSERT_NEAR(1.05536e+10, C(1,1), eps_C);
    ASSERT_NEAR(-8.94897e+09, C(1,2), eps_C);
    ASSERT_NEAR(2.92195e+09, C(2,0), eps_C);
    ASSERT_NEAR(-2.92195e+09, C(2,1), eps_C);
    ASSERT_NEAR(4.72319e+10, C(2,2), eps_C);
    ASSERT_NEAR(1.06608E+5, state->getShearYieldFunctionValue(), eps_sigma);
}

ext: cpp

max_stars_repo_name: wolf-pf/ogs

max_issues_repo_name: wolf-pf/ogs

max_forks_repo_name: wolf-pf/ogs


text:
// Licensed to the Apache Software Foundation (ASF) under one
// or more contributor license agreements.  See the NOTICE file
// distributed with this work for additional information
// regarding copyright ownership.  The ASF licenses this file
// to you 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.

#include <chrono>
#include <cstring>
#include <utility>

#ifdef _WIN32
#include "arrow/util/windows_compatibility.h"
#else
#include <errno.h>
#include <fcntl.h>
#include <stdio.h>
#include <sys/stat.h>
#endif

#include <boost/filesystem.hpp>

#include "arrow/filesystem/localfs.h"
#include "arrow/filesystem/util_internal.h"
#include "arrow/io/file.h"
#include "arrow/util/io_util.h"
#include "arrow/util/logging.h"

namespace arrow {
namespace fs {

namespace bfs = ::boost::filesystem;

using ::arrow::internal::NativePathString;
using ::arrow::internal::PlatformFilename;

namespace {

#define BOOST_FILESYSTEM_TRY try {
#define BOOST_FILESYSTEM_CATCH           \
  }                                      \
  catch (bfs::filesystem_error & _err) { \
    return ToStatus(_err);               \
  }

// NOTE: catching filesystem_error gives more context than system::error_code
// (it includes the file path(s) in the error message)

Status ToStatus(const bfs::filesystem_error& err) { return Status::IOError(err.what()); }

template <typename... Args>
Status ErrnoToStatus(Args&&... args) {
  auto err_string = ::arrow::internal::ErrnoMessage(errno);
  return Status::IOError(std::forward<Args>(args)..., err_string);
}

#ifdef _WIN32

std::string NativeToString(const NativePathString& ns) {
  PlatformFilename fn(ns);
  return fn.ToString();
}

template <typename... Args>
Status WinErrorToStatus(Args&&... args) {
  auto err_string = ::arrow::internal::WinErrorMessage(GetLastError());
  return Status::IOError(std::forward<Args>(args)..., err_string);
}

TimePoint ToTimePoint(FILETIME ft) {
  // Hundreds of nanoseconds between January 1, 1601 (UTC) and the Unix epoch.
  static constexpr int64_t kFileTimeEpoch = 11644473600LL * 10000000;

  int64_t hundreds = (static_cast<int64_t>(ft.dwHighDateTime) << 32) + ft.dwLowDateTime -
                     kFileTimeEpoch;  // hundreds of ns since Unix epoch
  std::chrono::nanoseconds ns_count(100 * hundreds);
  return TimePoint(std::chrono::duration_cast<TimePoint::duration>(ns_count));
}

FileStats FileInformationToFileStat(const BY_HANDLE_FILE_INFORMATION& info) {
  FileStats st;
  if (info.dwFileAttributes & FILE_ATTRIBUTE_DIRECTORY) {
    st.set_type(FileType::Directory);
    st.set_size(kNoSize);
  } else {
    // Regular file
    st.set_type(FileType::File);
    st.set_size((static_cast<int64_t>(info.nFileSizeHigh) << 32) + info.nFileSizeLow);
  }
  st.set_mtime(ToTimePoint(info.ftLastWriteTime));
  return st;
}

Status StatFile(const std::wstring& path, FileStats* out) {
  HANDLE h;
  std::string bytes_path = NativeToString(path);

  /* Inspired by CPython, see Modules/posixmodule.c */
  h = CreateFileW(path.c_str(), FILE_READ_ATTRIBUTES, /* desired access */
                  0,                                  /* share mode */
                  NULL,                               /* security attributes */
                  OPEN_EXISTING,
                  /* FILE_FLAG_BACKUP_SEMANTICS is required to open a directory */
                  FILE_ATTRIBUTE_NORMAL | FILE_FLAG_BACKUP_SEMANTICS, NULL);

  if (h == INVALID_HANDLE_VALUE) {
    DWORD err = GetLastError();
    if (err == ERROR_FILE_NOT_FOUND || err == ERROR_PATH_NOT_FOUND) {
      out->set_path(bytes_path);
      out->set_type(FileType::NonExistent);
      out->set_mtime(kNoTime);
      out->set_size(kNoSize);
      return Status::OK();
    } else {
      return WinErrorToStatus("Failed querying information for path '", bytes_path, "'");
    }
  }
  BY_HANDLE_FILE_INFORMATION info;
  if (!GetFileInformationByHandle(h, &info)) {
    Status st =
        WinErrorToStatus("Failed querying information for path '", bytes_path, "'");
    CloseHandle(h);
    return st;
  }
  CloseHandle(h);
  *out = FileInformationToFileStat(info);
  out->set_path(bytes_path);
  return Status::OK();
}

#else  // POSIX systems

TimePoint ToTimePoint(const struct timespec& s) {
  std::chrono::nanoseconds ns_count(static_cast<int64_t>(s.tv_sec) * 1000000000 +
                                    static_cast<int64_t>(s.tv_nsec));
  return TimePoint(std::chrono::duration_cast<TimePoint::duration>(ns_count));
}

FileStats StatToFileStat(const struct stat& s) {
  FileStats st;
  if (S_ISREG(s.st_mode)) {
    st.set_type(FileType::File);
    st.set_size(static_cast<int64_t>(s.st_size));
  } else if (S_ISDIR(s.st_mode)) {
    st.set_type(FileType::Directory);
    st.set_size(kNoSize);
  } else {
    st.set_type(FileType::Unknown);
    st.set_size(kNoSize);
  }
#ifdef __APPLE__
  // macOS doesn't use the POSIX-compliant spelling
  st.set_mtime(ToTimePoint(s.st_mtimespec));
#else
  st.set_mtime(ToTimePoint(s.st_mtim));
#endif
  return st;
}

Status StatFile(const std::string& path, FileStats* out) {
  struct stat s;
  int r = stat(path.c_str(), &s);
  if (r == -1) {
    if (errno == ENOENT || errno == ENOTDIR || errno == ELOOP) {
      out->set_type(FileType::NonExistent);
      out->set_mtime(kNoTime);
      out->set_size(kNoSize);
    } else {
      return ErrnoToStatus("Failed stat()ing path '", path, "'");
    }
  } else {
    *out = StatToFileStat(s);
  }
  out->set_path(path);
  return Status::OK();
}

#endif

Status StatSelector(const NativePathString& path, const Selector& select,
                    std::vector<FileStats>* out) {
  bfs::path p(path);

  if (select.allow_non_existent) {
    bfs::file_status st;
    BOOST_FILESYSTEM_TRY
    st = bfs::status(p);
    BOOST_FILESYSTEM_CATCH
    if (st.type() == bfs::file_not_found) {
      return Status::OK();
    }
  }

  BOOST_FILESYSTEM_TRY
  for (const auto& entry : bfs::directory_iterator(p)) {
    FileStats st;
    NativePathString ns = entry.path().native();
    RETURN_NOT_OK(StatFile(ns, &st));
    if (st.type() != FileType::NonExistent) {
      out->push_back(std::move(st));
    }
    if (select.recursive && st.type() == FileType::Directory) {
      RETURN_NOT_OK(StatSelector(ns, select, out));
    }
  }
  BOOST_FILESYSTEM_CATCH

  return Status::OK();
}

}  // namespace

LocalFileSystem::LocalFileSystem() {}

LocalFileSystem::~LocalFileSystem() {}

Status LocalFileSystem::GetTargetStats(const std::string& path, FileStats* out) {
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  return StatFile(fn.ToNative(), out);
}

Status LocalFileSystem::GetTargetStats(const Selector& select,
                                       std::vector<FileStats>* out) {
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(select.base_dir, &fn));
  out->clear();
  return StatSelector(fn.ToNative(), select, out);
}

Status LocalFileSystem::CreateDir(const std::string& path, bool recursive) {
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  if (recursive) {
    return ::arrow::internal::CreateDirTree(fn);
  } else {
    return ::arrow::internal::CreateDir(fn);
  }
}

Status LocalFileSystem::DeleteDir(const std::string& path) {
  bool deleted = false;
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  RETURN_NOT_OK(::arrow::internal::DeleteDirTree(fn, &deleted));
  if (deleted) {
    return Status::OK();
  } else {
    return Status::IOError("Directory does not exist: '", path, "'");
  }
}

Status LocalFileSystem::DeleteDirContents(const std::string& path) {
  bool deleted = false;
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  RETURN_NOT_OK(::arrow::internal::DeleteDirContents(fn, &deleted));
  if (deleted) {
    return Status::OK();
  } else {
    return Status::IOError("Directory does not exist: '", path, "'");
  }
}

Status LocalFileSystem::DeleteFile(const std::string& path) {
  bool deleted = false;
  PlatformFilename fn;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  RETURN_NOT_OK(::arrow::internal::DeleteFile(fn, &deleted));
  if (deleted) {
    return Status::OK();
  } else {
    return Status::IOError("File does not exist: '", path, "'");
  }
}

Status LocalFileSystem::Move(const std::string& src, const std::string& dest) {
  PlatformFilename sfn, dfn;
  RETURN_NOT_OK(PlatformFilename::FromString(src, &sfn));
  RETURN_NOT_OK(PlatformFilename::FromString(dest, &dfn));

#ifdef _WIN32
  if (!MoveFileExW(sfn.ToNative().c_str(), dfn.ToNative().c_str(),
                   MOVEFILE_REPLACE_EXISTING)) {
    return WinErrorToStatus("Failed renaming '", sfn.ToString(), "' to '", dfn.ToString(),
                            "': ");
  }
#else
  if (rename(sfn.ToNative().c_str(), dfn.ToNative().c_str()) == -1) {
    return ErrnoToStatus("Failed renaming '", sfn.ToString(), "' to '", dfn.ToString(),
                         "': ");
  }
#endif
  return Status::OK();
}

Status LocalFileSystem::CopyFile(const std::string& src, const std::string& dest) {
  PlatformFilename sfn, dfn;
  RETURN_NOT_OK(PlatformFilename::FromString(src, &sfn));
  RETURN_NOT_OK(PlatformFilename::FromString(dest, &dfn));
  // XXX should we use fstat() to compare inodes?
  if (sfn.ToNative() == dfn.ToNative()) {
    return Status::OK();
  }

#ifdef _WIN32
  if (!CopyFileW(sfn.ToNative().c_str(), dfn.ToNative().c_str(),
                 FALSE /* bFailIfExists */)) {
    return WinErrorToStatus("Failed copying '", sfn.ToString(), "' to '", dfn.ToString(),
                            "': ");
  }
  return Status::OK();
#else
  std::shared_ptr<io::InputStream> is;
  std::shared_ptr<io::OutputStream> os;
  RETURN_NOT_OK(OpenInputStream(src, &is));
  RETURN_NOT_OK(OpenOutputStream(dest, &os));
  RETURN_NOT_OK(internal::CopyStream(is, os, 1024 * 1024 /* chunk_size */));
  RETURN_NOT_OK(os->Close());
  return is->Close();
#endif
}

Status LocalFileSystem::OpenInputStream(const std::string& path,
                                        std::shared_ptr<io::InputStream>* out) {
  std::shared_ptr<io::ReadableFile> file;
  RETURN_NOT_OK(io::ReadableFile::Open(path, &file));
  *out = std::move(file);
  return Status::OK();
}

Status LocalFileSystem::OpenInputFile(const std::string& path,
                                      std::shared_ptr<io::RandomAccessFile>* out) {
  std::shared_ptr<io::ReadableFile> file;
  RETURN_NOT_OK(io::ReadableFile::Open(path, &file));
  *out = std::move(file);
  return Status::OK();
}

namespace {

Status OpenOutputStreamGeneric(const std::string& path, bool truncate, bool append,
                               std::shared_ptr<io::OutputStream>* out) {
  PlatformFilename fn;
  int fd;
  bool write_only = true;
  RETURN_NOT_OK(PlatformFilename::FromString(path, &fn));
  RETURN_NOT_OK(
      ::arrow::internal::FileOpenWritable(fn, write_only, truncate, append, &fd));
  Status st = io::FileOutputStream::Open(fd, out);
  if (!st.ok()) {
    ARROW_UNUSED(::arrow::internal::FileClose(fd));
  }
  return st;
}

}  // namespace

Status LocalFileSystem::OpenOutputStream(const std::string& path,
                                         std::shared_ptr<io::OutputStream>* out) {
  bool truncate = true;
  bool append = false;
  return OpenOutputStreamGeneric(path, truncate, append, out);
}

Status LocalFileSystem::OpenAppendStream(const std::string& path,
                                         std::shared_ptr<io::OutputStream>* out) {
  bool truncate = false;
  bool append = true;
  return OpenOutputStreamGeneric(path, truncate, append, out);
}

}  // namespace fs
}  // namespace arrow

ext: cc

max_stars_repo_name: yongbozhang/arrow

max_issues_repo_name: dewott-technologies/arrow

max_forks_repo_name: zettaflowio/arrow


text:
#include "client.h"

#include <clc_if.h>
#include <clc_if_types.h>
#include <utils/cwd.h>
#include <utils/logger.h>
#include <utils/paths.h>

#include <boost/make_shared.hpp>
#include <thrift/protocol/TCompactProtocol.h>
#include <thrift/transport/TBufferTransports.h>
#include <thrift/transport/TSocket.h>
#include <vector>

namespace tprotocol = apache::thrift::protocol;
namespace ttransport = apache::thrift::transport;

namespace clc { namespace client {

static clc::rpc::clc_ifClient* connect_to_server(void)
{
  auto socket = boost::make_shared<ttransport::TSocket>(utils::sock_path());
  auto transport = boost::make_shared<ttransport::TBufferedTransport>(socket);
  auto protocol = boost::make_shared<tprotocol::TCompactProtocol>(transport);
  auto client = new clc::rpc::clc_ifClient(protocol);

  do {
    try {
      /*
       * When the server is being launched, we might end-up here before
       * the server has created the unix socket. Let's loop until it has
       * been, this shouldn't be long.
       */
      transport->open();
    } catch (ttransport::TTransportException& e) {
      /* I wonder if that is a bug in thrift. */
      transport->close();
    }
  } while (!transport->isOpen());

  return client;
}

void register_compilation(int argc, char** argv)
{
  auto client = connect_to_server();
  auto transport = client->getOutputProtocol()->getTransport();

  std::vector<std::string> cmd(argv, argv + argc);
  client->register_compilation(cmd, utils::cwd());
  transport->close();

  delete client;
}

void find_definition(const char* usr)
{
  auto client = connect_to_server();
  auto transport = client->getOutputProtocol()->getTransport();

  rpc::cursor ret;
  std::string arg(usr);
  client->find_definition_from_usr(ret, arg);
  transport->close();

  delete client;

  LOG_INFO() << "symbol " << arg << " is defined in " << ret.file << ":" << ret.line;
}

}} // namespace clc::client

max_stars_repo_name: sas/clang-cache

max_issues_repo_name: sas/clang-cache

max_forks_repo_name: sas/clang-cache


text:
//
// GraphTools library
// Copyright 2017-2019 Illumina, Inc.
// All rights reserved.
//
// Author: Felix Schlesinger <fschlesinger@illumina.com>
//
// 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.
//

#include "graphio/GraphJson.hh"

#include <boost/optional.hpp>
#include <fstream>

using std::string;

namespace graphtools
{

Graph loadGraph(string const& jsonPath)
{
    std::ifstream jsonFile(jsonPath);
    assert(jsonFile.good());
    Json json;
    jsonFile >> json;
    return parseGraph((json.count("graph") == 1) ? json["graph"] : json);
}

Graph parseGraph(Json const& jsonGraph)
{
    // boost::optional<RefGenome> genome;
    // const auto refFasta = jsonGraph.find("reference_genome");
    // if (refFasta != jsonGraph.end())
    //{
    //    genome.emplace(refFasta->get<string>());
    //}

    Json::array_t nodes = jsonGraph.value("nodes", Json::array());
    auto const nNodes = nodes.size();
    string const& graphId = jsonGraph.value("graph_id", "");
    Graph graph(nNodes, graphId);

    std::unordered_map<string, NodeId> nodeIds; // NodeName -> NameID
    NodeId nodeIndex = 0;
    for (auto const& jsonNode : nodes)
    {
        string const name = jsonNode.at("name");
        graph.setNodeName(nodeIndex, name);
        auto const seq = jsonNode.find("sequence");
        if (seq != jsonNode.end())
        {
            graph.setNodeSeq(nodeIndex, *seq);
        }
        else
        {
            throw std::runtime_error("Node has an invalid sequence: " + graph.nodeName(nodeIndex));
            // auto const refRegion = jsonNode.find("reference");
            // if (refRegion == jsonNode.end())
            //{
            //    throw std::runtime_error("Node has no sequence: " + graph.nodeName(nodeIndex));
            //}
            // if (!genome)
            //{
            //    throw std::runtime_error("Need 'referenceGenome' FASTA file to use reference nodes");
            //}
            // auto const interval = ReferenceInterval::parseRegion(*refRegion);
            // graph.setNodeSeq(nodeIndex, genome->extractSeq(interval));
        }
        assert(nodeIds.count(name) == 0);
        nodeIds[name] = nodeIndex;
        nodeIndex++;
    }
    for (auto const& jsonEdge : jsonGraph.value("edges", Json::array()))
    {
        NodeId const sourceNode = nodeIds.at(jsonEdge.at("from"));
        NodeId const sinkNode = nodeIds.at(jsonEdge.at("to"));
        graph.addEdge(sourceNode, sinkNode);
        const auto labels = jsonEdge.find("labels");
        if (labels != jsonEdge.end())
        {
            for (auto const& label : *labels)
            {
                graph.addLabelToEdge(sourceNode, sinkNode, label.get<string>());
            }
        }
    }
    return graph;
}

Json graphToJson(Graph const& graph)
{
    Json json;
    if (!graph.graphId.empty())
    {
        json["graph_id"] = graph.graphId;
    }
    json["nodes"] = Json::array();
    for (size_t i = 0; i != graph.numNodes(); ++i)
    {
        json["nodes"].push_back({ { "name", graph.nodeName(i) }, { "sequence", graph.nodeSeq(i) } });
    }
    json["edges"] = Json::array();
    for (size_t n1 = 0; n1 < graph.numNodes(); ++n1)
    {
        for (NodeId n2 : graph.successors(n1))
        {
            Json edge = { { "from", graph.nodeName(n1) }, { "to", graph.nodeName(n2) } };
            auto const& labels = graph.edgeLabels(n1, n2);
            if (labels.size() > 0)
            {
                edge["labels"] = labels;
            }

            json["edges"].push_back(std::move(edge));
        }
    }
    return json;
}

GraphReferenceMapping parseReferenceMapping(Json const& jRefmap, Graph const& graph)
{
    Json::array_t nodes = jRefmap.value("nodes", Json::array());

    GraphReferenceMapping refmap(&graph);
    NodeId nodeIndex = 0;
    for (auto const& jNode : nodes)
    {
        const auto refInterval = jNode.find("reference");
        if (refInterval != jNode.end())
        {
            auto const region = ReferenceInterval::parseRegion(*refInterval);
            refmap.addMapping(nodeIndex, region);
        }
        ++nodeIndex;
    }
    return refmap;
}
}


text:
#include "Ancient_Code_Processor.h"
#include <functional>
//#include "textual/regex_additions.h"
//#include <boost/regex.hpp>

//using namespace textual;

namespace shading {

//  string Ancient_Code_Processor::process(Shader_Type type, const string source) {
//    auto included = process_includes(source, type);
//    return string("precision highp float;\n\n") + olden(included, type);
//  }
}

max_stars_repo_name: silentorb/mythic-cpp

max_issues_repo_name: silentorb/mythic-cpp

max_forks_repo_name: silentorb/mythic-cpp


text:
// Copyright 2018 The Simons Foundation, Inc. - All Rights Reserved.
//
// 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.

#include <Eigen/Dense>
#include <algorithm>
#include <cmath>
#include <iostream>
#include "Utils/random_utils.hpp"
#include <vector>
#include "Utils/json_utils.hpp"
#include "abstract_hilbert.hpp"

#ifndef NETKET_QUBITS_HPP
#define NETKET_QUBITS_HPP

namespace netket {

/**
  Hilbert space for qubits.
*/

class Qubit : public AbstractHilbert {
  std::vector<double> local_;

  int nqubits_;

 public:
  explicit Qubit(const json &pars) {
    int nqubits;

    if (!FieldExists(pars["Hilbert"], "Nqubits")) {
      std::cerr << "Nqubits is not defined" << std::endl;
    }

    nqubits = pars["Hilbert"]["Nqubits"];

    Init(nqubits);
  }

  void Init(int nqubits) {
    nqubits_ = nqubits;

    local_.resize(2);

    local_[0] = 0;
    local_[1] = 1;
  }

  bool IsDiscrete() const override { return true; }

  int LocalSize() const override { return 2; }

  int Size() const override { return nqubits_; }

  std::vector<double> LocalStates() const override { return local_; }

  void RandomVals(Eigen::VectorXd &state,
                  netket::default_random_engine &rgen) const override {
    std::uniform_int_distribution<int> distribution(0, 1);

    assert(state.size() == nqubits_);

    // unconstrained random
    for (int i = 0; i < state.size(); i++) {
      state(i) = distribution(rgen);
    }
  }

  void UpdateConf(Eigen::VectorXd &v, const std::vector<int> &tochange,
                  const std::vector<double> &newconf) const override {
    assert(v.size() == nqubits_);

    int i = 0;
    for (auto sf : tochange) {
      v(sf) = newconf[i];
      i++;
    }
  }
};

}  // namespace netket
#endif

max_stars_repo_name: artemborin/netket

max_issues_repo_name: artemborin/netket

max_forks_repo_name: artemborin/netket


text:
/*
 * Copyright (c) 2019 Opticks Team. All Rights Reserved.
 *
 * This file is part of Opticks
 * (see https://bitbucket.org/simoncblyth/opticks).
 *
 * 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.
 */

#include <cstring>
#include <iostream>
#include <sstream>

#include <boost/filesystem.hpp>
#include <boost/lexical_cast.hpp>

#include "BMap.hh"

#include "NLookup.hpp"
#include "PLOG.hh"

const plog::Severity NLookup::LEVEL = debug ;

namespace fs = boost::filesystem;

NLookup::NLookup()
   :
   m_alabel(NULL),
   m_blabel(NULL),
   m_closed(false)
{
}

const std::map<std::string, unsigned>& NLookup::getA() const
{
    return m_A ;
}
const std::map<std::string, unsigned>& NLookup::getB() const
{
    return m_B ;
}

void NLookup::mockup(const char* dir, const char* aname, const char* bname)
{
    mockA(dir, aname);
    mockB(dir, bname);
}
void NLookup::mockA(const char* adir, const char* aname)
{
    MSU mp ;
    mp["/dd/Materials/red"] = 10 ;
    mp["/dd/Materials/green"] = 20 ;
    mp["/dd/Materials/blue"] = 30 ;
    BMap<std::string, unsigned>::save(&mp, adir, aname );
}
void NLookup::mockB(const char* bdir, const char* bname)
{
    MSU mp;
    mp["red"] = 1 ;
    mp["green"] = 2 ;
    mp["blue"] = 3 ;
    BMap<std::string, unsigned>::save(&mp, bdir, bname );
}

void NLookup::setA( const char* json )
{
    BMap<std::string, unsigned>::LoadJSONString(&m_A, json, 0 );
    BMap<std::string, unsigned>::dump(&m_A, "NLookup::setA" );
}

void NLookup::dumpA(const char* msg)
{
    BMap<std::string, unsigned>::dump(&m_A, msg );
}
void NLookup::dumpB(const char* msg)
{
    BMap<std::string, unsigned>::dump(&m_B, msg );
}

void NLookup::setA( const std::map<std::string, unsigned>& A, const char* aprefix, const char* alabel )
{
    assert(!m_closed);

    m_alabel = strdup(alabel);
    LOG(debug) << "NLookup::setA " << alabel  ;

    for(MSU::const_iterator it=A.begin() ; it != A.end() ; it++)
    {
        const char* name = it->first.c_str() ;

        if(aprefix)
        {
            if(strncmp(name, aprefix, strlen(aprefix)) == 0)
            {
                std::string shortname = name + strlen(aprefix) ;
                m_A[shortname] = it->second ;
            }
        }
        else
        {
            std::string shortname = name  ;
            m_A[shortname] = it->second ;
        }
    }
}

void NLookup::setB( const std::map<std::string, unsigned>& B, const char* bprefix, const char* blabel)
{
    assert(!m_closed);

    m_blabel = strdup(blabel);
    LOG(debug) << "NLookup::setB " << blabel  ;
    for(MSU::const_iterator it=B.begin() ; it != B.end() ; it++)
    {
        const char* name = it->first.c_str() ;

        if(bprefix)
        {
            if(strncmp(name, bprefix, strlen(bprefix)) == 0)
            {
                std::string shortname = name + strlen(bprefix) ;
                m_B[shortname] = it->second ;
            }
        }
        else
        {
            std::string shortname = name  ;
            m_B[shortname] = it->second ;
        }
    }
}

std::string NLookup::brief() const
{
    std::stringstream ss ;

    ss
        << " A: " << ( m_alabel ? m_alabel : "NULL" ) << " " <<  m_A.size()
        << " B: " << ( m_blabel ? m_blabel : "NULL" ) << " " <<  m_B.size()
        ;

    return ss.str();
}

void NLookup::loadA(const char* adir, const char* aname, const char* aprefix)
{
    MSU A ;
    BMap<std::string, unsigned int>::load(&A, adir, aname );
    setA(A, aprefix, "NLookup::loadA");
}

void NLookup::loadB(const char* bdir, const char* bname, const char* /*bprefix*/)
{
    MSU B ;
    BMap<std::string, unsigned int>::load(&B, bdir, bname);
    setB(B, "", "NLookup::loadB");
}

bool NLookup::isClosed() const
{
    return m_closed ;
}

void NLookup::close(const char* msg)
{
    if(isClosed())
    {
        LOG(info)
                  << msg
                  << " CLOSED ALREADY : SKIPPING "
                  ;
        return ;

    }

    LOG(LEVEL) << msg << brief() ;

    bool complete = m_alabel && m_blabel ;

    if(!complete)
    {
        LOG(fatal)
           << " lookup missing A or B, cannot crossReference and close "
           << " msg : " << msg
           << " brief : " << brief()
           ;

        if(m_alabel) dumpA(m_alabel);
        if(m_blabel) dumpB(m_blabel);
    }

    assert(m_alabel) ; // have to setA and setB before close
    assert(m_blabel) ; // have to setA and setB before close

    crossReference();

    //dump(msg);
}

void NLookup::crossReference()
{
    assert(!m_closed);
    m_closed = true ;
    m_A2B = create(m_A, m_B);
    m_B2A = create(m_B, m_A);
}

std::string NLookup::acode2name(unsigned acode) const
{
    return find(m_A, acode);
}
std::string NLookup::bcode2name(unsigned bcode) const
{
    return find(m_B, bcode);
}

void NLookup::dump(const char* msg) const
{
    std::cout << msg
              << " A entries " <<  m_A.size()
              << " B entries " <<  m_B.size()
              << " A2B entries " <<  m_A2B.size()
              << std::endl
              ;

    for(NLookup_t::const_iterator it=m_A2B.begin() ; it != m_A2B.end() ; it++)
    {
        unsigned int acode = it->first ;
        unsigned int bcode = it->second ;
        std::string aname = acode2name(acode);
        std::string bname = bcode2name(bcode);
        printf("  A %4u : %25s  B %4u : %25s \n", acode, aname.c_str(), bcode, bname.c_str() );
    }
}

std::map<unsigned, unsigned> NLookup::create(MSU& a, MSU&b) const
{
    assert(m_closed);
    // cross referencing codes which correspond to the same names
    NLookup_t a2b_ ;
    for(MSU::const_iterator ia=a.begin() ; ia != a.end() ; ia++)
    {
        std::string aname = ia->first ;
        unsigned acode = ia->second ;
        if(b.find(aname) != b.end())
        {
            unsigned bcode = b.at(aname);
            a2b_[acode] = bcode ;
        }
    }
    return a2b_ ;
}

void NLookup::dumpMap(const char* msg, MSU& map) const
{
    printf("NLookup::dumpMap %s \n", msg);
    for(MSU::const_iterator it=map.begin() ; it != map.end() ; it++)
    {
        std::string name = it->first ;
        unsigned code = it->second ;
        printf("   %25s : %u \n", name.c_str(), code);
    }
}

std::string NLookup::find(const MSU& m, unsigned code) const
{
    std::string name ;
    for(MSU::const_iterator im=m.begin() ; im != m.end() ; im++)
    {
        unsigned mcode = im->second ;
        if(mcode == code)
        {
           name = im->first ;
           break;
        }
    }
    return name ;
}

int NLookup::a2b(unsigned a) const
{
    return lookup(m_A2B, a );
}
int NLookup::b2a(unsigned b) const
{
    return lookup(m_B2A, b );
}
int NLookup::lookup(const NLookup_t& lkup, unsigned x) const
{
    assert(m_closed);
    return lkup.find(x) == lkup.end() ? -1 : lkup.at(x) ;
}

max_stars_repo_name: hanswenzel/opticks

max_issues_repo_name: hanswenzel/opticks

max_forks_repo_name: hanswenzel/opticks


text:
/******************************************************************************
**
** Copyright (C) 2009-2011 Kyle Lutz <kyle.r.lutz@gmail.com>
** All rights reserved.
**
** This file is a part of the chemkit project. For more information
** see <http://www.chemkit.org>.
**
** Redistribution and use in source and binary forms, with or without
** modification, are permitted provided that the following conditions
** are met:
**
**   * Redistributions of source code must retain the above copyright
**     notice, this list of conditions and the following disclaimer.
**   * Redistributions in binary form must reproduce the above copyright
**     notice, this list of conditions and the following disclaimer in the
**     documentation and/or other materials provided with the distribution.
**   * Neither the name of the chemkit project nor the names of its
**     contributors may be used to endorse or promote products derived
**     from this software without specific prior written permission.
**
** THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
** "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
** LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
** A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
** OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
** SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
** LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
** DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
** THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
** (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
** OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
**
******************************************************************************/

#include "inchilineformat.h"

#include <string>
#include <vector>

#include <boost/thread/mutex.hpp>

#include <boost/algorithm/string.hpp>

#include "../../3rdparty/inchi/inchi_api.h"

#include <chemkit/atom.h>
#include <chemkit/bond.h>
#include <chemkit/foreach.h>

// BIL - the big inchi lock
//
// The InChI library is not reentrant. In order to safely use the
// library from multiple threads we introduce a single mutex which
// must be acquired to use the InChI library functions.
static boost::mutex BIL;

InchiLineFormat::InchiLineFormat()
    : chemkit::LineFormat("inchi")
{
}

chemkit::Molecule* InchiLineFormat::read(const std::string &formula)
{
    // verify formula
    if(formula.empty()){
        setErrorString("Formula is empty.");
        return 0;
    }

    // create molecule
    chemkit::Molecule *molecule = new chemkit::Molecule;

    std::string formulaString = formula;

    // add `InChI=` to the start if it is not there
    if(!boost::starts_with(formula,  "InChI=")){
        formulaString = "InChI=" + formulaString;
    }

    // setup input struct
    inchi_InputINCHI input;

    input.szInChI = const_cast<char *>(formulaString.c_str());
    input.szOptions = 0;

    // acquire inchi lock
    BIL.lock();

    // get inchi output
    inchi_OutputStruct output;
    int ret = GetStructFromStdINCHI(&input, &output);
    CHEMKIT_UNUSED(ret);

    // add atoms
    for(int i = 0; i < output.num_atoms; i++){
        molecule->addAtom(output.atom[i].elname);
    }

    // add bonds
    for(int i = 0; i < output.num_atoms; i++){
        inchi_Atom *inchiAtom = &output.atom[i];

        for(int j = 0; j < inchiAtom->num_bonds; j++){
            molecule->addBond(i, inchiAtom->neighbor[j], inchiAtom->bond_type[j]);
        }
    }

    // add implicit hydrogens (if enabled)
    bool addHydrogens = option("add-implicit-hydrogens").toBool();

    if(addHydrogens){
        for(int i = 0; i < output.num_atoms; i++){
            chemkit::Atom *atom = molecule->atom(i);
            inchi_Atom *inchiAtom = &output.atom[i];

            for(int j = 0; j < inchiAtom->num_iso_H[0]; j++){
                chemkit::Atom *hydrogen = molecule->addAtom(chemkit::Atom::Hydrogen);
                molecule->addBond(atom, hydrogen);
            }
        }
    }

    // free output structure
    FreeStructFromStdINCHI(&output);

    // release inchi lock
    BIL.unlock();

    return molecule;
}

std::string InchiLineFormat::write(const chemkit::Molecule *molecule)
{
    // check for valid molecule
    if(molecule->atomCount() > 1024){
        setErrorString("InChI does not support molecules with more that 1024 atoms.");
        return std::string();
    }

    // setup inchi input structure
    inchi_Input input;

    input.atom = new inchi_Atom[molecule->atomCount()];
    input.stereo0D = 0;
    input.szOptions = 0;
    input.num_atoms = molecule->atomCount();
    input.num_stereo0D = 0;

    inchi_Atom *inputAtom = &input.atom[0];

    std::vector<const chemkit::Atom *> chiralAtoms;

    foreach(const chemkit::Atom *atom, molecule->atoms()){

        // coordinates
        inputAtom->x = 0;
        inputAtom->y = 0;
        inputAtom->z = 0;

        // bonds and neighbors
        int neighborCount = 0;
        foreach(const chemkit::Bond *bond, atom->bonds()){
            const chemkit::Atom *neighbor = bond->otherAtom(atom);

            if(neighbor->index() < atom->index())
                continue;

            inputAtom->neighbor[neighborCount] = neighbor->index();
            inputAtom->bond_type[neighborCount] = bond->order();

            inputAtom->bond_stereo[neighborCount] = INCHI_BOND_STEREO_NONE;

            neighborCount++;
        }
        inputAtom->num_bonds = neighborCount;

        // element symbol
        strncpy(inputAtom->elname, atom->symbol().c_str(), ATOM_EL_LEN);

        // isotopic hydrogens
        inputAtom->num_iso_H[0] = -1;
        inputAtom->num_iso_H[1] = 0;
        inputAtom->num_iso_H[2] = 0;
        inputAtom->num_iso_H[3] = 0;

        // misc data
        inputAtom->isotopic_mass = 0;
        inputAtom->radical = 0;
        inputAtom->charge = 0;

        // chiral atoms
        if(atom->isChiral()){
            chiralAtoms.push_back(atom);
        }

        // move pointer to the next position in the atom array
        inputAtom++;
    }

    // count double bonds with stereochemistry
    std::vector<const chemkit::Bond *> stereogenicBonds;
    foreach(const chemkit::Bond *bond, molecule->bonds()){
        if(bond->order() == chemkit::Bond::Double &&
           bond->stereochemistry() != chemkit::Stereochemistry::None){
            stereogenicBonds.push_back(bond);
        }
    }

    // add stereochemistry if enabled
    bool stereochemistry = option("stereochemistry").toBool();

    if(stereochemistry){
        input.num_stereo0D = chiralAtoms.size() + stereogenicBonds.size();
        input.stereo0D = new inchi_Stereo0D[input.num_stereo0D];

        int chiralIndex = 0;

        foreach(const chemkit::Atom *atom, chiralAtoms){
            inchi_Stereo0D *stereo = &input.stereo0D[chiralIndex];
            memset(stereo, 0, sizeof(*stereo));
            stereo->central_atom = atom->index();

            int neighborIndex = 0;
            foreach(const chemkit::Atom *neighbor, atom->neighbors()){
                stereo->neighbor[neighborIndex++] = neighbor->index();
            }

            stereo->type = INCHI_StereoType_Tetrahedral;

            if(atom->chirality() == chemkit::Stereochemistry::R){
                stereo->parity = INCHI_PARITY_ODD;
            }
            else if(atom->chirality() == chemkit::Stereochemistry::S){
                stereo->parity = INCHI_PARITY_EVEN;
            }
            else if(atom->chirality() == chemkit::Stereochemistry::Unspecified){
                stereo->parity = INCHI_PARITY_UNDEFINED;
            }
            else{
                stereo->parity = INCHI_PARITY_UNKNOWN;
            }

            chiralIndex++;
        }

        foreach(const chemkit::Bond *bond, stereogenicBonds){
            inchi_Stereo0D *stereo = &input.stereo0D[chiralIndex];
            memset(stereo, 0, sizeof(*stereo));

            stereo->central_atom = NO_ATOM;
            stereo->type = INCHI_StereoType_DoubleBond;

            if(bond->stereochemistry() == chemkit::Stereochemistry::E){
                stereo->parity = INCHI_PARITY_EVEN;
            }
            else if(bond->stereochemistry() == chemkit::Stereochemistry::Z){
                stereo->parity = INCHI_PARITY_ODD;
            }
            else{
                stereo->parity = INCHI_PARITY_UNKNOWN;
            }

            // bond atoms
            stereo->neighbor[1] = bond->atom1()->index();
            stereo->neighbor[2] = bond->atom2()->index();

            // neighbor atoms
            const chemkit::Atom *highestPriorityNeighbor = 0;

            foreach(const chemkit::Atom *neighbor, bond->atom1()->neighbors()){
                if(neighbor == bond->atom2()){
                    continue;
                }

                if(!highestPriorityNeighbor){
                    highestPriorityNeighbor = neighbor;
                }
                else if(neighbor->atomicNumber() > highestPriorityNeighbor->atomicNumber()){
                    highestPriorityNeighbor = neighbor;
                }
            }

            stereo->neighbor[0] = highestPriorityNeighbor->index();

            highestPriorityNeighbor = 0;
            foreach(const chemkit::Atom *neighbor, bond->atom2()->neighbors()){
                if(neighbor == bond->atom1()){
                    continue;
                }

                if(!highestPriorityNeighbor){
                    highestPriorityNeighbor = neighbor;
                }
                else if(neighbor->atomicNumber() > highestPriorityNeighbor->atomicNumber()){
                    highestPriorityNeighbor = neighbor;
                }
            }

            stereo->neighbor[3] = highestPriorityNeighbor->index();

            chiralIndex++;
        }
    }
    else{
        input.stereo0D = 0;
        input.num_stereo0D = 0;
    }

    // acquire inchi lock
    BIL.lock();

    // create inchi generator object
    INCHIGEN_HANDLE generator = STDINCHIGEN_Create();

    // setup generator object
    INCHIGEN_DATA generatorData;
    STDINCHIGEN_Setup(generator, &generatorData, &input);

    // perform structure normalization
    STDINCHIGEN_DoNormalization(generator, &generatorData);

    // perform structure canonicalization
    STDINCHIGEN_DoCanonicalization(generator, &generatorData);

    // write inchi output structure
    inchi_Output output;
    STDINCHIGEN_DoSerialization(generator, &generatorData, &output);

    // get inchi string from output
    std::string inchiString;
    if(output.szInChI){
        inchiString = output.szInChI;
    }

    // destroy inchi input structure
    delete [] input.atom;
    delete [] input.stereo0D;

    // destroy inchi generator object
    STDINCHIGEN_Destroy(generator);

    // release inchi lock
    BIL.unlock();

    return inchiString;
}

chemkit::Variant InchiLineFormat::defaultOption(const std::string &name) const
{
    if(name == "stereochemistry")
        return true;
    else if(name == "add-implicit-hydrogens")
        return true;
    else
        return chemkit::Variant();
}

max_stars_repo_name: quizzmaster/chemkit

max_issues_repo_name: soplwang/chemkit

max_forks_repo_name: soplwang/chemkit


text:
// Copyright 2010-2015 RethinkDB, all rights reserved.
#ifndef CLUSTERING_QUERY_ROUTING_METADATA_HPP_
#define CLUSTERING_QUERY_ROUTING_METADATA_HPP_

#include <string>

#include "errors.hpp"
#include <boost/variant.hpp>

#include "clustering/generic/multi_client_metadata.hpp"
#include "clustering/generic/registration_metadata.hpp"
#include "concurrency/fifo_checker.hpp"
#include "concurrency/fifo_enforcer.hpp"
#include "rdb_protocol/protocol.hpp"
#include "rpc/mailbox/typed.hpp"

/* There is one `primary_query_bcard_t` per branch. It's created by the primary; parsers
use it to find the master. */

class primary_query_bcard_t {
public:
    class read_request_t {
    public:
        read_request_t() { }
        read_request_t(
                const read_t &r,
                order_token_t ot,
                fifo_enforcer_read_token_t ft,
                const mailbox_addr_t< void(boost::variant<read_response_t, cannot_perform_query_exc_t>)> &ca) :
            read(r), order_token(ot), fifo_token(ft), cont_addr(ca) { }
        read_t read;
        order_token_t order_token;
        fifo_enforcer_read_token_t fifo_token;
        mailbox_addr_t< void(boost::variant<read_response_t, cannot_perform_query_exc_t>)> cont_addr;
    };

    class write_request_t {
    public:
        write_request_t() { }
        write_request_t(
                const write_t &w,
                order_token_t ot,
                fifo_enforcer_write_token_t ft,
                const mailbox_addr_t< void(boost::variant<write_response_t, cannot_perform_query_exc_t>)> &ca) :
            write(w), order_token(ot), fifo_token(ft), cont_addr(ca) { }
        write_t write;
        order_token_t order_token;
        fifo_enforcer_write_token_t fifo_token;
        mailbox_addr_t< void(boost::variant<write_response_t, cannot_perform_query_exc_t>)> cont_addr;
    };

    typedef boost::variant< read_request_t, write_request_t > request_t;

    primary_query_bcard_t() { }
    primary_query_bcard_t(
            const region_t &r,
            const multi_client_business_card_t<request_t> &mc)
        : region(r), multi_client(mc) { }

    /* The region that this master covers */
    region_t region;

    /* Contact info for the master itself */
    multi_client_business_card_t<request_t> multi_client;
};

RDB_DECLARE_SERIALIZABLE_FOR_CLUSTER(primary_query_bcard_t::read_request_t);

RDB_DECLARE_SERIALIZABLE_FOR_CLUSTER(primary_query_bcard_t::write_request_t);

RDB_DECLARE_SERIALIZABLE(primary_query_bcard_t);
RDB_DECLARE_EQUALITY_COMPARABLE(primary_query_bcard_t);

/* Each replica exposes a `direct_query_bcard_t` for each shard that it is a primary or
secondary replica for. */

class direct_query_bcard_t {
public:
    typedef mailbox_t< void(
            read_t,
            mailbox_addr_t< void(read_response_t)>
            )> read_mailbox_t;

    direct_query_bcard_t() { }
    explicit direct_query_bcard_t(const read_mailbox_t::address_t &rm) : read_mailbox(rm) { }

    read_mailbox_t::address_t read_mailbox;
};

RDB_DECLARE_SERIALIZABLE(direct_query_bcard_t);
RDB_DECLARE_EQUALITY_COMPARABLE(direct_query_bcard_t);

/* `table_query_bcard_t` wraps the `primary_query_bcard_t` and/or `direct_query_bcard_t`
into a single object for convenience. */

class table_query_bcard_t {
public:
    region_t region;
    boost::optional<primary_query_bcard_t> primary;
    boost::optional<direct_query_bcard_t> direct;
};

RDB_DECLARE_SERIALIZABLE(table_query_bcard_t);
RDB_DECLARE_EQUALITY_COMPARABLE(table_query_bcard_t);

#endif /* CLUSTERING_QUERY_ROUTING_METADATA_HPP_ */

max_stars_repo_name: sauter-hq/rethinkdb

max_issues_repo_name: sauter-hq/rethinkdb

max_forks_repo_name: sauter-hq/rethinkdb


text:
/*
 * Copyright (C) 2015 Open Source Robotics Foundation
 *
 * 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.
 *
*/

#include <boost/bind.hpp>

#include "gazebo/gazebo_config.h"
#include "gazebo/common/Console.hh"

#include "gazebo/physics/Model.hh"
#include "gazebo/physics/Link.hh"
#include "gazebo/physics/ode/ODEFixedJoint.hh"

using namespace gazebo;
using namespace physics;

//////////////////////////////////////////////////
ODEFixedJoint::ODEFixedJoint(dWorldID _worldId, BasePtr _parent)
    : FixedJoint<ODEJoint>(_parent)
{
  this->jointId = dJointCreateFixed(_worldId, NULL);
}

//////////////////////////////////////////////////
ODEFixedJoint::~ODEFixedJoint()
{
}

//////////////////////////////////////////////////
void ODEFixedJoint::Load(sdf::ElementPtr _sdf)
{
  FixedJoint<ODEJoint>::Load(_sdf);
}

//////////////////////////////////////////////////
math::Vector3 ODEFixedJoint::GetAnchor(unsigned int /*index*/) const
{
  gzwarn << "SimbodyFixedJoint: called method "
         << "GetAnchor that is not valid for joints of type fixed.\n";
  return math::Vector3();
}

//////////////////////////////////////////////////
void ODEFixedJoint::SetAnchor(unsigned int /*index*/,
    const math::Vector3 &/*_anchor*/)
{
  if (this->childLink)
    this->childLink->SetEnabled(true);
  if (this->parentLink)
    this->parentLink->SetEnabled(true);
}

//////////////////////////////////////////////////
math::Vector3 ODEFixedJoint::GetGlobalAxis(unsigned int /*_index*/) const
{
  gzwarn << "SimbodyFixedJoint: called method "
         << "GetGlobalAxis that is not valid for joints of type fixed.\n";
  return math::Vector3();
}

//////////////////////////////////////////////////
void ODEFixedJoint::SetAxis(unsigned int /*_index*/,
                            const math::Vector3 &/*_axis*/)
{
  gzwarn << "ODEFixedJoint: called method "
         << "SetAxis that is not valid for joints of type fixed.\n";
}

//////////////////////////////////////////////////
math::Angle ODEFixedJoint::GetAngleImpl(unsigned int /*index*/) const
{
  gzwarn << "ODEFixedJoint: called method "
         << "GetAngleImpl that is not valid for joints of type fixed.\n";
  return math::Angle();
}

//////////////////////////////////////////////////
double ODEFixedJoint::GetVelocity(unsigned int /*index*/) const
{
  gzwarn << "ODEFixedJoint: called method "
         << "GetVelocity that is not valid for joints of type fixed.\n";
  return 0.0;
}

//////////////////////////////////////////////////
void ODEFixedJoint::SetVelocity(unsigned int /*_index*/, double /*_angle*/)
{
  gzwarn << "ODEFixedJoint: called method "
         << "SetVelocity that is not valid for joints of type fixed.\n";
}

//////////////////////////////////////////////////
void ODEFixedJoint::SetMaxForce(unsigned int /*index*/, double /*_t*/)
{
  gzwarn << "ODEFixedJoint: called method "
         << "SetMaxForce that is not valid for joints of type fixed.\n";
}

//////////////////////////////////////////////////
double ODEFixedJoint::GetMaxForce(unsigned int /*index*/)
{
  gzwarn << "ODEFixedJoint: called method "
         << "GetMaxForce that is not valid for joints of type fixed.\n";
  return 0.0;
}

//////////////////////////////////////////////////
void ODEFixedJoint::SetForceImpl(unsigned int /*_index*/, double /*_effort*/)
{
  gzwarn << "ODEFixedJoint: called method "
         << "SetForceImpl that is not valid for joints of type fixed.\n";
}

//////////////////////////////////////////////////
void ODEFixedJoint::Attach(LinkPtr _parent, LinkPtr _child)
{
    gazebo::physics::ODEJoint::Attach(_parent, _child);

    // Once we attach the links to the fixed joint,
    // we also call the dJointSetFixed method to
    // the current desired relative offset and desired
    // relative rotation between the bodies.

    dJointSetFixed(this->jointId);
}

max_stars_repo_name: harderthan/gazebo

max_issues_repo_name: harderthan/gazebo

max_forks_repo_name: harderthan/gazebo


text:
//  Copyright (c) 2007-2013 Hartmut Kaiser
//
//  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)

#if !defined(HPX_LCOS_SERVER_BARRIER_MAR_10_2010_0310PM)
#define HPX_LCOS_SERVER_BARRIER_MAR_10_2010_0310PM

#include <hpx/hpx_fwd.hpp>

#include <boost/version.hpp>
#include <boost/intrusive/slist.hpp>

#include <hpx/lcos/local/mutex.hpp>
#include <hpx/util/scoped_unlock.hpp>
#include <hpx/util/stringstream.hpp>
#include <hpx/runtime/threads/thread_data.hpp>
#include <hpx/runtime/threads/thread_helpers.hpp>
#include <hpx/runtime/components/component_type.hpp>
#include <hpx/runtime/components/server/managed_component_base.hpp>

///////////////////////////////////////////////////////////////////////////////
namespace hpx { namespace lcos { namespace server
{
    /// A barrier can be used to synchronize a specific number of threads,
    /// blocking all of the entering threads until all of the threads have
    /// entered the barrier.
    class barrier
      : public lcos::base_lco
      , public components::managed_component_base<barrier>
    {
    public:
        typedef lcos::base_lco base_type_holder;

    private:
        typedef components::managed_component_base<barrier> base_type;

        typedef hpx::lcos::local::mutex mutex_type;
        mutex_type mtx_;

        // define data structures needed for intrusive slist container used for
        // the queues
        struct barrier_queue_entry
        {
            typedef boost::intrusive::slist_member_hook<
                boost::intrusive::link_mode<boost::intrusive::normal_link>
            > hook_type;

            barrier_queue_entry(threads::thread_id_type id)
              : id_(id)
            {}

            threads::thread_id_type id_;
            hook_type slist_hook_;
        };

        typedef boost::intrusive::member_hook<
            barrier_queue_entry, barrier_queue_entry::hook_type,
            &barrier_queue_entry::slist_hook_
        > slist_option_type;

        typedef boost::intrusive::slist<
            barrier_queue_entry, slist_option_type,
            boost::intrusive::cache_last<true>,
            boost::intrusive::constant_time_size<false>
        > queue_type;

        struct reset_queue_entry
        {
            reset_queue_entry(barrier_queue_entry& e, queue_type& q)
              : e_(e), q_(q), last_(q.last())
            {}

            ~reset_queue_entry()
            {
                if (e_.id_)
                    q_.erase(last_);     // remove entry from queue
            }

            barrier_queue_entry& e_;
            queue_type& q_;
            queue_type::const_iterator last_;
        };

    public:
        // This is the component id. Every component needs to have an embedded
        // enumerator 'value' which is used by the generic action implementation
        // to associate this component with a given action.
        enum { value = components::component_barrier };

        barrier()
          : number_of_threads_(1)
        {}

        barrier(std::size_t number_of_threads)
          : number_of_threads_(number_of_threads)
        {}

        ~barrier()
        {
            if (!queue_.empty()) {
                LERR_(fatal) << "~barrier: thread_queue is not empty, aborting threads";

                mutex_type::scoped_lock l(mtx_);

                while (!queue_.empty()) {
                    threads::thread_id_type id = queue_.front().id_;
                    queue_.front().id_ = 0;
                    queue_.pop_front();

                    // we know that the id is actually the pointer to the thread
                    threads::thread_data_base* thrd =
                        static_cast<threads::thread_data_base*>(id);
                    LERR_(fatal) << "~barrier: pending thread: "
                            << get_thread_state_name(thrd->get_state())
                            << "(" << id << "): " << thrd->get_description();

                    // forcefully abort thread, do not throw
                    error_code ec(lightweight);
                    threads::set_thread_state(id, threads::pending,
                        threads::wait_abort, threads::thread_priority_default, ec);
                    if (ec) {
                        LERR_(fatal) << "~barrier: could not abort thread"
                            << get_thread_state_name(thrd->get_state())
                            << "(" << id << "): " << thrd->get_description();
                    }
                }
            }
        }

        // disambiguate base classes
        using base_type::finalize;
        typedef base_type::wrapping_type wrapping_type;

        static components::component_type get_component_type()
        {
            return components::component_barrier;
        }
        static void set_component_type(components::component_type) {}

        // standard LCO action implementations

        /// The function \a set_event will block the number of entering
        /// \a threads (as given by the constructor parameter \a number_of_threads),
        /// releasing all waiting threads as soon as the last \a thread
        /// entered this function.
        void set_event()
        {
            threads::thread_self& self = threads::get_self();

            mutex_type::scoped_lock l(mtx_);

            if (queue_.size() < number_of_threads_-1) {
                barrier_queue_entry e(self.get_thread_id());
                queue_.push_back(e);

                reset_queue_entry r(e, queue_);
                {
                    util::scoped_unlock<mutex_type::scoped_lock> ul(l);
                    this_thread::suspend(threads::suspended,
                        "barrier::set_event");
                }
            }
            else {
            // slist::swap has a bug in Boost 1.35.0
#if BOOST_VERSION < 103600
                // release the threads
                while (!queue_.empty()) {
                    threads::thread_id_type id = queue_.front().id_;
                    queue_.front().id_ = 0;
                    queue_.pop_front();
                    set_thread_state(id, threads::pending);
                }
#else
                // swap the list
                queue_type queue;
                queue.swap(queue_);
                l.unlock();

                // release the threads
                while (!queue.empty()) {
                    threads::thread_id_type id = queue.front().id_;
                    queue.front().id_ = 0;
                    queue.pop_front();
                    set_thread_state(id, threads::pending);
                }
#endif
            }
        }

        /// The \a function set_exception is called whenever a
        /// \a set_exception_action is applied on an instance of a LCO. This
        /// function just forwards to the virtual function \a set_exception, which
        /// is overloaded by the derived concrete LCO.
        ///
        /// \param e      [in] The exception encapsulating the error to report
        ///               to this LCO instance.
        void set_exception(boost::exception_ptr const& e)
        {
            try {
                mutex_type::scoped_lock l(mtx_);

                while (!queue_.empty()) {
                    threads::thread_id_type id = queue_.front().id_;
                    queue_.front().id_ = 0;
                    queue_.pop_front();

                    threads::set_thread_state(id, threads::pending,
                        threads::wait_abort);
                }

                boost::rethrow_exception(e);
            }
            catch (boost::exception const& be) {
                // rethrow again, but this time using the native hpx mechanics
                HPX_THROW_EXCEPTION(hpx::no_success, "barrier::set_exception",
                    boost::diagnostic_information(be));
            }
        }

        typedef
            hpx::components::server::create_component_action1<
                barrier
              , std::size_t
            >
            create_component_action;

    private:
        std::size_t const number_of_threads_;
        queue_type queue_;
    };
}}}

HPX_REGISTER_ACTION_DECLARATION(
    hpx::lcos::server::barrier::create_component_action
  , hpx_lcos_server_barrier_create_component_action
)

#endif

max_stars_repo_name: andreasbuhr/hpx

max_issues_repo_name: andreasbuhr/hpx

max_forks_repo_name: andreasbuhr/hpx


text:
#include "net/rpc/server.h"

#include <algorithm>
#include <cstddef>
#include <cstdint>
#include <functional>
#include <optional>
#include <system_error>
#include <utility>
#include <boost/icl/interval_set.hpp>
#include <boost/smart_ptr/intrusive_ptr.hpp>
#include <boost/smart_ptr/intrusive_ref_counter.hpp>
#include "base/logging.h"
#include "net/rpc/flatbuffers-handler.h"
#include "net/rpc/flatbuffers/list_generated.h"
#include "net/rpc/wire-structs.h"
#include "security/keys/zero.key.h"

namespace net {
namespace rpc {

class Server::PingHandler : public Handler {
public:
    void handle(
        std::vector<uint8_t> request,
        const security::Key &key,
        std::function<void(std::vector<uint8_t>)> callback) final;
};

class Server::ListHandler : public FlatbuffersHandler<
    ListHandler, fbs::ListRequest, fbs::ListResponse> {
public:
    ListHandler(const Server &server) : server_(server) {}

    void handle(
        const fbs::ListRequestT &,
        const security::Key &,
        std::function<void(fbs::ListResponseT)> callback);

private:
    const Server &server_;
};

class Server::Operation : public boost::intrusive_ref_counter<
    Operation, boost::thread_unsafe_counter> {
public:
    Operation(
        Server &server,
        const udp::endpoint &client_endpoint,
        uint64_t key_fingerprint,
        uint64_t request_id);
    Operation(const Operation &) = delete;
    Operation &operator=(const Operation &) = delete;
    ~Operation();

    void unpack();

    void set_client_endpoint(const udp::endpoint &client_endpoint) {
        client_endpoint_ = client_endpoint;
    }

private:
    void handle();
    void pack();
    void respond();

    Server &server_;
    steady_timer keep_alive_timer_;
    udp::endpoint client_endpoint_;
    uint64_t key_fingerprint_;
    uint64_t request_id_;
    std::optional<security::Key> key_;
    std::vector<uint8_t> request_;
    boost::icl::interval_set<size_t> request_intervals_;
    size_t request_size_ = std::numeric_limits<size_t>::max();
    std::vector<uint8_t> response_;
    std::vector<std::vector<uint8_t>> response_buffers_;
};

Server::Server(
    const any_io_executor &executor,
    const udp::endpoint &endpoint,
    const Options &options)
    : executor_(executor),
      options_(options),
      socket_(executor_, endpoint),
      receive_buffer_(
          std::make_unique<uint8_t[]>(options_.receive_buffer_size)) {
    handle("ping", std::make_unique<PingHandler>());
    handle("list", std::make_unique<ListHandler>(*this));
    add_key(security::keys::zero);
}

void Server::handle(std::string_view name, std::unique_ptr<Handler> handler) {
    handlers_.emplace(name, std::move(handler));
}

void Server::add_key(const security::Key &key) {
    keystore_.add(key);
}

void Server::receive() {
    socket_.async_receive_from(
        buffer(receive_buffer_.get(), options_.receive_buffer_size),
        receive_endpoint_,
        [this](std::error_code ec, size_t size) {
            if (ec) {
                LOG(error) << "async_receive_from failed: " << ec;
                if (ec == std::errc::operation_canceled) {
                    return;
                }
                receive();
                return;
            }
            receive_size_ = size;
            dispatch();
        });
}

void Server::dispatch() {
    if (receive_size_ < sizeof(wire::request::Header)) {
        receive();
        return;
    }
    const auto *request_header = reinterpret_cast<wire::request::Header *>(
        &receive_buffer_[0]);
    auto iter = operations_.find(
        {request_header->key_fingerprint, request_header->request_id});
    if (iter != operations_.end()) {
        iter->second->set_client_endpoint(receive_endpoint_);
        iter->second->unpack();
    } else {
        boost::intrusive_ptr<Operation>(new Operation(
            *this,
            receive_endpoint_,
            request_header->key_fingerprint,
            request_header->request_id))->unpack();
    }
    receive();
}

void Server::PingHandler::handle(
    std::vector<uint8_t> request,
    const security::Key &,
    std::function<void(std::vector<uint8_t>)> callback) {
    callback(std::move(request));
}

void Server::ListHandler::handle(
    const fbs::ListRequestT &,
    const security::Key &,
    std::function<void(fbs::ListResponseT)> callback) {
    fbs::ListResponseT response;
    for (const auto &pair : server_.handlers_) {
        response.methods.push_back(pair.first);
    }
    std::sort(response.methods.begin(), response.methods.end());
    callback(std::move(response));
}

Server::Operation::Operation(
    Server &server,
    const udp::endpoint &client_endpoint,
    uint64_t key_fingerprint,
    uint64_t request_id)
    : server_(server),
      keep_alive_timer_(server_.executor_),
      client_endpoint_(client_endpoint),
      key_fingerprint_(key_fingerprint),
      request_id_(request_id) {
    server_.operations_.emplace(
        OperationKey{key_fingerprint_, request_id_}, this);
}

Server::Operation::~Operation() {
    server_.operations_.erase({key_fingerprint_, request_id_});
}

void Server::Operation::unpack() {
    const auto *request_header = reinterpret_cast<wire::request::Header *>(
        &server_.receive_buffer_[0]);
    size_t offset = request_header->offset;
    size_t request_body_size =
        server_.receive_size_ - sizeof(wire::request::Header);
    if (request_.size() < offset + request_body_size) {
        request_.resize(offset + request_body_size);
    }
    security::NonceArray nonce{};
    *reinterpret_cast<uint64_t *>(&nonce[0]) = request_header->request_id;
    *reinterpret_cast<uint16_t *>(&nonce[8]) = request_header->offset;
    size_t decrypt_result;
    if (key_) {
        decrypt_result = key_->decrypt(
            {
                &server_.receive_buffer_[sizeof(wire::request::Header)],
                request_body_size,
            },
            &nonce[0],
            &request_[offset]);
        if (decrypt_result == std::numeric_limits<size_t>::max()) {
            return;
        }
    } else {
        std::vector<security::Key> keys;
        server_.keystore_.find(request_header->key_fingerprint, keys);
        decrypt_result = std::numeric_limits<size_t>::max();
        for (const security::Key &key : keys) {
            decrypt_result = key.decrypt(
                {
                    &server_.receive_buffer_[sizeof(wire::request::Header)],
                    request_body_size,
                },
                &nonce[0],
                &request_[offset]);
            if (decrypt_result != std::numeric_limits<size_t>::max()) {
                key_.emplace(key);
                break;
            }
        }
        if (decrypt_result == std::numeric_limits<size_t>::max()) {
            return;
        }
        keep_alive_timer_.expires_after(std::chrono::minutes(1));
        keep_alive_timer_.async_wait(
            [_ = boost::intrusive_ptr<Operation>(this)](std::error_code) {});
    }
    if (contains(request_intervals_,
                 boost::icl::interval<size_t>::right_open(0, request_size_))) {
        if (!(request_header->flags & wire::request::flags::partial)) {
            respond();
        }
        return;
    }
    request_intervals_ += boost::icl::interval<size_t>::right_open(
        offset, offset + decrypt_result);
    if (!(request_header->flags & wire::request::flags::partial)) {
        request_size_ = offset + decrypt_result;
    }
    if (contains(request_intervals_,
                 boost::icl::interval<size_t>::right_open(0, request_size_))) {
        handle();
    }
}

void Server::Operation::handle() {
    const uint8_t *request_body_end = &request_[request_size_];
    const uint8_t *request_end = request_body_end;
    while (request_end != &request_[0] && request_end[-1]) {
        --request_end;
    }
    std::string_view method(
        reinterpret_cast<const char *>(request_end),
        request_body_end - request_end);
    if (request_end != &request_[0]) {
        --request_end;
    }
    auto handler_iter = server_.handlers_.find(method);
    if (handler_iter == server_.handlers_.end()) {
        // TODO(iceboy): Invalid request?
        return;
    }
    request_.resize(request_end - &request_[0]);
    handler_iter->second->handle(
        std::move(request_),
        *key_,
        [this, _ = boost::intrusive_ptr<Operation>(this)](
            std::vector<uint8_t> response) {
            response_ = std::move(response);
            pack();
        });
}

void Server::Operation::pack() {
    security::NonceArray nonce{};
    *reinterpret_cast<uint64_t *>(&nonce[0]) = request_id_;
    nonce[11] = 1;
    size_t offset = 0;
    do {
        size_t size = std::min(response_.size() - offset,
                               server_.options_.send_fragment_size);
        response_buffers_.emplace_back(sizeof(wire::response::Header) + size +
                                       security::Key::encrypt_overhead);
        auto *response_header = reinterpret_cast<wire::response::Header *>(
            &response_buffers_.back()[0]);
        response_header->key_fingerprint = key_fingerprint_;
        response_header->request_id = request_id_;
        response_header->offset = static_cast<uint16_t>(offset);
        response_header->flags = 0;
        if (offset + size < response_.size()) {
            response_header->flags |= wire::response::flags::partial;
        }
        *reinterpret_cast<uint16_t *>(&nonce[8]) =
            static_cast<uint16_t>(offset);
        size_t encrypt_result = key_->encrypt(
            {&response_[offset], size},
            &nonce[0],
            &response_buffers_.back()[sizeof(wire::response::Header)]);
        response_buffers_.back().resize(
            sizeof(wire::response::Header) + encrypt_result);
        offset += size;
    } while (offset < response_.size());
    respond();
}

void Server::Operation::respond() {
    for (const std::vector<uint8_t> &response_buffer : response_buffers_) {
        server_.socket_.async_send_to(
            buffer(response_buffer),
            client_endpoint_,
            [_ = boost::intrusive_ptr<Operation>(this)](
                std::error_code, size_t) {});
    }
}

}  // namespace rpc
}  // namespace net