線程安全函數指針gsl monte carlo集成

我想使用gsl monte carlo集成庫，我目前正在運行線程安全問題，目前我還無法解決這個問題。線程安全函數指針gsl monte carlo集成

我有一個靜態函數gsl_func_wrapper，它包含具有gsl（非成員函數指針）所需的正確類型。這是傳遞給gsl蒙特卡洛積分器內部的積分成員函數。

class IntegralStrategyGSL2D: public IntegralStrategy2D { 
    static Model2D *current_model; 

    const gsl_rng_type *T; 
    gsl_rng *r; 
    gsl_monte_vegas_state *s; 
    size_t calls; 
    size_t maxcalls; 

    static double gsl_func_wrapper(double *x, size_t dim, void *params) { 
    return current_model->eval(x); 
    } 
... 
} 

IntegralStrategyGSL2D::IntegralStrategyGSL2D() : 
    T(gsl_rng_default) { 
    calls = 500; // keep calls low at first 
    maxcalls = 100000; 
    gsl_rng_env_setup(); 
    r = gsl_rng_alloc(T); 
    s = gsl_monte_vegas_alloc(2); 
} 

double IntegralStrategyGSL2D::Integral(Model2D *model2d, double xlow, 
    double xhigh, double ylow, double yhigh, double precision) { 
    current_model = model2d; 
    double result, error; 

    double xl[2] = { xlow, ylow }; 
    double xu[2] = { xhigh, yhigh } 

    gsl_monte_function G = { &gsl_func_wrapper, 2, 0 }; 

    gsl_monte_vegas_init(s); 

    while (true) { 
    gsl_monte_vegas_integrate(&G, xl, xu, 2, calls, r, s, &result, &error); 
    ... 
    } 
    ... 
}

所有這些工作，如果我運行1線程，並且如果我使用多個線程，我得到nan積分一次。因爲gsl需要這個非成員函數指針，所以我現在並沒有真正看到如何讓這個線程安全......有人能指出我正確的方向嗎？ Thx

編輯：我添加了將由gsl_func_wrapper調用的eval函數。

double GaussianModel2D::eval(const double *x) const { 
    // see wikipedia definition 

    double normalization = 1.0 
     /(2.0 * M_PI * gauss_sigma_var1->getValue() 
       * gauss_sigma_var2->getValue() 
       * sqrt(1 - pow(gauss_rho->getValue(), 2.0))); 

    double exp_value = exp(
     -(pow(x[0] - gauss_mean_var1->getValue(), 2.0) 
       /pow(gauss_sigma_var1->getValue(), 2.0) 
       + pow(x[1] - gauss_mean_var2->getValue(), 2.0) 
         /pow(gauss_sigma_var2->getValue(), 2.0) 
       + 2.0 * gauss_rho->getValue() 
         * (x[0] - gauss_mean_var1->getValue()) 
         * (x[1] - gauss_mean_var2->getValue()) 
         /(gauss_sigma_var1->getValue() 
           * gauss_sigma_var2->getValue())) 
       /(2.0 * (1 - pow(gauss_rho->getValue(), 2.0)))); 

    return gauss_amplitude->getValue() * normalization * exp_value; 
}

EDIT2： 我加入了呼籲的整體功能由Mike的要求。所以下面是進行Integral調用的兩個相關代碼塊。我真的不明白爲什麼current_model變量的thread_local關鍵字修復了這個問題。事件雖然代碼有點難看，但是current_model變量不應該總是指向相同的model2d指針嗎？所以，即使兩個線程可能有衝突的寫入，寫入的值始終是相同的，所以它應該不重要或？

double smearing_probability = divergence_model->Integral(int_range, 1e-4); 

double Model2D::Integral(const std::vector<DataStructs::DimensionRange &ranges, 
    double precision) { 
    return integral_strategy->Integral(this, ranges[0].range_low, 
     ranges[0].range_high, ranges[1].range_low, ranges[1].range_high, 
     precision); 
}

來源

2016-02-19 steve

顯示'Model2D :: eval'的代碼。 –

你可以使用類作爲函數指針來避免這個問題，但不是很好。 –

@Pete：Model2D :: eval是一個純虛擬函數，但對於我來說它是一個非常複雜的函數，所以我不確定是否應該在這裏發佈所有這些...是否有關聯？ – steve

假設Integral被稱爲在同一個線程，執行gsl_func_wrapper，您可以使用thread_local代替static作爲current_model存儲類。

這會給你一個單獨的current_model變量爲每個線程。

來源

2016-02-19 16:31:15 MikeMB

我測試了它，它似乎工作，我會明天測試一下更徹底，並給予更新。非常感謝邁克。 – steve

好的，所以它只是偶爾比平時跑得長一些。它仍然給我nan的積分一次偶爾的結果... – steve

@steve：那麼問題很可能不是（只）在你告訴我們的代碼中。請提供[mcve]（/ help/mcve）。你正在討論多線程的問題，但是現在它完全不清楚，你在哪裏產生這些線程以及每個線程訪問什麼數據 – MikeMB

好吧，我設法弄清楚了這個問題。這是這兩個變量

const gsl_rng_type *T; gsl_rng *r; gsl_monte_vegas_state *s;

都是實例變量，因此以同樣的方式使用的所有線程。如果我在每個線程的Integral（）函數內分別創建這些函數，那麼每個線程都會運行良好...

來源

2016-02-21 12:03:10 steve

線程安全函數指針gsl monte carlo集成

回答

相關問題