3D Gaussian CUDA Source

The below code is meant to accompany the 3D Gaussian Convolution post. Note that the code is not generic and only calculates a 5x5x5 Gaussian on a set of 256x256 planes. It is very easy to modify to support other sizes and could even be templated. The Z convolution is not included in the code, but remember in this step you would need to do the final division.

As blogger doesn't seem to like keeping the format of code, please find it on GitHub https://github.com/Repmov/3D-Gaussian

CUDA 9

NVidia announced CUDA 9 a few weeks ago. I've been using CUDA since v1.1 and compute capability 1, and things have matured significantly over the years.

The new CUDA adds support for the new Volta architecture, C++14, faster libraries and Tensor core matrix multiply, which is clearly targeting deep learning applications. But, for me, there is one stand out feature:  Cooperative Groups.

The release says that it is a new programming model for managing groups of communicating threads. What does that really mean?

Previously you could synchronize threads across a thread block with the __syncthreads() function. Cooperative groups allow you to define groups of threads at the sub-block and multi block levels and synchronization across the entire grid.

The grid sync means you now longer have to have multiple kernels operating in successive launches in order to complete a complex task on a data set. A single kernel can now operate on the data and using something like:

thread_group group_grid = this_grid();

//do something here

grid.sync();

//do something else here

//etc

You also get a this_multi_grid() variant which will synchronize the kernel across all GPU's its been launched on!

You no longer have to wait to the sync at the end of a kernel launch and launch another kernel from the cpu code.  Presumably you will still be restricted by the timeout on the driver for your primary device.

This coupled with the pinned memory / zero copy means you can have long running kernels running all sorts of operations on memory that can be streamed into the device.

Exciting stuff!

New blog home

BV2 has a new home on blogger (google) as I just couldn't be bothered to run my own backups, email hosting etc etc anymore and the expense of hosting my own server was getting a bit excessive.

The blog has been renamed to ComputeCube as the idea is to merge the two seperate blogs into one, both links should point here.

The automated wordpress->blogger import tool did a good job but some of the image links have broken. I will be working on fixing these. Just to be clear my move to blogger has nothing to do with wordpress. I still think its a really good blogging system, just at the moment my needs are better served here.

Google now pretty much hosts everything of mine. Email, blog and google+ for the social networking thingies. In my push to move with the times I even have a twitter account now (ComputeCube) which contains personal, blog and GPU related stuff.

On the GPU front: I'm still running 2x 8800GT's and a Tesla C1060 but am considering getting a 560Ti as they seem to offer the best price/performance ratio and allow me to use the new CUDA features.

The server crash did knock the wind out of my blogging sails somewhat but since my last post I've not been idle and have implemented a lot of my ideas in CUDA, watch this space for more details.

Catch up

No apologies for the long delays between posts, or even checking the blog. It just has to fit in with life at the moment and there is so much going on!

Those of you who have left comments on the blog and emails for me should have got an answer last night or this morning. A bit of a delay for some of you.  I do apologise for those of you who sent an email and haven't got a reply back. My computer that used to handle all my email died and I haven't had a chance to fix it yet. I think its just the power supply.... hopefully! so will get the emails back in a few weeks when I eventually get around to fixing it.

I've subsequently upgraded my email server - so any forthcoming mails will get to me. 

In development work I've been working on my SPH simulations, and some GP stuff whenever I get a chance. GP stuff is traditionally recursive - well the equation trees anyway and have needed a substantial amount of reworking to get working efficiently on the GPU.

Speaking of recursive.... in order to be Turing Complete (assuming infinite memory for now) do you need to support / include recursion? Some posters on certain forums seem to think it is needed, but personally I can't see why?  Most recursion with a bit of effort can be iterative - although possibly not very pretty or efficient.

For example. A GPU doesn't really support recursion*, but I would consider cuda / GPU combination as Turing complete. Admittedly not very efficient in certain cases - single thread for example. And again ignoring the infinite memory issue.  *You can if you implement your own stack type system in global memory...

I'd be interested in knowing others views on this - email the normal place or comment here :)

To all the regular readers of the blog - anyone else being amazed by the absolute explosion of GPU / CUDA related code / products / hardware.  Very exciting indeed!

SPH Screenshot

Finally the promised screenshot :)

SPH with symmetry

It's not all that impressive to look at as I've restricted all the particles to 2d although it does use 3d calculations. I do this to help look for any issues in the code as I find it hard to spot errors in a 3d particle rendering.

This particular screenshot has 64000 particles that have been dropped into the box in a column formation and are now starting to slosh around at the bottom.

