추력 및 부스트 조합 odeint; 잘못된 포인터 전달 segfault를 일으킬 것으로 생각했습니다

Question

많은 초기 조건에 대한 초기 값 문제를 해결하기 위해 lorenz_parameters.cu 예제를 boost::odeint에서 수정하려고합니다. 이 프로그램은 segfaults, 나는 내 우연히 포인터를 integrate_adaptive()의 호출에 전달하여 파생물을 저장하기위한 벡터의 크기를 인쇄하려고 할 때 발생한다고 믿습니다.

참고로 여기에 documentation for the example 및 code for the same examples이 있습니다 (어떤 이유로 원본 문서의 코드 링크가 작동하지 않음).

아래에는 프로그램의 전체 소스 코드가 첨부되어 있습니다. UNEXPECTED OUTPUT!으로 표시된 줄이 있는데, 여기에서 파생 벡터의 크기를 인쇄하려고 시도하고 상태 벡터 (512)와 동일한 크기 여야하는 0으로 표시됩니다. LINE I DO NOT UNDERSTAND이라고 표시된 줄이 있는데, 나는 이것이 내가 잘못 생각하고있는 부분 일 것이라고 생각합니다. 이 줄은 intergrate_adaptive()에 대한 호출의 일부이며이 인수는 parallel_initial_condition_problem의 operator() 함수에 전달됩니다 (결국/odeint 패키지를 통해 전달되는 것입니다). 그러나 나는 그 예에서와 똑같이하려고 노력하고 있지만, 아직 효과가 없습니다. 어떤 도움이라도 대단히 감사합니다!

#include <iostream> 
#include <vector> 
#include <thrust/device_vector.h> 
#include <boost/numeric/odeint.hpp> 
#include <stdio.h> 
#include <iostream> 
#include <fstream> 
#include <typeinfo> 

typedef double value_type; 

//change this to device_vector<...> of you want to run on GPU 
typedef thrust::host_vector<value_type> state_type; 

using namespace std; 
using namespace boost::numeric::odeint; 

const size_t D = 16; // how many values of each initial condition to try 
const size_t N_ICS = D * D; // number of initial conditions 
const size_t N_VARS = 2; // number of variables in the system 

// time constants 
const value_type t_0 = 0.0; 
const value_type t_final = 0.1;  
const value_type dt = 0.01; 

struct parallel_initial_condition_problem { 
    struct functor 
    { 
    template< class T > 
    __host__ __device__ 
    void operator()(T t) const 
    { 
     const value_type k_b = 0.45/6.0; 
     const value_type k_f = 1.0; 
     const value_type k_d = 1.0; 

     // unpack the variables (and the parameters we want to vary, if any) 
     value_type f = thrust::get<0>(t); 
     value_type a = thrust::get<1>(t); 
     // value_type df = thrust::get<2>(t); // this also causes a segfault 
     // value_type da = thrust::get<3>(t); 
     // calculate the derivative 
     thrust::get<2>(t) = -0.01; // minimal value for debugging 
     thrust::get<3>(t) = -0.01; 
    } 
    }; 

    parallel_initial_condition_problem(size_t N) 
    : m_N(N) { } 

    template< class State , class Deriv > 
    void operator()(const State &x , Deriv &dxdt , value_type t) const 
    { 
    //// Printing out debugging info (here is the error..but why?) 
    cout << "INSIDE operator()" << endl; 
    cout << "t: " << t << endl;  
    cout << "size of x: \t" << boost::end(x) - boost::begin(x) << endl; 
    cout << "size of dxdt: \t" << boost::end(dxdt) - boost::begin(dxdt) << endl; // UNEXPECTED OUTPUT!! 
    cout << endl; 

    thrust::for_each(
     thrust::make_zip_iterator( 
     thrust::make_tuple(
       boost::begin(x) + 0*m_N, 
       boost::begin(x) + 1*m_N , 
       boost::begin(dxdt) + 0*m_N, 
       boost::begin(dxdt) + 1*m_N 
       ) 
       ) , 
     thrust::make_zip_iterator( 
     thrust::make_tuple(
       boost::begin(x) + 1*m_N, 
       boost::begin(x) + 2*m_N , 
       boost::begin(dxdt) + 3*m_N, 
       boost::begin(dxdt) + 4*m_N 
       ) 
       ) , 
     functor()); 
    } 

    size_t m_N; 
}; 

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

    // x stores both (initial) state (x) [0:N_VARS*N_ICS] and dxdt [N_VARS*N_ICS:END] 
    state_type x(2 * N_VARS * N_ICS); 

    // DUMMY INITIAL CONDITIONS FOR DEBUGGING PURPOSES 
    // fill each value with its index 
    for(int i=0;i<x.size();i++) { 
    x[i] = 0.01+i; 
    } 

    // system function 
    parallel_initial_condition_problem init_con_solver(N_ICS); 

    //////////////////////////////// Debugging info: (all seems correct here) 
    cout << "POINTERS" << endl; 
    cout << "size of x: " << boost::end(x) - boost::begin(x) << endl; 
    cout << "x :" << &(*x.begin()) << endl; 
    cout << "dxdt :" << &(*x.begin()) + N_VARS*N_ICS << endl; 
    cout << "end :" << &(*x.end()) << endl; 
    cout << endl; 

    /////////////////////////// STATE_T  VALUE_T  DERIV_T  TIME_T 
    typedef runge_kutta_dopri5< state_type , value_type , state_type , value_type > stepper_type; 
    const value_type abs_err = 1.0e-6; 
    const value_type rel_err = 1.0e-6; 
    integrate_adaptive(make_controlled(abs_err , rel_err , stepper_type()), 
       init_con_solver , 
       std::make_pair(x.begin() , x.begin() + N_VARS*N_ICS ), /// LINE I DO NOT UNDERSTAND 
       t_0, t_final, dt); 

} 

이 실행의 출력은 다음과 같습니다

$ g++ main.c && ./a.out 
POINTERS 
size of x: 1024 
x :0x9a2010 
dxdt :0x9a3010 
end :0x9a4010 

INSIDE operator() 
t: 0 
size of x: 512 
size of dxdt: 0 

Segmentation fault (core dumped) 

Answer 1

나는 일부가 누락 포함되어 있다고 생각합니다. 사용해보기

#include <boost/numeric/odeint/external/thrust/thrust.hpp> 

사실, 왜 코드가 컴파일되는지 궁금합니다. 당신이 부스트 (1.55)에서 현재 odeint 버전을 사용 당신은 또한 대수를 지정해야하고, 작업을 입력하면 :

typedef runge_kutta_dopri5< 
    state_type , 
    value_type , 
    state_type , 
    value_type , 
    thrust_algebra , 
    thrust_operations > stepper_type; 

이 다음 부스트 릴리스에서 변경해야합니다. oihint의 gihub 버전에서는 코드가 잘되어야합니다.