text:
//  (C) Copyright Jeremy Siek 2006
//  Distributed under the Boost Software License, Version 1.0. (See
//  accompanying file LICENSE_1_0.txt or copy at
//  http://www.boost.org/LICENSE_1_0.txt)

#ifndef BOOST_PROPERTY_SERIALIZE_HPP
#define BOOST_PROPERTY_SERIALIZE_HPP

#include <boost/pending/property.hpp>
#include <boost/serialization/is_bitwise_serializable.hpp>
#include <boost/serialization/base_object.hpp>
#include <boost/serialization/nvp.hpp>

namespace boost {
  template<class Archive>
  inline void serialize(Archive&, no_property&, const unsigned int) { }

  template<class Archive, class Tag, class T, class Base>
  void
  serialize(Archive& ar, property<Tag, T, Base>& prop,
            const unsigned int /*version*/)
  {
    ar & serialization::make_nvp( "property_value" , prop.m_value );
    ar & serialization::make_nvp( "property_base" , prop.m_base );
  }

#ifdef BOOST_GRAPH_USE_MPI

// Setting the serialization properties of boost::property<> and
// boost::no_property to is_bitwise_serializable, object_serializable,
// track_never only when BOOST_GRAPH_USE_MPI is defined is dubious.
//
// This changes the serialization format of these classes, and hence
// of boost::adjacency_list, depending on whether BOOST_GRAPH_USE_MPI
// is defined.
//
// These serialization properties should probably be set in either case.
//
// Unfortunately, doing that now will change the serialization format
// of boost::adjacency_list in the non-MPI case, and could potentially
// break software that reads files serialized with an older release.

  namespace mpi {

    // forward declaration, to avoid including mpi
    template<typename T> struct is_mpi_datatype;

    template<typename Tag, typename T, typename Base>
    struct is_mpi_datatype<property<Tag, T, Base> >
      : mpl::and_<is_mpi_datatype<T>,
                  is_mpi_datatype<Base> > { };
  }

  namespace serialization {
    template<typename Tag, typename T, typename Base>
    struct is_bitwise_serializable<property<Tag, T, Base> >
      : mpl::and_<is_bitwise_serializable<T>,
                  is_bitwise_serializable<Base> > { };

  template<typename Tag, typename T, typename Base>
  struct implementation_level<property<Tag, T, Base>  >
   : mpl::int_<object_serializable> {} ;

  template<typename Tag, typename T, typename Base>
  struct tracking_level<property<Tag, T, Base>  >
   : mpl::int_<track_never> {} ;

  }
#endif // BOOST_GRAPH_USE_MPI

} // end namespace boost

#ifdef BOOST_GRAPH_USE_MPI
namespace boost { namespace mpi {
    template<>
    struct is_mpi_datatype<boost::no_property> : mpl::true_ { };

} } // end namespace boost::mpi

BOOST_IS_BITWISE_SERIALIZABLE(boost::no_property)
BOOST_CLASS_IMPLEMENTATION(boost::no_property,object_serializable)
BOOST_CLASS_TRACKING(boost::no_property,track_never)
#endif // BOOST_GRAPH_USE_MPI

#endif // BOOST_PROPERTY_SERIALIZE_HPP

ext: hpp

max_stars_repo_name: cpp-pm/boost

max_issues_repo_name: cpp-pm/boost

max_forks_repo_name: cpp-pm/boost


text:
#pragma once
#ifndef TABLESTORE_CORE_IMPL_BATCH_WRITER_HPP
#define TABLESTORE_CORE_IMPL_BATCH_WRITER_HPP
/*
BSD 3-Clause License

Copyright (c) 2017, Alibaba Cloud
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

* Redistributions of source code must retain the above copyright notice, this
  list of conditions and the following disclaimer.

* Redistributions in binary form must reproduce the above copyright notice,
  this list of conditions and the following disclaimer in the documentation
  and/or other materials provided with the distribution.

* Neither the name of the copyright holder 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 HOLDER 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 "tablestore/core/batch_writer.hpp"
#include "tablestore/util/threading.hpp"
#include "tablestore/util/random.hpp"
#include "tablestore/util/optional.hpp"
#include "tablestore/util/logger.hpp"
#include <boost/atomic.hpp>
#include <boost/noncopyable.hpp>
#include <tr1/functional>
#include <deque>

namespace aliyun {
namespace tablestore {
namespace core {
namespace impl {

enum SingleRowType
{
    kPutRow,
    kUpdateRow,
    kDeleteRow,
};

template<typename T>
struct WriteTraits
{};

template<>
struct WriteTraits<BatchWriteRowRequest::Put>
{
    typedef BatchWriteRowRequest::Put TypeInBatchWriteRequest;
    typedef PutRowRequest SingleRowRequest;
    typedef PutRowResponse SingleRowResponse;
    typedef RowPutChange SingleRowChange;
    typedef std::tr1::function<
        void(SingleRowRequest&, util::Optional<OTSError>&, SingleRowResponse&)> Callback;
    static const SingleRowType kType = kPutRow;
};
template<>
struct WriteTraits<PutRowRequest> : public WriteTraits<BatchWriteRowRequest::Put>
{};
template<>
struct WriteTraits<PutRowResponse> : public WriteTraits<BatchWriteRowRequest::Put>
{};
template<>
struct WriteTraits<BatchWriteRowRequest::Update>
{
    typedef BatchWriteRowRequest::Update TypeInBatchWriteRequest;
    typedef UpdateRowRequest SingleRowRequest;
    typedef UpdateRowResponse SingleRowResponse;
    typedef RowUpdateChange SingleRowChange;
    typedef std::tr1::function<
        void(SingleRowRequest&, util::Optional<OTSError>&, SingleRowResponse&)> Callback;
    static const SingleRowType kType = kUpdateRow;
};
template<>
struct WriteTraits<UpdateRowRequest> : public WriteTraits<BatchWriteRowRequest::Update>
{};
template<>
struct WriteTraits<UpdateRowResponse> : public WriteTraits<BatchWriteRowRequest::Update>
{};
template<>
struct WriteTraits<BatchWriteRowRequest::Delete>
{
    typedef BatchWriteRowRequest::Delete TypeInBatchWriteRequest;
    typedef DeleteRowRequest SingleRowRequest;
    typedef DeleteRowResponse SingleRowResponse;
    typedef RowDeleteChange SingleRowChange;
    typedef std::tr1::function<
        void(SingleRowRequest&, util::Optional<OTSError>&, SingleRowResponse&)> Callback;
    static const SingleRowType kType = kDeleteRow;
};
template<>
struct WriteTraits<DeleteRowRequest> : public WriteTraits<BatchWriteRowRequest::Delete>
{};
template<>
struct WriteTraits<DeleteRowResponse> : public WriteTraits<BatchWriteRowRequest::Delete>
{};

class AsyncBatchWriter : public core::AsyncBatchWriter, private boost::noncopyable
{
public:
    typedef std::tr1::function<void(
        PutRowRequest&, util::Optional<OTSError>&, PutRowResponse&)> PutRowCallback;
    typedef std::tr1::function<void(
        UpdateRowRequest&, util::Optional<OTSError>&, UpdateRowResponse&)> UpdateRowCallback;
    typedef std::tr1::function<void(
        DeleteRowRequest&, util::Optional<OTSError>&, DeleteRowResponse&)> DeleteRowCallback;

    struct CallbackCarrier
    {
        std::deque<PutRowCallback> mPutCallbacks;
        std::deque<UpdateRowCallback> mUpdateCallbacks;
        std::deque<DeleteRowCallback> mDeleteCallbacks;
    };

    struct Item
    {
        enum ItemType
        {
            kInvalid,
            kPutRow,
            kUpdateRow,
            kDeleteRow,
        };

        ItemType mType;
        RowPutChange mRowPutChange;
        PutRowCallback mPutRowCallback;
        RowUpdateChange mRowUpdateChange;
        UpdateRowCallback mUpdateRowCallback;
        RowDeleteChange mRowDeleteChange;
        DeleteRowCallback mDeleteRowCallback;

        explicit Item()
          : mType(kInvalid)
        {}

        explicit Item(RowPutChange&, PutRowCallback&);
        explicit Item(RowUpdateChange&, UpdateRowCallback&);
        explicit Item(RowDeleteChange&, DeleteRowCallback&);

        explicit Item(const util::MoveHolder<Item>& a)
        {
            *this = a;
        }

        Item& operator=(const util::MoveHolder<Item>& a)
        {
            if (this == &*a) {
                return *this;
            }
            mType = a->mType;
            a->mType = kInvalid;
            util::moveAssign(mRowPutChange, util::move(a->mRowPutChange));
            util::moveAssign(mPutRowCallback, util::move(a->mPutRowCallback));
            util::moveAssign(mRowUpdateChange, util::move(a->mRowUpdateChange));
            util::moveAssign(mUpdateRowCallback, util::move(a->mUpdateRowCallback));
            util::moveAssign(mRowDeleteChange, util::move(a->mRowDeleteChange));
            util::moveAssign(mDeleteRowCallback, util::move(a->mDeleteRowCallback));
            return *this;
        }
    };

    class IContext
    {
    public:
        explicit IContext(SingleRowType tp)
          : mType(tp)
        {}

        virtual ~IContext()
        {}

        SingleRowType type() const;
        SingleRowType& mutableType();

    private:
        SingleRowType mType;
    };

    template<class Request>
    struct Context : public IContext
    {
        typedef typename WriteTraits<Request>::SingleRowResponse Response;
        typedef typename WriteTraits<Request>::Callback Callback;

        explicit Context();

        const Callback& callback() const;
        Callback& mutableCallback();
        const Request& request() const;
        Request& mutableRequest();
        const util::Optional<OTSError>& error() const;
        util::Optional<OTSError>& mutableError();
        const Response& response() const;
        Response& mutableResponse();

    private:
        Callback mCallback;
        Request mRequest;
        util::Optional<OTSError> mError;
        Response mResponse;
    };

public:
    explicit AsyncBatchWriter(
        AsyncClient& client,
        const BatchWriterConfig& cfg);
    virtual ~AsyncBatchWriter();

    void putRow(
        PutRowRequest&,
        const std::tr1::function<void(
            PutRowRequest&, util::Optional<OTSError>&, PutRowResponse&)>&);
    void updateRow(
        UpdateRowRequest&,
        const std::tr1::function<void(
            UpdateRowRequest&, util::Optional<OTSError>&, UpdateRowResponse&)>&);
    void deleteRow(
        DeleteRowRequest&,
        const std::tr1::function<void(
            DeleteRowRequest&, util::Optional<OTSError>&, DeleteRowResponse&)>&);

    // internal use only. tests require them to be exposed to public.
    void flush();
    std::tr1::tuple<util::Duration, int64_t> nextNapAndConcurrency(
        boost::atomic<bool>& backoff,
        int64_t maxConcurrency,
        util::Duration);

private:
    void callbackOnBatch(
        CallbackCarrier*,
        BatchWriteRowRequest&,
        util::Optional<OTSError>&,
        BatchWriteRowResponse&);
    void aggregator();
    void takeSomeNap(util::Duration);
    void batch(
        BatchWriteRowRequest&,
        CallbackCarrier&,
        int64_t& remains,
        int64_t& num);
    void send(int64_t concurrency);
    void waitAgain(CallbackCarrier&, BatchWriteRowRequest&);
    void prependWaitingList(std::deque<Item>&);
    void triggerCallback(const std::tr1::function<void()>& cb);

    template<class Request>
    void feedbackFromBatchReq(
        std::deque<Item>&,
        CallbackCarrier&,
        BatchWriteRowRequest&,
        BatchWriteRowResponse&);
    template<class Request>
    void feedbackAllError(
        CallbackCarrier&,
        BatchWriteRowRequest&,
        const OTSError&);
    template<class Request>
    void feedbackOkRequest(
        typename WriteTraits<Request>::SingleRowChange&,
        typename WriteTraits<Request>::Callback&,
        util::Optional<Row>&,
        const std::string& requestId,
        const std::string& traceId);
    template<class Request>
    void feedbackErrRequest(
        std::deque<Item>&,
        typename WriteTraits<Request>::SingleRowChange&,
        typename WriteTraits<Request>::Callback&,
        OTSError&,
        const std::string& requestId,
        const std::string& traceId);

    template<class Request>
    void issue(
        Request&,
        const typename WriteTraits<Request>::Callback&);

public:
    static int64_t sConcurrencyIncStep;
    static int64_t sDefaultActors;

private:
    std::auto_ptr<util::Logger> mLogger;
    core::AsyncClient& mClient;

    int64_t mMaxConcurrency;
    int64_t mMaxBatchSize;
    util::Duration mRegularNap;
    util::Duration mMaxNap;
    util::Duration mNapShrinkStep;
    std::deque<std::tr1::shared_ptr<util::Actor> > mActors;

    util::Thread mAggregateThread;
    util::Semaphore mAggregateSem;
    boost::atomic<bool> mExit;
    boost::atomic<int64_t> mOngoingRequests;
    std::auto_ptr<util::Random> mRng;
    util::Mutex mMutex;
    std::deque<Item> mWaitingList;
    boost::atomic<bool> mShouldBackoff;
    boost::atomic<uint64_t> mActorSelector;
};

} // namespace impl
} // namespace core
} // namespace tablestore
} // namespace aliyun

#endif

max_stars_repo_name: TimeExceed/aliyun-tablestore-cpp-sdk

max_issues_repo_name: TimeExceed/aliyun-tablestore-cpp-sdk

max_forks_repo_name: TimeExceed/aliyun-tablestore-cpp-sdk


text:
/*
 * Utility functions for opencv-stereo
 *
 * Copyright 2013, Andrew Barry <abarry@csail.mit.edu>
 *
 */

