Поиск локальных минимумов выборочной функции CUDA Thrust

Question

Я хочу написать программу с использованием Thrust, которая должна вычислять локальные минимумы данных функций, f.i. sin(x). Я сделал это, аппроксимируя производную функции конечными разностями, а затем искал те абсциссы, где производная меняет знак. Теперь я хочу собрать местные минимумы. Я отметил локальные минимумы с «1» , а остальные точки с «0». Я сделал inclusive_scan (для расчета мест на новой вкладке). Моя проблема заключается в сборе локальных минимумов с gather_if (условие трафарета, минимумы карты), , но код не компилируется, и я не знаю почему.

Может кто-нибудь объяснить, почему?

/** 
* Copyright 1993-2012 NVIDIA Corporation. All rights reserved. 
* 
* Please refer to the NVIDIA end user license agreement (EULA) associated 
* with this source code for terms and conditions that govern your use of 
* this software. Any use, reproduction, disclosure, or distribution of 
* this software and related documentation outside the terms of the EULA 
* is strictly prohibited. 
*/ 

#include <stdio.h> 
#include <thrust/device_vector.h> 
#include <thrust/gather.h> 
#include <thrust/host_vector.h> 
#include <thrust/reduce.h> 
#include <thrust/copy.h>          
#include <thrust/remove.h> 
#include <thrust/functional.h> 
#include <thrust/iterator/constant_iterator.h> 
#include <thrust/iterator/counting_iterator.h> 
#include <thrust/scan.h> 
#include <sys/time.h> 

__host__ __device__ unsigned int bitreverse(unsigned int number) { 
    number = ((0xf0f0f0f0 & number) >> 4) | ((0x0f0f0f0f & number) << 4); 
    number = ((0xcccccccc & number) >> 2) | ((0x33333333 & number) << 2); 
    number = ((0xaaaaaaaa & number) >> 1) | ((0x55555555 & number) << 1); 
    return number; 
} 

struct is_even 
{ 
    __host__ __device__ 
    bool operator()(const int x) { 
     return (x % 2) == 0; 
    } 
}; 

struct select_mine 
{ 
    __host__ __device__ 
    float operator()(const float x) { 
     return (x < 0) ? 1.0f : 0.0f; 
    } 
}; 

struct bitreverse_functor 
{ 
    __host__ __device__ unsigned int operator()(const unsigned int &x) { 
     return bitreverse(x); 
    } 
}; 

struct sign 
{ 
    __host__ __device__ float operator()(const float x) { 
     if (x > 0.0f) 
      return 1.0f; 
     if (x < 0.0f) 
      return -1.0f; 
     return 0.0f; 
    } 
}; 

struct sine: public thrust::unary_function<float, float> 
{ 
    __host__ __device__ 
    float operator()(float x) { 
     return sinf(x); 
    } 
}; 

struct absolute: public thrust::unary_function<float, float> 
{ 
    __host__ __device__ 
    float operator()(float x) { 
     if (x < 0.0f) 
      x = -x; 
     return x; 
    } 
}; 

struct lokalne_minimum : public thrust::binary_function<float,float,float> 
{ 
    __host__ __device__ 
    float operator()(float x, float y) 
    { 
     if (x > 0 && y < 0) 
      return 1.0f; 
     return 0.0f; 
    } 
}; 

struct conv : public thrust::unary_function<float,int> 
{ 
    __host__ __device__ 
    int operator()(float x) 
    { 
     return (int)(x); 
    } 
}; 

using namespace thrust; 