The unusual thing with regards to a CUDA implementation is that it is using symmetry in the interactions thereby decreasing the memory/processing load. I've still got more work to do but its showing a lot of promise in running superfast particle interaction simulations.

I've aso been doing a bit of work on my second version of my raytracer. I've once again stepped away from KD-trees and Octrees and am using a type of BVH, ray marching system. Screenshots once I have a decent scene rendered :)

In other news I'm now compiling all my new C++/CUDA code in 64bit with the CUDA 3.0 beta. Although I think putting in c++ object support into CUDA was a mistake the new version does produce decent code.

Poor neglected blog...

Nearly 3 months since my last post :(

Work has been exceptionally busy: In the last two months on top of my normal product maintenance and improvement duties I have prepared and filed a patent application, architected and largely completed a distributed, resilient document processing framework and found a bit of time to eat and sleep!

I've noticed other blogs in the raytracing / graphics / visualization space have been very quiet lately - maybe everyone else is also working like crazy?

Not a huge amount has happened in my raytracer and SPH projects although got some interesting effects running with a non-uniform mass particle system when I had time over Christmas. Screenshots soon.

I do have the beta release of Nexus (the NVidia Visual Studio plugin)  but sadly it only runs on Windows Vista or Windows 7 which leads nicely on to my next point:

I am a bit irritated with Microsoft for two reasons:  Even though I purchased a 64 bit Windows XP professional about 6 or 8 months ago there is no upgrade path to Windows 7...  Secondly even though visual studio 2008 standard has a switch for openMP it doesnt contain the openmp headers. Only the more expensive professional version does. Not something that was immediately obvious from the documentation before I purchased...

Although I also run Linux (centos) I prefer to develop on a Windows GUI - less buggy and more responsive than gnome / kde in my opinion. For running code the Linux os does usually win though! I would really like to run Nexus so am a bit stuck about what to do....  Succumb and buy Windows 7 and get Nexus on Visual Studio? or just forget entirely about Windows development / environment and use Linux / gcc / Intel compilers instead?  While the Intel compilers are great (if a bit expensive) for an IDE I really do like Visual Studio.  Most of my code is cross platform and for graphics I mostly use openGL so could switch without too much trouble...    But direct compute is so tempting.....

Arrrrgh what to do!

Fermi

I've just completed reading the white paper released by nvidia which you can find here.

Rather interestingly no mention of graphics performance has been made which, in a way, is really exciting. This has clearly been aimed at the high performance or throughput computing markets with the notable inclusion of ECC memory and increased double precision throughput along with the updated IEEE 754-2008 floating point support.

Concurrent kernel execution and faster context switching will allow, with the use of DAG's, the optimization of execution on the devices rather than just working out the most efficient order of kernels to execute sequentially.

Also tucked away in the white paper is the mention of have predication at the instruction level which should give greater control of divergent paths in your kernels.

The inclusion of C++ support will appeal to a lot of people but am I rather unconvinced this is the correct way to go for throughput computing as it will encourage the use of all the old patterns that may work well in serial cases but are often rather poor for enabling maximum throughput.

There is a lot more in the paper and already an announcement by Oak Ridge that they will be using it in a new supercomputer.

All in all its a wonderful development and I can't help feeling that computing took a substantial leap forward today.

3D Gaussian Convolution

There hasn't been much in the way of posts here lately as I've been really busy at work getting some new components built into the systems I work on. Not really hard but it's frustrating things like trying to get various components and libraries written in different languages to work together. So lately I've not had the energy to do much work on the computer once I get home...

There has been a bit of interest in my 3D Gaussian convolution kernels. Although I explained the technique mathematically in an earlier post I never actually posted the code. As it is rather quick and quite a novel way of calculating the convolution for a xy plane I decided to post it so everyone can benefit from / improve upon the technique.  As always comments / bug reports etc are always welcome :)

Glass half full or half empty?

Or as this is a hpc/cuda/parallel processing site:

Gustafson's Law or Amdahl's law?

Personally I prefer Gustafson's Law .... it seems more logical to me or is this just because I'm inherently an optimist?

I would be quite interested on hearing your views on this - so comments /forum posts most welcome.

In other news:  The thermal monitor downloads have gone over 80! :)  The updated version (v0.2) is ready after mucking around with subclassing a control.... I will release it soon...

