CUDA function calling, what's the best approach?

Compile it and look at the assembled function. My guess is that the second method would be slower, but only the good Lord knows what the compiler will do with it.

What?! Kernel2 still needs to upload the entire array into device memory to be able to access value[0], value[1], .. value[N]. It doesn't fetch them from host memory, that would be horribly slow. All data needed by the kernel is always cached into device memory. The transfer may be done before the kernel call in the second case, though (upon device memory allocation and data upload) but it would still be there, however Kernel2 is better in the sense that it can be reused more easily, and I believe there is a limit on the number of kernel arguments you can have. Also, I believe the arguments in Kernel1 are in a different memory address space than in Kernel2.

Although this is based on my knowledge of OpenCL. CUDA should be largely the same at this level, though.

What?! Kernel2 still needs to upload the entire array into device memory to be able to access value[0], value[1], .. value[N]. It doesn't fetch them from host memory, that would be horribly slow. All data needed by the kernel is always cached into device memory. The transfer may be done before the kernel call in the second case, though (upon device memory allocation and data upload) but it would still be there, however Kernel2 is better in the sense that it can be reused more easily, and I believe there is a limit on the number of kernel arguments you can have. Also, I believe the arguments in Kernel1 are in a different memory address space than in Kernel2.

Although this is based on my knowledge of OpenCL. CUDA should be largely the same at this level, though.

Those values aren't actually global memory, they are private. I'm very sorry for the confusion and will try to provide some context:

I have a 3D mesh with each position representing a float value and I use GPGPU to update this mesh T times.

//This is the function called from host to device.
__global__ void Update3DMesh(float *mesh, int time)
{
    //Retrieve the index from the thread ID that represents the position being updated.

    //Retrieve the center, top, bottom, left, right, front, back position values relative to the index.

    //Compute N factors using these surrounding position values according to time.
    float factor1 = ComputeFactor1(parameters);
    float factor2 = ComputeFactor2(parameters);
    ...
    float factorN = ComputeFactorN(parameters);

   //Update the position.
   float[index] = factor1+factor2 ... +factorN;
}

I'm struggling to figure the best way for ComputeFactor to receive the parameters. I thought that maybe packing all these surrounding positions into a (float *) vector and passing the pointer (sizeof(float *) bytes) instead of passing all positions separatelly (sizeof(float)*7 bytes) would result in a faster program. I'll just have to compare both methods.

Appreciate the responses, Khatharr and Bacterius.

On Compute Capability < 2.0 devices, all function parameters are stored in shared memory and immediately loaded into registers. That means, that each parameters has a very high cost, because it limits the number of threads that can run in parallel. I had a lot of cases, where it was beneficial to store parameters in constant memory.

For Compute Capability >= 2.0, the spec says that parameters are stored in constant memory, but I haven't done any testing, if they still block an entire register for each parameter.

I would expect the single pointer to be either equally fast or faster.

PS: If you have every thread read a different value from the array and cache them in shared memory, you can significantly reduce the number of global memory accesses for the second (pointer based) variant.

I had a lot of cases, where it was beneficial to store parameters in constant memory.

That's what I'm doing, prior to start executing the program I load static settings (i.e mesh size, mesh length...) as constants using cudaMemcpySymbol and just send the mesh and time parameters with the kernel. I'm trying to minimize transfers from host to device as much as possible.

PS: If you have every thread read a different value from the array and cache them in shared memory, you can significantly reduce the number of global memory accesses for the second (pointer based) variant.

I didn't know about shared memory (I'm basically a starter in CUDA) and, by doing some research over the internet, I learned the maximum cache size for each thread block is around 16-32KB. My kernels run 1024 threads per block, needing 60KB, therefore, to make it work with my kernel I'd have to reduce the number of threads per block from 1024 to 512, is this generally a good tradeoff?

My kernels run 1024 threads per block, needing 60KB, therefore, to make it work with my kernel I'd have to reduce the number of threads per block from 1024 to 512, is this generally a good tradeoff?

