GPU coder support for GPU arrayfun in parallel computation toolbox

I wonder if the GPU coder support the arrayfun that is commonly used in the GPU parallel computation toolbox ?

In the following testFun, the variable I is seen by both the parent function testFun and sub-function myFun (It works just like the Stencil Operations). The GPU coder app generates the following error message when calling myFun function handle within the arrayfun :

??? Function handles cannot be passed to extrinsic functions.

function output = testFun(I,d) %#codegen

    I = gpuArray(I);
    Row = (1:size(I,1))'; Col = (1:size(I,2);
    StepSize = d;
    output = arrayfun(@myFun,Row,Col,StepSize);
    
    function  Iavg = myFun(X,Y,S)
       Iavg = (I(X,Y)/S + I(X-S,Y-S))/2;
    end
end

I don't know how to write such a GPU code without using the arrayfun combined with the extrinic funcion handle to share the variable between GPU workers ... Is there anyway to make the mex of this function?

[Chaowen Liang]

This question has been posted for a while but no response. I appreciate a lot if there is someone there help me with this issue?... Even just a hint..

Answer 1

[Infinite_king]

Hi Chaowen Liang,

I understand that you are trying to generate a GPU MEX for a function in which the 'arrayfun' function has been used.

The 'Parallel Computing Toolbox' (PCT) is typically employed to transition computations to the GPU with minimal adjustments to a standard MATLAB script. By employing 'gpuArray,' data can be transferred to the GPU, allowing operations to be performed on the GPU itself. PCT supports most standard MATLAB functions and efficiently utilizes the GPU.

GPU Coder was employed to generate CUDA code from the MATLAB script. The resulting code is adaptable and can be seamlessly integrated into existing code bases. Furthermore, GPU Coder can be utilized to generate GPU MEX files, offering significantly faster execution compared to the original MATLAB script.

It's important to note that 'GPU Coder' does not currently support the 'gpuArray' and 'arrayfun' functions. This limitation arises because 'GPU Coder' expects a standard MATLAB script and automatically generates optimized CUDA code, including necessary CPU to GPU and GPU to CPU memory transfers.

So, the code can be simply written as,

function output = testFun(I,d) %#codegen

    Row = (1:size(I,1))'; Col = (1:size(I,2);
    
    S = d;
    X = Row;
    Y = Col;
    
    Iavg = (I(X,Y)/S + I(X-S,Y-S))/2;

    output = Iavg;
end

Note :- A quick observation: Ensure that the indices are always positive and within range. Based on the code, it appears that the indices for the expression I(X-S,Y-S)) may include 0 or negative values.

For more information on how to use GPU Coder, refer the following MATLAB Documentation, https://www.mathworks.com/help/gpucoder/getting-started-with-gpu-coder.html

For more information on how to PCT to shift computation from CPU to GPU, refer the following documentation

https://www.mathworks.com/products/parallel-computing.html#:~:text=Accelerate%20MATLAB%20with%20GPUs

Hope this is helpful.

Answer 2

[Ram Kokku]

@Infinite_king - Thank you for your response contrasting PCT and GPU Coder approaches. GPU Coder automatically identifies parallel regions, creates kernels, and inserts necessary CPU and GPU copies. Hence, marking variables gpuArray will not work. This is different from how Parallel Computing Toolbox works.

@Chaowen Liang - apologies for the late response.

Stecil like operations can be expressed using GPU Coder's stencilfun idiom. Here is the documentation for this https://www.mathworks.com/help/gpucoder/ref/stencilfun.html. However, for your example, we just need to remove gpuArray constructs from the testFcn. This generates kernels for arrayfun as expected. This is the modified code that works with GPU Coder

function output = testFun(I,d) %#codegen
%% I = gpuArray(I); commented out this line of code 
Row = (1:size(I,1))'; Col = (1:size(I,2))';
StepSize = d;
output = arrayfun(@myFun,Row,Col,repmat(StepSize, 1000,1));

    function  Iavg = myFun(X,Y,S)
        Iavg = (I(X,Y)/S + I(X-S,Y-S))/2;
    end
end