#ifndef OPENCV_STEREO_UTIL
#define OPENCV_STEREO_UTIL

#include <cv.h>
#include <highgui.h>
#include "opencv2/legacy/legacy.hpp"
#include "opencv2/opencv.hpp"

#include <sys/time.h>

#include <boost/filesystem.hpp>

#include <string>
#include <glib.h> // for configuration files

#include "lcmtypes/bot_core_image_t.h" // from libbot for images over LCM

#include "../../LCM/lcmt_stereo.h"
#include "../../utils/utils/RealtimeUtils.hpp"

extern "C"
{
    #include <stdio.h>
    #include <stdint.h>
    #include <stdlib.h>
    #include <inttypes.h>

    #include <dc1394/dc1394.h>

    #include <camera.h>
    #include <utils.h>

    #include "../../externals/jpeg-utils/jpeg-utils.h"
}

using namespace cv;

struct OpenCvStereoConfig
{
    uint64 guidLeft;
    uint64 guidRight;

    string lcmUrl;
    string stereoControlChannel;
    string calibrationDir;
    float calibrationUnitConversion;
    int lastValidPixelRow;
    string videoSaveDir;
    string fourcc;

    bool usePGM;

    string stereo_replay_channel;
    string baro_airspeed_channel;
    string pose_channel;
    string gps_channel;
    string battery_status_channel;
    string servo_out_channel;
    string optotrak_channel;
    string cpu_info_channel1;
    string cpu_info_channel2;
    string cpu_info_channel3;

    string log_size_channel1;
    string log_size_channel2;
    string log_size_channel3;

    int disparity;
    int infiniteDisparity;
    int interestOperatorLimit;
    int blockSize;
    int sadThreshold;
    float horizontalInvarianceMultiplier;

    int displayOffsetX;
    int displayOffsetY;

};

struct OpenCvStereoCalibration
{
    Mat mx1fp;
    Mat mx2fp;
    Mat qMat;

    Mat M1;
    Mat D1;
    Mat R1;
    Mat P1;

    Mat M2;
    Mat D2;
    Mat R2;
    Mat P2;
};

Mat GetFrameFormat7(dc1394camera_t *camera);

void FlushCameraBuffer(dc1394camera_t *camera);

int64_t getTimestampNow();

bool ParseConfigFile(string configFile, OpenCvStereoConfig *configStruct);

bool LoadCalibration(string calibrationDir, OpenCvStereoCalibration *stereoCalibration);

void StopCapture(dc1394_t *dcContext, dc1394camera_t *camera);

void InitBrightnessSettings(dc1394camera_t *camera1, dc1394camera_t *camera2, bool enable_gamma = false);

void MatchBrightnessSettings(dc1394camera_t *camera1, dc1394camera_t *camera2, bool complete_set = false, int force_brightness = -1, int force_exposure = -1);

void SendImageOverLcm(lcm_t* lcm, string channel, Mat image, int compression_quality = 80);

void Get3DPointsFromStereoMsg(const lcmt_stereo *msg, vector<Point3f> *points_out);

void Draw3DPointsOnImage(Mat camera_image, vector<Point3f> *points_list_in, Mat cam_mat_m, Mat cam_mat_d, Mat cam_mat_r, Scalar outline_color = 128, Scalar inside_color = 255, Point2d box_top = Point2d(-1, -1), Point2d box_bottom = Point2d(-1, -1), vector<int> *points_in_box = NULL, float min_z = 0, float max_z = 0, int box_size = 4);

int GetDisparityForDistance(double distance, const OpenCvStereoCalibration &calibration, int *inf_disparity = NULL);

#endif

max_stars_repo_name: titanomachy/flight

max_issues_repo_name: qiuwenhui/flight

max_forks_repo_name: qiuwenhui/flight


text:
#include "Common/Common.h"
#include "Demos/Visualization/MiniGL.h"
#include "Demos/Visualization/Selection.h"
#include "GL/glut.h"
#include "Demos/Simulation/TimeManager.h"
#include <Eigen/Dense>
#include "Demos/Simulation/SimulationModel.h"
#include "Demos/Simulation/TimeStepController.h"
#include <iostream>
#include "Demos/Utils/OBJLoader.h"
#include "Demos/Visualization/Visualization.h"
#include "Demos/Simulation/DistanceFieldCollisionDetection.h"
#include "Demos/Utils/Logger.h"
#include "Demos/Utils/Timing.h"
#include "Demos/Utils/FileSystem.h"

#define _USE_MATH_DEFINES
#include "math.h"

// Enable memory leak detection
#if defined(_DEBUG) && !defined(EIGEN_ALIGN)
	#define new DEBUG_NEW
#endif

INIT_TIMING
INIT_LOGGING

using namespace PBD;
using namespace Eigen;
using namespace std;

void timeStep ();
void buildModel ();
void createBodyModel();
void render ();
void reset();
void cleanup();
void initShader();
void selection(const Eigen::Vector2i &start, const Eigen::Vector2i &end);
void TW_CALL setTimeStep(const void *value, void *clientData);
void TW_CALL getTimeStep(void *value, void *clientData);
void TW_CALL setVelocityUpdateMethod(const void *value, void *clientData);
void TW_CALL getVelocityUpdateMethod(void *value, void *clientData);
void TW_CALL setMaxIterations(const void *value, void *clientData);
void TW_CALL getMaxIterations(void *value, void *clientData);
void TW_CALL setContactTolerance(const void *value, void *clientData);
void TW_CALL getContactTolerance(void *value, void *clientData);
void TW_CALL setContactStiffnessRigidBody(const void *value, void *clientData);
void TW_CALL getContactStiffnessRigidBody(void *value, void *clientData);
void TW_CALL setContactStiffnessParticleRigidBody(const void *value, void *clientData);
void TW_CALL getContactStiffnessParticleRigidBody(void *value, void *clientData);

SimulationModel model;
DistanceFieldCollisionDetection cd;
TimeStepController sim;

bool doPause = true;
std::vector<unsigned int> selectedBodies;
Vector3r oldMousePos;
Shader *shader;
Shader *shaderTex;
bool renderContacts = false;
string exePath;
string dataPath;

// main
int main( int argc, char **argv )
{
	REPORT_MEMORY_LEAKS

	std::string logPath = FileSystem::normalizePath(FileSystem::getProgramPath() + "/log");
	FileSystem::makeDirs(logPath);
	logger.addSink(unique_ptr<ConsoleSink>(new ConsoleSink(LogLevel::INFO)));
	logger.addSink(unique_ptr<FileSink>(new FileSink(LogLevel::DEBUG, logPath + "/PBD.log")));

	exePath = FileSystem::getProgramPath();
	dataPath = exePath + "/" + std::string(PBD_DATA_PATH);

	// OpenGL
	MiniGL::init (argc, argv, 1024, 768, 0, 0, "Rigid body collision demo");
	MiniGL::initLights ();
	MiniGL::initTexture();
	MiniGL::setClientIdleFunc (50, timeStep);
	MiniGL::setKeyFunc(0, 'r', reset);
	MiniGL::setSelectionFunc(selection);
	initShader();

	buildModel ();

	MiniGL::setClientSceneFunc(render);
	MiniGL::setViewport (40.0f, 0.1f, 500.0, Vector3r (5.0, 30.0, 70.0), Vector3r (5.0, 0.0, 0.0));

	TwAddVarRW(MiniGL::getTweakBar(), "Pause", TW_TYPE_BOOLCPP, &doPause, " label='Pause' group=Simulation key=SPACE ");
	TwAddVarRW(MiniGL::getTweakBar(), "RenderContacts", TW_TYPE_BOOLCPP, &renderContacts, " label='Render contacts' group=Simulation ");
	TwAddVarCB(MiniGL::getTweakBar(), "TimeStepSize", TW_TYPE_REAL, setTimeStep, getTimeStep, &model, " label='Time step size'  min=0.0 max = 0.1 step=0.001 precision=4 group=Simulation ");
	TwType enumType = TwDefineEnum("VelocityUpdateMethodType", NULL, 0);
	TwAddVarCB(MiniGL::getTweakBar(), "VelocityUpdateMethod", enumType, setVelocityUpdateMethod, getVelocityUpdateMethod, &sim, " label='Velocity update method' enum='0 {First Order Update}, 1 {Second Order Update}' group=Simulation");
	TwAddVarCB(MiniGL::getTweakBar(), "MaxIter", TW_TYPE_UINT32, setMaxIterations, getMaxIterations, &sim, " label='Max. iterations'  min=1 step=1 group=Simulation ");
	TwAddVarCB(MiniGL::getTweakBar(), "ContactTolerance", TW_TYPE_REAL, setContactTolerance, getContactTolerance, &cd, " label='Contact tolerance'  min=0.0 step=0.001 precision=3 group=Simulation ");
	TwAddVarCB(MiniGL::getTweakBar(), "ContactStiffnessRigidBody", TW_TYPE_REAL, setContactStiffnessRigidBody, getContactStiffnessRigidBody, &model, " label='Contact stiffness RB'  min=0.0 step=0.1 precision=2 group=Simulation ");
	TwAddVarCB(MiniGL::getTweakBar(), "ContactStiffnessParticleRigidBody", TW_TYPE_REAL, setContactStiffnessParticleRigidBody, getContactStiffnessParticleRigidBody, &model, " label='Contact stiffness Particle-RB'  min=0.0 step=0.1 precision=2 group=Simulation ");

	glutMainLoop ();

	cleanup ();

	Timing::printAverageTimes();

	return 0;
}