void help(char *arg) { 
    fprintf(stderr, 
      "Nieprawidlowe uzycie: %s [x1] [x2] [n]\nx1 - zakres od\nx2 - zakres do\nn - liczba podzialow zakresu\n", 
      arg); 
} 

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

    if (argc != 4) { 
     help(argv[0]); 
     return 1; 
    } 

    int n = atoi(argv[3]); 
    float x1 = (float) atof(argv[1]); 
    float x2 = (float) atof(argv[2]); 

    if (n < 0 || x2 < x1) { 
     help(argv[0]); 
     return 1; 
    } 

    float step = (x2 - x1)/n; 

    fprintf(stderr, "Step: %f\n", step); 

    thrust::device_vector<float> oxdata(n); 
    thrust::device_vector<float> oydata(n); 
    thrust::device_vector<float> diff(n); 
    thrust::host_vector<float> ixdata(n); 

    // FIXME change it 
    for (int i = 0; i < n; i++) 
     ixdata[i] = x1 + i * step; 

    thrust::copy(ixdata.begin(), ixdata.end(), oxdata.begin()); 

    thrust::transform(oxdata.begin(),oxdata.end(),oydata.begin(),sine()); 

    thrust::transform(oydata.begin() + 1, oydata.end(), oydata.begin(), 
      diff.begin()+1, thrust::minus<float>()); 

    thrust::copy(diff.begin(), diff.end(), ixdata.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%f, ", ixdata[i]); 
    printf ("\n"); 

    thrust::transform(diff.begin()+1,diff.end(), diff.begin(),diff.begin(),lokalne_minimum()); 


    for (int i = 0; i < n; i++) 
     printf ("%f, ", ixdata[i]); 
    printf ("\n"); 

    thrust::copy(oydata.begin(), oydata.end(), ixdata.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%f, ", ixdata[i]); 
    printf ("\n"); 

    thrust::copy(diff.begin(), diff.end(), ixdata.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%f, ", ixdata[i]); 
    printf ("\n"); 

    //thrust::inclusive_scan(diff.begin(),diff.end(),diff.begin()); 

    thrust::copy(diff.begin(), diff.end(), ixdata.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%f, ", ixdata[i]); 
    printf ("\n"); 

    thrust::device_vector<int> minima(n); 
    thrust::device_vector<int> stencil(n); 
    thrust::host_vector<int> hminima(n); 
    thrust::transform(diff.begin(),diff.end(),minima.begin(),conv()); 

    thrust::copy(minima.begin(),minima.end(),hminima.begin()); 
    thrust::copy(minima.begin(),minima.end(),stencil.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%d, ", hminima[i]); 
    printf ("\n"); 

    thrust::inclusive_scan(minima.begin(), minima.end(),minima.begin()); 

    thrust::copy(minima.begin(),minima.end(),hminima.begin()); 

    for (int i = 0; i < n; i++) 
     printf ("%d, ", hminima[i]); 
    printf ("\n"); 

    //thrust::gather_if(minima.begin(),minima.end(),stencil.begin(),ixdata.begin(),ixdata.begin()); 

    return 0; 
} 

Answer 1

Это очень поздний ответ, предоставленный для устранения этого вопроса из списка без ответа. Я получаю удовольствие от комментариев Роберта Кровеллы и показывая ниже полный рабочий код, чтобы найти локальные минимумы выборочной функции с CUDA Thrust. Обоснование кода выглядит следующим образом:

Производная функции аппроксимируется центральными различиями как приложение thrust::transform;

Функции точек отбора проб отмечены знаком «1» как приложение thrust::transform путем поиска знаковых изменений производной по предикату local_minima_check();

Число локальных минимумов засчитывается как приложение thrust::count;

Местные минимумы изолированы как приложение thrust::copy_if.

#include <stdio.h> 
#include <thrust/count.h> 
#include <thrust/device_vector.h> 
#include <thrust/host_vector.h> 
#include <thrust/copy.h>          
#include <thrust/iterator/constant_iterator.h> 
#include <thrust/iterator/counting_iterator.h> 

/**************/ 
/* COS STRUCT */ 
/**************/ 
struct cosine: public thrust::unary_function<float, float> 
{ 
    __host__ __device__ float operator()(float h_x) { return cosf(h_x); } 
}; 

/******************************************/ 
/* SECOND ORDER CENTRAL DIFFERENCE STRUCT */ 
/******************************************/ 
struct second_order_central_difference 
{ 
    __host__ __device__ float operator()(thrust::tuple<float,float,float> t) 
    { 
     float f_1, f0, f1; 

     thrust::tie(f_1,f0,f1) = t; 

     return f_1 - 2.0f * f0 + f1; 
    } 
}; 

/******************************/ 
/* LOCAL MINIMUM CHECK STRUCT */ 
/******************************/ 
struct local_minimum_check:public thrust::binary_function<float,float,float> 
{ 
    __host__ __device__ float operator()(float x, float y) 
    { 
     if (x < 0 && y > 0) return 1.0f; 
     return 0.0f; 
    } 
}; 

/****************************************/ 
/* LOCAL MINIMUM PREDICATE CHECK STRUCT */ 
/****************************************/ 
struct pred 
{ 
    __host__ __device__ bool operator()(const int d_x) { return (d_x == 1.f); } 
}; 

void main() { 

    // --- Input parameters 
    int n = 100;   // Number of sampling points 
    float x1 = 3.14f/2.f;  // (x1,x2) is the sampling interval 
    float x2 = 1.5f * 3.14f; 

    // --- Calculating the sampling points x 
    thrust::host_vector<float> h_x(n); 
    float step = (x2 - x1)/n; 
    for (int i = 0; i < n; i++) h_x[i] = x1 + (float)i * step; 
    thrust::device_vector<float> d_x = h_x; 

    // --- Evaluating the function values y = f(x) 
    thrust::device_vector<float> d_y(n); 
    thrust::transform(d_x.begin(),d_x.end(),d_y.begin(),cosine()); 

    // --- Computing first order central finite differences 
    //  In Matlab's notation, it calculates d_diff1(1:n-2) = d_y(3:n,:) - d_y(1:n-2,:); 
    thrust::device_vector<float> d_diff1(n-2); 
    thrust::transform(d_y.begin() + 2, d_y.end(), d_y.begin(), d_diff1.begin(), thrust::minus<float>()); 

    // --- Computing second order central finite differences 
    //  In Matlab's notation, it calculates d_diff2(1:n-2) = d_y(3:n) - 2. * d_y(2:n-1) + d_y(1:n-2); 
    thrust::device_vector<float> d_diff2(n-2); 
    thrust::transform(thrust::make_zip_iterator(
      thrust::make_tuple(d_y.begin(), d_y.begin() + 1, d_y.begin() + 2)), 
         thrust::make_zip_iterator(
      thrust::make_tuple(d_y.end() - 2, d_y.end() - 1, d_y.end())), 
         d_diff2.begin(),second_order_central_difference()); 

    // --- Setting a flag for all those points for which the derivative changes sign from negative to positive 
    thrust::device_vector<float> d_fo_derivative(n-3); 
    thrust::transform(d_diff1.begin(), d_diff1.end() - 1, d_diff1.begin() + 1, d_fo_derivative.begin(), local_minimum_check()); 

    // --- Counting the number of local minima and selecting the local minima coordinates 
    int min_number = thrust::count(d_fo_derivative.begin(), d_fo_derivative.end(), 1.f); 
    thrust::device_vector<float> d_x_minima(min_number); 
    thrust::copy_if(d_x.begin() + 1, d_x.end() - 1, d_fo_derivative.begin(), d_x_minima.begin(), pred()); 

    for (int i = 0; i < d_x_minima.size(); i++) { 
     printf ("Local minimum # %i = %f\n ", i+1, (float)d_x_minima[i]); 
    } 

    getchar(); 

}