C++ Library for Linear Algebra on Supercomputers (libelemental.org)
Elemental is an open-source C++ library meant for performing dense linear algebra on tens of thousands of processes (via the Message Passing Interface). I have been the primary developer for the past four years, but a community is starting to emerge. In addition to the website, the project is also hosted on github: http://github.com/elemental/Elemental
38 comments
[ 4.1 ms ] story [ 86.7 ms ] threadI'm also happy to answer general questions about parallel linear algebra and to point people to the appropriate literature.
With that said, the main advantage of the library is its high-level of software-engineering, which tends to encourage the rapid development of new features. It is actively used within a large number of research projects.
From my perspective, Elemental is primarily about human performance. That it scales to thousands of nodes of Blue Gene and achieves performance on par with (and in some cases significantly better than) ScaLAPACK is really icing on the cake. Granted, I really like scalability and performance, but I will not ruin my code in the name of these things.
That's a shame, because having collections of useful tools like that makes it a lot easier to keep track of options.
I imagine the authors could make a single node of the cluster use a GPU, if they wanted.
Distributed memory and GPUs are not mutually exclusive. Multi-GPU clusters are extremely common. In fact the latest devices (e.g. Tesla K10) have multiple GPU processors packaged in a single card, so it is necessary for applications to target multiple GPUs. There is explicit support for distributed-memory applications in GPUs through the "GPUDirect" technology that allows peer-to-peer DMA and RDMA transfers between GPUs.
Given that reports of 30-50x GPU performance gains (versus CPUs) are common, the issue is important because it means solving a problem with (say) $10,000 of kit instead of $500,000.
For examples, look at the MKL benchmarks for Xeon Phi (http://software.intel.com/en-us/intel-mkl/) and normalize as 3 Sandy Bridge sockets per Xeon Phi (common configurations pair one SB socket with one Xeon Phi, the Phi has more than twice the TDP and costs more than twice as much). Don't forget to look at QR and Cholesky, for which the Phi at best breaks even, but only for enormous matrix sizes.
Let me explain: One of the most common issues with x86 HPC applications coming from the scientific crowd is a lack of vector optimization such as loop unrolling. Even having the right compiler flags is rather difficult for this kind of thing. Another reason is a lack of understanding on how to program for memory bandwidth optimization. GPU programming on the other hand, especially with CUDA, is hard to get into at first, but once you have the right formula you can apply it pretty easily to most common tasks. Getting to, say, 70% of model performance on GPU is much easier than on x86. One reason is the implicitly bandwidth optimized idiomatic way of writing CUDA/OpenCL programs as a set of scalar kernels applied over a whole data region - this allows the programmer to think of block dimensions in an abstract way - no need to fiddle around manually with loops to achieve this. There is also no need to use any intrinsics, just plain C in idiomatic CUDA is enough.
So, to wrap it up, there is more to GPU programming than just the hardware itself, the software model actually makes a lot more sense than traditional OOP/procedural programming for HPC - often resulting in higher than expected speedups when going from idiomatic x86 to idiomatic GPGPU (since there is no such thing as easy to learn idioms for HPC x86 programming).
And btw. Xeon Phi is the result of Intel not understanding exactly this interrelation, since it doesn't even have OpenCL support as of now.
My guess is that you've spent a lot more time trying to understand CUDA and NVIDIA hardware. Or maybe all of your arrays are dimensions that are multiples of 16.
While it's true that there are kernels that are easier to optimize for a GPU, I still think most claims of huge speedup are due to neglecting the CPU implementation. See "Debunking the 100x GPU vs CPU Myth" and other papers on this. And amidst all of that, there are a lot of applications that utterly fail on GPUs, despite having lots of parallelism, just not at quite the right granularity.
Also, not all supercomputers have accelerators (consider Blue Gene/Q), and often simply having access to more memory is more of a concern than solving the problem at the absolute fastest rate.
The reoccurring bifurcation of talent and resources in the open source community is really disheartening. Can't we focus on one or two libraries and make them actually good? Or at least fork off of something that already exists and add your own features. I look at benchmarks of the existing tools and one library will do one operation very efficiently, while another will work well with something. Often the differences in speed are huge (more than a factor of 10). So I end up having to flip a coin in choosing which library to use.
the C++ linear algebra libraries are basically syntactic sugar for interacting with optimized Fortran libraries like BLAS, LAPACK, etc. (I'm sure I'm overlooking some complexity here, because the C++ libraries, while linking to the same Fortran libraries have very different run times)
Yours interacts with ScaLAPACK and other distributed memory libraries.
Why would it not be possible to simply extend say armadillo or eigen to interact with distributed memory libraries? If you need more syntax, then extend the interface.
The other major difference is that sequential libraries tend to get away with letting users not have to worry about where data resides. This is of fundamental importance in distributed libraries, and, for this reason, it is usually a bad idea to think of simply modifying sequential APIs.
Looking at ScaLAPACK, it's been developed since 1995. I've never touched it, but it's probably many many lines of code (and maybe a few PhD thesis) with all sorts of kinks worked out that will take you decades to iron out yourself. To throw out all that knowledge/work/man-hours and to start from scratch seems like a waste.
1. What granularity to distribute the entries in the matrix is a long and subtle argument which doesn't provide a clear winner for every operation (the current conclusions are different for LU with partial pivoting vs. reduction to tridiagonal form). I would by no means say that the approach used by ScaLAPACK is wrong, but only that it is unnecessarily complex and only one operation purposefully targets the finest granularity case.
2. Again, I appreciate the complement, but I don't think that arguments from authority are valid, nor do I think that one can be the "best". I have a large number of colleagues doing wonderful work, much of which I find extremely impressive.
Also, it would be good to disclose that you're affiliated with the project you were promoting.
With that said, if you were familiar with the source code and APIs of both libraries, I think that you would see a clear benefit. This is supported by the way that Elemental is growing and arguably has more functionality than ScaLAPACK (with the notable exception being a parallel Schur decomposition). Over the past couple of years, the library has primarily been developed to support my research goals, but a large number of research groups are actively using it now.
I think you have unintentionally violated the last provision. Please don't post about things you have not taken the time to investigate at all, which is true by your own admission: "I've never touched [ScaLAPACK]".