Otherwise I have been extremely busy on debugging a sparse matrix solver - bugs in huge datasets can be hard to find! Even with the aid of Gold / Sivler / Bronze kernels mentioned in the last post they have been proving remarkably tricky to isolate. Rather surprising to me is the fact that the long long data type doesn't consume a lot more processing time than a normal unsigned int - so I have been using that wherever there is a risk of exceeding 2^32.

CFD code coming soon too - although I have unwound a lot of the optimizations in order to make it easier to understand and possibly be a good foundation for your own optimizations.

Right .... back to the grindstone!

Gold, Silver and Bronze

Kernels of course! :)

Most of the readers of this blog should be familiar with a "Gold" kernel in which your data is processed on the CPU (usually) and the output is carefully checked. This kernel and its associated outputs form the basis of the regression testing of subsequent implementations on the GPU including algorithmic optimizations.

Personally I like most of my gold kernels to be naive implementations of an algorithm. This causes them to be  easily verifiable and usually easy to debug if there is a problem.

If you currently don't implement a Gold kernel before writing your CUDA implementations and/or adapting you algorithm I strongly suggest you do.

The purpose of this post is to suggest two other debugging techniques I use when needed and where possible. I call them my Silver and Bronze kernels.

A Silver kernel is implemented on the GPU without any optimizations or algorithmic enhancements. The grid / block structure is as simple as possible making sure we don't vary from the Gold kernels implementation too much - only unwinding the loops into grid/blocks is allowed where possible. This type of kernel I use when I am writing something that depends on numerical precision. Once written and verified within acceptable numerical limits against the Gold kernel it becomes the new baseline kernel before later optimizations. This allows exact matching of later kernel outputs rather than using an "acceptable deviation" approach.

CUDA Emulator Output

As many people have noticed the same code executed in Emulator mode gives different floating point results from the kernels run in Debug or Release mode.

Although I know what causes this I have never bothered to investigate the actual differences as most of the stuff I write runs entirely on the GPU. Recently I have had to compare results on the CPU<->GPU and wrote some code to change the FPU settings. Firstly a quick explanation:

By default the CPU (FPU) is set to use 80 bit floating point internally. This means that when you load in an integer (fild) or a single / double float (fld) it gets converted to a 80 bit number inside the FPU stack. All operations are performed internally at 80 bits and when storing the result it converts back to the correct floating point width (single / double)  (fst / fstp). 

This method of operation is desirable as it reduces the effect of rounding / truncating on the intermediate results.  Of course while very useful for computing on the CPU this is not how the CUDA devices operate.

GPU Temperature Monitor

I was hoping to release the GPU temperature monitor to the downloads section sometime during this last bank holiday weekend. I had also planned to sort out my somewhat ailing / overheating computer and perform some upgrades. Unfortunately the repair / upgrades took almost 2 full days. My pc is now running a lot cooler and a bit quieter and I am now something of an expert with heat spreader / heat sink cleaning and thermal grease application.  By the way artic silver is really good! Even before the 200 hour break in period mentioned on their site I'm already seeing more than 10 degrees lower temperature on the CPU.

As to the upgrades:  more posts on this later :) But it did include me purchasing Windows XP Professional 64bit.

The longer than expected pc maintenance time has impacted the GPU Thermal Monitor application and it won't be ready for a few more days. A bit of good news about it though:  it works perfectly on Windows XP 64bit without any modifications. Don't you just love asm :)  Some of my other C/C++ applications were not quite so happy on the new OS.

Keep watching this space for release date :)

GPU Thermal Monitor

(update 27/5/09 - I've now released v0.1 - see this post )

It's been a while since my last post as I've been rather busy with work and some other projects - so time for a quick update. 

While working on some long running kernels I wanted to keep track of the GPU's temperature and ran the monitor app that came on their installation disks. The problem with the monitoring apps that I have is that they are quite large and take up a lot of screen space and for some unknown reason they crash / stop working when viewed over VNC (or similar) remote control app.  Now as I mostly use my CUDA machines over the network this is not a good situation.

So decided to write my own :)  Here is a screenshot of the initial version. It still needs config and some cleaning up (can you spot the glitch by the minimize button?) and possibly some remote reporting over the network functionality. Currently it can display the temperature from two GPU's and updates every 500ms. When minimized it sits in the system tray and updates its tooltip with the reported temperatures. I will be releasing it to the downloads section of the website sometime over the weekend. The current .exe is 146k (mostly the skin) and doesn't require any installer.

[caption id="attachment_823" align="aligncenter" width="252" caption="GPU Thermal Monitoring Application"][/caption]

Nothing wrong with a bit of self promotion sometimes :)

Tesla C1060 memory performance #2

