c++ - 当使用超过 1 个逻辑核心时，TBB 减少会导致段错误

Question

我使用 TBB 实现了 Jacobi 算法，它工作得很好。然后我使用缩减并行化收敛计算，但出于某种原因，如果我使用超过 1 个逻辑核心，我会遇到段错误，我无法弄清楚原因。

我可以在只有 1 个逻辑核心的系统上使用 1 个以上的线程。

使用 OpenMP 的相同实现可以毫不费力地工作

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <time.h>
#include <string.h>
#include <tbb/parallel_for.h>
#include <tbb/parallel_reduce.h>
#include <tbb/blocked_range.h>
#include <tbb/task_scheduler_init.h>
#include <tbb/tick_count.h>

// ----------------------------------------------------------------

#define SIZE 1024
#define RESIDUO 0.0009f*SIZE
#define THREADS 2

using namespace tbb;

// ----------------------------------------------------------------

struct Sum {
    float ret;
    float (*a)[SIZE];
    float (*x);
    float (*b);
    Sum(float A[SIZE][SIZE], float X[SIZE], float B[SIZE]) : ret(0), a(A), x(X), b(B) {}
    Sum( Sum&, split ) {ret = 0;}
    void operator()( const blocked_range<int>& r ) {
        float temp = ret;
        for( int i = r.begin(); i != r.end(); i++ ) {
            float sum = 0.0f;
            for(int j = 0; j < SIZE; j++){
                sum += a[i][j] * x[j];
            }
            temp += pow(b[i] - sum, 2);
        }
        ret = temp;
    }
    void join( Sum& rhs ) {ret += rhs.ret;}
};


/*
    // || b - Ax ||
*/

int converge(float a[SIZE][SIZE], float x[SIZE], float b[SIZE]){

    Sum total(a, x, b);
    parallel_reduce( blocked_range<int>(0, SIZE), total );

    float norm = sqrt(total.ret);
    printf("Ret: %f | Residuo: %f\n", total.ret, norm);

    return (norm <= RESIDUO);
}

// ----------------------------------------------------------------

float randomFloat()
{
    float r = (float)rand()/(float)RAND_MAX;
    return r;
}

// ----------------------------------------------------------------

int check_ddm(float (*a)[SIZE]){
    float sum = 0.0f;
    int i = 0, j = 0;

    for(i = 0; i < SIZE; i++){
        sum = 0.0f;
        for(j = 0; j < SIZE; j++){
            if(i != j){
                sum += a[i][j];
            }
        }
        if(a[i][i] < sum){
            printf("line: %d, sum: %f, a[i][i]: %f  \n", i, sum, a[i][i]);
            for(j = 0; j < SIZE; j++){
                if(i != j) printf("%f ", a[i][j]);
                else printf("(%f) ", a[i][j]);
            }
            printf("\n");
            return 0;
        }
    }
    return 1;
}

// ----------------------------------------------------------------

int generate_ddm(float (*a)[SIZE], float *b)
{
    int i = 0, j = 0;
    float line = 0.0f;
    for(i = 0; i < SIZE; i++){
        line = 0.0f;
        for(j = 0; j < SIZE; j++){
            if(i != j){
                a[i][j] = randomFloat();
            }
            line += a[i][j];
        }
        a[i][i] = SIZE;
        b[i] = line + SIZE;
    }

    return check_ddm(a);
}

// ----------------------------------------------------------------

int main(  )
{
    float (*x)[SIZE] = (float(*)[SIZE])malloc(sizeof *x * 2);
    float (*a)[SIZE] = (float(*)[SIZE])malloc(sizeof *a * SIZE);
    float (*b) = (float*)malloc(sizeof(float) * SIZE);
    int i = 0, j = 0;
    float delta = 0.0f;
    int read = 0;
    int write = 1;

    srand(time(NULL));
    tbb::task_scheduler_init init(THREADS);

    // set up initial solution
    for(i = 0; i < SIZE; i++){
        x[0][i] = i;
        x[1][i] = i;
    }

    // generate a diagonal dominant matrix
    if(!generate_ddm(a, b)){
        printf("Array generated is not ddm!\n");
        return 1;
    }

    tick_count startTime = tick_count::now();        
    while(!converge(a, x[write], b)){
        read = !read;
        write = !write;

        parallel_for(blocked_range<int>(0,SIZE), 
            [&] (const blocked_range<int>& r) {
            for (int i = r.begin(); i < r.end(); i++) {
                float delta = 0.0f;
                for(int j = 0; j < SIZE; j++){
                    if(j != i){
                        delta += a[i][j] * x[read][j];
                    }
                }
                x[write][i] = (b[i] - delta) / a[i][i];
            }
        });
    }


    tick_count lastTime = tick_count::now();
    float walltime = (lastTime - startTime).seconds();

    printf("tbb %f\n", walltime);
    converge(a, x[write], b);
    printf("x0: %f | x%d: %f\n", x[write][0], SIZE-1, x[write][SIZE-1]);
    free(a);
    free(b);
    free(x);

    return 0;
}

段错误发生在 Sum 类中的以下行:

sum += a[i][j] * x[j];

如果我将该行更改为

float tmpa = a[i][j];
float tmpx = x[j];
sum += tmpa * tmpx;

错误一直在

sum += tmpa * tmpx;

Answer 1

在原始版本中，“拆分构造函数”未定义 a、x 和 b。它们需要从传入的 Sum& 参数中复制。例如，将拆分构造函数更改为:

Sum( Sum& s, split ) {a=s.a; b=s.b; x=s.x; ret = 0;}