void initShader()
{
	std::string vertFile = dataPath + "/shaders/vs_smooth.glsl";
	std::string fragFile = dataPath + "/shaders/fs_smooth.glsl";
	shader = MiniGL::createShader(vertFile, "", fragFile);

	if (shader == NULL)
		return;

	shader->begin();
	shader->addUniform("modelview_matrix");
	shader->addUniform("projection_matrix");
	shader->addUniform("surface_color");
	shader->addUniform("shininess");
	shader->addUniform("specular_factor");
	shader->end();

	std::string vertFileTex = dataPath + "/shaders/vs_smoothTex.glsl";
	std::string fragFileTex = dataPath + "/shaders/fs_smoothTex.glsl";
	shaderTex = MiniGL::createShader(vertFileTex, "", fragFileTex);

	if (shaderTex == NULL)
		return;

	shaderTex->begin();
	shaderTex->addUniform("modelview_matrix");
	shaderTex->addUniform("projection_matrix");
	shaderTex->addUniform("surface_color");
	shaderTex->addUniform("shininess");
	shaderTex->addUniform("specular_factor");
	shaderTex->end();
}

void cleanup()
{
	delete TimeManager::getCurrent();
	delete shader;
	delete shaderTex;
}

void reset()
{
	Timing::printAverageTimes();
	Timing::reset();

	model.reset();
	sim.reset();
	TimeManager::getCurrent()->setTime(0.0);
}

void mouseMove(int x, int y)
{
	Vector3r mousePos;
	MiniGL::unproject(x, y, mousePos);
	const Vector3r diff = mousePos - oldMousePos;

	TimeManager *tm = TimeManager::getCurrent();
	const Real h = tm->getTimeStepSize();

	SimulationModel::RigidBodyVector &rb = model.getRigidBodies();
	for (size_t j = 0; j < selectedBodies.size(); j++)
	{
		if (rb[selectedBodies[j]]->getMass() != 0.0)
			rb[selectedBodies[j]]->getVelocity() += 1.0 / h * diff;
	}
	oldMousePos = mousePos;
}

void selection(const Eigen::Vector2i &start, const Eigen::Vector2i &end)
{
	std::vector<unsigned int> hits;
	selectedBodies.clear();

	SimulationModel::RigidBodyVector &rb = model.getRigidBodies();
	std::vector<Vector3r, Eigen::aligned_allocator<Vector3r> > x;
	x.resize(rb.size());
	for (unsigned int i = 0; i < rb.size(); i++)
	{
		x[i] = rb[i]->getPosition();
	}

	Selection::selectRect(start, end, &x[0], &x[rb.size() - 1], selectedBodies);
	if (selectedBodies.size() > 0)
		MiniGL::setMouseMoveFunc(GLUT_MIDDLE_BUTTON, mouseMove);
	else
		MiniGL::setMouseMoveFunc(-1, NULL);

	MiniGL::unproject(end[0], end[1], oldMousePos);
}

void timeStep ()
{
	if (doPause)
		return;

	// Simulation code
	for (unsigned int i = 0; i < 8; i++)
		sim.step(model);
}

void buildModel ()
{
	TimeManager::getCurrent ()->setTimeStepSize (0.005);

	sim.setCollisionDetection(model, &cd);

	createBodyModel();
}

void renderBallJoint(BallJoint &bj)
{
	float jointColor[4] = { 0.0, 0.6f, 0.2f, 1 };
	MiniGL::drawSphere(bj.m_jointInfo.col(2), 0.15f, jointColor);
}

void renderContact(RigidBodyContactConstraint &cc)
{
	float col1[4] = { 0.0, 0.6f, 0.2f, 1 };
	float col2[4] = { 0.6f, 0.0, 0.2f, 1 };
	MiniGL::drawPoint(cc.m_constraintInfo.col(0), 5.0, col1);
	MiniGL::drawPoint(cc.m_constraintInfo.col(1), 5.0, col2);
	MiniGL::drawVector(cc.m_constraintInfo.col(1), cc.m_constraintInfo.col(1) + cc.m_constraintInfo.col(2), 1.0, col2);
}

void render ()
{
	MiniGL::coordinateSystem();

	// Draw sim model

	SimulationModel::RigidBodyVector &rb = model.getRigidBodies();
	SimulationModel::ConstraintVector &constraints = model.getConstraints();
	SimulationModel::RigidBodyContactConstraintVector &contacts = model.getRigidBodyContactConstraints();

	float selectionColor[4] = { 0.8f, 0.0f, 0.0f, 1 };
	float surfaceColor[4] = { 0.3f, 0.5f, 0.8f, 1 };
	float staticColor[4] = { 0.5f, 0.5f, 0.5f, 1 };

	if (renderContacts)
	{
		for (unsigned int i = 0; i < contacts.size(); i++)
			renderContact(contacts[i]);
	}

	if (shader)
	{
		shader->begin();
		glUniform1f(shader->getUniform("shininess"), 5.0);
		glUniform1f(shader->getUniform("specular_factor"), 0.2f);

		GLfloat matrix[16];
		glGetFloatv(GL_MODELVIEW_MATRIX, matrix);
		glUniformMatrix4fv(shader->getUniform("modelview_matrix"), 1, GL_FALSE, matrix);
		GLfloat pmatrix[16];
		glGetFloatv(GL_PROJECTION_MATRIX, pmatrix);
		glUniformMatrix4fv(shader->getUniform("projection_matrix"), 1, GL_FALSE, pmatrix);
	}

	for (size_t i = 0; i < rb.size(); i++)
	{
		bool selected = false;
		for (unsigned int j = 0; j < selectedBodies.size(); j++)
		{
			if (selectedBodies[j] == i)
				selected = true;
		}

		const VertexData &vd = rb[i]->getGeometry().getVertexData();
		const IndexedFaceMesh &mesh = rb[i]->getGeometry().getMesh();
		if (rb[i]->getMass() != 0.0)
		{
			if (!selected)
			{
				if (shader)
					glUniform3fv(shader->getUniform("surface_color"), 1, surfaceColor);
				Visualization::drawMesh(vd, mesh, 0, surfaceColor);
			}
			else
			{
				if (shader)
					glUniform3fv(shader->getUniform("surface_color"), 1, selectionColor);
				Visualization::drawMesh(vd, mesh, 0, selectionColor);
			}
		}
	}

	if (shader)
		shader->end();

	if (shaderTex)
	{
		shaderTex->begin();
		glUniform1f(shaderTex->getUniform("shininess"), 5.0);
		glUniform1f(shaderTex->getUniform("specular_factor"), 0.2f);

		GLfloat matrix[16];
		glGetFloatv(GL_MODELVIEW_MATRIX, matrix);
		glUniformMatrix4fv(shaderTex->getUniform("modelview_matrix"), 1, GL_FALSE, matrix);
		GLfloat pmatrix[16];
		glGetFloatv(GL_PROJECTION_MATRIX, pmatrix);
		glUniformMatrix4fv(shaderTex->getUniform("projection_matrix"), 1, GL_FALSE, pmatrix);
	}

	for (size_t i = 0; i < rb.size(); i++)
	{
		bool selected = false;
		for (unsigned int j = 0; j < selectedBodies.size(); j++)
		{
			if (selectedBodies[j] == i)
				selected = true;
		}

		const VertexData &vd = rb[i]->getGeometry().getVertexData();
		const IndexedFaceMesh &mesh = rb[i]->getGeometry().getMesh();
		if (rb[i]->getMass() == 0.0)
		{
			if (!selected)
			{
				if (shaderTex)
					glUniform3fv(shaderTex->getUniform("surface_color"), 1, staticColor);
				Visualization::drawTexturedMesh(vd, mesh, 0, staticColor);
			}
			else
			{
				if (shaderTex)
					glUniform3fv(shaderTex->getUniform("surface_color"), 1, selectionColor);
				Visualization::drawTexturedMesh(vd, mesh, 0, selectionColor);
			}
		}
	}

	if (shaderTex)
		shaderTex->end();

	for (size_t i = 0; i < constraints.size(); i++)
	{
		if (constraints[i]->getTypeId() == BallJoint::TYPE_ID)
		{
			renderBallJoint(*(BallJoint*)constraints[i]);
		}
	}

	MiniGL::drawTime( TimeManager::getCurrent ()->getTime ());
}