In my post a few days ago I mentioned that some of the numbers reported by the Visual Profiler needed further investigation.

In particular the device to device memory bandwith reported by the profiler differed from the value reported by the bandwith test sample. This was easily tracked down:  according to the documentation, the profiler divides the bytes read/written by 10^9 whereas the bandwith test sample divides by ( time in milliseconds * 2^20).    So: 1000*2^20/10^9 = 1.048576    i.e. the bandwith sample reports a lower number.

D3Q19 Lattice

I managed to get some time last night to convert my LBM implementation to CUDA. Its far from optimal at the moment. Here is a screen shot showing the lid-driven-cavity with two obstacles, one horizontal and one vertical in the box. The visualization is from the middle plane. The obstacles are not shaded but are easily seen in the picture as an area of black. This image was taken after 16000 timeslices and some nice stable vortices have developed.

Lid driven cavity with obstacles

Once again the importance of making a gold kernel cannot be understated as I had quite a number of bugs in my initial CUDA implementation. I use a "multi-tap" type approach when debugging where the kernels write intermediate results to device memory as they go along. This can be easily compared to data coming from the various stages of the gold kernel and makes it much easier to identify the source of the error. Keep in mind the CUDA floats will never be the same as the CPU's floats (CPU uses 80 bits to compute intermediate results).

Tesla C1060 memory performance

According to the CUDA programming guide the memory coalescing rules have been relaxed in devices with compute capability 1.2 or greater.  Chapter 5 has a subsection entitled "Coalescing on Devices with Compute Capability 1.2 and Higher" which gives more information.

As I now have a C1060 which has compute capability of 1.3 I thought I'd run my old coalescing tests on it to see how it has improved.  I modified the launch configuration to 32760 blocks in order to maximize utilization of all the Multiprocessors and increased the thread count to 256. These changes cause the kernels to report 100% utilization in profiler.  I expected the uncoalesced part of the old tests to be quite kind to the Tesla as they fall within the parameters given in the programming guide and indeed the memory transfer rates were similar to a pure device to device copy. I then modified the uncoalesced kernels to avoid the auto coalescing of memory accesses in a half-warp.

D3Q19 Lattice

After playing with a D2Q9 lattice mentioned in a previous post I felt I'd learned enough to progress to the wonderful 3D world :) The arrival of my Tesla has also given me more processing power to move into three dimensions.

So far, in order to get data I can compare the CUDA kernels against, I have come up with a cpu gold kernel. For the first test have constructed a Box with 5 sides and constantly moving flow in the top layer of the box. I am only simulating incompressible at the moment.

I have made a very simple OpenGL viewer that can either show me the entire box or a plane through it. The direction of the velocity vectors is indicated both by their colour and direction of line while the magnitude of the velocity is indicated by the length of the line (thats why some spill over the edges of the box). The below image shows a section through the middle of the box.

[caption id="attachment_763" align="aligncenter" width="292" caption="Visualization from the gold kernel"][/caption]

The next step is to construct the CUDA kernels and compare their outputs. I'm hoping for a massive increase in lattice update speed whilst maintaining numerical stability. The gold kernel uses doubles for reference purposes.

Visual Studio 2008 and CUDA

As promised I received my Visual Studio 2008 on Friday from polyhedron and set about installing it on Saturday. This post will describe how to get VS 2008 to work with CUDA.  I use CUDA 2.2 which is still under NDA so I wont be making any comments about its performance improvements today. Rather I will describe how to set up syntax highlighting, building and Intellisense.

3D Gaussian - in sections

After installing VS2008 it took me a little while to get all my projects to compile nicely again. Quite a big jump from 2003->2008 :) and a lot of my paths were a bit out. I'll post later about getting CUDA 2.2 to work in VS2008 and the problems I had making it cross compile.

Just time this evening for a quick update: I finished my 3D Gaussian convolution on my generated data set - the attached image has clearly visible joins in it, these are to make sure my segmented approach of processing big data sets works correctly. The are easily removed by overlapping the non-boundary segments. The colours are based on my transfer function and don't reflect the underlying data accurately yet. It renders the volume at around 40fps.

3D Gaussian convolution (segmented) output

CUDA Processing of Large Data Sets

As those who follow this blog will know, I have been playing with large generated volumetric data sets (around 16GB). The first operation I performed on them was a 3D Gaussian convolution in order to smooth a lot of the randomness out. My first implementation of the kernel was just the naive one - take a point and look at all the surrounding points. As I'm using a 5x5x5 kernel this means 125 memory reads per point - very poor!