#include <CudaMatrixPrivate.cuh>
#include <math.h>
#include <thrust/transform.h>
#include <thrust/functional.h>
#include <thrust/execution_policy.h>
using namespace thrust::placeholders;

struct PowOp: public thrust::unary_function<float, float>{
    float exponent;
    PowOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return powf(x, exponent);
    }
};

struct CompareGOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareGOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent<x?1.0:.0);
    }
};

struct CompareGEOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareGEOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent<=x?1.0:.0);
    }
};

struct CompareEOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareEOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent==x?1.0:.0);
    }
};

struct CompareNEOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareNEOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent!=x?1.0:.0);
    }
};

struct CompareLOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareLOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent>x?1.0:.0);
    }
};

struct CompareLEOp: public thrust::unary_function<float, float>{
    float exponent;
    CompareLEOp(float v):exponent(v) {}
    void setExponent(float v){
        exponent = v;
    }
    
    __host__ __device__
    float operator()(const float& x)  {
        return (exponent>=x?1.0:.0);
    }
};


struct CompareAGOp{
    __host__ __device__
    float operator()(const float& x,const float& y)  {
        return x>y?1:0;
    }
};

struct CompareAGEOp{
    __host__ __device__
    float operator()(const float& x,const float& y)  {
        return x>=y?1:0;
    }
};

struct CompareAEOp{
    __host__ __device__
    float operator()(const float& x,const float& y)  {
        return x==y?1:0;
    }
};

struct CompareANEOp{
    __host__ __device__
    float operator()(const float& x,const float& y)  {
        return x!=y?1:0;
    }
};


void unaryAdd(float* in1, float* in2, float* out, unsigned long  length)
{
    thrust::plus<float> op;
    thrust::transform(thrust::device,in1,in1+length,in2,out,op);
}

void unaryAdd(float* in1, const float& in2, float* out, unsigned long  length)
{
    thrust::transform(thrust::device,in1,in1+length,out,in2 + _1);
}

void unaryMul(float* in1, float* in2, float* out, unsigned long  length)
{
    thrust::multiplies<float> op;
    thrust::transform(thrust::device,in1,in1+length,in2,out,op);
}

void unaryMul(float* in1, const float& in2, float* out, unsigned long  length)
{
    thrust::transform(thrust::device,in1, in1+length, out, in2 * _1);
}

void unaryNeg(float* in1, float* out, unsigned long length){
    thrust::negate<float> op;
    thrust::transform(thrust::device,in1,in1+length,out,op);
}

void unarySub(float* in1, float* in2, float* out, unsigned long length){
    thrust::minus<float> op;
    thrust::transform(thrust::device,in1,in1+length,in2,out,op);
}

void unaryDiv(float* in1, float* in2, float* out, unsigned long length){
    thrust::divides<float> op;
    thrust::transform(thrust::device,in1,in1+length,in2,out,op);
}

void unarySub(const float& in1, float* in2, float* out, unsigned long length){
    thrust::transform(thrust::device,in2,in2+length,out,in1-_1);
}

void unaryDiv(const float& in1, float* in2, float* out, unsigned long length){
    thrust::transform(thrust::device,in2,in2+length,out,in1/_1);

}

void unarySub(float* in1, const float& in2, float* out, unsigned long length){
    thrust::transform(thrust::device,in1,in1+length,out,_1-in2);
}

void unaryDiv(float* in1, const float& in2, float* out, unsigned long length){
    thrust::transform(thrust::device,in1,in1+length,out,_1/in2);
}

void unaryPow(float* in1, float N,float* out, unsigned long length){
    if (N == 0.0f)
    {
            thrust::fill(thrust::device,out,out+length,1);
            return;
    }
    if (N == 1.0f)
    {
        thrust::copy(thrust::device,in1,in1+length,out);
        return;
    }
    if (N == 2.0f){
        thrust::square<float> op;
        thrust::transform(thrust::device,in1,in1+length,out,op);
        return;
    }
    thrust::transform(thrust::device,in1,in1+length,out,PowOp(N));
}

void unaryCompare(float* in1, const float& in2, float* out, unsigned long length, int type){
    switch (type)
    {
        case G:
            thrust::transform(thrust::device,in1,in1+length,out,CompareGOp(in2));
            break;
        case GE:
            thrust::transform(thrust::device,in1,in1+length,out,CompareGEOp(in2));
            break;
        case E:
            thrust::transform(thrust::device,in1,in1+length,out,CompareEOp(in2));
            break;
        case NE:
            thrust::transform(thrust::device,in1,in1+length,out,CompareNEOp(in2));
            break;
        case LE:
            thrust::transform(thrust::device,in1,in1+length,out,CompareLEOp(in2));
            break;
        case L:
            thrust::transform(thrust::device,in1,in1+length,out,CompareLOp(in2));
            break;
        default:
            break;
    }
}
void unaryCompare(const float& in1, float* in2, float* out, unsigned long length, int type){
    switch (type)
    {
        case G:
            thrust::transform(thrust::device,in2,in2+length,out,CompareLOp(in1));
            break;
        case GE:
            thrust::transform(thrust::device,in2,in2+length,out,CompareLEOp(in1));
            break;
        case E:
            thrust::transform(thrust::device,in2,in2+length,out,CompareEOp(in1));
            break;
        case NE:
            thrust::transform(thrust::device,in2,in2+length,out,CompareNEOp(in1));
            break;
        case LE:
            thrust::transform(thrust::device,in2,in2+length,out,CompareGEOp(in1));
            break;
        case L:
            thrust::transform(thrust::device,in2,in2+length,out,CompareGOp(in1));
            break;
        default:
            break;
    }
}
void unaryCompare(float* in1, float* in2, float* out, unsigned long length, int type){
    switch (type)
    {
        case G:
            thrust::transform(thrust::device,in1,in1+length,in2,out,CompareAGOp());
            break;
        case GE:
            thrust::transform(thrust::device,in1,in1+length,in2,out,CompareAGEOp());
            break;
        case E:
            thrust::transform(thrust::device,in1,in1+length,in2,out,CompareAEOp());
            break;
        case NE:
            thrust::transform(thrust::device,in1,in1+length,in2,out,CompareANEOp());
            break;
        case LE:
            thrust::transform(thrust::device,in2,in2+length,in1,out,CompareAGEOp());
            break;
        case L:
            thrust::transform(thrust::device,in2,in2+length,in1,out,CompareAGOp());
            break;
        default:
            break;
    }
}

void thrustFill(float* aBegin, float* aEnd, float aValue)
{
    thrust::fill(thrust::device, aBegin, aEnd, aValue);
}

void thrustFill(float* aBegin, float* aEnd, std::complex<float> aValue)
{
    thrust::fill(thrust::device, (std::complex<float>*)aBegin, (std::complex<float>*)aEnd, aValue);
}