/** Create the rigid body model
*/
void createBodyModel()
{
	SimulationModel::RigidBodyVector &rb = model.getRigidBodies();
	SimulationModel::ConstraintVector &constraints = model.getConstraints();

	string fileNameBox = FileSystem::normalizePath(dataPath + "/models/cube.obj");
	IndexedFaceMesh meshBox;
	VertexData vdBox;
	OBJLoader::loadObj(fileNameBox, vdBox, meshBox);

	string fileNameCylinder = FileSystem::normalizePath(dataPath + "/models/cylinder.obj");
	IndexedFaceMesh meshCylinder;
	VertexData vdCylinder;
	OBJLoader::loadObj(fileNameCylinder, vdCylinder, meshCylinder);

	string fileNameTorus = FileSystem::normalizePath(dataPath + "/models/torus.obj");
	IndexedFaceMesh meshTorus;
	VertexData vdTorus;
	OBJLoader::loadObj(fileNameTorus, vdTorus, meshTorus);

	string fileNameCube = FileSystem::normalizePath(dataPath + "/models/cube_5.obj");
	IndexedFaceMesh meshCube;
	VertexData vdCube;
	OBJLoader::loadObj(fileNameCube, vdCube, meshCube);

	string fileNameSphere = FileSystem::normalizePath(dataPath + "/models/sphere.obj");
	IndexedFaceMesh meshSphere;
	VertexData vdSphere;
	OBJLoader::loadObj(fileNameSphere, vdSphere, meshSphere);

	const unsigned int num_piles_x = 5;
	const unsigned int num_piles_z = 5;
	const Real dx_piles = 4.0;
	const Real dz_piles = 4.0;
	const Real startx_piles = -0.5 * (Real)(num_piles_x - 1)*dx_piles;
	const Real startz_piles = -0.5 * (Real)(num_piles_z - 1)*dz_piles;
	const unsigned int num_piles = num_piles_x * num_piles_z;
	const unsigned int num_bodies_x = 3;
	const unsigned int num_bodies_y = 5;
	const unsigned int num_bodies_z = 3;
	const Real dx_bodies = 6.0;
	const Real dy_bodies = 6.0;
	const Real dz_bodies = 6.0;
	const Real startx_bodies = -0.5 * (Real)(num_bodies_x - 1)*dx_bodies;
	const Real starty_bodies = 14.0;
	const Real startz_bodies = -0.5 * (Real)(num_bodies_z - 1)*dz_bodies;
	const unsigned int num_bodies = num_bodies_x * num_bodies_y * num_bodies_z;
	rb.resize(num_piles + num_bodies + 1);
	unsigned int rbIndex = 0;

	// floor
	rb[rbIndex] = new RigidBody();
	rb[rbIndex]->initBody(1.0,
		Vector3r(0.0, -0.5, 0.0),
		Quaternionr(1.0, 0.0, 0.0, 0.0),
		vdBox, meshBox, Vector3r(100.0, 1.0, 100.0));
	rb[rbIndex]->setMass(0.0);

	const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
	const unsigned int nVert = static_cast<unsigned int>(vertices->size());

	cd.addCollisionBox(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, Vector3r(100.0, 1.0, 100.0));
	rbIndex++;

	Real current_z = startz_piles;
	for (unsigned int i = 0; i < num_piles_z; i++)
	{
		Real current_x = startx_piles;
		for (unsigned int j = 0; j < num_piles_x; j++)
		{
			rb[rbIndex] = new RigidBody();
			rb[rbIndex]->initBody(100.0,
				Vector3r(current_x, 5.0, current_z),
				Quaternionr(1.0, 0.0, 0.0, 0.0),
				vdCylinder, meshCylinder,
				Vector3r(0.5, 10.0, 0.5));
			rb[rbIndex]->setMass(0.0);

			const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
			const unsigned int nVert = static_cast<unsigned int>(vertices->size());
			cd.addCollisionCylinder(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, Vector2r(0.5, 10.0));
			current_x += dx_piles;
			rbIndex++;
		}
		current_z += dz_piles;
	}

	Real current_y = starty_bodies;
	unsigned int currentType = 0;
	for (unsigned int i = 0; i < num_bodies_y; i++)
	{
		Real current_x = startx_bodies;
		for (unsigned int j = 0; j < num_bodies_x; j++)
		{
			Real current_z = startz_bodies;
			for (unsigned int k = 0; k < num_bodies_z; k++)
			{
				rb[rbIndex] = new RigidBody();

				Real ax = static_cast <Real> (rand()) / static_cast <Real> (RAND_MAX);
				Real ay = static_cast <Real> (rand()) / static_cast <Real> (RAND_MAX);
				Real az = static_cast <Real> (rand()) / static_cast <Real> (RAND_MAX);
				Real w = static_cast <Real> (rand()) / static_cast <Real> (RAND_MAX);
				Quaternionr q(w, ax, ay, az);
				q.normalize();

				currentType = rand() % 4;
				if (currentType == 0)
				{
					rb[rbIndex]->initBody(100.0,
						Vector3r(current_x, current_y, current_z),
						q, //Quaternionr(1.0, 0.0, 0.0, 0.0),
						vdTorus, meshTorus,
						Vector3r(2.0, 2.0, 2.0));

					const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
					const unsigned int nVert = static_cast<unsigned int>(vertices->size());
					cd.addCollisionTorus(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, Vector2r(2.0, 1.0));
				}
				else if (currentType == 1)
				{
					rb[rbIndex]->initBody(100.0,
						Vector3r(current_x, current_y, current_z),
						q, //Quaternionr(1.0, 0.0, 0.0, 0.0),
						vdCube, meshCube,
						Vector3r(4.0, 1.0, 1.0));

					const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
					const unsigned int nVert = static_cast<unsigned int>(vertices->size());
					cd.addCollisionBox(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, Vector3r(4.0, 1.0, 1.0));
				}
				else if (currentType == 2)
				{
					rb[rbIndex]->initBody(100.0,
						Vector3r(current_x, current_y, current_z),
						q, //Quaternionr(1.0, 0.0, 0.0, 0.0),
						vdSphere, meshSphere);

					const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
					const unsigned int nVert = static_cast<unsigned int>(vertices->size());
					cd.addCollisionSphere(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, 1.0);
				}
				else if (currentType == 3)
				{
					rb[rbIndex]->initBody(100.0,
						Vector3r(current_x, current_y, current_z),
						q, //Quaternionr(1.0, 0.0, 0.0, 0.0),
						vdCylinder, meshCylinder,
						Vector3r(0.75, 5.0, 0.75));

					const std::vector<Vector3r> *vertices = rb[rbIndex]->getGeometry().getVertexDataLocal().getVertices();
					const unsigned int nVert = static_cast<unsigned int>(vertices->size());
					cd.addCollisionCylinder(rbIndex, CollisionDetection::CollisionObject::RigidBodyCollisionObjectType, &(*vertices)[0], nVert, Vector2r(0.75, 5.0));
			}
				currentType = (currentType + 1) % 4;
				current_z += dz_bodies;
				rbIndex++;
			}
			current_x += dx_bodies;
		}
		current_y += dy_bodies;
	}
}

void TW_CALL setTimeStep(const void *value, void *clientData)
{
	const Real val = *(const Real *)(value);
	TimeManager::getCurrent()->setTimeStepSize(val);
}

void TW_CALL getTimeStep(void *value, void *clientData)
{
	*(Real *)(value) = TimeManager::getCurrent()->getTimeStepSize();
}

void TW_CALL setVelocityUpdateMethod(const void *value, void *clientData)
{
	const short val = *(const short *)(value);
	((TimeStepController*)clientData)->setVelocityUpdateMethod((unsigned int)val);
}

void TW_CALL getVelocityUpdateMethod(void *value, void *clientData)
{
	*(short *)(value) = (short)((TimeStepController*)clientData)->getVelocityUpdateMethod();
}

void TW_CALL setMaxIterations(const void *value, void *clientData)
{
	const unsigned int val = *(const unsigned int *)(value);
	((TimeStepController*)clientData)->setMaxIterations(val);
}

void TW_CALL getMaxIterations(void *value, void *clientData)
{
	*(unsigned int *)(value) = ((TimeStepController*)clientData)->getMaxIterations();
}

void TW_CALL setContactStiffnessRigidBody(const void *value, void *clientData)
{
	const Real val = *(const Real *)(value);
	((SimulationModel*)clientData)->setContactStiffnessRigidBody(val);
}

void TW_CALL getContactStiffnessRigidBody(void *value, void *clientData)
{
	*(Real *)(value) = ((SimulationModel*)clientData)->getContactStiffnessRigidBody();
}

void TW_CALL setContactStiffnessParticleRigidBody(const void *value, void *clientData)
{
	const Real val = *(const Real *)(value);
	((SimulationModel*)clientData)->setContactStiffnessParticleRigidBody(val);
}

void TW_CALL getContactStiffnessParticleRigidBody(void *value, void *clientData)
{
	*(Real *)(value) = ((SimulationModel*)clientData)->getContactStiffnessParticleRigidBody();
}

void TW_CALL setContactTolerance(const void *value, void *clientData)
{
	const Real val = *(const Real *)(value);
	((DistanceFieldCollisionDetection*)clientData)->setTolerance(val);
}

void TW_CALL getContactTolerance(void *value, void *clientData)
{
	*(Real *)(value) = ((DistanceFieldCollisionDetection*)clientData)->getTolerance();
}

ext: cpp


text:
/**
 * Includes
 */
#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <iomanip>
#include <complex>
#include <typeinfo>
#include <boost/math/tools/polynomial.hpp>
#include <boost/program_options/options_description.hpp>
#include <boost/program_options/parsers.hpp>
#include <boost/program_options/variables_map.hpp>
#include <boost/lexical_cast.hpp>
#include "string.h"
/**
 * Compiler const
 */
#define deflationMethod 1
#define deflationMethod2 2
#define mullerMethod 3
#define laguerreMethod 4
#define debug

namespace po = boost::program_options;
namespace anpi {
	/**
	 * Deflation class
	 */
	template<typename T> class deflation {
		public :
			/**
			 * Deflation class operator () overload
			 * @param poly: the polynomial array it can be almost any numerical type.
			 * @param root: a known root of the polynomial, has to be the same type as poly.
			 * @param order: the size of the polynomial -1, has to be a int.
			 * @return
			 */
			inline std::vector<T> operator () ( std::vector<T>& poly, T& root, T& res, int order) {
				res =poly[order];
				poly[order]=0;
				#ifdef debug
					std::cout << "Deflation method started." << std::endl;
				#endif
				for (int i=order-1;i>=0;i--){
					T s = poly[i];
					poly[i]=res;
					res=s+res*root;
					std::cout << poly[i] << "x^" << i;
					if(i!=0) std::cout << " + ";
				}
				std::cout << "\n";
				#ifdef debug
					std::cout << "Deflation method ended." << std::endl;
				#endif
				return poly;
			}
	};
	/**
	 * Deflation2 class, for the case that both complex conjugate roots have to be deleted
	 */
	template<typename T> class deflation2 {
		public :
			/**
			 * Deflation2 class operator () overload
			 * @param poly: the polynomial array, has to be a complex<T> type.
			 * @param root: a known root of the polynomial, has to be  the same type as poly.
			 * @param order: the size of the polynomial -1, has to be a int.
			 * @return
			 */
			inline std::vector<T> operator () (std::vector<T>& poly, T& root, T& res, int order) {
				res =poly[order];
				poly[order]=0;
				#ifdef debug
					std::cout << "Deflation method 2 started." << std::endl;
				#endif
				for (int i=order-1;i>=0;i--){
					T s = poly[i];
					poly[i]=res;
					res=s+res*root;
				}
				res =poly[--order];
				root=std::conj(root);
				poly[order]=0;
				for (int i=order-1;i>=0;i--){
					T s = poly[i];
					poly[i]=res;
					res=s+res*root;
					std::cout << poly[i] << "x^" << i;
					if(i!=0) std::cout << " + ";
				}
				std::cout << "\n";
				#ifdef debug
					std::cout << "Deflation method 2 ended." << std::endl;
				#endif
				return poly;
			}
	};
	/**
	 * Muller class
	 */
	template <typename T>
	class muller{
		public:
			/**
			 * Muller class operator () overload
			 * @param x0:
			 * @param x1:
			 * @param x2:
			 * @param coefs:
			 * @param max:
			 * @param pTolerance:
			 * @return
			 */
			std::vector<std::complex<T> > operator()(T x0, T x1, T x2,const std::vector< std::complex<T> >& coefs,unsigned int max, double pTolerance) {
			std::vector<std::complex<T> > res;
			std::vector<std::complex<T> > q = coefs;
			std::complex<T> z;
			while (q.size() > 2) {
			    z = polinom_zero(x0, x1, x2, q, max, pTolerance);
			    //z = polinom_zero(p, z);
			    q = horner(q, z);
			    res.push_back(z);
			  }
			  res.push_back(-q[0] / q[1]);
			  return res;
			}

			std::complex<T> evaluatePol(std::vector<std::complex<T> >& coef, std::complex<T> * x){
				std::complex<T> result = 0;
				for (int ite = 0; ite<coef.size(); ite++){
					result+= coef[ite]*std::pow(*x, ite);
				}
				return result;
			}

			std::vector<std::complex<T> > horner(const std::vector<std::complex<T> > &a, std::complex<T> x0) {
			  int n = a.size();
			  std::vector<std::complex<T> >  b(std::max(1, n - 1),0);

			  for(int i = n - 1; i > 0; i--)
			    b[i - 1] = a[i] + (i < n - 1 ? b[i] * x0 : 0);
			  return b;
			}

