수치 적 통합 - 어떻게 병렬화할까요?

Question

OpenCL을 시작하면서 벡터 추가 예제를보고 이해할 수있었습니다. 그러나 나는 사지의 방법을 생각하고있었습니다. 이것은 [a, b]의 x^2에 대한 적분 계산을위한 코드 (C)입니다.

double f(double x) 
{ 
    return x*x; 
} 

double Simple_Trap(double a, double b) 
{ 
    double fA, fB; 
    fA = f(a); 
    fB = f(b); 
    return ((fA + fB) * (b-a))/2; 
} 

double Comp_Trap(double a, double b) 
{ 
    double Suma = 0; 
    double i = 0; 
    i = a + INC; 
    Suma += Simple_Trap(a,i); 
    while(i < b) 
    { 
     i+=INC; 
     Suma += Simple_Trap(i,i + INC); 
    } 
    return Suma; 
} 

질문은 사다리꼴 방법을 사용한 적분 계산을위한 커널을 얻는 방법입니까?

---

그래서, 난 아이디어에 대해 생각했다 : 파셜은 [I] = 언급 Patrick87로 파셜의 합을 계산하는 커널을 다음 (A, A + 오프셋) 통합합니다.

하지만 이것이 최선의 방법입니까?

Answer 1

다음은 내가 생각해 낸 것입니다. 필자는이 커널에 대한 엔드 투 엔드 테스트를 수행하지 않았습니다. 나는 조금 더 시간이있을 때 나는 갱신을 할 것이다.

comp_trap은 위의 코드를 기반으로 기본 나누기 & 정복 방법입니다. comp_trap_multi는 각 작업 항목이 하위 섹션을 나누어 가짐으로써 정확도를 높입니다.

각 작업 그룹이 결과를 반환 할 수 있도록 두 배의 배열 만 호스트에 할당하면됩니다. 이것은 피해야 할 벡터 할당을 줄이는 데 도움이됩니다.

문제가 있으면 알려주십시오.

업데이트 :

1) 부유 모든 이중 참조 변경, 오픈 CL의 선택 사항입니다

2) 하드 코딩이 값이 내 시스템에 최적 (64)에 작업 그룹 크기, 더블 때문에해야 실험적으로 결정되어야한다. 호스트 프로그램이 궁극적으로 대상 시스템의 최적 값만 사용하기 때문에이 값을 하드 코딩하여 로컬 float 배열을 사용하는 것보다 선호합니다.

3) 잘못된 계산 (A1 잘못했다,

/* 
numerical-integration.cl 
*/ 

float f(float x) 
{ 
    return x*x; 
} 

float simple_trap(float a, float b) 
{ 
    float fA, fB; 
    fA = f(a); 
    fB = f(b); 
    return ((fA + fB) * (b-a))/2; 
} 

__kernel void comp_trap(
    float a, 
    float b, 
    __global float* sums) 
{ 
/* 
- assumes 1D global and local work dimensions 
- each work unit will calculate 1/get_global_size of the total sum 
- the 0th work unit of each group then accumulates the sum for the 
group and stores it in __global * sums 
- memory allocation: sizeof(sums) = get_num_groups(0) * sizeof(float) 
- assumes local scratchpad size is at lease 8 bytes per work unit in the group 
ie sizeof(wiSums) = get_local_size(0) * sizeof(float) 
*/ 
    __local float wiSums[64]; 
    int l_id = get_local_id(0); 

    //cumpute range for this work item is: a1, b1 
    float a1 = a+((b-a)/get_global_size(0))*get_global_id(0); 
    float b1 = a1+(b-a)/get_global_size(0); 

    wiSums[l_id] = simple_trap(a1,b1); 

    barrier(CLK_LOCAL_MEM_FENCE); 

    int i; 
    if(l_id == 0){ 
     for(i=1;i<get_local_size(0);i++){ 
      wiSums[0] += wiSums[i]; 
     } 
     sums[get_group_id(0)] = wiSums[0]; 
    } 
} 

__kernel void comp_trap_multi(
    float a, 
    float b, 
    __global float* sums, 
    int divisions) 
{ 
/* 
- same as above, but each work unit further divides its range into 
'divisions' equal parts, yielding a more accurate result 
- work units still store only one sum in the local array, which is 
used later for the final group accumulation 
*/ 
    __local float wiSums[64]; 
    int l_id = get_local_id(0); 

    float a1 = a+((b-a)/get_global_size(0))*get_global_id(0); 
    float b1 = a1+(b-a)/get_global_size(0); 
    float range; 
    if(divisions > 0){ 
     range = (b1-a1)/divisions; 
    }else{ 
     range = (b1-a1); 
    } 

    int i; 
    wiSums[l_id] = 0; 
    for(i=0;i<divisions;i++){ 
     wiSums[l_id] += simple_trap(a1+range*i,a1+range*(i+1)); 
    } 

    barrier(CLK_LOCAL_MEM_FENCE); 

    if(l_id == 0){ 
     for(i=1;i<get_local_size(0);i++){ 
      wiSums[0] += wiSums[i]; 
     } 
     sums[get_group_id(0)] = wiSums[0]; 
    } 
} 

Answer 2

사다리꼴 방법은 리만 합계를 수행하는 것의 약간의 미세 조정입니다. 이 작업을 병렬로 수행하려면 간격을 스레드만큼 원하는만큼의 하위 간격으로 나누고 싶을 것입니다. 그런 다음 각 스레드가 해당 하위 구간에서 함수를 통합하게하십시오. 마지막으로, 이전 단계에서 계산 된 모든 적분에 대한 전체 합계 감소를 수행합니다. 각 스테이지에 사용할 스레드 수를 실험 할 수 있습니다.