I don't think there is an alternative to run-time profiling of the kernels. Which for me means exhausting all possible block sizes and comparing times.

__device__ float SimpleKernel1(float value1, float value2, float value3 ...float valueN)
{
return value1 + value2 + value3 ... + valueN;
}

I don't know what happens when pass by value, I'm curious though. I'll look at ptx when I get to work. Questions such as these should be answered on the clock.

The pointer version will simply access the thread's local memory, which I'm guessing is better/best. It never dawned on me to pass, say, 6-7 values like that instead of using an array.

I don't see right off how shared memory would provide a benefit unless you are passing data between threads within the block. I'll admit to a lot of voodoo optimization of CUDA code -- buffering in shared memory never worked for me. At best it will mimic the use of private local memory. The use of local memory and registers will force you to vary the block size also.

Oh, I see, I misread your original post. You are not talking about the kernel launch, but about a device function.

The actual block size is not as relevant because a MultiProcessor can run several blocks in parallel. The maximal amount of shared memory is actually 16KB for ComputeCapability < 2.0 and 48KB otherwise. Note, that cache and shared memory are two distinct things, although they are implemented by the same hardware.

Which hardware level are you actually targeting? How many parameters do you have, and how many "factors" are computed?

The way I interpret your code, the "parameters" are sampled from a 3D grid. Each of these grid cells are sampled 6 times by 6 different threads so without storing them in shared memory, you have 6 times as many global loads, as you would need in theory. In addition, if I interpret your snippets correctly, you have not one parameter per grid cell but many and you read them one by one. That is, up to 6 threads may each read the same value, and you do this up to N times, where N is the number of parameters per grid cell. So you end up with 6 * N global memory loads per grid cell.

If you use shared memory, you can share those grid cells to some degree among the threads. Also, all the threads in a block can join forces to read the cells parameters sequentially from memory, thus further reducing the number of global loads.

The actual block size is not as relevant because a MultiProcessor can run several blocks in parallel. The maximal amount of shared memory is actually 16KB for ComputeCapability < 2.0 and 48KB otherwise.

Good, didn't know that. If many blocks are allowed to run in parallel, does it mean each block will have it's unique shared memory?

Which hardware level are you actually targeting? How many parameters do you have, and how many "factors" are computed?

The program output says it a Tesla 2090, it's compute capability is 2.0. The parameters on the kernel are float *mesh3D and int time and the "factor functions" are two actually; Laplace and Chemotaxis; their parameters are the surrouding values (top, bottom, left, right, front, back and middle).

Your guess is correct with the exception that there are 7 samplings. Also, each point which belongs to the grid has 7 cells (the grid is actually a microscopic portion of human skin tissue), so it's 7*7*N redundant acesses.

Laplace and Chemotaxis are calculated on each cell on a mesh position, so there are lots of acesses to global values. I want to move these values to a faster memory before Laplace and Chemotaxis start acessing them.

This is pretty much what cuda's so-called texture memory is for, you should look into it. If I recall correctly, the introduction to texture memory in some of the NVIDIA docs is a similar problem. The short of it, you'll get these neighbors into cache for fast access. Requires a bit of tuning though...

CUDA by Example

Chapter 7. Texture Memory - Simulating Heat Transfer

This is pretty much what cuda's so-called texture memory is for, you should look into it. If I recall correctly, the introduction to texture memory in some of the NVIDIA docs is a similar problem. The short of it, you'll get these neighbors into cache for fast access. Requires a bit of tuning though...

Thanks for the suggestion but, correct me if I'm wrong, texture memory is a read-only (my mesh has to change it's value at each iteration) and special type of global memory. So shared memory acess is still going to be faster than texture memory.

Also, from looking at the simpleCubemapTexture sample, texture memory would seem to suit better for surfaces, not volumes as in my case and there's also memory size limitations I'd have to consider for a single texture, the meshes in my code can easily reach 100MB.