			std::complex<T>  polinom_zero (T x0, T x1, T x2,std::vector< std::complex<T> >& coefs,unsigned int max, double pTolerance){
				std::complex<T> xi0, xi1, xi2,xi3, h0,h1,q0,q1,A, B, C,result, tmp1, tmp2;
				long double err;
				xi0=x0;
				xi1=x1;
				xi2=x2;
				err=2;
				tmp1=0;
				tmp2=0;
				int ite =0;
				while(ite<max && (long double)err >= pTolerance){

					h0=xi1-xi0;
					h1=xi2-xi1;
					q0=(evaluatePol(coefs, &xi1)-evaluatePol(coefs, &xi0))/h0;
					q1=(evaluatePol(coefs, &xi2)-evaluatePol(coefs, &xi1))/h1;

					A=(q1-q0)/(h1+h0);
					B=A*h1+q1;
					C=(evaluatePol(coefs, &xi2));

					tmp1=xi2+((std::complex<T> )-2*C)/(B+std::sqrt(B*B-(std::complex<T> )4*A*C));
					tmp2=xi2+((std::complex<T> )-2*C)/(B-std::sqrt(B*B-(std::complex<T> )4*A*C));

					if(B.real()>=0)
						xi3=tmp1;
					else
						xi3=tmp2;
					err=std::abs((((xi3-xi2)/xi3)).real())*100;

					xi0=xi1;
					xi1=xi2;
					xi2=xi3;
					ite++;

				}
				return xi3;
			};

	};
	typedef std::complex<double> cdouble;
	typedef std::vector<cdouble> polyn;
	/**
	 * Laguerre class
	 */
	std::pair<polyn, cdouble> horner(const polyn &a, cdouble x0) {
	  int n = a.size();
	  polyn b = polyn(std::max(1, n - 1));

	  for(int i = n - 1; i > 0; i--)
	    b[i - 1] = a[i] + (i < n - 1 ? b[i] * x0 : 0);
	  return make_pair(b, a[0] + b[0] * x0);
	}

	cdouble eval(const polyn &p, cdouble x) {
	  return horner(p, x).second;
	}

	polyn derivative(const polyn &p) {
	  int n = p.size();
	  polyn r = polyn(std::max(1, n - 1));
	  for(int i = 1; i < n; i++)
	    r[i - 1] = p[i] * cdouble(i);
	  return r;
	}

	const double EPS = 1e-9;

	int cmp(cdouble x, cdouble y) {
	  double diff = abs(x) - abs(y);
	  return diff < -EPS ? -1 : (diff > EPS ? 1 : 0);
	}

	cdouble find_one_root(const polyn &p0, cdouble x) {
	  int n = p0.size() - 1;
	  polyn p1 = derivative(p0);
	  polyn p2 = derivative(p1);
	  for (int step = 0; step < 10000; step++) {
	    cdouble y0 = eval(p0, x);
	    if (cmp(y0, 0) == 0) break;
	    cdouble G = eval(p1, x) / y0;
	    cdouble H = G * G - eval(p2, x) - y0;
	    cdouble R = sqrt(cdouble(n - 1) * (H * cdouble(n) - G * G));
	    cdouble D1 = G + R;
	    cdouble D2 = G - R;
	    cdouble a = cdouble(n) / (cmp(D1, D2) > 0 ? D1 : D2);
	    x -= a;
	    if (cmp(a, 0) == 0) break;
	  }
	  return x;
	}

	std::vector<cdouble> find_all_roots(const polyn &p) {
	  std::vector<cdouble> res;
	  polyn q = p;
	  while (q.size() > 2) {
	    cdouble z(rand() / double(RAND_MAX), rand() / double(RAND_MAX));
	    z = find_one_root(q, z);
	    z = find_one_root(p, z);
	    q = horner(q, z).first;
	    res.push_back(z);
	  }
	  res.push_back(-q[0] / q[1]);
	  return res;
	}
}

/**
 * Auxiliar method for the main
 * @param method          Which method will be run (number from 1 to 4 or none for default).
 * @param poly            The polynomial that the method will be run in.
 * @param order           The orther of the polynomial.
 * @param root            The root of the polinomial (use for deflation).
 * @param X               The positions where the method will be started.
 * @param maxIt           The max times the method will be run(muller).
 * @param tolerance       The tolerance the method will have(muller).
 * @param singlePrecision Single precision flag.
 * @param polish          Polish flag.
 */
void mainAux(unsigned short method, std::vector<double>& poly, short order, double root, std::vector<double> X, int maxIt=50, double tolerance=0.00001, bool singlePrecision=0,bool polish=0){
	if(poly.size()==0 and method!=0){
		std::cout << "To use a method you need to enter a polynomial with the option --polynomial and the coeficients."<< std::endl;
		return;
	}
	switch(method){
		case deflationMethod: //Deflation
		{
			if(singlePrecision){
				std::complex<float> res =*(new std::complex<float>(0,0));
				std::complex<float> roots=*(new std::complex<float>(root,0));
				std::vector<std::complex<float> > polys;
				anpi::deflation<std::complex<float> > defla;
				for(int i=0; i<poly.size(); i++){
					std::complex<float> tmp=*(new std::complex<float>(poly[i],0));
					polys.push_back(tmp);
				}
				defla(polys,roots, res, (int)order);
			}else{
				std::complex<double> res =*(new std::complex<double>(0,0));
				std::complex<double> rootd=*(new std::complex<double>(root,0));
				std::vector<std::complex<double> > polyd;
				anpi::deflation<std::complex<double> > defla;
				for(int i=0; i<poly.size(); i++){
					std::complex<float> tmp=*(new std::complex<float>(poly[i],0));
					polyd.push_back(tmp);
				}
				defla(polyd,rootd, res, (int)order);
			}
			break;
		}
		case deflationMethod2: //Deflation2
		{
			if(singlePrecision){
				std::complex<float> res =*(new std::complex<float>(0,0));
				std::complex<float> roots=*(new std::complex<float>(root,0));
				std::vector<std::complex<float> > polys;
				anpi::deflation2<std::complex<float> > defla;
				for(int i=0; i<poly.size(); i++){
					std::complex<float> tmp=*(new std::complex<float>(poly[i],0));
					polys.push_back(tmp);
				}
				defla(polys,roots, res, (int)order);
			}else{
				std::complex<double> res =*(new std::complex<double>(0,0));
				std::complex<double> rootd=*(new std::complex<double>(root,0));
				std::vector<std::complex<double> > polyd;
				anpi::deflation2<std::complex<double> > defla;
				for(int i=0; i<poly.size(); i++){
					std::complex<float> tmp=*(new std::complex<float>(poly[i],0));
					polyd.push_back(tmp);
				}
				defla(polyd,rootd, res, (int)order);
			}
			break;
		}
		case mullerMethod: //Muller
		{

			if(X.size()!=3){
				std::cout << "Three positions has to be set to use this method. Use the option: --position and the positions space separated." << std::endl;
				return ;
			}
			if(singlePrecision){
				anpi::muller<float> mull;
				std::vector<std::complex<float> > polyd;
				for(int i=0; i<poly.size(); i++){
					std::complex<float> tmp=*(new std::complex<float>(poly[i],0));
					polyd.push_back(tmp);
				}
				std::vector<std::complex<float> > res;
				res = mull((float)X[0],
								(float)X[1],
								(float)X[2],
								polyd,
								(unsigned int)maxIt,
								(float)tolerance);
				for(int i = 0; i < res.size(); i++){
					std::cout<< "x" << i <<": " << res[i] << std::endl;
				}
				break;
			}else{
				anpi::muller<double> mull;
				std::vector<std::complex<double> > polyd;
				for(int i=0; i<poly.size(); i++){
					std::complex<double> tmp=*(new std::complex<double>(poly[i],0));
					polyd.push_back(tmp);
				}
				std::vector<std::complex<double> > res;
				res = mull((double)X[0],
								(double)X[1],
								(double)X[2],
								polyd,
								(unsigned int)maxIt,
								(double)tolerance);
				for(int i = 0; i < res.size(); i++){
					std::cout<< "x" << i <<": " << res[i] << std::endl;
				}
				break;
			}
		}
		case laguerreMethod:
		{
			anpi::polyn p;
			// x^3 - 8x^2 - 13x + 140 = (x+4)(x-5)(x-7)
			//X**3 + 2x**2 -109x -110 = (x+11)(x+8)(x-4)
			for(int i=0; i<poly.size(); i++){
				p.push_back(poly[i]);
			}
			std::vector<anpi::cdouble> roots = anpi::find_all_roots(p);

			for(size_t i = 0; i < roots.size(); i++) {
				if (abs(roots[i].real()) < anpi::EPS) roots[i] -= anpi::cdouble(roots[i].real(), 0);
					if (abs(roots[i].imag()) < anpi::EPS) roots[i] -= anpi::cdouble(0, roots[i].imag());
						std::cout << std::setprecision(3) << roots[i] << std::endl;
			}
			break;
		}
		default:
		{
			//Deflation method
			std::cout << "Deflation of 1x^3+3x^2+4x+12, with 2i as root." <<std::endl;
			std::vector<std::complex<double> >  polyd;
			std::complex<double> tmp= *(new std::complex<double>(12,0));
			polyd.push_back(tmp);
			tmp= *(new std::complex<double>(4,0));
			polyd.push_back(tmp);
			tmp= *(new std::complex<double>(3,0));
			polyd.push_back(tmp);
			tmp= *(new std::complex<double>(1,0));
			polyd.push_back(tmp);
			std::complex<double> rootd =*(new std::complex<double>(0,2));
			std::complex<double> res =*(new std::complex<double>(0,0));
			anpi::deflation<std::complex<double> > defla;
			defla(polyd,rootd, res, (int)order);
			//Deflation method 2
			std::cout << "Deflation method 2 of 1x^3+3x^2+4x+12, with 2i as root." <<std::endl;
			rootd =*(new std::complex<double>(0,2));
			res =*(new std::complex<double>(0,0));
			anpi::deflation2<std::complex<double> > defla2;
			defla2(polyd,rootd, res, (int)order);
			//Muller method
			std::cout << "Muller method roots of 1x^3+2x^2+-5x+-6." <<std::endl;
			anpi::muller<double> mull;
			std::vector<std::complex<double> >coefs;
			tmp=-6;
			coefs.push_back(tmp);
			tmp=-5;
			coefs.push_back(tmp);
			tmp=2;
			coefs.push_back(tmp);
			tmp=1;
			coefs.push_back(tmp);
			std::vector<std::complex<double> > result;
			result = mull((double)-4,
							(double)-2,
							(double)0,
							coefs,
							(unsigned int)50,
							(long double)0.000001);
			for(int i = 0; i < result.size(); i++){
				std::cout<< "x" << i <<": " << result[i] << std::endl;
			}

			//Laguerre method
			std::cout << "Laguerre roots of 1x^3-2x^2+-109x+-110." <<std::endl;
			anpi::polyn p;
			// x^3 - 8x^2 - 13x + 140 = (x+4)(x-5)(x-7)
			//X**3 + 2x**2 -109x -110 = (x+11)(x+8)(x-4)
			p.push_back(-110);
			p.push_back(-109);
			p.push_back(-2);
			p.push_back(1);

			std::vector<anpi::cdouble> roots = anpi::find_all_roots(p);

			for(size_t i = 0; i < roots.size(); i++) {
				if (abs(roots[i].real()) < anpi::EPS) roots[i] -= anpi::cdouble(roots[i].real(), 0);
					if (abs(roots[i].imag()) < anpi::EPS) roots[i] -= anpi::cdouble(0, roots[i].imag());
						std::cout << std::setprecision(3) << roots[i] << std::endl;
			}
		}

	}
}

