将openmp与odeint和自适应步长配合使用 [英] Using openmp with odeint and adaptative step sizes
问题描述
我正在尝试使用openmp来并行化我的代码.当我使用恒定的步长时，一切工作都很好，但是当我使用自适应步进器运行相同的代码时，会出现我不理解的错误.
I am trying to use openmp to parallelize my code. Everything works just fine when I use constant step sizes, however when I run the same code using an adaptative stepper I get errors that I don't understand.
以下是代码的重要部分:
Here are the essential parts of the code :
using namespace std;
using namespace boost::numeric::odeint;
const int jmax = 10;
typedef double value_type;
typedef boost::array<value_type ,2*(jmax+1) > state_type;
//The step function
void rhs( const state_type A , state_type &dAdt , const value_type t )
{
value_type RHStemp0;
value_type RHStemp1;
//We will write the RHS of the equations a big sum
#pragma omp parallel for schedule(runtime)
for(int j = 0; j < jmax+1 ; j++ ) //Real part
{
RHStemp0 = value_type(0.0);
RHStemp1 = value_type(0.0);
for (int k = 0; k< jmax+1 ;k++)
{
for (int l = max(0,j+k-jmax); l < 1 + min(jmax,j+k);l++)
{
RHStemp0 = RHStemp0 + S[j*SIZE_S*SIZE_S + k*SIZE_S + l]*(-A[k+jmax+1]*A[l]*A[j+k-l] + A[k]*A[l+jmax+1]*A[j+k-l]
+ A[k]*A[l]*A[j+k-l+jmax+1] +A[k+jmax+1]*A[l+jmax+1]*A[j+k-l+jmax+1]);
RHStemp1 = RHStemp1 + S[j*SIZE_S*SIZE_S + k*SIZE_S + l]*(A[k]*A[l]*A[j+k-l] - A[k]*A[l+jmax+1]*A[j+k-l+jmax+1]
+ A[k+jmax+1]*A[l]*A[j+k-l+jmax+1] +A[k+jmax+1]*A[l+jmax+1]*A[j+k-l]);
}
}
dAdt[j] = (-1/(value_type((2*(2*j+3)))))*RHStemp0;
dAdt[j+jmax+1] = (1/(value_type((2*(2*j+3)))))*RHStemp1;
}
}
int main()
{
const state_type initial = loadInitialData(); //Initial condition
omp_set_num_threads(jmax+1);
int chunk_size = jmax/omp_get_max_threads();
omp_set_schedule( omp_sched_dynamic, chunk_size );
//I define my controlled error steppers
typedef runge_kutta_fehlberg78< state_type , value_type ,
state_type , value_type,openmp_range_algebra> error_stepper_type;
typedef controlled_runge_kutta< error_stepper_type > controlled_stepper_type;
controlled_stepper_type controlled_stepper;
int steps = integrate_adaptive( controlled_stepper ,rhs ,
initial, TINITIAL , TFINAL,INITIAL_STEP , push_back_state_and_time( A_vec , times ) );
}
我没有显示所有变量的定义，但是我怀疑它们是问题所在，因为如果我只是从error_stepper_type的定义中删除openmp_range_algebra选项，这将很好地工作.如果我使用具有恒定步进大小的openmp_range_algebra(例如4阶的Runge Kutta)，这也可以正常工作.
I do not show the definition of all variables but I doubt that they are the problem, since this works fine if I just remove the openmp_range_algebra option from the definition of error_stepper_type. This works also fine if I use openmp_range_algebra with a constant stepper size, like the Runge Kutta of order 4.
但是，使用此代码，我得到以下错误:
However, with this code I get the following errors :
invalid conversion from 'boost::range_iterator<const boost::array<double, 22ull>, void>::type {aka const double*}' to 'boost::range_iterator<boost::array<double, 22ull>, void>::type {aka double*}' [-fpermissive]|
因此，看来我尝试分配常量.此错误出现在文件以下文件的openmp_range_algebra.hpp中:
So it seems that I try to allocate sonmething that is constant. This errors appears in the file openmp_range_algebra.hpp, in the following code :
template< class S >
static typename norm_result_type< S >::type norm_inf( const S &s )
{
using std::max;
using std::abs;
typedef typename norm_result_type< S >::type result_type;
result_type init = static_cast< result_type >( 0 );
const size_t len = boost::size(s);
typename boost::range_iterator<S>::type beg = boost::begin(s);
#pragma omp parallel for reduction(max: init) schedule(dynamic)
for( size_t i = 0 ; i < len ; ++i )
init = max( init , abs( beg[i] ) );
return init;
}
我希望我已经足够清楚了，我只希望能够对我的并行化代码使用自适应步进器.
I hope I have been clear enough, I just would like to be able to use adaptative stepper with my parallelized code.
非常感谢您的帮助.
推荐答案
This is a bug in odeint, I've filed it on github: https://github.com/headmyshoulder/odeint-v2/issues/166 and I will try to fix it asap. Thanks for posting.
已修复，您的程序现在应该可以编译.另外，我更改了一个使用自适应步进器的示例: https://github.com/headmyshoulder/odeint-v2 /blob/master/examples/openmp/lorenz_ensemble_simple.cpp
edit: fixed, your program should now compile. Also, I've changed an example to use adaptive stepper: https://github.com/headmyshoulder/odeint-v2/blob/master/examples/openmp/lorenz_ensemble_simple.cpp
这篇关于将openmp与odeint和自适应步长配合使用的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
