I have a problem with automatic GPU memory management when using CUDA. The situation is the following:

I create a simple CUDA function that operates on a 2D array (image), which does not take much memory space. This 2D array is passed as input then modified and output again. The function is loaded into mathematica using CUDAFuncitonLoad (pseudocode):

CUDATestFunction = CUDAFunctionLoad[...,{"Float", 2, "InputOutput"},{16,16}];

The function compiles fine and works. Then I create some sort of loop that involves calling the CUDA-function mentioned before, assigning the output to a new variable and then doing the whole thing again till some number of iterations is reached (again only pseudocode):

Array2Dold = imageData; (*assign raw image data*)
    Array2Dnew = CUDATestFunction[Array2Dold];
    Array2Dold = Array2Dnew;

The problem I have now is that when the number of iterations gets high I am running out of GPU memory and I get the error-message:"CUDAFunction::outmem: CUDALink ran out of available memory, possibly due to not freeing memory using the memory manager".

Does anybody know why that happens?

My guess is that mathematica for some reason keeps the Array2Dold or some modified version of it after calculation in CUDA in GPU memory and does not free the memory when the CUDA-funtion (CUDATestFunction) is finished. I know that there is a CUDALink Memory Manager, but from what I read in the mathematica documentation mathematica should at least be able to handle the above situation automatically in a way that does not lead to a memory overflow.

Real Code for the actual problem with some modification to simplify the input (which does not change the occurring problem):

(*data initialization*)
imageData = Table[RandomReal[], {i, 200}, {j, 300}]; (*this would usually be some GrayScale ImageData*)
f = Table[RandomReal[], {i, 200}, {j, 300}]; (*this would be some function evaluated on the imageData*)
maxiter = 10000;
testImageData = Table[{RandomReal[], RandomReal[]}, {i, 200}, {j, 300}]; (*would be the gradient of the imageData*)

(*the for loop that leads to the memory overflow*)
For[i = 1, i < maxiter, i++,
 test = First@CUDATestFunction[testImageData, 0.1, f, imageGradientNormalized, 1, Sequence @@ Dimensions[imageData]];
 testImageData = test;

CUDATestFunction = CUDAFunctionLoad["
  __device__ float length(const float2& a)
    return sqrtf(a.x*a.x+a.y*a.y);
  __device__ float2 operator+(const float2& a, const float2& b)
    return make_float2(a.x + b.x, a.y + b.y);
  __device__ float2 operator-(const float2& a, const float2& b)
    return make_float2(a.x - b.x, a.y - b.y);
  __device__ float2 operator*(const float& a, const float2& b)
    return make_float2(a * b.x, a * b.y);
  __device__ float2 operator/(const float2& a, const float& b)
    return make_float2(a.x / b, a.y / b);
  __global__ void resolvFs(float* p, float sigma, float* f, float* imageGradientNormalized, float lambda1, mint width, mint height)
    int xIndex = threadIdx.x + blockIdx.x * blockDim.x;
      int yIndex = threadIdx.y + blockIdx.y * blockDim.y;
      int index = 2*(xIndex + yIndex * width);
    if(xIndex < width && yIndex < height)
        float2 vecP = make_float2(p[index], p[index+1]);
        float2 vecN = make_float2(imageGradientNormalized[index], imageGradientNormalized[index+1]);

        vecP = vecP + 2 * sqrtf(lambda1) * vecN;
        vecP = vecP/fmaxf(1, length(vecP)/(2*sqrtf(f[index/2]+lambda1)));
        vecP = vecP - 2 * sqrtf(lambda1) * vecN;

        p[index] = vecP.x;
        p[index+1] = vecP.y;
  }", "resolvFs", {{"Float", 3, "InputOutput"}, 
   "Float", {"Float", 2, "Input"}, {"Float", 3, "Input"}, 
   "Float", _Integer, _Integer}, {16, 16}]
  • $\begingroup$ It would be interesting to have a little bit more than psedocode. I would think that declaring as "Input" an argument (instead of what you described which seems to me should be "InputOutput") wouldn't let it be returned. $\endgroup$
    – Rojo
    Commented Aug 30, 2013 at 16:25
  • $\begingroup$ I modified the question above and appended some real code, that should run in mathematica with loading CUDALink first. The Code uses different input data than I do, but that should not have any influence on the problem at hand. $\endgroup$
    – Wizard
    Commented Aug 30, 2013 at 20:54
  • $\begingroup$ I would guess the memory manager is buggy. I'm sending Wolfram an email to report something that might be part of your issue here (but probably not all of it). Anyway, you can in the meantime use CUDAMemory objects explicitly for your Output and InputOutput arguments. $\endgroup$
    – Rojo
    Commented Sep 13, 2013 at 2:37
  • $\begingroup$ Thanks for the comment. Too bad the memory manager has some issues there, but for the final version of the algorithm I intended to use explicit memory management anyway. A working memory manager would of course help in testing phase, when the main purpose is to test if the CUDA Code is working. $\endgroup$
    – Wizard
    Commented Sep 16, 2013 at 13:41

2 Answers 2


There's no fixing for that memory leak problem so far (ver 10.4). Here is the workaround I post in another place https://stackoverflow.com/a/35922539/6046035


I believe that CUDALink uses cuBLAS from NVIDIA. cuBLAS has a memory leak when a handle is free'd. When running some code that I wrote which calls cuBLAS directly, I need to do a device reset about every 5 to 10 times that I free a handle to force the memory to be free'd. This problem has been documented on the NVIDIA forums.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.