texture memory is a read-only (my mesh has to change it's value at each iteration) and special type of global memory

It's more a caching-scheme than it is a special type of memory. I'm pretty sure your Tesla is capable of writing to it, though I've not done that in any of my solutions.

I wouldn't get too caught up in the real world objects you are modeling in considering this, it is really a matter of controlling your memory access pattern so that caching is beneficial, ie threads in the same block access data that are 'close' to each other in global memory. You would transpose your input data to accommodate.

there's also memory size limitations I'd have to consider for a single texture, the meshes in my code can easily reach 100MB

If you transpose your data to a 1d array then a texture can cover ~130MB. If that isn't enough you can use more than 1, but (at least with the last version of CUDA I used) this led to some awkward code and branching because there was no way to represent an array of textures.

Textures are actually really good in cuda, so I agree with jms bc, you should use them if possible. They are a lot easier to handle then shared memory and you get caching, interpolation and handling of out of bounds sampling for free.

Textures aren't really read only. You can not write to a texture using a global store or surface store, and expect to read that data using a texture fetch within the same kernel call. If you read from texture1 and write to the memory region of texture2 and in the next iteration/kernel call read from texture2 and write to texture1, then you are fine. Also consider using half floats for storing the data in textures. They can cut your bandwidth needs in half, if the precision is sufficient.

About the blocks and occupancy thing:

As a rule of thumb, you should watch out for two things in cuda: Occupancy and global loads/stores. For the global loads/stores, textures are a good optimization you should try first. If that doesn't work out, you can still try manual caching in shared memory.

Occupancy is the amount of threads, that can be held in the MultiProcessor in relation to it's maximum. A CC 2.x MP can hold up to 1536 threads. These are not actually all running in parallel. Rather think of it the way, SMT/Hyper threading works on Intel CPUs: The thread contexts are held in the chip and the chip tries to keep all of it's execution units busy by fetching instructions from those threads. For example, if a couple of threads are performing a global load/store, which will take multiple cycles, all the other threads can use the ALUs to perform computations.

This means that using only 60% of those 1536 threads does not necessarily mean a decrease in performance. But it does decrease the flexibility, that the warp schedulers have and might make global loads/stores and texture reads more expensive.

The amount of threads, that you can put on a MP, is not only limited by that magical number 1536, but also on the amount of registers and shared memory, that the kernel needs.

A CC2.0 MP has a total of 32000 registers. So if your kernel (and thus each thread) needs 32 registers, you can not have more than 1000 threads running on the MP, an occupancy of 65%.

Also, a CC2.0 MP has a total of 48KB of shared memory, that can be split among the blocks. If each block needs 12KB of shared memory (including whatever overhead CUDA produces), then you can have a maximum of 4 blocks present in your MP. If each block only consisty of 16 threads, then you have the horrible situation of only 64 threads being active (Occupancy 1%), which isn't even enough to hide the latency of the arithmetic instructions.

All of these numbers are nicely structured in tables on wikipedia:

http://en.wikipedia.org/wiki/CUDA

IIRC, NVidia offers a spread sheet for computing the occupancy (I never used it) and the NVidia Visual Profiler will also tell you everything you need to know.

Thanks for the tips Ohforf sake and jms bc. By doing some research I've learned more about the textures and there's a CUDA binding called surface which allows write acess to texture memory, it's seems like a good solution to my problem.

My only issue right know is that I've already started implementing the shared cache version of my code, so changing yet again to another approach would cost me a lot of time (I'd have to implement a work around those 130MB limit for a single 1D texture), so, for now, I'll test shared memory and check the speed-up, if it's not good enough I'll switch to texture memory.

Also, didn't know CUDA supported half float variables (thought that was only an OpenCL feature), so I'll switch to half since it could effectivelly double the reach of a single CUDA block to the mesh values.

Thanks again, learned a lot.