/**
 * Main method, runs all the methods with a default or you can enter options
 * @param Uses program options, call the program plus --help to see detailed description
 * @return the return is depending which method is called
 */
int main(int argc, char** argv){
	using boost::lexical_cast;
	using boost::bad_lexical_cast;
	//Default values for the methods:
	short order=3;
	unsigned short method=0;
	double tolerance=0.000001;
	int maxIt=50;
	std::vector<double> poly;
	double root;
	std::vector<double> position;
	bool polish=0;
	bool singlePrecision=0;
	//Program options definitions
	po::options_description description("Project1 Usage");
	description.add_options()
	("help", "Display this help message")
	("polynomial", po::value<std::vector<double> >(&poly)->multitoken(), "The polynomial that the program will use, just add the coeficients space separated.")
	("position", po::value<std::vector<double> >(&position)->multitoken(), "The positions to start the method space separated.")
	("root", po::value<double>(&root), "The root to find.")
	("order", po::value<short>(&order), "The order of the polynomial.")
	("muller", "Use this option to use the muller method.")
	("deflation", "Use this option to use deflation method.")
	("deflation2", "Use this option to use deflation2 method.")
	("laguerre", "Use this option to use laguerre method.")
	("polish", "Use this option to polish the roots.")
	("tolerance", po::value<double>(&tolerance), "The order of the polynomial.")
	("maxIterations", po::value<int>(&maxIt), "The order of the polynomial.")
	("singlePrecision", "Use this to set the methods to double singlePrecision, otherwise double precision will be used.");
	po::positional_options_description p;
	p.add("file", -1);
	po::variables_map vm;
	//Try to use the entered options
	try {
		po::store(po::parse_command_line(argc, argv, description, po::command_line_style::unix_style ^ po::command_line_style::allow_short), vm);		po::notify(vm);
		if(vm.count("help")){
			std::cout << description;
			return 0;
		}
		if(vm.count("order")){
			#ifdef debug
				std::cout << "Polynomial order changed to: " << order << std::endl;
			#endif
			if(order<0){
				std::cout << "The order has to be a positive number." << std::endl;
				return -1;
			}
		}
		if(vm.count("tolerance")){
			#ifdef debug
				std::cout << "The tolerance for the method will be: " << tolerance << std::endl;
			#endif
			if(tolerance<0){
				std::cout << "The tolerance has to be a positive number." << std::endl;
				return -1;
			}

		}
		if(vm.count("maxIterations")){
			#ifdef debug
				std::cout << "Max iterations will be changed to: " << maxIt << std::endl;
			#endif
			if(maxIt<0){
				std::cout << "The iterations has to be a positive number." << std::endl;
				return -1;
			}

		}
		if(vm.count("singlePrecision")){
			#ifdef debug
				std::cout << "Polynomial order changed to: " << singlePrecision << std::endl;
			#endif
			singlePrecision = 1;
		}
		if(vm.count("deflation")){
			#ifdef debug
				std::cout << "Deflation method will be used." << std::endl;
			#endif
			method = deflationMethod;
		}
		if(vm.count("deflation2")){
			#ifdef debug
				std::cout << "Deflation method 2 will be used." << std::endl;
			#endif
			method = deflationMethod2;
		}
		if(vm.count("muller")){
			#ifdef debug
				std::cout << "Muller method will be used." << std::endl;
			#endif
			method = mullerMethod;
		}
		if(vm.count("laguerre")){
			#ifdef debug
				std::cout << "Laguerre method will be used." << std::endl;
			#endif
			method = laguerreMethod;
		}
		if(vm.count("poly")){
			#ifdef debug
				std::cout << "A custom polynomial will be used. With size: " <<poly.size() << std::endl;
			#endif
			if(poly.size()<=0){
				std::cout << "The polynomial coeficients must be more than 0." << std::endl;
				return -1;
			}
		}
		if(vm.count("root")){
			#ifdef debug
				std::cout << "A custom root will be used: " << root << std::endl;
			#endif
		}
		if(vm.count("position")){
			#ifdef debug
				std::cout << "Custom positions will be used." << std::endl;
			#endif
			if(poly.size()<=0){
				std::cout << "The positions must be more than 0." << std::endl;
				return -1;
			}
		}
		mainAux(method, poly, order, root, position,maxIt, tolerance, singlePrecision,polish);
		return 0;
	}
	catch(const po::error &ex){
		std::cerr << ex.what() << '\n';
	}

}


text:
/**************************************************************************\
|
|    Copyright (C) 2009 Marc Stevens
|
|    This program is free software: you can redistribute it and/or modify
|    it under the terms of the GNU General Public License as published by
|    the Free Software Foundation, either version 3 of the License, or
|    (at your option) any later version.
|
|    This program is distributed in the hope that it will be useful,
|    but WITHOUT ANY WARRANTY; without even the implied warranty of
|    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
|    GNU General Public License for more details.
|
|    You should have received a copy of the GNU General Public License
|    along with this program.  If not, see <http://www.gnu.org/licenses/>.
|
\**************************************************************************/

#include <stdexcept>
#include <iostream>
#include <fstream>
#include <vector>
#include <string>

#include <boost/program_options.hpp>

#include "main.hpp"

#include <hashclash/md5detail.hpp>
#include <hashclash/timer.hpp>
#include <hashclash/rng.hpp>

using namespace hashclash;
using namespace std;
namespace po = boost::program_options;
namespace fs = boost::filesystem;

unsigned load_block_fillrandom(istream& i, uint32 block[])
{
	unsigned bytes = load_block(i, block);
	unsigned len = 0;
	for (unsigned k = 0; k < 16; ++k)
		for (unsigned c = 0; c < 4; ++c)
		{
			if (len >= bytes)
				block[k] ^= (xrng128()&0xFF) << (c*8);
			++len;
		}
	return bytes;
}

unsigned load_block(istream& i, uint32 block[])
{
	for (unsigned k = 0; k < 16; ++k)
		block[k] = 0;

	unsigned len = 0;
	char uc;
	for (unsigned k = 0; k < 16; ++k)
		for (unsigned c = 0; c < 4; ++c)
		{
			i.get(uc);
			if (i) {
				++len;
				block[k] += uint32((unsigned char)(uc)) << (c*8);
			} else {
				i.putback(uc);
				i.setstate(ios::failbit);
				return len;
			}
		}
	return len;
}

void save_block(ostream& o, uint32 block[])
{
	for (unsigned k = 0; k < 16; ++k)
		for (unsigned c = 0; c < 4; ++c)
			o << (unsigned char)((block[k]>>(c*8))&0xFF);
}

int dostep(birthday_parameters& parameters)
{
	uint32 ihv1[4];
	uint32 ihv2[4];
	uint32 msg1[16];
	uint32 msg2[16];
	{
		ifstream if1(parameters.inputfile1.c_str(), ios::binary);
		if (!if1) {
			cerr << "Error: cannot open inputfile 1 '" << parameters.inputfile1 << "'!" << endl;
			return 1;
		}
		ifstream if2(parameters.inputfile2.c_str(), ios::binary);
		if (!if2) {
			cerr << "Error: cannot open inputfile 2 '" << parameters.inputfile2 << "'!" << endl;
			return 1;
		}
		ofstream* of1 = 0;
		ofstream* of2 = 0;
		if (parameters.modi == 0) {
			of1 = new ofstream(parameters.outputfile1.c_str(), ios::binary);
			if (!(*of1)) {
				cerr << "Error: cannot open outputfile 1 '" << parameters.outputfile1 << "'!" << endl;
				return 1;
			}
			of2 = new ofstream(parameters.outputfile2.c_str(), ios::binary);
			if (!(*of2)) {
				delete of1;
				cerr << "Error: cannot open outputfile 2 '" << parameters.outputfile1 << "'!" << endl;
				return 1;
			}
		}

		for (unsigned k = 0; k < 4; ++k)
			ihv1[k] = ihv2[k] = md5_iv[k];

		seed(0x12345678);
		// load, md5 and save inputfile1
		unsigned file1blocks = 0;
		while (load_block_fillrandom(if1, msg1) > 64-12) { // stop when at least 12 bytes are not used
			md5compress(ihv1, msg1);
			if (of1) save_block(*of1, msg1);
			++file1blocks;
			addseed(ihv1[0]);
		}
		// load, md5 and save inputfile2
		unsigned file2blocks = 0;
		while (load_block_fillrandom(if2, msg2) > 64-12) { // stop when at least 12 bytes are not used
			md5compress(ihv2, msg2);
			if (of2) save_block(*of2, msg2);
			++file2blocks;
			addseed(ihv2[0]);
		}
		// flush last partial message block when padding blocks are required
		if (file1blocks < file2blocks) {
			md5compress(ihv1, msg1);
			if (of1) save_block(*of1, msg1);
			++file1blocks;
		}
		if (file2blocks < file1blocks) {
			md5compress(ihv2, msg2);
			if (of2) save_block(*of2, msg2);
			++file2blocks;
		}
		// expand length outputfile1 to length outputfile2, by saving and md5'ing
		while (file1blocks < file2blocks) {
			for (unsigned k = 0; k < 16; ++k)
				msg1[k] = xrng128()+xrng128();
			md5compress(ihv1, msg1);
			if (of1) save_block(*of1, msg1);
			++file1blocks;
		}
		// expand length outputfile2 to length outputfile1, by saving and md5'ing
		while (file2blocks < file1blocks) {
			for (unsigned k = 0; k < 16; ++k)
				msg2[k] = xrng128()+xrng128();
			md5compress(ihv2, msg2);
			if (of2) save_block(*of2, msg2);
			++file2blocks;
		}

		if (of1) delete of1;
		if (of2) delete of2;
	} // close if1, if2, of1, of2
	addseed(uint32(time(0)));
	addseed(parameters.modi);

	for (unsigned k = 0; k < 4; ++k)
	{
		parameters.ihv1[k] = ihv1[k];
		parameters.ihv2[k] = ihv2[k];
	}
	for (unsigned k = 0; k < 16; ++k)
	{
		parameters.msg1[k] = msg1[k];
		parameters.msg2[k] = msg2[k];
	}

	cout << "IHV1 = {" << ihv1[0] << "," << ihv1[1] << "," << ihv1[2] << "," << ihv1[3] << "}" << endl;
	cout << "IHV1 = " << hex;
	for (unsigned k = 0; k < 4; ++k)
		for (unsigned c = 0; c < 4; ++c)
		{
			cout.width(2); cout.fill('0');
			cout << ((ihv1[k]>>(c*8))&0xFF);
		}
	cout << dec << endl << endl;

	cout << "IHV2 = {" << ihv2[0] << "," << ihv2[1] << "," << ihv2[2] << "," << ihv2[3] << "}" << endl;
	cout << "IHV2 = " << hex;
	for (unsigned k = 0; k < 4; ++k)
		for (unsigned c = 0; c < 4; ++c)
		{
			cout.width(2); cout.fill('0');
			cout << ((ihv2[k]>>(c*8))&0xFF);
		}
	cout << dec << endl << endl;

	// Start job with given parameters
	birthday(parameters);

	ofstream* of1 = 0;
	ofstream* of2 = 0;
	if (parameters.modi == 0) {
		of1 = new ofstream(parameters.outputfile1.c_str(), ios::binary | ios::app);
		if (!(*of1)) {
			cerr << "Error: cannot open outputfile 1 '" << parameters.outputfile1 << "'!" << endl;
			return 1;
		}
		of2 = new ofstream(parameters.outputfile2.c_str(), ios::binary | ios::app);
		if (!(*of2)) {
			delete of1;
			cerr << "Error: cannot open outputfile 2 '" << parameters.outputfile1 << "'!" << endl;
			return 1;
		}
		msg1[13] = collc1;
		msg1[14] = colla1;
		msg1[15] = collb1;
		msg2[13] = collc2;
		msg2[14] = colla2;
		msg2[15] = collb2;
		if (of1) save_block(*of1, msg1);
		if (of2) save_block(*of2, msg2);
	}
	if (of1) delete of1;
	if (of2) delete of2;

	return 0;
}


text:
#include "CommandBuffer.hpp"
#include "Command.hpp"
#include "DrawCommand.hpp"
#include "ClearDrawingCommand.hpp"
#include "RemoveDrawingCommand.hpp"
#include "ClearPlotCommand.hpp"
#include "vizkit3d_debug_drawings/DrawingManager.hpp"

//need to be included. otherwise BOOST_CLASS_EXPORT wont work (http://www.boost.org/doc/libs/1_46_1/libs/serialization/doc/special.html)
#include <boost/archive/binary_oarchive.hpp>
#include <boost/archive/binary_iarchive.hpp>

#include "BoostSerializationExports.hpp"
#include <iostream>

using namespace vizkit3dDebugDrawings;

CommandBuffer::CommandBuffer(size_t maxBufferSize) : maxBufferSize(maxBufferSize)
{
}

CommandBuffer::CommandBuffer() : maxBufferSize(50000)
{
}

void CommandBuffer::addCommand(boost::shared_ptr<Command> cmd)
{
    if(dynamic_cast<ClearDrawingCommand*>(cmd.get()) ||
       dynamic_cast<RemoveDrawingCommand*>(cmd.get()) ||
       dynamic_cast<ClearPlotCommand*>(cmd.get()))
    {
        commands.clear();
    }

    commands.push_back(cmd);
    if(commands.size() > maxBufferSize)
    {
        commands.pop_front();
        std::cerr << "Commandbuffer is full. Dropping old commands." << std::endl;
    }
}

void CommandBuffer::executeAll(vizkit3dDebugDrawings::DrawingManager* manager) const
{

    for(const auto& cmd : commands)
    {
        cmd->execute(manager);
    }
}

const std::deque<boost::shared_ptr<Command>>& CommandBuffer::getCommands() const
{
    return commands;
}

BOOST_CLASS_EXPORT(vizkit3dDebugDrawings::CommandBuffer);


text:
#include "common.hpp"

#include <MetaObject/runtime_reflection/StructTraits.hpp>
#include <MetaObject/runtime_reflection/visitor_traits/filesystem.hpp>
#include <MetaObject/runtime_reflection/visitor_traits/string.hpp>
#include <MetaObject/runtime_reflection/visitor_traits/vector.hpp>

#include "MetaObject/core/metaobject_config.hpp"

#include "MetaObject/types/file_types.hpp"
#include <MetaObject/params/AccessToken.hpp>

#include <boost/lexical_cast.hpp>

#include <cereal/types/string.hpp>
#include <cereal/types/vector.hpp>

#ifdef MO_EXPORTS
#undef MO_EXPORTS
#endif
#if (defined WIN32 || defined _WIN32 || defined WINCE || defined __CYGWIN__)
#define MO_EXPORTS __declspec(dllexport)
#elif defined __GNUC__ && __GNUC__ >= 4
#define MO_EXPORTS __attribute__((visibility("default")))
#else
#define MO_EXPORTS
#endif

using namespace mo;

#ifdef MO_HAVE_PYTHON
#include <boost/python.hpp>
namespace ct
{

    template <>
    inline bool convertFromPython(const boost::python::object& obj, EnumParam& param)
    {
        boost::python::extract<std::string> str_ext(obj);
        if (str_ext.check())
        {
            auto string = str_ext();
            auto itr = std::find(param.enumerations.begin(), param.enumerations.end(), string);
            if (itr != param.enumerations.end())
            {
                param.current_selection = itr - param.enumerations.begin();
                return true;
            }
            return false;
        }
        boost::python::extract<int> int_ext(obj);

        if (int_ext.check())
        {
            int val = int_ext();
            param.current_selection = static_cast<size_t>(val);
            return true;
        }

        return false;
    }

    template <>
    inline bool convertFromPython(const boost::python::object& obj, ReadFile& result)
    {
        boost::python::extract<std::string> extractor(obj);
        result = extractor();
        return true;
    }

    template <>
    inline bool convertFromPython(const boost::python::object& obj, WriteFile& result)
    {
        boost::python::extract<std::string> extractor(obj);
        result = extractor();
        return true;
    }

    template <>
    inline bool convertFromPython(const boost::python::object& obj, WriteDirectory& result)
    {
        boost::python::extract<std::string> extractor(obj);
        result = extractor();
        return true;
    }

    template <>
    inline bool convertFromPython(const boost::python::object& obj, ReadDirectory& result)
    {
        boost::python::extract<std::string> extractor(obj);
        result = extractor();
        return true;
    }

    template <>
    inline boost::python::object convertToPython(const ReadFile& file)
    {
        return boost::python::object(file.string());
    }

    template <>
    inline boost::python::object convertToPython(const WriteFile& file)
    {
        return boost::python::object(file.string());
    }

    template <>
    inline boost::python::object convertToPython(const ReadDirectory& file)
    {
        return boost::python::object(file.string());
    }

    template <>
    inline boost::python::object convertToPython(const WriteDirectory& file)
    {
        return boost::python::object(file.string());
    }
} // namespace ct

#endif

namespace mo
{
    namespace metaparams
    {
        void instMOTypes(SystemTable* table)
        {
            registerTrait<ReadFile>();
            registerTrait<WriteFile>();
            registerTrait<ReadDirectory>();
            registerTrait<WriteDirectory>();
            registerTrait<EnumParam>();
        }
    } // namespace metaparams
} // namespace mo

max_stars_repo_name: dtmoodie/MetaObject

max_issues_repo_name: dtmoodie/MetaObject

max_forks_repo_name: dtmoodie/MetaObject


text:
/**
 * Copyright (c) 2016, Adrien Devresse <adrien.devresse@epfl.ch>
 *
 * Boost Software License - Version 1.0
 *
 * Permission is hereby granted, free of charge, to any person or organization
 * obtaining a copy of the software and accompanying documentation covered by
 * this license (the "Software") to use, reproduce, display, distribute,
 * execute, and transmit the Software, and to prepare derivative works of the
 * Software, and to permit third-parties to whom the Software is furnished to
 * do so, all subject to the following:
 *
 * The copyright notices in the Software and this entire statement, including
 * the above license grant, this restriction and the following disclaimer,
 * must be included in all copies of the Software, in whole or in part, and
 * all derivative works of the Software, unless such copies or derivative
 * works are solely in the form of machine-executable object code generated by
 * a source language processor.
 *
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 * FITNESS FOR A PARTICULAR PURPOSE, TITLE AND NON-INFRINGEMENT. IN NO EVENT
 * SHALL THE COPYRIGHT HOLDERS OR ANYONE DISTRIBUTING THE SOFTWARE BE LIABLE
 * FOR ANY DAMAGES OR OTHER LIABILITY, WHETHER IN CONTRACT, TORT OR OTHERWISE,
 * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
 * DEALINGS IN THE SOFTWARE.
 *
 */

#define BOOST_TEST_MODULE parallelTests
#define BOOST_TEST_MAIN

#include <algorithm>
#include <future>
#include <iostream>
#include <numeric>
#include <random>
#include <sstream>
#include <stdexcept>

#include <chrono>

#include <boost/test/unit_test.hpp>

#include <hadoken/parallel/algorithm.hpp>

//#include <parallel/algorithm>

using namespace hadoken;
using cl = std::chrono::system_clock;

BOOST_AUTO_TEST_CASE(parallel_for_each_simple_test) {

    std::vector<double> values;
    parallel::for_each(parallel::seq, values.begin(), values.end(), [](double& v) { v += 42; });

    parallel::for_each(parallel::parallel_execution_policy(), values.begin(), values.end(), [](double& v) { v += 422; });

    parallel::for_each(parallel::parallel_vector_execution_policy(), values.begin(), values.end(),
                       [](double& v) { v += 4222; });

    values = std::vector<double>(12800, 128);

    auto t1 = cl::now();

    auto fops = [](double& v) {
        double res = std::pow(v, 4);
        res = std::sqrt(res);
        v = res;
    };

    parallel::for_each(parallel::seq, values.begin(), values.end(), fops);

    auto t2 = cl::now();

    for (auto i : values) {
        BOOST_CHECK_CLOSE(i, 128 * 128, 0.01);
    }

    std::cout << "sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;

    values = std::vector<double>(12800, 128);

    t1 = cl::now();

    parallel::for_each(parallel::parallel_execution_policy(), values.begin(), values.end(), fops);
    //__gnu_parallel::for_each(values.begin(), values.end(), [](int & v){ v += 100; });

    t2 = cl::now();

    for (auto i : values) {
        BOOST_CHECK_CLOSE(i, 128 * 128, 0.01);
    }

    std::cout << "parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;

    values = std::vector<double>(12800, 128);

    parallel::for_each(parallel::parallel_vector_execution_policy(), values.begin(), values.end(), fops);

    for (auto i : values) {
        BOOST_CHECK_CLOSE(i, 128 * 128, 0.01);
    }
}

BOOST_AUTO_TEST_CASE(parallel_for_each_nested) {

    std::vector<std::vector<double>> vec_vec_values(1024);

    std::size_t index = 0;
    for (auto& vec_values : vec_vec_values) {
        vec_values.resize(2048);
        std::fill(vec_values.begin(), vec_values.end(), index);
        index++;
    }

    parallel::for_each(parallel::parallel_execution_policy(), vec_vec_values.begin(), vec_vec_values.end(),
                       [](std::vector<double>& vec_values) {
                           parallel::for_each(parallel::parallel_execution_policy(), vec_values.begin(), vec_values.end(),
                                              [](double& v) { v += 42; });
                       });

    index = 0;
    for (auto& vec_values : vec_vec_values) {
        for (auto i : vec_values) {
            BOOST_CHECK_CLOSE(index + 42, i, 0.001);
        }
        index++;
    }
}

BOOST_AUTO_TEST_CASE(parallel_fill_test) {

    using namespace hadoken;

    std::size_t n = 5000;

    std::size_t val = 42;

    std::vector<double> values(n), values_seq(n), values_par(n), values_par_n(n);

    std::fill(values.begin(), values.end(), val);

    parallel::fill(parallel::seq, values_seq.begin(), values_seq.end(), val);

    parallel::fill(parallel::parallel_execution_policy(), values_par.begin(), values_par.end(), val);

    parallel::fill_n(parallel::parallel_vector_execution_policy(), values_par_n.begin(), n, val);

    for (auto& v : values) {
        BOOST_CHECK_EQUAL(v, val);
    }

    BOOST_CHECK_EQUAL_COLLECTIONS(values.begin(), values.end(), values_seq.begin(), values_seq.end());
    BOOST_CHECK_EQUAL_COLLECTIONS(values.begin(), values.end(), values_par.begin(), values_par.end());
    BOOST_CHECK_EQUAL_COLLECTIONS(values.begin(), values.end(), values_par_n.begin(), values_par_n.end());
}

BOOST_AUTO_TEST_CASE(parallel_count_test) {

    using namespace hadoken;

    std::size_t n = 50000;

    std::size_t num_4, num_2, p_num_4, p_num_2;

    std::vector<std::size_t> values(n, 0);

    std::size_t counter = 0;
    for (auto& v : values) {
        v = (counter++) % 5;
    }

    {
        auto t1 = cl::now();

        num_4 = std::count(values.begin(), values.end(), 4);
        num_2 = std::count_if(values.begin(), values.end(), [&](const std::size_t& v) { return v == 1; });

        auto t2 = cl::now();
        std::cout << " std::count " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;
    }

    {
        auto t1 = cl::now();

        p_num_4 = parallel::count(parallel::par, values.begin(), values.end(), 4);
        p_num_2 = parallel::count_if(parallel::par, values.begin(), values.end(), [&](const std::size_t& v) { return v == 1; });

        auto t2 = cl::now();
        std::cout << " parallel::count " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;
    }

    std::size_t s_num_4 = parallel::count(parallel::seq, values.begin(), values.end(), 4);
    std::size_t s_num_2 =
        parallel::count_if(parallel::seq, values.begin(), values.end(), [&](const std::size_t& v) { return v == 1; });

    BOOST_CHECK_EQUAL(num_4, p_num_4);
    BOOST_CHECK_EQUAL(num_4, s_num_4);

    BOOST_CHECK_EQUAL(num_2, p_num_2);
    BOOST_CHECK_EQUAL(num_2, s_num_2);

    std::cout << " n " << s_num_4 << std::endl;
}

BOOST_AUTO_TEST_CASE(parallel_transform_test) {

    using namespace hadoken;

    std::size_t n = 80000;

    auto dummy_ops = [](std::size_t i, std::size_t j) { return (i * 100 + j) / 7; };

    std::vector<std::size_t> v1(n), v2(n), res(n), res2(n), res3(n), res4(n);

    std::size_t counter = 0;
    for (std::size_t i = 0; i < n; ++i) {
        v1[i] = counter++;
        v2[i] = counter++;
    }

    {
        auto t1 = cl::now();

        hadoken::parallel::transform(hadoken::parallel::par, v1.begin(), v1.end(), v2.begin(), res.begin(), dummy_ops);

        auto t2 = cl::now();

        std::cout << " transform parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;

        hadoken::parallel::transform(hadoken::parallel::seq, v1.begin(), v1.end(), res2.begin(),
                                     [](const std::size_t& v) { return v + 100; });
    }

    {
        auto t1 = cl::now();

        std::transform(v1.begin(), v1.end(), v2.begin(), res3.begin(), dummy_ops);

        auto t2 = cl::now();

        std::cout << " transform sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;

        std::transform(v1.begin(), v1.end(), res4.begin(), [](const std::size_t& v) { return v + 100; });
    }

    BOOST_CHECK_EQUAL_COLLECTIONS(res.begin(), res.end(), res3.begin(), res3.end());

    BOOST_CHECK_EQUAL_COLLECTIONS(res2.begin(), res2.end(), res4.begin(), res4.end());
}

BOOST_AUTO_TEST_CASE(parallel_all_of_test) {

    using namespace hadoken;

    std::size_t n = 50000;

    auto is_positive = [](int v) { return (v >= 0); };

    std::vector<int> values;

    bool res_seq = parallel::all_of(parallel::seq, values.begin(), values.end(), is_positive);

    bool res_par = parallel::all_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_positive);

    bool res_std = std::all_of(values.begin(), values.end(), is_positive);

    BOOST_CHECK_EQUAL(res_seq, true);
    BOOST_CHECK_EQUAL(res_par, true);
    BOOST_CHECK_EQUAL(res_std, true);

    // fullfil vector with positive values excepted one
    for (std::size_t i = 0; i < n; ++i) {
        values.push_back(i);
    }

    values.push_back(-1);

    res_seq = parallel::all_of(parallel::seq, values.begin(), values.end(), is_positive);

    {
        auto t1 = cl::now();
        res_par = parallel::all_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_positive);
        auto t2 = cl::now();

        std::cout << " all_of parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;
    }

    {
        auto t1 = cl::now();
        res_std = std::all_of(values.begin(), values.end(), is_positive);
        auto t2 = cl::now();

        std::cout << " all_of sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    BOOST_CHECK_EQUAL(res_seq, false);
    BOOST_CHECK_EQUAL(res_par, false);
    BOOST_CHECK_EQUAL(res_std, false);

    // remove the positive values
    values.pop_back();

    res_seq = parallel::all_of(parallel::seq, values.begin(), values.end(), is_positive);

    res_par = parallel::all_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_positive);

    res_std = std::all_of(values.begin(), values.end(), is_positive);

    BOOST_CHECK_EQUAL(res_seq, true);
    BOOST_CHECK_EQUAL(res_par, true);
    BOOST_CHECK_EQUAL(res_std, true);
}

BOOST_AUTO_TEST_CASE(parallel_none_of_test) {

    using namespace hadoken;

    std::size_t n = 50000;

    auto is_negative = [](int v) { return (v < 0); };

    std::vector<int> values;

    bool res_seq = parallel::none_of(parallel::seq, values.begin(), values.end(), is_negative);

    bool res_par = parallel::none_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);

    bool res_std = std::none_of(values.begin(), values.end(), is_negative);

    BOOST_CHECK_EQUAL(res_seq, true);
    BOOST_CHECK_EQUAL(res_par, true);
    BOOST_CHECK_EQUAL(res_std, true);

    // fullfil vector with positive values excepted one
    for (std::size_t i = 0; i < n; ++i) {
        values.push_back(i);
    }

    values.push_back(-1);

    res_seq = parallel::none_of(parallel::seq, values.begin(), values.end(), is_negative);

    {
        auto t1 = cl::now();
        res_par = parallel::none_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);
        auto t2 = cl::now();

        std::cout << " none_of parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    {
        auto t1 = cl::now();
        res_std = std::none_of(values.begin(), values.end(), is_negative);
        auto t2 = cl::now();

        std::cout << " none_of sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    BOOST_CHECK_EQUAL(res_seq, false);
    BOOST_CHECK_EQUAL(res_par, false);
    BOOST_CHECK_EQUAL(res_std, false);

    // remove the positive values
    values.pop_back();

    res_seq = parallel::none_of(parallel::seq, values.begin(), values.end(), is_negative);

    res_par = parallel::none_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);

    res_std = std::none_of(values.begin(), values.end(), is_negative);

    BOOST_CHECK_EQUAL(res_seq, true);
    BOOST_CHECK_EQUAL(res_par, true);
    BOOST_CHECK_EQUAL(res_std, true);
}

BOOST_AUTO_TEST_CASE(parallel_any_of_test) {

    using namespace hadoken;

    std::size_t n = 50000;

    auto is_negative = [](int v) { return (v < 0); };

    std::vector<int> values;

    bool res_seq = parallel::any_of(parallel::seq, values.begin(), values.end(), is_negative);

    bool res_par = parallel::any_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);

    bool res_std = std::any_of(values.begin(), values.end(), is_negative);

    BOOST_CHECK_EQUAL(res_seq, false);
    BOOST_CHECK_EQUAL(res_par, false);
    BOOST_CHECK_EQUAL(res_std, false);

    // fullfil vector with positive values excepted one
    for (std::size_t i = 0; i < n; ++i) {
        values.push_back(i);
    }

    values.push_back(-1);

    res_seq = parallel::any_of(parallel::seq, values.begin(), values.end(), is_negative);

    {
        auto t1 = cl::now();
        res_par = parallel::any_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);
        auto t2 = cl::now();

        std::cout << " any_of parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;
    }

    {
        auto t1 = cl::now();
        res_std = std::any_of(values.begin(), values.end(), is_negative);
        auto t2 = cl::now();

        std::cout << " any_of sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    BOOST_CHECK_EQUAL(res_seq, true);
    BOOST_CHECK_EQUAL(res_par, true);
    BOOST_CHECK_EQUAL(res_std, true);

    // remove the positive values
    values.pop_back();

    res_seq = parallel::any_of(parallel::seq, values.begin(), values.end(), is_negative);

    res_par = parallel::any_of(parallel::parallel_execution_policy(), values.begin(), values.end(), is_negative);

    res_std = std::any_of(values.begin(), values.end(), is_negative);

    BOOST_CHECK_EQUAL(res_seq, false);
    BOOST_CHECK_EQUAL(res_par, false);
    BOOST_CHECK_EQUAL(res_std, false);
}

template <typename Container>
bool is_ordered(Container& c) {
    for (std::size_t i = 0; i < c.size() - 1; ++i) {
        if (c[i] > c[i + 1]) {
            return false;
        }
    }
    return true;
}

BOOST_AUTO_TEST_CASE(parallel_sort) {

    using namespace hadoken;

    std::size_t n = 50000;

    std::vector<int> values(n);

    std::mt19937_64 mt;
    std::uniform_int_distribution<int> dist;

    std::generate(values.begin(), values.end(), [&]() { return dist(mt); });

    BOOST_CHECK(is_ordered(values) == false);

    {
        auto v1 = values;
        parallel::sort(parallel::seq, v1.begin(), v1.end());
        BOOST_CHECK(is_ordered(v1) == true);
    }

    {
        auto v2 = values;

        auto t1 = cl::now();
        parallel::sort(parallel::parallel_execution_policy(), v2.begin(), v2.end());

        auto t2 = cl::now();

        std::cout << " sort parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;

        BOOST_CHECK(is_ordered(v2) == true);
    }

    {

        auto v3 = values;

        auto t1 = cl::now();

        std::sort(v3.begin(), v3.end());

        auto t2 = cl::now();

        std::cout << " sort sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() << std::endl;

        BOOST_CHECK(is_ordered(v3) == true);
    }
}

BOOST_AUTO_TEST_CASE(parallel_inclusive_scan) {

    using namespace hadoken;

    std::size_t n = 50000;

    std::vector<std::uint64_t> values(n);

    decltype(values) v1(n), v2(n), v3(n);

    std::iota(values.begin(), values.end(), 1);

    { parallel::inclusive_scan(parallel::seq, values.begin(), values.end(), v1.begin()); }

    {

        auto t1 = cl::now();

        parallel::inclusive_scan(parallel::par, values.begin(), values.end(), v2.begin());

        auto t2 = cl::now();

        std::cout << " inclusive scan parallel " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    {

        auto t1 = cl::now();

        std::partial_sum(values.begin(), values.end(), v3.begin());

        auto t2 = cl::now();

        std::cout << " inclusive scan sequential " << std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count()
                  << std::endl;
    }

    BOOST_CHECK_EQUAL(v1.size(), n);
    BOOST_CHECK_EQUAL(v2.size(), n);
    BOOST_CHECK_EQUAL(v3.size(), n);

    decltype(values)::value_type res = 0;
    for (std::size_t i = 0; i < n; ++i) {
        BOOST_CHECK_EQUAL(v1[i], v2[i]);
        BOOST_CHECK_EQUAL(v2[i], v3[i]);

        BOOST_CHECK_LE(res, v1[i]);
        res = v1[i];
    }
}

max_stars_repo_name: adevress/hadoken

max_issues_repo_name: mrexodia/hadoken

max_forks_repo_name: mrexodia/hadoken