Merge pull request #262 from CNugteren/CLBlast-237-tuning-api

CLBlast #237: Tuning API

34 files changed, 1972 insertions, 1179 deletions

diff --git a/CHANGELOG b/CHANGELOG
index c23bf4c0..7c22009c 100644
--- a/CHANGELOG
+++ b/CHANGELOG
@@ -2,6 +2,7 @@
 Development (next version)
 - Added Python interface to CLBlast 'PyCLBlast'
 - Added CLBlast to Ubuntu PPA and macOS Homebrew package managers
+- Added an API to run the tuners programmatically without any I/O
 - Updated and reorganised the CLBlast documentation
 - Various minor fixes and enhancements
 - Added non-BLAS level-1 routines:
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 2e685d76..eb04287e 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -210,7 +210,7 @@ set(PRECISIONS 32 64 3232 6464 16)
 
 # Sample programs
 if(OPENCL)
-  set(SAMPLE_PROGRAMS_CPP sgemm sgemm_batched dtrsm)
+  set(SAMPLE_PROGRAMS_CPP sgemm sgemm_batched dtrsm tuning_api)
   set(SAMPLE_PROGRAMS_C sasum dgemv sgemm haxpy cache)
   if(NETLIB)
     set(SAMPLE_PROGRAMS_C ${SAMPLE_PROGRAMS_C} sgemm_netlib)
@@ -235,6 +235,8 @@ set(SOURCES
   src/kernel_preprocessor.cpp
   src/routine.cpp
   src/routines/levelx/xinvert.cpp  # only source, don't include it as a test
+  src/tuning/configurations.cpp
+  src/tuning/tuning_api.cpp
 )
 set(HEADERS  # such that they can be discovered by IDEs such as CLion and Visual Studio
   include/clblast_half.h
@@ -258,6 +260,9 @@ set(HEADERS  # such that they can be discovered by IDEs such as CLion and Visual
   src/kernel_preprocessor.hpp
   src/cxpp11_common.hpp
   src/routine.hpp
+  src/tuning/configurations.hpp
+  src/tuning/tuning.hpp
+  src/tuning/routines/routine_tuner.hpp
 )
 if(OPENCL)
   set(SOURCES ${SOURCES} src/clblast.cpp src/clblast_c.cpp)
@@ -295,6 +300,9 @@ foreach(DATABASE ${DATABASES})
   set(HEADERS ${HEADERS} src/database/kernels/${DATABASE}/${DATABASE}_3232.hpp)
   set(HEADERS ${HEADERS} src/database/kernels/${DATABASE}/${DATABASE}_6464.hpp)
 endforeach()
+foreach(KERNEL ${KERNELS})
+  set(HEADERS ${HEADERS} src/tuning/kernels/${KERNEL}.hpp)
+endforeach()
 
 # Creates and links the library
 if(BUILD_SHARED_LIBS)
diff --git a/ROADMAP.md b/ROADMAP.md
index bcae8bc4..28a4afa2 100644
--- a/ROADMAP.md
+++ b/ROADMAP.md
@@ -17,6 +17,6 @@ This file gives an overview of the main features planned for addition to CLBlast
 | [#233](https://github.com/CNugteren/CLBlast/issues/233)        | Feb '18     | CNugteren | ✔      | Add CLBlast to common package managers |
 | [#223](https://github.com/CNugteren/CLBlast/issues/223)        | Feb '18     | CNugteren | ✔      | Python OpenCL interface |
 | [#228](https://github.com/CNugteren/CLBlast/issues/228)        | Mar '18     | CNugteren |        | Improving performance for Qualcomm Adreno GPUs |
-| [#237](https://github.com/CNugteren/CLBlast/issues/237)        | Mar '18     | CNugteren |        | Making tuning possible from the CLBlast API |
+| [#237](https://github.com/CNugteren/CLBlast/issues/237)        | Mar '18     | CNugteren | ✔      | Making tuning possible from the CLBlast API |
 | [#136](https://github.com/CNugteren/CLBlast/issues/136)        | Apr '18     | CNugteren |        | Implement xAXPBY and xSET |
 | [#169](https://github.com/CNugteren/CLBlast/issues/169)        | ??          | dividiti  |        | Problem-specific tuning parameter selection |
diff --git a/doc/api.md b/doc/api.md
index 0fbdeaa0..a60e16ce 100644
--- a/doc/api.md
+++ b/doc/api.md
@@ -3497,3 +3497,77 @@ Arguments to OverrideParameters (C++ version):
 * `const std::string &kernel_name`: The target kernel name. This has to be one of the existing CLBlast kernels (Xaxpy, Xdot, Xgemv, XgemvFast, XgemvFastRot, Xgemv, Xger, Copy, Pad, Transpose, Padtranspose, Xgemm, or XgemmDirect). If this argument is incorrect, this function will return with the `clblast::kInvalidOverrideKernel` status-code.
 * `const Precision precision`: The CLBlast precision enum to set the new parameters for.
 * `const std::unordered_map<std::string,size_t> &parameters`: An unordered map of strings to integers. This has to contain all the tuning parameters for a specific kernel as reported by the included tuners (e.g. `{ {"COPY_DIMX",8}, {"COPY_DIMY",32}, {"COPY_VW",4}, {"COPY_WPT",8} }` for the `Copy` kernel). If this argument is incorrect, this function will return with the `clblast::kMissingOverrideParameter` status-code.
+
+
+
+Tune<kernel_name>: Run the tuner for a particular kernel (advanced usage)
+-------------
+
+The CLBlast kernels can be tuned using the tuning binaries, but also programmatically through an API. This is only recommended for advanced usage, see for more information [the tuning docs](tuning.md).
+
+C++ API:
+```
+// Tunes the "Xaxpy" kernel, used for many level-1 routines such as XAXPY, XCOPY, and XSWAP
+template <typename T>
+StatusCode PUBLIC_API TuneXaxpy(cl_command_queue* queue, const size_t n,
+                                const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xdot" kernel, used for level-1 reduction routines such as XDOT, XMAX, and XSUM
+template <typename T>
+StatusCode PUBLIC_API TuneXdot(cl_command_queue* queue, const size_t n,
+                               const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemv" kernel, used for matrix-vector level-2 routines such as XGEMV, XGBMV, and XHEMV
+template <typename T>
+StatusCode PUBLIC_API TuneXgemv(cl_command_queue* queue, const size_t m, const size_t n,
+                                const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xger" kernel, used for matrix update level-2 routines such as XGER, XHER, and XSYR2
+template <typename T>
+StatusCode PUBLIC_API TuneXger(cl_command_queue* queue, const size_t m, const size_t n,
+                               const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemm" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TuneXgemm(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                               const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "XgemmDiret" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TuneXgemmDirect(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                                      const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Copy" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TuneCopy(cl_command_queue* queue, const size_t m, const size_t n,
+                               const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Pad" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TunePad(cl_command_queue* queue, const size_t m, const size_t n,
+                              const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Transpose" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TuneTranspose(cl_command_queue* queue, const size_t m, const size_t n,
+                                    const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Padtranspose" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode PUBLIC_API TunePadtranspose(cl_command_queue* queue, const size_t m, const size_t n,
+                                       const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemm" kernel, used for the level-3 routine XTRSM
+template <typename T>
+StatusCode PUBLIC_API TuneInvert(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                                 const double fraction, std::unordered_map<std::string,size_t> &parameters);
+```
+
+Arguments to Tune<kernel_name> (C++ version):
+
+* `cl_command_queue* queue`: Pointer to an OpenCL command queue associated with a context and device to tune the kernel for.
+* `const size_t m`: The routine argument `m` to tune for (not applicable for all kernels)
+* `const size_t n`: The routine argument `n` to tune for
+* `const size_t k`: The routine argument `k` to tune for (not applicable for all kernels)
+* `const double fraction`: A value between 0.0 and 1.0 which determines the fraction of the tuning search space to explore.
+* `std::unordered_map<std::string,size_t> &parameters`: An unordered map of strings to integers. This will return the best found tuning parameters.
diff --git a/doc/tuning.md b/doc/tuning.md
index 88c4fc4c..ebf3cb0c 100644
--- a/doc/tuning.md
+++ b/doc/tuning.md
@@ -100,6 +100,14 @@ In summary, tuning the entire library for your device can be done as follows (st
 After the kernels are tuned, you can run the `clblast_tuner_routine_xgemm` tuner to optimize the high-level GEMM routine, i.e. selecting which method to use: the direct kernel or the in-direct kernel.
 
 
+Tuning using the API (advanced users only)
+-------------
+
+Apart from running the tuning binaries, it is also possible to run the tuners programmatically through the CLBlast API. This could be useful if you want to tune for non-standard arguments (e.g. a rectangular or very small matrix). The tuning results can then also be set programmatically using `OverrideParameters`.
+
+The tuning API does not perform any disk or stdout I/O, thus it is not possible to track progress. Running the regular tuner binaries should give an idea of the amount of configurations to explore for a particular device, thus giving an indication of a good value for the `fraction` argument (see the [API documentation](api.md) for more details).
+
+
 Inspecting and changing tuning parameters at run-time
 -------------
 
@@ -120,3 +128,19 @@ Tuning OpenCL compiler options
 -------------
 
 For all of CLBlast's APIs, it is possible to optionally set an OS environmental variable `CLBLAST_BUILD_OPTIONS` to pass specific build options to the OpenCL compiler. Also make sure this is set in the same way when running the tuners.
+
+
+Which kernels are used for which routines?
+-------------
+
+To find out which tuners to run for which routines, you can use the table below. The kernel names correspond to the tuner binaries, the tuner API, and to the arguments for `OverrideParameters` and `RetrieveParameters`.
+
+| Routines                                                                 | Kernel(s) / Tuner(s)            |
+| -------------------------------------------------------------------------|---------------------------------|
+| AXPY COPY SCAL SWAP OMATCOPY AXPYBATCHED                                 | Xaxpy                           |
+| AMAX ASUM DOT DOTC DOTU NRM2 SUM MAX MIN AMIN                            | Xdot                            |
+| GBMV GEMV HBMV HEMV HPMV SBMV SPMV SYMV TMBV TPMV TRMV TRSV              | Xgemv                           |
+| GER GERC GERU HER HER2 HPR HPR2 SPR SPR2 SYR SYR2                        | Xger                            |
+| GEMM HEMM HER2K HERK SYMM SYR2K SYRK TRMM GEMMBATCHED GEMMSTRIDEDBATCHED | Xgemm XgemmDirect Copy Pad Transpose Padtranspose |
+| TRSM                                                                     | Xgemm XgemmDirect Copy Pad Transpose Padtranspose Invert |
+| IM2COL                                                                   | Copy                            |
diff --git a/include/clblast.h b/include/clblast.h
index 9d3b9ea0..ce64b37a 100644
--- a/include/clblast.h
+++ b/include/clblast.h
@@ -705,6 +705,63 @@ StatusCode PUBLIC_API OverrideParameters(const cl_device_id device, const std::s
 
 // =================================================================================================
 
+// Tunes the "Xaxpy" kernel, used for many level-1 routines such as XAXPY, XCOPY, and XSWAP
+template <typename T>
+StatusCode TuneXaxpy(cl_command_queue* queue, const size_t n,
+                     const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xdot" kernel, used for level-1 reduction routines such as XDOT, XMAX, and XSUM
+template <typename T>
+StatusCode TuneXdot(cl_command_queue* queue, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemv" kernel, used for matrix-vector level-2 routines such as XGEMV, XGBMV, and XHEMV
+template <typename T>
+StatusCode TuneXgemv(cl_command_queue* queue, const size_t m, const size_t n,
+                     const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xger" kernel, used for matrix update level-2 routines such as XGER, XHER, and XSYR2
+template <typename T>
+StatusCode TuneXger(cl_command_queue* queue, const size_t m, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemm" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TuneXgemm(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "XgemmDiret" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TuneXgemmDirect(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                           const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Copy" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TuneCopy(cl_command_queue* queue, const size_t m, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Pad" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TunePad(cl_command_queue* queue, const size_t m, const size_t n,
+                   const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Transpose" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TuneTranspose(cl_command_queue* queue, const size_t m, const size_t n,
+                         const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Padtranspose" kernel, used for most level-3 routines such as XGEMM, XSYMM, and XHER2K
+template <typename T>
+StatusCode TunePadtranspose(cl_command_queue* queue, const size_t m, const size_t n,
+                            const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// Tunes the "Xgemm" kernel, used for the level-3 routine XTRSM
+template <typename T>
+StatusCode TuneInvert(cl_command_queue* queue, const size_t m, const size_t n, const size_t k,
+                      const double fraction, std::unordered_map<std::string,size_t> &parameters);
+
+// =================================================================================================
+
 } // namespace clblast
 
 // CLBLAST_CLBLAST_H_
diff --git a/samples/tuning_api.cpp b/samples/tuning_api.cpp
new file mode 100644
index 00000000..f92b6909
--- /dev/null
+++ b/samples/tuning_api.cpp
@@ -0,0 +1,77 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file demonstrates the use of the runtime tuning API. It is a stand-alone example, but it
+// does require the Khronos C++ OpenCL API header file (downloaded by CMake).
+//
+// =================================================================================================
+
+#include <cstdio>
+#include <chrono>
+#include <vector>
+
+#define CL_USE_DEPRECATED_OPENCL_1_1_APIS // to disable deprecation warnings
+#define CL_USE_DEPRECATED_OPENCL_1_2_APIS // to disable deprecation warnings
+
+// Includes the C++ OpenCL API. If not yet available, it can be found here:
+// https://www.khronos.org/registry/cl/api/1.1/cl.hpp
+#include "cl.hpp"
+
+// Includes the CLBlast library
+#include <clblast.h>
+
+// =================================================================================================
+
+int main() {
+
+  // OpenCL platform/device settings
+  const auto platform_id = 0;
+  const auto device_id = 0;
+
+  // Example arguments
+  const size_t m = 128;
+  const size_t n = 64;
+  const auto fraction = 1.0; // between 0.0 and 1.0
+
+  // Initializes the OpenCL platform
+  auto platforms = std::vector<cl::Platform>();
+  cl::Platform::get(&platforms);
+  if (platforms.size() == 0 || platform_id >= platforms.size()) { return 1; }
+  auto platform = platforms[platform_id];
+
+  // Initializes the OpenCL device
+  auto devices = std::vector<cl::Device>();
+  platform.getDevices(CL_DEVICE_TYPE_ALL, &devices);
+  if (devices.size() == 0 || device_id >= devices.size()) { return 1; }
+  auto device = devices[device_id];
+
+  // Creates the OpenCL context, queue, and an event
+  auto device_as_vector = std::vector<cl::Device>{device};
+  auto context = cl::Context(device_as_vector);
+  auto queue = cl::CommandQueue(context, device);
+
+  // Performs the tuning
+  printf("Starting the tuning...\n");
+  std::unordered_map<std::string,size_t> parameters;
+  auto queue_plain = queue();
+  auto status = clblast::TuneCopy<float>(&queue_plain, m, n, fraction, parameters);
+
+  // Tuning completed. See "clblast.h" for status codes (0 -> success).
+  printf("Completed TuneCopy with status %d (0 == OK), found parameters:\n", static_cast<int>(status));
+  for (const auto parameter: parameters) {
+    printf(">  %s = %zu\n", parameter.first.c_str(), parameter.second);
+  }
+
+  // Set the new parameters
+  status = clblast::OverrideParameters(device(), "Copy", clblast::Precision::kSingle, parameters);
+  printf("Completed OverrideParameters with status %d (0 == OK)\n", static_cast<int>(status));
+  return 0;
+}
+
+// =================================================================================================
diff --git a/scripts/generator/generator.py b/scripts/generator/generator.py
index 5ead9b73..c4db2e0c 100755
--- a/scripts/generator/generator.py
+++ b/scripts/generator/generator.py
@@ -50,9 +50,9 @@ FILES = [
     "/src/pyclblast/src/pyclblast.pyx"
 ]
 HEADER_LINES = [123, 21, 127, 24, 29, 41, 29, 65, 32, 95, 21, 288]
-FOOTER_LINES = [41, 56, 112, 275, 6, 6, 6, 9, 2, 41, 55, 1]
+FOOTER_LINES = [98, 56, 112, 275, 6, 6, 6, 9, 2, 41, 55, 1]
 HEADER_LINES_DOC = 0
-FOOTER_LINES_DOC = 158
+FOOTER_LINES_DOC = 232
 
 # Different possibilities for requirements
 ald_m = "The value of `a_ld` must be at least `m`."
diff --git a/src/tuning/kernels/copy_fast.cpp b/src/tuning/kernels/copy_fast.cpp
index d046c9e9..0314113c 100644
--- a/src/tuning/kernels/copy_fast.cpp
+++ b/src/tuning/kernels/copy_fast.cpp
@@ -11,86 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "copy";
-  settings.kernel_name = "CopyMatrixFast";
-  settings.sources =
-#include "../src/kernels/level3/level3.opencl"
-#include "../src/kernels/level3/copy_fast.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.n;
-  settings.size_b = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3};
-  settings.outputs = {3};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"COPY_DIMX", "COPY_DIMY"}};
-  settings.div_global = {{"COPY_VW", "COPY_WPT"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"COPY_DIMX", {8, 16, 32}},
-    {"COPY_DIMY", {8, 16, 32}},
-    {"COPY_WPT", {1, 2, 4, 8}},
-    {"COPY_VW", {1, 2, 4, 8}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(2, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(3, GetRealArg(args.alpha));
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/copy_fast.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -101,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::CopyGetTunerDefaults, clblast::CopyGetTunerSettings<half>, clblast::CopyTestValidArguments<half>, clblast::CopySetConstraints, clblast::CopySetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::CopyGetTunerDefaults, clblast::CopyGetTunerSettings<float>, clblast::CopyTestValidArguments<float>, clblast::CopySetConstraints, clblast::CopySetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::CopyGetTunerDefaults, clblast::CopyGetTunerSettings<double>, clblast::CopyTestValidArguments<double>, clblast::CopySetConstraints, clblast::CopySetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::CopyGetTunerDefaults, clblast::CopyGetTunerSettings<float2>, clblast::CopyTestValidArguments<float2>, clblast::CopySetConstraints, clblast::CopySetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::CopyGetTunerDefaults, clblast::CopyGetTunerSettings<double2>, clblast::CopyTestValidArguments<double2>, clblast::CopySetConstraints, clblast::CopySetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/copy_fast.hpp b/src/tuning/kernels/copy_fast.hpp
new file mode 100644
index 00000000..f9a58bc7
--- /dev/null
+++ b/src/tuning/kernels/copy_fast.hpp
@@ -0,0 +1,93 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the copy OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults CopyGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings CopyGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "copy";
+  settings.kernel_name = "CopyMatrixFast";
+  settings.sources =
+#include "../src/kernels/level3/level3.opencl"
+#include "../src/kernels/level3/copy_fast.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.n;
+  settings.size_b = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3};
+  settings.outputs = {3};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"COPY_DIMX", "COPY_DIMY"}};
+  settings.div_global = {{"COPY_VW", "COPY_WPT"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"COPY_DIMX", {8, 16, 32}},
+    {"COPY_DIMY", {8, 16, 32}},
+    {"COPY_WPT", {1, 2, 4, 8}},
+    {"COPY_VW", {1, 2, 4, 8}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void CopyTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> CopySetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void CopySetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(2, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(3, GetRealArg(args.alpha));
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/copy_pad.cpp b/src/tuning/kernels/copy_pad.cpp
index 1b483e86..909a71c8 100644
--- a/src/tuning/kernels/copy_pad.cpp
+++ b/src/tuning/kernels/copy_pad.cpp
@@ -11,94 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "pad";
-  settings.kernel_name = "CopyPadMatrix";
-  settings.sources =
-#include "../src/kernels/level3/level3.opencl"
-#include "../src/kernels/level3/copy_pad.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.n;
-  settings.size_b = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3};
-  settings.outputs = {3};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"PAD_DIMX", "PAD_DIMY"}};
-  settings.div_global = {{"PAD_WPTX", "PAD_WPTY"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"PAD_DIMX", {8, 16, 32}},
-    {"PAD_DIMY", {8, 16, 32}},
-    {"PAD_WPTX", {1, 2, 4}},
-    {"PAD_WPTY", {1, 2, 4}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, static_cast<int>(args.m));
-  kernel.SetArgument(3, 0);
-  kernel.SetArgument(4, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(5, static_cast<int>(args.m));
-  kernel.SetArgument(6, static_cast<int>(args.n));
-  kernel.SetArgument(7, static_cast<int>(args.m));
-  kernel.SetArgument(8, 0);
-  kernel.SetArgument(9, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(10, GetRealArg(args.alpha));
-  kernel.SetArgument(11, 0);
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/copy_pad.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -109,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::PadGetTunerDefaults, clblast::PadGetTunerSettings<half>, clblast::PadTestValidArguments<half>, clblast::PadSetConstraints, clblast::PadSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::PadGetTunerDefaults, clblast::PadGetTunerSettings<float>, clblast::PadTestValidArguments<float>, clblast::PadSetConstraints, clblast::PadSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::PadGetTunerDefaults, clblast::PadGetTunerSettings<double>, clblast::PadTestValidArguments<double>, clblast::PadSetConstraints, clblast::PadSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::PadGetTunerDefaults, clblast::PadGetTunerSettings<float2>, clblast::PadTestValidArguments<float2>, clblast::PadSetConstraints, clblast::PadSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::PadGetTunerDefaults, clblast::PadGetTunerSettings<double2>, clblast::PadTestValidArguments<double2>, clblast::PadSetConstraints, clblast::PadSetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/copy_pad.hpp b/src/tuning/kernels/copy_pad.hpp
new file mode 100644
index 00000000..e612ca9e
--- /dev/null
+++ b/src/tuning/kernels/copy_pad.hpp
@@ -0,0 +1,101 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the pad OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults PadGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings PadGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "pad";
+  settings.kernel_name = "CopyPadMatrix";
+  settings.sources =
+#include "../src/kernels/level3/level3.opencl"
+#include "../src/kernels/level3/copy_pad.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.n;
+  settings.size_b = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3};
+  settings.outputs = {3};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"PAD_DIMX", "PAD_DIMY"}};
+  settings.div_global = {{"PAD_WPTX", "PAD_WPTY"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"PAD_DIMX", {8, 16, 32}},
+    {"PAD_DIMY", {8, 16, 32}},
+    {"PAD_WPTX", {1, 2, 4}},
+    {"PAD_WPTY", {1, 2, 4}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void PadTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> PadSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void PadSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, static_cast<int>(args.m));
+  kernel.SetArgument(3, 0);
+  kernel.SetArgument(4, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(5, static_cast<int>(args.m));
+  kernel.SetArgument(6, static_cast<int>(args.n));
+  kernel.SetArgument(7, static_cast<int>(args.m));
+  kernel.SetArgument(8, 0);
+  kernel.SetArgument(9, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(10, GetRealArg(args.alpha));
+  kernel.SetArgument(11, 0);
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/invert.cpp b/src/tuning/kernels/invert.cpp
index d846fdf7..3dfeb508 100644
--- a/src/tuning/kernels/invert.cpp
+++ b/src/tuning/kernels/invert.cpp
@@ -11,99 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgN, kArgM, kArgK};
-  settings.default_n = 128; // dimension of input matrix 'n'
-  settings.default_m = 64; // block size
-  settings.default_k = 16; // current size
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "invert";
-  settings.kernel_name = "TripleMatMul16Part1Lower";
-  settings.sources =
-"#define ROUTINE_INVERT"
-#include "../src/kernels/level3/invert_diagonal_blocks_part1.opencl"
-#include "../src/kernels/level3/invert_diagonal_blocks_part2.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.n * args.n + args.a_offset;
-  settings.size_b = Ceil(args.n, args.m) * args.m; // Ceil(n, block_size) * block_size
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3};
-  settings.outputs = {3};
-
-  // Sets the base thread configuration
-  const auto num_pages = CeilDiv(args.n, args.k * 2); // CeilDiv(n, current_size*2)
-  settings.global_size = {args.k / 4, num_pages * (args.k / 16) * 4};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {4, 4};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"TMMWGSX", "TMMWGSY"}};
-  settings.div_global = {{}};
-
-  // Sets the tuning parameters and their possible values
-  // TODO: Make these actually tunable, apart from LOCALPAD
-  settings.parameters = {
-    {"INTERNAL_BLOCK_SIZE", {16}},
-    {"LOCALPAD", {0, 1}},
-    {"TMMWGSX", {4}},
-    {"TMMWGSY", {4}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 1 * GetBytes(args.precision);
-  settings.performance_unit = "N/A";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &args) {
-  if (!(args.k == 16)) {
-    throw std::runtime_error("'TripleMatMul16Part1Lower' requires 'k' to be 16");
-  }
-}
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  const auto num_pages = CeilDiv(args.n, args.k * 2); // CeilDiv(n, current_size*2)
-  kernel.SetArgument(0, static_cast<int>(args.n)); // n
-  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(2, 0); // a_offset
-  kernel.SetArgument(3, static_cast<int>(args.n)); // a_ld
-  kernel.SetArgument(4, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(5, static_cast<int>(args.k)); // current_size
-  kernel.SetArgument(6, static_cast<int>(num_pages)); // num_pages
-  kernel.SetArgument(7, static_cast<int>(args.m)); // block_size
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/invert.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -114,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::InvertGetTunerDefaults, clblast::InvertGetTunerSettings<half>, clblast::InvertTestValidArguments<half>, clblast::InvertSetConstraints, clblast::InvertSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::InvertGetTunerDefaults, clblast::InvertGetTunerSettings<float>, clblast::InvertTestValidArguments<float>, clblast::InvertSetConstraints, clblast::InvertSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::InvertGetTunerDefaults, clblast::InvertGetTunerSettings<double>, clblast::InvertTestValidArguments<double>, clblast::InvertSetConstraints, clblast::InvertSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::InvertGetTunerDefaults, clblast::InvertGetTunerSettings<float2>, clblast::InvertTestValidArguments<float2>, clblast::InvertSetConstraints, clblast::InvertSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::InvertGetTunerDefaults, clblast::InvertGetTunerSettings<double2>, clblast::InvertTestValidArguments<double2>, clblast::InvertSetConstraints, clblast::InvertSetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/invert.hpp b/src/tuning/kernels/invert.hpp
new file mode 100644
index 00000000..0a0c9ce2
--- /dev/null
+++ b/src/tuning/kernels/invert.hpp
@@ -0,0 +1,106 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the invert OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults InvertGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgN, kArgM, kArgK};
+  settings.default_n = 128; // dimension of input matrix 'n'
+  settings.default_m = 64; // block size
+  settings.default_k = 16; // current size
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings InvertGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "invert";
+  settings.kernel_name = "TripleMatMul16Part1Lower";
+  settings.sources =
+"#define ROUTINE_INVERT"
+#include "../src/kernels/level3/invert_diagonal_blocks_part1.opencl"
+#include "../src/kernels/level3/invert_diagonal_blocks_part2.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.n * args.n + args.a_offset;
+  settings.size_b = Ceil(args.n, args.m) * args.m; // Ceil(n, block_size) * block_size
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3};
+  settings.outputs = {3};
+
+  // Sets the base thread configuration
+  const auto num_pages = CeilDiv(args.n, args.k * 2); // CeilDiv(n, current_size*2)
+  settings.global_size = {args.k / 4, num_pages * (args.k / 16) * 4};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {4, 4};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"TMMWGSX", "TMMWGSY"}};
+  settings.div_global = {{}};
+
+  // Sets the tuning parameters and their possible values
+  // TODO: Make these actually tunable, apart from LOCALPAD
+  settings.parameters = {
+    {"INTERNAL_BLOCK_SIZE", {16}},
+    {"LOCALPAD", {0, 1}},
+    {"TMMWGSX", {4}},
+    {"TMMWGSY", {4}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 1 * GetBytes(args.precision);
+  settings.performance_unit = "N/A";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void InvertTestValidArguments(const int, const Arguments<T> &args) {
+  if (!(args.k == 16)) {
+    throw std::runtime_error("'TripleMatMul16Part1Lower' requires 'k' to be 16");
+  }
+}
+std::vector<Constraint> InvertSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void InvertSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  const auto num_pages = CeilDiv(args.n, args.k * 2); // CeilDiv(n, current_size*2)
+  kernel.SetArgument(0, static_cast<int>(args.n)); // n
+  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(2, 0); // a_offset
+  kernel.SetArgument(3, static_cast<int>(args.n)); // a_ld
+  kernel.SetArgument(4, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(5, static_cast<int>(args.k)); // current_size
+  kernel.SetArgument(6, static_cast<int>(num_pages)); // num_pages
+  kernel.SetArgument(7, static_cast<int>(args.m)); // block_size
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/transpose_fast.cpp b/src/tuning/kernels/transpose_fast.cpp
index 5b701a5b..6b37a31d 100644
--- a/src/tuning/kernels/transpose_fast.cpp
+++ b/src/tuning/kernels/transpose_fast.cpp
@@ -11,86 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "transpose";
-  settings.kernel_name = "TransposeMatrixFast";
-  settings.sources =
-#include "../src/kernels/level3/level3.opencl"
-#include "../src/kernels/level3/transpose_fast.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.n;
-  settings.size_b = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3};
-  settings.outputs = {3};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"TRA_DIM", "TRA_DIM"}};
-  settings.div_global = {{"TRA_WPT", "TRA_WPT"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"TRA_DIM", {4, 8, 16, 32, 64}},
-    {"TRA_WPT", {1, 2, 4, 8, 16}},
-    {"TRA_PAD", {0, 1}},
-    {"TRA_SHUFFLE", {0, 1}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(2, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(3, GetRealArg(args.alpha));
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/transpose_fast.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -101,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::TransposeGetTunerDefaults, clblast::TransposeGetTunerSettings<half>, clblast::TransposeTestValidArguments<half>, clblast::TransposeSetConstraints, clblast::TransposeSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::TransposeGetTunerDefaults, clblast::TransposeGetTunerSettings<float>, clblast::TransposeTestValidArguments<float>, clblast::TransposeSetConstraints, clblast::TransposeSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::TransposeGetTunerDefaults, clblast::TransposeGetTunerSettings<double>, clblast::TransposeTestValidArguments<double>, clblast::TransposeSetConstraints, clblast::TransposeSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::TransposeGetTunerDefaults, clblast::TransposeGetTunerSettings<float2>, clblast::TransposeTestValidArguments<float2>, clblast::TransposeSetConstraints, clblast::TransposeSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::TransposeGetTunerDefaults, clblast::TransposeGetTunerSettings<double2>, clblast::TransposeTestValidArguments<double2>, clblast::TransposeSetConstraints, clblast::TransposeSetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/transpose_fast.hpp b/src/tuning/kernels/transpose_fast.hpp
new file mode 100644
index 00000000..e8917ad2
--- /dev/null
+++ b/src/tuning/kernels/transpose_fast.hpp
@@ -0,0 +1,93 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the transpose OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults TransposeGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings TransposeGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "transpose";
+  settings.kernel_name = "TransposeMatrixFast";
+  settings.sources =
+#include "../src/kernels/level3/level3.opencl"
+#include "../src/kernels/level3/transpose_fast.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.n;
+  settings.size_b = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3};
+  settings.outputs = {3};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"TRA_DIM", "TRA_DIM"}};
+  settings.div_global = {{"TRA_WPT", "TRA_WPT"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"TRA_DIM", {4, 8, 16, 32, 64}},
+    {"TRA_WPT", {1, 2, 4, 8, 16}},
+    {"TRA_PAD", {0, 1}},
+    {"TRA_SHUFFLE", {0, 1}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void TransposeTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> TransposeSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void TransposeSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(2, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(3, GetRealArg(args.alpha));
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/transpose_pad.cpp b/src/tuning/kernels/transpose_pad.cpp
index ed24fb04..fc7244f6 100644
--- a/src/tuning/kernels/transpose_pad.cpp
+++ b/src/tuning/kernels/transpose_pad.cpp
@@ -11,93 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "padtranspose";
-  settings.kernel_name = "TransposePadMatrix";
-  settings.sources =
-#include "../src/kernels/level3/level3.opencl"
-#include "../src/kernels/level3/transpose_pad.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.n;
-  settings.size_b = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3};
-  settings.outputs = {3};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"PADTRA_TILE", "PADTRA_TILE"}};
-  settings.div_global = {{"PADTRA_WPT", "PADTRA_WPT"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"PADTRA_TILE", {8, 16, 32, 64}},
-    {"PADTRA_WPT", {1, 2, 4, 8, 16}},
-    {"PADTRA_PAD", {0, 1}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, static_cast<int>(args.m));
-  kernel.SetArgument(3, 0);
-  kernel.SetArgument(4, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(5, static_cast<int>(args.n));
-  kernel.SetArgument(6, static_cast<int>(args.m));
-  kernel.SetArgument(7, static_cast<int>(args.n));
-  kernel.SetArgument(8, 0);
-  kernel.SetArgument(9, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(10, GetRealArg(args.alpha));
-  kernel.SetArgument(11, 0);
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/transpose_pad.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -108,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::PadtransposeGetTunerDefaults, clblast::PadtransposeGetTunerSettings<half>, clblast::PadtransposeTestValidArguments<half>, clblast::PadtransposeSetConstraints, clblast::PadtransposeSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::PadtransposeGetTunerDefaults, clblast::PadtransposeGetTunerSettings<float>, clblast::PadtransposeTestValidArguments<float>, clblast::PadtransposeSetConstraints, clblast::PadtransposeSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::PadtransposeGetTunerDefaults, clblast::PadtransposeGetTunerSettings<double>, clblast::PadtransposeTestValidArguments<double>, clblast::PadtransposeSetConstraints, clblast::PadtransposeSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::PadtransposeGetTunerDefaults, clblast::PadtransposeGetTunerSettings<float2>, clblast::PadtransposeTestValidArguments<float2>, clblast::PadtransposeSetConstraints, clblast::PadtransposeSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::PadtransposeGetTunerDefaults, clblast::PadtransposeGetTunerSettings<double2>, clblast::PadtransposeTestValidArguments<double2>, clblast::PadtransposeSetConstraints, clblast::PadtransposeSetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/transpose_pad.hpp b/src/tuning/kernels/transpose_pad.hpp
new file mode 100644
index 00000000..8d24a0dc
--- /dev/null
+++ b/src/tuning/kernels/transpose_pad.hpp
@@ -0,0 +1,100 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the pad-transpose OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults PadtransposeGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings PadtransposeGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "padtranspose";
+  settings.kernel_name = "TransposePadMatrix";
+  settings.sources =
+#include "../src/kernels/level3/level3.opencl"
+#include "../src/kernels/level3/transpose_pad.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.n;
+  settings.size_b = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3};
+  settings.outputs = {3};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"PADTRA_TILE", "PADTRA_TILE"}};
+  settings.div_global = {{"PADTRA_WPT", "PADTRA_WPT"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"PADTRA_TILE", {8, 16, 32, 64}},
+    {"PADTRA_WPT", {1, 2, 4, 8, 16}},
+    {"PADTRA_PAD", {0, 1}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void PadtransposeTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> PadtransposeSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void PadtransposeSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, static_cast<int>(args.m));
+  kernel.SetArgument(3, 0);
+  kernel.SetArgument(4, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(5, static_cast<int>(args.n));
+  kernel.SetArgument(6, static_cast<int>(args.m));
+  kernel.SetArgument(7, static_cast<int>(args.n));
+  kernel.SetArgument(8, 0);
+  kernel.SetArgument(9, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(10, GetRealArg(args.alpha));
+  kernel.SetArgument(11, 0);
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xaxpy.cpp b/src/tuning/kernels/xaxpy.cpp
index dd44018c..6a95600d 100644
--- a/src/tuning/kernels/xaxpy.cpp
+++ b/src/tuning/kernels/xaxpy.cpp
@@ -11,88 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgN, kArgAlpha};
-  settings.default_n = 4096*1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "xaxpy";
-  settings.kernel_name = "XaxpyFastest";
-  settings.sources =
-#include "../src/kernels/level1/level1.opencl"
-#include "../src/kernels/level1/xaxpy.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_x = args.n;
-  settings.size_y = args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {0, 1};
-  settings.outputs = {1};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1};
-  settings.local_size_ref = {64};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"WGS"}};
-  settings.div_global = {{"WPT"},{"VW"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"WGS", {64, 128, 256, 512, 1024, 2048}},
-    {"WPT", {1, 2, 4, 8}},
-    {"VW", {1, 2, 4, 8}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 3 * args.n * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &args) {
-  if (!IsMultiple(args.n, 64)) {
-    throw std::runtime_error("'XaxpyFastest' requires 'n' to be a multiple of WGS*WPT*VW");
-  }
-}
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.n));
-  kernel.SetArgument(1, GetRealArg(args.alpha));
-  kernel.SetArgument(2, buffers[0]()); // 0 == X vector
-  kernel.SetArgument(3, buffers[1]()); // 1 == Y vector
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xaxpy.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -103,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::XaxpyGetTunerDefaults, clblast::XaxpyGetTunerSettings<half>, clblast::XaxpyTestValidArguments<half>, clblast::XaxpySetConstraints, clblast::XaxpySetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::XaxpyGetTunerDefaults, clblast::XaxpyGetTunerSettings<float>, clblast::XaxpyTestValidArguments<float>, clblast::XaxpySetConstraints, clblast::XaxpySetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::XaxpyGetTunerDefaults, clblast::XaxpyGetTunerSettings<double>, clblast::XaxpyTestValidArguments<double>, clblast::XaxpySetConstraints, clblast::XaxpySetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::XaxpyGetTunerDefaults, clblast::XaxpyGetTunerSettings<float2>, clblast::XaxpyTestValidArguments<float2>, clblast::XaxpySetConstraints, clblast::XaxpySetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::XaxpyGetTunerDefaults, clblast::XaxpyGetTunerSettings<double2>, clblast::XaxpyTestValidArguments<double2>, clblast::XaxpySetConstraints, clblast::XaxpySetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/xaxpy.hpp b/src/tuning/kernels/xaxpy.hpp
new file mode 100644
index 00000000..24550ed9
--- /dev/null
+++ b/src/tuning/kernels/xaxpy.hpp
@@ -0,0 +1,95 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the xaxpy OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XaxpyGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgN, kArgAlpha};
+  settings.default_n = 4096*1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XaxpyGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "xaxpy";
+  settings.kernel_name = "XaxpyFastest";
+  settings.sources =
+#include "../src/kernels/level1/level1.opencl"
+#include "../src/kernels/level1/xaxpy.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_x = args.n;
+  settings.size_y = args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {0, 1};
+  settings.outputs = {1};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1};
+  settings.local_size_ref = {64};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"WGS"}};
+  settings.div_global = {{"WPT"},{"VW"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"WGS", {64, 128, 256, 512, 1024, 2048}},
+    {"WPT", {1, 2, 4, 8}},
+    {"VW", {1, 2, 4, 8}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 3 * args.n * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XaxpyTestValidArguments(const int, const Arguments<T> &args) {
+  if (!IsMultiple(args.n, 64)) {
+    throw std::runtime_error("'XaxpyFastest' requires 'n' to be a multiple of WGS*WPT*VW");
+  }
+}
+std::vector<Constraint> XaxpySetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void XaxpySetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.n));
+  kernel.SetArgument(1, GetRealArg(args.alpha));
+  kernel.SetArgument(2, buffers[0]()); // 0 == X vector
+  kernel.SetArgument(3, buffers[1]()); // 1 == Y vector
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xdot.cpp b/src/tuning/kernels/xdot.cpp
index 635d012a..6d10c4d8 100644
--- a/src/tuning/kernels/xdot.cpp
+++ b/src/tuning/kernels/xdot.cpp
@@ -12,94 +12,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgN};
-  settings.default_n = 2*1024*1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int V, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "xdot_"+std::to_string(V);
-  settings.kernel_name = (V==1) ? "Xdot" : "XdotEpilogue";
-  settings.sources =
-#include "../src/kernels/level1/xdot.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_x = args.n;
-  settings.size_y = args.n;
-  settings.size_temp = args.n; // Worst case
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {0, 1, 5};
-  settings.outputs = {}; // no output checking
-
-  // Sets the base thread configuration
-  settings.global_size = (V==1) ? std::vector<size_t>{2*64} : std::vector<size_t>{1};
-  settings.global_size_ref = (V==1) ? std::vector<size_t>{2*64*64} : std::vector<size_t>{64};
-  settings.local_size = {1};
-  settings.local_size_ref = {64};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = (V==1) ? TransformVector{{"WGS1"}} : TransformVector{{"WGS2"}};
-  settings.mul_global = (V==1) ? TransformVector{{"WGS1"}} : TransformVector{{"WGS2"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"WGS"+std::to_string(V), {32, 64, 128, 256, 512, 1024}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = (V==1) ? (2*args.n + 1) * GetBytes(args.precision) : 1 * GetBytes(args.precision);
-  settings.performance_unit = (V==1) ? "GB/s" : "N/A";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int V, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  if (V == 1) {
-    kernel.SetArgument(0, static_cast<int>(args.n));
-    kernel.SetArgument(1, buffers[0]()); // 0 == X vector
-    kernel.SetArgument(2, 0);
-    kernel.SetArgument(3, 1);
-    kernel.SetArgument(4, buffers[1]()); // 1 == Y vector
-    kernel.SetArgument(5, 0);
-    kernel.SetArgument(6, 1);
-    kernel.SetArgument(7, buffers[5]()); // 5 == temp; no output checking - size varies
-    kernel.SetArgument(8, static_cast<int>(false));
-  }
-  else {
-    kernel.SetArgument(0, buffers[5]()); // 5 == temp
-    kernel.SetArgument(1, buffers[0]()); // 0 == X vector; no output checking - size varies
-    kernel.SetArgument(2, 0);
-  }
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xdot.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -111,11 +24,11 @@ template <int V>
 void StartVariation(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::XdotGetTunerDefaults, clblast::XdotGetTunerSettings<half>, clblast::XdotTestValidArguments<half>, clblast::XdotSetConstraints, clblast::XdotSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::XdotGetTunerDefaults, clblast::XdotGetTunerSettings<float>, clblast::XdotTestValidArguments<float>, clblast::XdotSetConstraints, clblast::XdotSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::XdotGetTunerDefaults, clblast::XdotGetTunerSettings<double>, clblast::XdotTestValidArguments<double>, clblast::XdotSetConstraints, clblast::XdotSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::XdotGetTunerDefaults, clblast::XdotGetTunerSettings<float2>, clblast::XdotTestValidArguments<float2>, clblast::XdotSetConstraints, clblast::XdotSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::XdotGetTunerDefaults, clblast::XdotGetTunerSettings<double2>, clblast::XdotTestValidArguments<double2>, clblast::XdotSetConstraints, clblast::XdotSetArguments<double2>); break;
   }
 }
 
diff --git a/src/tuning/kernels/xdot.hpp b/src/tuning/kernels/xdot.hpp
new file mode 100644
index 00000000..15673c79
--- /dev/null
+++ b/src/tuning/kernels/xdot.hpp
@@ -0,0 +1,102 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the xdot OpenCL kernels. Note that the results are
+// not verified, since the result is not final and depends on the WGS2 parameter.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XdotGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgN};
+  settings.default_n = 2*1024*1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XdotGetTunerSettings(const int V, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "xdot_"+std::to_string(V);
+  settings.kernel_name = (V==1) ? "Xdot" : "XdotEpilogue";
+  settings.sources =
+#include "../src/kernels/level1/xdot.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_x = args.n;
+  settings.size_y = args.n;
+  settings.size_temp = args.n; // Worst case
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {0, 1, 5};
+  settings.outputs = {}; // no output checking
+
+  // Sets the base thread configuration
+  settings.global_size = (V==1) ? std::vector<size_t>{2*64} : std::vector<size_t>{1};
+  settings.global_size_ref = (V==1) ? std::vector<size_t>{2*64*64} : std::vector<size_t>{64};
+  settings.local_size = {1};
+  settings.local_size_ref = {64};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = (V==1) ? TransformVector{{"WGS1"}} : TransformVector{{"WGS2"}};
+  settings.mul_global = (V==1) ? TransformVector{{"WGS1"}} : TransformVector{{"WGS2"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"WGS"+std::to_string(V), {32, 64, 128, 256, 512, 1024}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = (V==1) ? (2*args.n + 1) * GetBytes(args.precision) : 1 * GetBytes(args.precision);
+  settings.performance_unit = (V==1) ? "GB/s" : "N/A";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XdotTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> XdotSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void XdotSetArguments(const int V, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  if (V == 1) {
+    kernel.SetArgument(0, static_cast<int>(args.n));
+    kernel.SetArgument(1, buffers[0]()); // 0 == X vector
+    kernel.SetArgument(2, 0);
+    kernel.SetArgument(3, 1);
+    kernel.SetArgument(4, buffers[1]()); // 1 == Y vector
+    kernel.SetArgument(5, 0);
+    kernel.SetArgument(6, 1);
+    kernel.SetArgument(7, buffers[5]()); // 5 == temp; no output checking - size varies
+    kernel.SetArgument(8, static_cast<int>(false));
+  }
+  else {
+    kernel.SetArgument(0, buffers[5]()); // 5 == temp
+    kernel.SetArgument(1, buffers[0]()); // 0 == X vector; no output checking - size varies
+    kernel.SetArgument(2, 0);
+  }
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xgemm.cpp b/src/tuning/kernels/xgemm.cpp
index 5c242757..d365ce6d 100644
--- a/src/tuning/kernels/xgemm.cpp
+++ b/src/tuning/kernels/xgemm.cpp
@@ -7,162 +7,11 @@
 // Author(s):
 //   Cedric Nugteren <www.cedricnugteren.nl>
 //
-// This file uses the auto-tuner to tune the xgemm OpenCL kernels. There are two variations:
-// - V==1: This tests some limited set of tuning parameters exhaustively.
-// - V==2: This tests a much larger set of tuning parameters by randomly sampling a subset.
+// This file uses the auto-tuner to tune the xgemm OpenCL kernels.
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int V) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgK, kArgAlpha, kArgBeta, kArgFraction,
-                      kArgHeuristicSelection, kArgPsoSwarmSize,
-                      kArgPsoInfGlobal, kArgPsoInfLocal, kArgPsoInfRandom};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  settings.default_k = 1024;
-  settings.default_fraction = (V==1) ? 1.0 : 512.0; // test all or sample randomly
-  settings.default_num_runs = 2;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int V, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = (V==1) ? "xgemm_1" : "xgemm_2";
-  settings.kernel_name = "Xgemm";
-  settings.sources =
-#include "../src/kernels/level3/xgemm_part1.opencl"
-#include "../src/kernels/level3/xgemm_part2.opencl"
-#include "../src/kernels/level3/xgemm_part3.opencl"
-#include "../src/kernels/level3/xgemm_part4.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.k;
-  settings.size_b = args.n * args.k;
-  settings.size_c = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3, 4};
-  settings.outputs = {4};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"MDIMC", "NDIMC"}};
-  settings.mul_global = {{"MDIMC", "NDIMC"}};
-  settings.div_global = {{"MWG", "NWG"}};
-
-  // Sets the tuning parameters and their possible values
-  if (V==1) { // limited subset of tuning parameters - but explorable exhaustively
-    settings.parameters = {
-      {"MWG", {16, 32, 64}},
-      {"NWG", {16, 32, 64}},
-      {"KWG", {32}},
-      {"MDIMC", {8, 16, 32}},
-      {"NDIMC", {8, 16, 32}},
-      {"MDIMA", {8, 16, 32}},
-      {"NDIMB", {8, 16, 32}},
-      {"KWI", {2}},
-      {"VWM", {1, 2, 4}},
-      {"VWN", {1, 2, 4}},
-      {"STRM", {0}},
-      {"STRN", {0}},
-      {"SA", {0, 1}},
-      {"SB", {0, 1}},
-    };
-  }
-  else { // a lot more tuning parameters - has to be sampled randomly, too much to test all
-    settings.parameters = {
-      {"MWG", {16, 32, 64, 128}},
-      {"NWG", {16, 32, 64, 128}},
-      {"KWG", {16, 32}},
-      {"MDIMC", {8, 16, 32}},
-      {"NDIMC", {8, 16, 32}},
-      {"MDIMA", {8, 16, 32}},
-      {"NDIMB", {8, 16, 32}},
-      {"KWI", {2}},
-      {"VWM", {1, 2, 4, 8}},
-      {"VWN", {1, 2, 4, 8}},
-      {"STRM", {0, 1}},
-      {"STRN", {0, 1}},
-      {"SA", {0, 1}},
-      {"SB", {0, 1}},
-    };
-  }
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * args.k;
-  settings.performance_unit = "GFLOPS";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int V) {
-  auto constraints = std::vector<Constraint>();
-  auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
-  auto MultipleOfXMulY = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]*v[2]); };
-  auto MultipleOfXMulYDivZ = [] (std::vector<size_t> v) { return IsMultiple(v[0], (v[1]*v[2])/v[3]); };
-  // Requirement for unrolling the KWG loop
-  constraints.push_back({MultipleOfX, {"KWG", "KWI"}});
-  // Required for integer MWI and NWI
-  constraints.push_back({MultipleOfXMulY, {"MWG", "MDIMC", "VWM"}});
-  constraints.push_back({MultipleOfXMulY, {"NWG", "NDIMC", "VWN"}});
-  // Required for integer MWIA and NWIB
-  constraints.push_back({MultipleOfXMulY, {"MWG", "MDIMA", "VWM"}});
-  constraints.push_back({MultipleOfXMulY, {"NWG", "NDIMB", "VWN"}});
-  // KWG has to be a multiple of KDIMA = ((MDIMC*NDIMC)/(MDIMA)) and KDIMB = (...)
-  constraints.push_back({MultipleOfXMulYDivZ, {"KWG", "MDIMC", "NDIMC", "MDIMA"}});
-  constraints.push_back({MultipleOfXMulYDivZ, {"KWG", "MDIMC", "NDIMC", "NDIMB"}});
-
-  // Extra constraints for variation 1 to limit the set of options significantly
-  if (V==1) {
-    auto IsEqual = [] (std::vector<size_t> v) { return v[0] == v[1]; };
-    constraints.push_back({IsEqual, {"MDIMC", "MDIMA"}});
-    constraints.push_back({IsEqual, {"NDIMC", "NDIMB"}});
-    constraints.push_back({IsEqual, {"SA", "SB"}});
-  }
-  return constraints;
-}
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, static_cast<int>(args.k));
-  kernel.SetArgument(3, GetRealArg(args.alpha));
-  kernel.SetArgument(4, GetRealArg(args.beta));
-  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(6, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(7, buffers[4]()); // 4 == C matrix
-  kernel.SetArgument(8, 0);
-  kernel.SetArgument(9, 0);
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xgemm.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -174,11 +23,11 @@ template <int V>
 void StartVariation(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::XgemmGetTunerDefaults, clblast::XgemmGetTunerSettings<half>, clblast::XgemmTestValidArguments<half>, clblast::XgemmSetConstraints, clblast::XgemmSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::XgemmGetTunerDefaults, clblast::XgemmGetTunerSettings<float>, clblast::XgemmTestValidArguments<float>, clblast::XgemmSetConstraints, clblast::XgemmSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::XgemmGetTunerDefaults, clblast::XgemmGetTunerSettings<double>, clblast::XgemmTestValidArguments<double>, clblast::XgemmSetConstraints, clblast::XgemmSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::XgemmGetTunerDefaults, clblast::XgemmGetTunerSettings<float2>, clblast::XgemmTestValidArguments<float2>, clblast::XgemmSetConstraints, clblast::XgemmSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::XgemmGetTunerDefaults, clblast::XgemmGetTunerSettings<double2>, clblast::XgemmTestValidArguments<double2>, clblast::XgemmSetConstraints, clblast::XgemmSetArguments<double2>); break;
   }
 }
 
diff --git a/src/tuning/kernels/xgemm.hpp b/src/tuning/kernels/xgemm.hpp
new file mode 100644
index 00000000..66e197e1
--- /dev/null
+++ b/src/tuning/kernels/xgemm.hpp
@@ -0,0 +1,165 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the xgemm OpenCL kernels. There are two variations:
+// - V==1: This tests some limited set of tuning parameters exhaustively.
+// - V==2: This tests a much larger set of tuning parameters by randomly sampling a subset.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XgemmGetTunerDefaults(const int V) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgK, kArgAlpha, kArgBeta, kArgFraction,
+                      kArgHeuristicSelection, kArgPsoSwarmSize,
+                      kArgPsoInfGlobal, kArgPsoInfLocal, kArgPsoInfRandom};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  settings.default_k = 1024;
+  settings.default_fraction = (V==1) ? 1.0 : 512.0; // test all or sample randomly
+  settings.default_num_runs = 2;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XgemmGetTunerSettings(const int V, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = (V==1) ? "xgemm_1" : "xgemm_2";
+  settings.kernel_name = "Xgemm";
+  settings.sources =
+#include "../src/kernels/level3/xgemm_part1.opencl"
+#include "../src/kernels/level3/xgemm_part2.opencl"
+#include "../src/kernels/level3/xgemm_part3.opencl"
+#include "../src/kernels/level3/xgemm_part4.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.k;
+  settings.size_b = args.n * args.k;
+  settings.size_c = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3, 4};
+  settings.outputs = {4};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"MDIMC", "NDIMC"}};
+  settings.mul_global = {{"MDIMC", "NDIMC"}};
+  settings.div_global = {{"MWG", "NWG"}};
+
+  // Sets the tuning parameters and their possible values
+  if (V==1) { // limited subset of tuning parameters - but explorable exhaustively
+    settings.parameters = {
+      {"MWG", {16, 32, 64}},
+      {"NWG", {16, 32, 64}},
+      {"KWG", {32}},
+      {"MDIMC", {8, 16, 32}},
+      {"NDIMC", {8, 16, 32}},
+      {"MDIMA", {8, 16, 32}},
+      {"NDIMB", {8, 16, 32}},
+      {"KWI", {2}},
+      {"VWM", {1, 2, 4}},
+      {"VWN", {1, 2, 4}},
+      {"STRM", {0}},
+      {"STRN", {0}},
+      {"SA", {0, 1}},
+      {"SB", {0, 1}},
+    };
+  }
+  else { // a lot more tuning parameters - has to be sampled randomly, too much to test all
+    settings.parameters = {
+      {"MWG", {16, 32, 64, 128}},
+      {"NWG", {16, 32, 64, 128}},
+      {"KWG", {16, 32}},
+      {"MDIMC", {8, 16, 32}},
+      {"NDIMC", {8, 16, 32}},
+      {"MDIMA", {8, 16, 32}},
+      {"NDIMB", {8, 16, 32}},
+      {"KWI", {2}},
+      {"VWM", {1, 2, 4, 8}},
+      {"VWN", {1, 2, 4, 8}},
+      {"STRM", {0, 1}},
+      {"STRN", {0, 1}},
+      {"SA", {0, 1}},
+      {"SB", {0, 1}},
+    };
+  }
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * args.k;
+  settings.performance_unit = "GFLOPS";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XgemmTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> XgemmSetConstraints(const int V) {
+  auto constraints = std::vector<Constraint>();
+  auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
+  auto MultipleOfXMulY = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]*v[2]); };
+  auto MultipleOfXMulYDivZ = [] (std::vector<size_t> v) { return IsMultiple(v[0], (v[1]*v[2])/v[3]); };
+  // Requirement for unrolling the KWG loop
+  constraints.push_back({MultipleOfX, {"KWG", "KWI"}});
+  // Required for integer MWI and NWI
+  constraints.push_back({MultipleOfXMulY, {"MWG", "MDIMC", "VWM"}});
+  constraints.push_back({MultipleOfXMulY, {"NWG", "NDIMC", "VWN"}});
+  // Required for integer MWIA and NWIB
+  constraints.push_back({MultipleOfXMulY, {"MWG", "MDIMA", "VWM"}});
+  constraints.push_back({MultipleOfXMulY, {"NWG", "NDIMB", "VWN"}});
+  // KWG has to be a multiple of KDIMA = ((MDIMC*NDIMC)/(MDIMA)) and KDIMB = (...)
+  constraints.push_back({MultipleOfXMulYDivZ, {"KWG", "MDIMC", "NDIMC", "MDIMA"}});
+  constraints.push_back({MultipleOfXMulYDivZ, {"KWG", "MDIMC", "NDIMC", "NDIMB"}});
+
+  // Extra constraints for variation 1 to limit the set of options significantly
+  if (V==1) {
+    auto IsEqual = [] (std::vector<size_t> v) { return v[0] == v[1]; };
+    constraints.push_back({IsEqual, {"MDIMC", "MDIMA"}});
+    constraints.push_back({IsEqual, {"NDIMC", "NDIMB"}});
+    constraints.push_back({IsEqual, {"SA", "SB"}});
+  }
+  return constraints;
+}
+
+// Sets the kernel's arguments
+template <typename T>
+void XgemmSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, static_cast<int>(args.k));
+  kernel.SetArgument(3, GetRealArg(args.alpha));
+  kernel.SetArgument(4, GetRealArg(args.beta));
+  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(6, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(7, buffers[4]()); // 4 == C matrix
+  kernel.SetArgument(8, 0);
+  kernel.SetArgument(9, 0);
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xgemm_direct.cpp b/src/tuning/kernels/xgemm_direct.cpp
index 0bd2e94d..7298a6c3 100644
--- a/src/tuning/kernels/xgemm_direct.cpp
+++ b/src/tuning/kernels/xgemm_direct.cpp
@@ -7,159 +7,11 @@
 // Author(s):
 //   Cedric Nugteren <www.cedricnugteren.nl>
 //
-// This file uses the auto-tuner to tune the direct xgemm kernels. There are two variations:
-// - V==1: This tests some limited set of tuning parameters exhaustively.
-// - V==2: This tests a much larger set of tuning parameters by randomly sampling a subset.
+// This file uses the auto-tuner to tune the direct xgemm kernels.
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int V) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgK, kArgAlpha, kArgBeta, kArgFraction,
-                      kArgHeuristicSelection, kArgPsoSwarmSize,
-                      kArgPsoInfGlobal, kArgPsoInfLocal, kArgPsoInfRandom};
-  settings.default_m = 256;
-  settings.default_n = 256;
-  settings.default_k = 256;
-  settings.default_fraction = (V==1) ? 1.0 : 64.0; // test all or sample randomly
-  settings.default_num_runs = 4;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int V, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = (V==1) ? "xgemm_direct_1" : "xgemm_direct_2";
-  settings.kernel_name = "XgemmDirectTN";
-  settings.sources =
-#include "../src/kernels/level3/xgemm_direct_part1.opencl"
-#include "../src/kernels/level3/xgemm_direct_part2.opencl"
-#include "../src/kernels/level3/xgemm_direct_part3.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_a = args.m * args.k;
-  settings.size_b = args.n * args.k;
-  settings.size_c = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {2, 3, 4};
-  settings.outputs = {4};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"MDIMCD", "NDIMCD"}};
-  settings.mul_global = {{"MDIMCD", "NDIMCD"}};
-  settings.div_global = {{"WGD", "WGD"}};
-
-  // Sets the tuning parameters and their possible values
-  if (V==1) { // limited subset of tuning parameters - but explorable exhaustively
-    settings.parameters = {
-      {"WGD", {8, 16, 32}},
-      {"MDIMCD", {8, 16, 32}},
-      {"NDIMCD", {8, 16, 32}},
-      {"MDIMAD", {8, 16, 32}},
-      {"NDIMBD", {8, 16, 32}},
-      {"KWID", {2}},
-      {"VWMD", {1, 2, 4, 8}},
-      {"VWND", {1, 2, 4, 8}},
-      {"PADA", {1}},
-      {"PADB", {1}},
-    };
-  }
-  else { // a lot more tuning parameters - has to be sampled randomly, too much to test all
-    settings.parameters = {
-      {"WGD", {8, 16, 32, 64}},
-      {"MDIMCD", {8, 16, 32}},
-      {"NDIMCD", {8, 16, 32}},
-      {"MDIMAD", {8, 16, 32}},
-      {"NDIMBD", {8, 16, 32}},
-      {"KWID", {2, 8, 16}},
-      {"VWMD", {1, 2, 4, 8}},
-      {"VWND", {1, 2, 4, 8}},
-      {"PADA", {0, 1}},
-      {"PADB", {0, 1}},
-    };
-  }
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = 2 * args.m * args.n * args.k;
-  settings.performance_unit = "GFLOPS";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int V) {
-  auto constraints = std::vector<Constraint>();
-  auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
-  auto MultipleOfXMulY = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]*v[2]); };
-  auto MultipleOfXMulYDivZ = [] (std::vector<size_t> v) { return IsMultiple(v[0], (v[1]*v[2])/v[3]); };
-  // Requirement for unrolling the WGD loop
-  constraints.push_back({MultipleOfX, {"WGD", "KWID"}});
-  // Required for integer MWID and NWID
-  constraints.push_back({MultipleOfXMulY, {"WGD", "MDIMCD", "VWMD"}});
-  constraints.push_back({MultipleOfXMulY, {"WGD", "NDIMCD", "VWND"}});
-  // Required for integer MWIAD and NWIBD
-  constraints.push_back({MultipleOfXMulY, {"WGD", "MDIMAD", "VWMD"}});
-  constraints.push_back({MultipleOfXMulY, {"WGD", "NDIMBD", "VWND"}});
-  // WGD has to be a multiple of KDIMAD = ((MDIMCD*NDIMCD)/(MDIMAD)) and KDIMBD = (...)
-  constraints.push_back({MultipleOfXMulYDivZ, {"WGD", "MDIMCD", "NDIMCD", "MDIMAD"}});
-  constraints.push_back({MultipleOfXMulYDivZ, {"WGD", "MDIMCD", "NDIMCD", "NDIMBD"}});
-
-  // Extra constraints for variation 1 to limit the set of options significantly
-  if (V==1) {
-    auto IsEqual = [] (std::vector<size_t> v) { return v[0] == v[1]; };
-    constraints.push_back({IsEqual, {"MDIMCD", "MDIMAD"}});
-    constraints.push_back({IsEqual, {"NDIMCD", "NDIMBD"}});
-  }
-  return constraints;
-}
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, static_cast<int>(args.k));
-  kernel.SetArgument(3, GetRealArg(args.alpha));
-  kernel.SetArgument(4, GetRealArg(args.beta));
-  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(6, 0); // a_offset
-  kernel.SetArgument(7, static_cast<int>(args.k)); // a_ld
-  kernel.SetArgument(8, buffers[3]()); // 3 == B matrix
-  kernel.SetArgument(9, 0); // b_offset
-  kernel.SetArgument(10, static_cast<int>(args.n)); // b_ld
-  kernel.SetArgument(11, buffers[4]()); // 4 == C matrix
-  kernel.SetArgument(12, 0); // c_offset
-  kernel.SetArgument(13, static_cast<int>(args.n)); // c_ld
-  kernel.SetArgument(14, 1); // c_do_transpose
-  kernel.SetArgument(15, 0); // a_conjugate
-  kernel.SetArgument(16, 0); // b_conjugate
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xgemm_direct.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -171,11 +23,11 @@ template <int V>
 void StartVariation(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::XgemmDirectGetTunerDefaults, clblast::XgemmDirectGetTunerSettings<half>, clblast::XgemmDirectTestValidArguments<half>, clblast::XgemmDirectSetConstraints, clblast::XgemmDirectSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::XgemmDirectGetTunerDefaults, clblast::XgemmDirectGetTunerSettings<float>, clblast::XgemmDirectTestValidArguments<float>, clblast::XgemmDirectSetConstraints, clblast::XgemmDirectSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::XgemmDirectGetTunerDefaults, clblast::XgemmDirectGetTunerSettings<double>, clblast::XgemmDirectTestValidArguments<double>, clblast::XgemmDirectSetConstraints, clblast::XgemmDirectSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::XgemmDirectGetTunerDefaults, clblast::XgemmDirectGetTunerSettings<float2>, clblast::XgemmDirectTestValidArguments<float2>, clblast::XgemmDirectSetConstraints, clblast::XgemmDirectSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::XgemmDirectGetTunerDefaults, clblast::XgemmDirectGetTunerSettings<double2>, clblast::XgemmDirectTestValidArguments<double2>, clblast::XgemmDirectSetConstraints, clblast::XgemmDirectSetArguments<double2>); break;
   }
 }
 
diff --git a/src/tuning/kernels/xgemm_direct.hpp b/src/tuning/kernels/xgemm_direct.hpp
new file mode 100644
index 00000000..ecb10bc6
--- /dev/null
+++ b/src/tuning/kernels/xgemm_direct.hpp
@@ -0,0 +1,162 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the direct xgemm kernels. There are two variations:
+// - V==1: This tests some limited set of tuning parameters exhaustively.
+// - V==2: This tests a much larger set of tuning parameters by randomly sampling a subset.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XgemmDirectGetTunerDefaults(const int V) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgK, kArgAlpha, kArgBeta, kArgFraction,
+                      kArgHeuristicSelection, kArgPsoSwarmSize,
+                      kArgPsoInfGlobal, kArgPsoInfLocal, kArgPsoInfRandom};
+  settings.default_m = 256;
+  settings.default_n = 256;
+  settings.default_k = 256;
+  settings.default_fraction = (V==1) ? 1.0 : 64.0; // test all or sample randomly
+  settings.default_num_runs = 4;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XgemmDirectGetTunerSettings(const int V, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = (V==1) ? "xgemm_direct_1" : "xgemm_direct_2";
+  settings.kernel_name = "XgemmDirectTN";
+  settings.sources =
+#include "../src/kernels/level3/xgemm_direct_part1.opencl"
+#include "../src/kernels/level3/xgemm_direct_part2.opencl"
+#include "../src/kernels/level3/xgemm_direct_part3.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_a = args.m * args.k;
+  settings.size_b = args.n * args.k;
+  settings.size_c = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {2, 3, 4};
+  settings.outputs = {4};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"MDIMCD", "NDIMCD"}};
+  settings.mul_global = {{"MDIMCD", "NDIMCD"}};
+  settings.div_global = {{"WGD", "WGD"}};
+
+  // Sets the tuning parameters and their possible values
+  if (V==1) { // limited subset of tuning parameters - but explorable exhaustively
+    settings.parameters = {
+      {"WGD", {8, 16, 32}},
+      {"MDIMCD", {8, 16, 32}},
+      {"NDIMCD", {8, 16, 32}},
+      {"MDIMAD", {8, 16, 32}},
+      {"NDIMBD", {8, 16, 32}},
+      {"KWID", {2}},
+      {"VWMD", {1, 2, 4, 8}},
+      {"VWND", {1, 2, 4, 8}},
+      {"PADA", {1}},
+      {"PADB", {1}},
+    };
+  }
+  else { // a lot more tuning parameters - has to be sampled randomly, too much to test all
+    settings.parameters = {
+      {"WGD", {8, 16, 32, 64}},
+      {"MDIMCD", {8, 16, 32}},
+      {"NDIMCD", {8, 16, 32}},
+      {"MDIMAD", {8, 16, 32}},
+      {"NDIMBD", {8, 16, 32}},
+      {"KWID", {2, 8, 16}},
+      {"VWMD", {1, 2, 4, 8}},
+      {"VWND", {1, 2, 4, 8}},
+      {"PADA", {0, 1}},
+      {"PADB", {0, 1}},
+    };
+  }
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = 2 * args.m * args.n * args.k;
+  settings.performance_unit = "GFLOPS";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XgemmDirectTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> XgemmDirectSetConstraints(const int V) {
+  auto constraints = std::vector<Constraint>();
+  auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
+  auto MultipleOfXMulY = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]*v[2]); };
+  auto MultipleOfXMulYDivZ = [] (std::vector<size_t> v) { return IsMultiple(v[0], (v[1]*v[2])/v[3]); };
+  // Requirement for unrolling the WGD loop
+  constraints.push_back({MultipleOfX, {"WGD", "KWID"}});
+  // Required for integer MWID and NWID
+  constraints.push_back({MultipleOfXMulY, {"WGD", "MDIMCD", "VWMD"}});
+  constraints.push_back({MultipleOfXMulY, {"WGD", "NDIMCD", "VWND"}});
+  // Required for integer MWIAD and NWIBD
+  constraints.push_back({MultipleOfXMulY, {"WGD", "MDIMAD", "VWMD"}});
+  constraints.push_back({MultipleOfXMulY, {"WGD", "NDIMBD", "VWND"}});
+  // WGD has to be a multiple of KDIMAD = ((MDIMCD*NDIMCD)/(MDIMAD)) and KDIMBD = (...)
+  constraints.push_back({MultipleOfXMulYDivZ, {"WGD", "MDIMCD", "NDIMCD", "MDIMAD"}});
+  constraints.push_back({MultipleOfXMulYDivZ, {"WGD", "MDIMCD", "NDIMCD", "NDIMBD"}});
+
+  // Extra constraints for variation 1 to limit the set of options significantly
+  if (V==1) {
+    auto IsEqual = [] (std::vector<size_t> v) { return v[0] == v[1]; };
+    constraints.push_back({IsEqual, {"MDIMCD", "MDIMAD"}});
+    constraints.push_back({IsEqual, {"NDIMCD", "NDIMBD"}});
+  }
+  return constraints;
+}
+
+// Sets the kernel's arguments
+template <typename T>
+void XgemmDirectSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, static_cast<int>(args.k));
+  kernel.SetArgument(3, GetRealArg(args.alpha));
+  kernel.SetArgument(4, GetRealArg(args.beta));
+  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(6, 0); // a_offset
+  kernel.SetArgument(7, static_cast<int>(args.k)); // a_ld
+  kernel.SetArgument(8, buffers[3]()); // 3 == B matrix
+  kernel.SetArgument(9, 0); // b_offset
+  kernel.SetArgument(10, static_cast<int>(args.n)); // b_ld
+  kernel.SetArgument(11, buffers[4]()); // 4 == C matrix
+  kernel.SetArgument(12, 0); // c_offset
+  kernel.SetArgument(13, static_cast<int>(args.n)); // c_ld
+  kernel.SetArgument(14, 1); // c_do_transpose
+  kernel.SetArgument(15, 0); // a_conjugate
+  kernel.SetArgument(16, 0); // b_conjugate
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xgemv.cpp b/src/tuning/kernels/xgemv.cpp
index 965fc53d..9e45d73f 100644
--- a/src/tuning/kernels/xgemv.cpp
+++ b/src/tuning/kernels/xgemv.cpp
@@ -7,135 +7,11 @@
 // Author(s):
 //   Cedric Nugteren <www.cedricnugteren.nl>
 //
-// This file uses the auto-tuner to tune the xgemv OpenCL kernels. Three variants are tuned:
-// 1: The full version of the kernel
-// 2: The fast version for non-transposed matrices
-// 3: The fast version for transposed matrices
+// This file uses the auto-tuner to tune the xgemv OpenCL kernels.
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha, kArgBeta};
-  settings.default_m = 2048;
-  settings.default_n = 2048;
-  settings.default_num_runs = 4;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int V, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = (V==1) ? "xgemv" : ((V==2) ? "xgemv_fast" : "xgemv_fast_rot");
-  settings.kernel_name = (V==1) ? "Xgemv" : ((V==2) ? "XgemvFast" : "XgemvFastRot");
-  settings.sources =
-#include "../src/kernels/level2/xgemv.opencl"
-#include "../src/kernels/level2/xgemv_fast.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_x = args.n;
-  settings.size_y = args.m;
-  settings.size_a = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {0, 1, 2};
-  settings.outputs = {1};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1};
-  settings.local_size_ref = {64};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"WGS"+std::to_string(V)}};
-  settings.div_global = (V==1 || V==2) ? TransformVector{{"WPT"+std::to_string(V)}} : TransformVector{};
-
-  // Sets the tuning parameters and their possible values
-  if (V==1) {
-    settings.parameters = {
-      {"WGS"+std::to_string(V), {32, 64, 128, 256}},
-      {"WPT"+std::to_string(V), {1, 2, 4}},
-    };
-  }
-  if (V==2) {
-    settings.parameters = {
-      {"WGS"+std::to_string(V), {16, 32, 64, 128, 256}},
-      {"WPT"+std::to_string(V), {1, 2, 4}},
-      {"VW"+std::to_string(V), {1, 2, 4, 8}},
-    };
-  }
-  if (V==3) {
-    settings.parameters = {
-      {"WGS"+std::to_string(V), {16, 32, 64, 128}},
-      {"WPT"+std::to_string(V), {1, 2, 4, 8, 16, 32}},
-      {"VW"+std::to_string(V), {1, 2, 4, 8}},
-    };
-  }
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = (args.m*args.n + 2*args.m + args.n) * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int V) {
-  auto constraints = std::vector<Constraint>();
-  if (V==2 || V==3) {
-    auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
-    constraints.push_back({MultipleOfX, {"WPT"+std::to_string(V), "VW"+std::to_string(V)}});
-  }
-  if (V==3) {
-    auto LargerOrEqual = [] (std::vector<size_t> v) { return v[0] >= v[1]; };
-    constraints.push_back({LargerOrEqual, {"WGS"+std::to_string(V), "WPT"+std::to_string(V)}});
-  }
-  return constraints;
-}
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int V, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  auto a_rotated = (V==3) ? 1 : 0;
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, GetRealArg(args.alpha));
-  kernel.SetArgument(3, GetRealArg(args.beta));
-  kernel.SetArgument(4, a_rotated);
-  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(6, 0);
-  kernel.SetArgument(7, static_cast<int>(args.m));
-  kernel.SetArgument(8, buffers[0]()); // 0 == X vector
-  kernel.SetArgument(9, 0);
-  kernel.SetArgument(10, 1);
-  kernel.SetArgument(11, buffers[1]()); // 1 == Y vector
-  kernel.SetArgument(12, 0);
-  kernel.SetArgument(13, 1);
-  kernel.SetArgument(14, 0); // Conjugate transpose
-  kernel.SetArgument(15, 0); // Additional parameter
-  kernel.SetArgument(16, 0); // Banded 'kl'
-  kernel.SetArgument(17, 0); // Banded 'ku'
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xgemv.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -147,11 +23,11 @@ template <int V>
 void StartVariation(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, V, clblast::XgemvGetTunerDefaults, clblast::XgemvGetTunerSettings<half>, clblast::XgemvTestValidArguments<half>, clblast::XgemvSetConstraints, clblast::XgemvSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, V, clblast::XgemvGetTunerDefaults, clblast::XgemvGetTunerSettings<float>, clblast::XgemvTestValidArguments<float>, clblast::XgemvSetConstraints, clblast::XgemvSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, V, clblast::XgemvGetTunerDefaults, clblast::XgemvGetTunerSettings<double>, clblast::XgemvTestValidArguments<double>, clblast::XgemvSetConstraints, clblast::XgemvSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, V, clblast::XgemvGetTunerDefaults, clblast::XgemvGetTunerSettings<float2>, clblast::XgemvTestValidArguments<float2>, clblast::XgemvSetConstraints, clblast::XgemvSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, V, clblast::XgemvGetTunerDefaults, clblast::XgemvGetTunerSettings<double2>, clblast::XgemvTestValidArguments<double2>, clblast::XgemvSetConstraints, clblast::XgemvSetArguments<double2>); break;
   }
 }
 
diff --git a/src/tuning/kernels/xgemv.hpp b/src/tuning/kernels/xgemv.hpp
new file mode 100644
index 00000000..e44efe32
--- /dev/null
+++ b/src/tuning/kernels/xgemv.hpp
@@ -0,0 +1,138 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the xgemv OpenCL kernels. Three variants are tuned:
+// 1: The full version of the kernel
+// 2: The fast version for non-transposed matrices
+// 3: The fast version for transposed matrices
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XgemvGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha, kArgBeta};
+  settings.default_m = 2048;
+  settings.default_n = 2048;
+  settings.default_num_runs = 4;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XgemvGetTunerSettings(const int V, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = (V==1) ? "xgemv" : ((V==2) ? "xgemv_fast" : "xgemv_fast_rot");
+  settings.kernel_name = (V==1) ? "Xgemv" : ((V==2) ? "XgemvFast" : "XgemvFastRot");
+  settings.sources =
+#include "../src/kernels/level2/xgemv.opencl"
+#include "../src/kernels/level2/xgemv_fast.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_x = args.n;
+  settings.size_y = args.m;
+  settings.size_a = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {0, 1, 2};
+  settings.outputs = {1};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1};
+  settings.local_size_ref = {64};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"WGS"+std::to_string(V)}};
+  settings.div_global = (V==1 || V==2) ? TransformVector{{"WPT"+std::to_string(V)}} : TransformVector{};
+
+  // Sets the tuning parameters and their possible values
+  if (V==1) {
+    settings.parameters = {
+      {"WGS"+std::to_string(V), {32, 64, 128, 256}},
+      {"WPT"+std::to_string(V), {1, 2, 4}},
+    };
+  }
+  if (V==2) {
+    settings.parameters = {
+      {"WGS"+std::to_string(V), {16, 32, 64, 128, 256}},
+      {"WPT"+std::to_string(V), {1, 2, 4}},
+      {"VW"+std::to_string(V), {1, 2, 4, 8}},
+    };
+  }
+  if (V==3) {
+    settings.parameters = {
+      {"WGS"+std::to_string(V), {16, 32, 64, 128}},
+      {"WPT"+std::to_string(V), {1, 2, 4, 8, 16, 32}},
+      {"VW"+std::to_string(V), {1, 2, 4, 8}},
+    };
+  }
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = (args.m*args.n + 2*args.m + args.n) * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XgemvTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> XgemvSetConstraints(const int V) {
+  auto constraints = std::vector<Constraint>();
+  if (V==2 || V==3) {
+    auto MultipleOfX = [] (std::vector<size_t> v) { return IsMultiple(v[0], v[1]); };
+    constraints.push_back({MultipleOfX, {"WPT"+std::to_string(V), "VW"+std::to_string(V)}});
+  }
+  if (V==3) {
+    auto LargerOrEqual = [] (std::vector<size_t> v) { return v[0] >= v[1]; };
+    constraints.push_back({LargerOrEqual, {"WGS"+std::to_string(V), "WPT"+std::to_string(V)}});
+  }
+  return constraints;
+}
+
+// Sets the kernel's arguments
+template <typename T>
+void XgemvSetArguments(const int V, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  auto a_rotated = (V==3) ? 1 : 0;
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, GetRealArg(args.alpha));
+  kernel.SetArgument(3, GetRealArg(args.beta));
+  kernel.SetArgument(4, a_rotated);
+  kernel.SetArgument(5, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(6, 0);
+  kernel.SetArgument(7, static_cast<int>(args.m));
+  kernel.SetArgument(8, buffers[0]()); // 0 == X vector
+  kernel.SetArgument(9, 0);
+  kernel.SetArgument(10, 1);
+  kernel.SetArgument(11, buffers[1]()); // 1 == Y vector
+  kernel.SetArgument(12, 0);
+  kernel.SetArgument(13, 1);
+  kernel.SetArgument(14, 0); // Conjugate transpose
+  kernel.SetArgument(15, 0); // Additional parameter
+  kernel.SetArgument(16, 0); // Banded 'kl'
+  kernel.SetArgument(17, 0); // Banded 'ku'
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/kernels/xger.cpp b/src/tuning/kernels/xger.cpp
index a88fb5d6..6dfc9ffa 100644
--- a/src/tuning/kernels/xger.cpp
+++ b/src/tuning/kernels/xger.cpp
@@ -11,95 +11,7 @@
 //
 // =================================================================================================
 
-#include <string>
-#include <vector>
-
-#include "utilities/utilities.hpp"
-#include "tuning/tuning.hpp"
-
-namespace clblast {
-// =================================================================================================
-
-// Settings for this kernel (default command-line arguments)
-TunerDefaults GetTunerDefaults(const int) {
-  auto settings = TunerDefaults();
-  settings.options = {kArgM, kArgN, kArgAlpha};
-  settings.default_m = 1024;
-  settings.default_n = 1024;
-  return settings;
-}
-
-// Settings for this kernel (general)
-template <typename T>
-TunerSettings GetTunerSettings(const int, const Arguments<T> &args) {
-  auto settings = TunerSettings();
-
-  // Identification of the kernel
-  settings.kernel_family = "xger";
-  settings.kernel_name = "Xger";
-  settings.sources =
-#include "../src/kernels/level2/level2.opencl"
-#include "../src/kernels/level2/xger.opencl"
-  ;
-
-  // Buffer sizes
-  settings.size_x = args.m;
-  settings.size_y = args.n;
-  settings.size_a = args.m * args.n;
-
-  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
-  settings.inputs = {0, 1, 2};
-  settings.outputs = {2};
-
-  // Sets the base thread configuration
-  settings.global_size = {args.m, args.n};
-  settings.global_size_ref = settings.global_size;
-  settings.local_size = {1, 1};
-  settings.local_size_ref = {8, 8};
-
-  // Transforms the thread configuration based on the parameters
-  settings.mul_local = {{"WGS1", "WGS2"}};
-  settings.div_global = {{"WPT", "WPT"}};
-
-  // Sets the tuning parameters and their possible values
-  settings.parameters = {
-    {"WGS1", {4, 8, 16, 32, 64, 128, 256, 512}},
-    {"WGS2", {1, 2, 4, 8, 16, 32, 64, 128, 256}},
-    {"WPT", {1, 2, 4}},
-  };
-
-  // Describes how to compute the performance metrics
-  settings.metric_amount = (2*args.m*args.n + args.m + args.n) * GetBytes(args.precision);
-  settings.performance_unit = "GB/s";
-
-  return settings;
-}
-
-// Tests for valid arguments
-template <typename T>
-void TestValidArguments(const int, const Arguments<T> &) { }
-std::vector<Constraint> SetConstraints(const int) { return {}; }
-
-// Sets the kernel's arguments
-template <typename T>
-void SetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
-  kernel.SetArgument(0, static_cast<int>(args.m));
-  kernel.SetArgument(1, static_cast<int>(args.n));
-  kernel.SetArgument(2, GetRealArg(args.alpha));
-  kernel.SetArgument(3, buffers[0]()); // 0 == X vector
-  kernel.SetArgument(4, 0); // x_offset
-  kernel.SetArgument(5, 1); // x_increment
-  kernel.SetArgument(6, buffers[1]()); // 1 == Y vector
-  kernel.SetArgument(7, 0); // y_offset
-  kernel.SetArgument(8, 1); // y_increment
-  kernel.SetArgument(9, buffers[2]()); // 2 == A matrix
-  kernel.SetArgument(10, 0); // a_offset
-  kernel.SetArgument(11, static_cast<int>(args.m)); // a_ld
-  kernel.SetArgument(12, 0); // a_is_rowmajor
-}
-
-// =================================================================================================
-} // namespace clblast
+#include "tuning/kernels/xger.hpp"
 
 // Shortcuts to the clblast namespace
 using half = clblast::half;
@@ -110,11 +22,11 @@ using double2 = clblast::double2;
 int main(int argc, char *argv[]) {
   const auto command_line_args = clblast::RetrieveCommandLineArguments(argc, argv);
   switch(clblast::GetPrecision(command_line_args)) {
-    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<half>, clblast::TestValidArguments<half>, clblast::SetConstraints, clblast::SetArguments<half>); break;
-    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float>, clblast::TestValidArguments<float>, clblast::SetConstraints, clblast::SetArguments<float>); break;
-    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double>, clblast::TestValidArguments<double>, clblast::SetConstraints, clblast::SetArguments<double>); break;
-    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<float2>, clblast::TestValidArguments<float2>, clblast::SetConstraints, clblast::SetArguments<float2>); break;
-    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::GetTunerDefaults, clblast::GetTunerSettings<double2>, clblast::TestValidArguments<double2>, clblast::SetConstraints, clblast::SetArguments<double2>); break;
+    case clblast::Precision::kHalf: clblast::Tuner<half>(argc, argv, 0, clblast::XgerGetTunerDefaults, clblast::XgerGetTunerSettings<half>, clblast::XgerTestValidArguments<half>, clblast::XgerSetConstraints, clblast::XgerSetArguments<half>); break;
+    case clblast::Precision::kSingle: clblast::Tuner<float>(argc, argv, 0, clblast::XgerGetTunerDefaults, clblast::XgerGetTunerSettings<float>, clblast::XgerTestValidArguments<float>, clblast::XgerSetConstraints, clblast::XgerSetArguments<float>); break;
+    case clblast::Precision::kDouble: clblast::Tuner<double>(argc, argv, 0, clblast::XgerGetTunerDefaults, clblast::XgerGetTunerSettings<double>, clblast::XgerTestValidArguments<double>, clblast::XgerSetConstraints, clblast::XgerSetArguments<double>); break;
+    case clblast::Precision::kComplexSingle: clblast::Tuner<float2>(argc, argv, 0, clblast::XgerGetTunerDefaults, clblast::XgerGetTunerSettings<float2>, clblast::XgerTestValidArguments<float2>, clblast::XgerSetConstraints, clblast::XgerSetArguments<float2>); break;
+    case clblast::Precision::kComplexDouble: clblast::Tuner<double2>(argc, argv, 0, clblast::XgerGetTunerDefaults, clblast::XgerGetTunerSettings<double2>, clblast::XgerTestValidArguments<double2>, clblast::XgerSetConstraints, clblast::XgerSetArguments<double2>); break;
   }
   return 0;
 }
diff --git a/src/tuning/kernels/xger.hpp b/src/tuning/kernels/xger.hpp
new file mode 100644
index 00000000..afd2f36e
--- /dev/null
+++ b/src/tuning/kernels/xger.hpp
@@ -0,0 +1,102 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file uses the auto-tuner to tune the xger OpenCL kernels.
+//
+// =================================================================================================
+
+#include <string>
+#include <vector>
+
+#include "utilities/utilities.hpp"
+#include "tuning/tuning.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+// Settings for this kernel (default command-line arguments)
+TunerDefaults XgerGetTunerDefaults(const int) {
+  auto settings = TunerDefaults();
+  settings.options = {kArgM, kArgN, kArgAlpha};
+  settings.default_m = 1024;
+  settings.default_n = 1024;
+  return settings;
+}
+
+// Settings for this kernel (general)
+template <typename T>
+TunerSettings XgerGetTunerSettings(const int, const Arguments<T> &args) {
+  auto settings = TunerSettings();
+
+  // Identification of the kernel
+  settings.kernel_family = "xger";
+  settings.kernel_name = "Xger";
+  settings.sources =
+#include "../src/kernels/level2/level2.opencl"
+#include "../src/kernels/level2/xger.opencl"
+  ;
+
+  // Buffer sizes
+  settings.size_x = args.m;
+  settings.size_y = args.n;
+  settings.size_a = args.m * args.n;
+
+  // Inputs and outputs IDs (X:0, Y:1, A:2, B:3, C:4, temp:5)
+  settings.inputs = {0, 1, 2};
+  settings.outputs = {2};
+
+  // Sets the base thread configuration
+  settings.global_size = {args.m, args.n};
+  settings.global_size_ref = settings.global_size;
+  settings.local_size = {1, 1};
+  settings.local_size_ref = {8, 8};
+
+  // Transforms the thread configuration based on the parameters
+  settings.mul_local = {{"WGS1", "WGS2"}};
+  settings.div_global = {{"WPT", "WPT"}};
+
+  // Sets the tuning parameters and their possible values
+  settings.parameters = {
+    {"WGS1", {4, 8, 16, 32, 64, 128, 256, 512}},
+    {"WGS2", {1, 2, 4, 8, 16, 32, 64, 128, 256}},
+    {"WPT", {1, 2, 4}},
+  };
+
+  // Describes how to compute the performance metrics
+  settings.metric_amount = (2*args.m*args.n + args.m + args.n) * GetBytes(args.precision);
+  settings.performance_unit = "GB/s";
+
+  return settings;
+}
+
+// Tests for valid arguments
+template <typename T>
+void XgerTestValidArguments(const int, const Arguments<T> &) { }
+std::vector<Constraint> XgerSetConstraints(const int) { return {}; }
+
+// Sets the kernel's arguments
+template <typename T>
+void XgerSetArguments(const int, Kernel &kernel, const Arguments<T> &args, std::vector<Buffer<T>>& buffers) {
+  kernel.SetArgument(0, static_cast<int>(args.m));
+  kernel.SetArgument(1, static_cast<int>(args.n));
+  kernel.SetArgument(2, GetRealArg(args.alpha));
+  kernel.SetArgument(3, buffers[0]()); // 0 == X vector
+  kernel.SetArgument(4, 0); // x_offset
+  kernel.SetArgument(5, 1); // x_increment
+  kernel.SetArgument(6, buffers[1]()); // 1 == Y vector
+  kernel.SetArgument(7, 0); // y_offset
+  kernel.SetArgument(8, 1); // y_increment
+  kernel.SetArgument(9, buffers[2]()); // 2 == A matrix
+  kernel.SetArgument(10, 0); // a_offset
+  kernel.SetArgument(11, static_cast<int>(args.m)); // a_ld
+  kernel.SetArgument(12, 0); // a_is_rowmajor
+}
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/tuning/tuning.hpp b/src/tuning/tuning.hpp
index ee7e0087..cbecc300 100644
--- a/src/tuning/tuning.hpp
+++ b/src/tuning/tuning.hpp
@@ -121,6 +121,16 @@ void Tuner(int argc, char* argv[], const int V,
            SetConstraintsFunc SetConstraints,
            SetArgumentsFunc<T> SetArguments);
 
+// Function to run the tuners through the CLBlast API, no I/O
+template <typename T>
+StatusCode TunerAPI(Queue &queue, const Arguments<T> &args, const int V,
+                    const GetTunerDefaultsFunc GetTunerDefaults,
+                    const GetTunerSettingsFunc<T> GetTunerSettings,
+                    const TestValidArgumentsFunc<T> TestValidArguments,
+                    const SetConstraintsFunc SetConstraints,
+                    const SetArgumentsFunc<T> SetArguments,
+                    std::unordered_map<std::string,size_t> &parameters);
+
 // =================================================================================================
 } // namespace clblast
 
diff --git a/src/tuning/tuning_api.cpp b/src/tuning/tuning_api.cpp
new file mode 100644
index 00000000..d03b428c
--- /dev/null
+++ b/src/tuning/tuning_api.cpp
@@ -0,0 +1,385 @@
+
+// =================================================================================================
+// This file is part of the CLBlast project. The project is licensed under Apache Version 2.0. This
+// project loosely follows the Google C++ styleguide and uses a tab-size of two spaces and a max-
+// width of 100 characters per line.
+//
+// Author(s):
+//   Cedric Nugteren <www.cedricnugteren.nl>
+//
+// This file implements the parameter configurations for the CLBlast auto-tuner (taken from CLTune).
+// This is only used for the optional tuner binaries and not part of the core of CLBlast.
+//
+// =================================================================================================
+
+#include <vector>
+#include <string>
+#include <random>
+#include <utility>
+#include <algorithm>
+
+#include "tuning/tuning.hpp"
+#include "tuning/kernels/xaxpy.hpp"
+#include "tuning/kernels/xdot.hpp"
+#include "tuning/kernels/xgemv.hpp"
+#include "tuning/kernels/xger.hpp"
+#include "tuning/kernels/xgemm.hpp"
+#include "tuning/kernels/xgemm_direct.hpp"
+#include "tuning/kernels/copy_fast.hpp"
+#include "tuning/kernels/copy_pad.hpp"
+#include "tuning/kernels/transpose_fast.hpp"
+#include "tuning/kernels/transpose_pad.hpp"
+#include "tuning/kernels/invert.hpp"
+
+namespace clblast {
+// =================================================================================================
+
+template <typename T>
+StatusCode TuneXaxpy(RawCommandQueue * queue, const size_t n,
+                     const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, XaxpyGetTunerDefaults, XaxpyGetTunerSettings<T>,
+                     XaxpyTestValidArguments<T>, XaxpySetConstraints, XaxpySetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXaxpy<half>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXaxpy<float>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXaxpy<double>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXaxpy<float2>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXaxpy<double2>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneXdot(RawCommandQueue * queue, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  auto status = TunerAPI<T>(queue_cpp, args, 1, XdotGetTunerDefaults, XdotGetTunerSettings<T>,
+                            XdotTestValidArguments<T>, XdotSetConstraints, XdotSetArguments<T>, parameters);
+  if (status != StatusCode::kSuccess) { return status; }
+  return TunerAPI<T>(queue_cpp, args, 2, XdotGetTunerDefaults, XdotGetTunerSettings<T>,
+                     XdotTestValidArguments<T>, XdotSetConstraints, XdotSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXdot<half>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXdot<float>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXdot<double>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXdot<float2>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXdot<double2>(RawCommandQueue*, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneXgemv(RawCommandQueue * queue, const size_t m, const size_t n,
+                     const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  auto status = TunerAPI<T>(queue_cpp, args, 1, XgemvGetTunerDefaults, XgemvGetTunerSettings<T>,
+                            XgemvTestValidArguments<T>, XgemvSetConstraints, XgemvSetArguments<T>, parameters);
+  if (status != StatusCode::kSuccess) { return status; }
+  status = TunerAPI<T>(queue_cpp, args, 2, XgemvGetTunerDefaults, XgemvGetTunerSettings<T>,
+                       XgemvTestValidArguments<T>, XgemvSetConstraints, XgemvSetArguments<T>, parameters);
+  if (status != StatusCode::kSuccess) { return status; }
+  return TunerAPI<T>(queue_cpp, args, 3, XgemvGetTunerDefaults, XgemvGetTunerSettings<T>,
+                     XgemvTestValidArguments<T>, XgemvSetConstraints, XgemvSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXgemv<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemv<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemv<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemv<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemv<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneXger(RawCommandQueue * queue, const size_t m, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, XgerGetTunerDefaults, XgerGetTunerSettings<T>,
+                     XgerTestValidArguments<T>, XgerSetConstraints, XgerSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXger<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXger<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXger<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXger<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXger<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneXgemm(RawCommandQueue * queue, const size_t m, const size_t n, const size_t k,
+                     const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n; args.k = k;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 2, XgemmGetTunerDefaults, XgemmGetTunerSettings<T>,
+                     XgemmTestValidArguments<T>, XgemmSetConstraints, XgemmSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXgemm<half>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemm<float>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemm<double>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemm<float2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemm<double2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneXgemmDirect(RawCommandQueue * queue, const size_t m, const size_t n, const size_t k,
+                           const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n; args.k = k;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 2, XgemmDirectGetTunerDefaults, XgemmDirectGetTunerSettings<T>,
+                     XgemmDirectTestValidArguments<T>, XgemmDirectSetConstraints, XgemmDirectSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneXgemmDirect<half>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemmDirect<float>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemmDirect<double>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemmDirect<float2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneXgemmDirect<double2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneCopy(RawCommandQueue * queue, const size_t m, const size_t n,
+                    const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, CopyGetTunerDefaults, CopyGetTunerSettings<T>,
+                     CopyTestValidArguments<T>, CopySetConstraints, CopySetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneCopy<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneCopy<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneCopy<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneCopy<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneCopy<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TunePad(RawCommandQueue * queue, const size_t m, const size_t n,
+                   const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, PadGetTunerDefaults, PadGetTunerSettings<T>,
+                     PadTestValidArguments<T>, PadSetConstraints, PadSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TunePad<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePad<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePad<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePad<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePad<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneTranspose(RawCommandQueue * queue, const size_t m, const size_t n,
+                         const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, TransposeGetTunerDefaults, TransposeGetTunerSettings<T>,
+                     TransposeTestValidArguments<T>, TransposeSetConstraints, TransposeSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneTranspose<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneTranspose<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneTranspose<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneTranspose<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneTranspose<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TunePadtranspose(RawCommandQueue * queue, const size_t m, const size_t n,
+                            const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, PadtransposeGetTunerDefaults, PadtransposeGetTunerSettings<T>,
+                     PadtransposeTestValidArguments<T>, PadtransposeSetConstraints, PadtransposeSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TunePadtranspose<half>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePadtranspose<float>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePadtranspose<double>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePadtranspose<float2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TunePadtranspose<double2>(RawCommandQueue*, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+template <typename T>
+StatusCode TuneInvert(RawCommandQueue * queue, const size_t m, const size_t n, const size_t k,
+                      const double fraction, std::unordered_map<std::string,size_t> &parameters) {
+  auto args = Arguments<T>(); args.fraction = fraction; args.m = m; args.n = n; args.k = k;
+  auto queue_cpp = Queue(*queue);
+  return TunerAPI<T>(queue_cpp, args, 0, InvertGetTunerDefaults, InvertGetTunerSettings<T>,
+                     InvertTestValidArguments<T>, InvertSetConstraints, InvertSetArguments<T>, parameters);
+}
+template StatusCode PUBLIC_API TuneInvert<half>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneInvert<float>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneInvert<double>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneInvert<float2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+template StatusCode PUBLIC_API TuneInvert<double2>(RawCommandQueue*, const size_t, const size_t, const size_t, const double, std::unordered_map<std::string,size_t>&);
+
+// =================================================================================================
+
+// The main tuner API, similar to the one in tuning.cpp, but without I/O
+template <typename T>
+StatusCode TunerAPI(Queue &queue, const Arguments<T> &args, const int V,
+                    const GetTunerDefaultsFunc GetTunerDefaults,
+                    const GetTunerSettingsFunc<T> GetTunerSettings,
+                    const TestValidArgumentsFunc<T> TestValidArguments,
+                    const SetConstraintsFunc SetConstraints,
+                    const SetArgumentsFunc<T> SetArguments,
+                    std::unordered_map<std::string,size_t> &parameters) {
+
+  // Sets the parameters and platform/device for which to tune (command-line options)
+  const TunerDefaults defaults = GetTunerDefaults(V);
+  const TunerSettings settings = GetTunerSettings(V, args);
+
+  // Tests validity of the given arguments
+  TestValidArguments(V, args);
+
+  // Retrieves OpenCL classes
+  const auto device = queue.GetDevice();
+  const auto context = queue.GetContext();
+
+  // Inspects whether or not FP64 is supported in case of double precision
+  if ((PrecisionValue<T>() == Precision::kDouble && !PrecisionSupported<double>(device)) ||
+      (PrecisionValue<T>() == Precision::kComplexDouble && !PrecisionSupported<double2>(device))) {
+    return StatusCode::kNoDoublePrecision;
+  }
+
+  // As above, but for FP16 (half precision)
+  if (PrecisionValue<T>() == Precision::kHalf && !PrecisionSupported<half>(device)) {
+    return StatusCode::kNoHalfPrecision;
+  }
+
+  // Retrieves properties
+  const auto device_type = GetDeviceType(device);
+  const auto device_vendor = GetDeviceVendor(device);
+  const auto device_architecture = GetDeviceArchitecture(device);
+  const auto device_name = GetDeviceName(device);
+
+  // Creates input buffers with random data
+  const auto buffer_sizes = std::vector<size_t>{
+      settings.size_x, settings.size_y,
+      settings.size_a, settings.size_b, settings.size_c,
+      settings.size_temp
+  };
+  const auto seed = static_cast<unsigned long>(time(nullptr));
+  std::mt19937 mt(seed);
+  std::uniform_real_distribution<double> dist(kTestDataLowerLimit, kTestDataUpperLimit);
+  auto source_buffers = std::vector<std::vector<T>>();
+  auto reference_buffers = std::vector<std::vector<T>>();
+  auto result_buffers = std::vector<std::vector<T>>();
+  auto device_buffers = std::vector<Buffer<T>>();
+  for (const auto size : buffer_sizes) {
+    auto host_buffer = std::vector<T>(size);
+    PopulateVector(host_buffer, mt, dist);
+    source_buffers.push_back(host_buffer);
+    reference_buffers.push_back(std::vector<T>(size));
+    result_buffers.push_back(std::vector<T>(size));
+    device_buffers.push_back(Buffer<T>(context, size));
+  }
+
+  // Sets the tunable parameters and their possible values
+  auto configurations = SetConfigurations(settings.parameters, SetConstraints(V));
+
+  // Select the search method (full search or a random fraction)
+  if (args.fraction != 0.0 && args.fraction != 1.0) {
+    const auto new_size = static_cast<size_t>(configurations.size() * args.fraction);
+    auto rng = std::default_random_engine{};
+    std::shuffle(std::begin(configurations), std::end(configurations), rng);
+    configurations.resize(new_size);
+  }
+
+  // First runs a reference example to compare against
+  try {
+
+    // Sets the input
+    for (const auto id : settings.inputs) {
+      device_buffers[id].Write(queue, buffer_sizes[id], source_buffers[id]);
+    }
+
+    // Compiles the kernel
+    auto compiler_options = std::vector<std::string>();
+    const auto program = CompileFromSource(settings.sources, args.precision, settings.kernel_name,
+                                           device, context, compiler_options, 0);
+    auto kernel = Kernel(program, settings.kernel_name);
+    SetArguments(V, kernel, args, device_buffers);
+
+    // Runs the kernel
+    const auto time_ms = TimeKernel(args.num_runs, kernel, queue, device,
+                                    settings.global_size_ref, settings.local_size_ref, true);
+    if (time_ms == -1.0) { throw std::runtime_error("Error in reference implementation"); }
+
+    // Saves the result
+    for (const auto id : settings.outputs) {
+      device_buffers[id].Read(queue, buffer_sizes[id], reference_buffers[id]);
+    }
+  }
+  catch (...) {
+    const auto status_code = DispatchExceptionCatchAll(true);
+    return status_code;
+  }
+
+  // Starts the tuning process
+  auto results = std::vector<TuningResult>();
+  for (auto config_id = size_t{0}; config_id < configurations.size(); ++config_id) {
+    try {
+      auto configuration = configurations[config_id];
+
+      // Sets the input
+      for (const auto id : settings.inputs) {
+        device_buffers[id].Write(queue, buffer_sizes[id], source_buffers[id]);
+      }
+
+      // Sets the thread configuration
+      const auto global = SetThreadConfiguration(configuration, settings.global_size,
+                                                 settings.mul_global, settings.div_global);
+      const auto local = SetThreadConfiguration(configuration, settings.local_size,
+                                                settings.mul_local, settings.div_local);
+
+      // Sets the parameters for this configuration
+      auto kernel_source = std::string{""};
+      for (const auto &parameter : configuration) {
+        kernel_source += "#define " + parameter.first + " " + ToString(parameter.second) + "\n";
+      }
+      kernel_source += settings.sources;
+
+      // Compiles the kernel
+      auto compiler_options = std::vector<std::string>();
+      const auto program = CompileFromSource(kernel_source, args.precision, settings.kernel_name,
+                                             device, context, compiler_options, 0, true);
+      auto kernel = Kernel(program, settings.kernel_name);
+
+      // Runs the kernel
+      SetArguments(V, kernel, args, device_buffers);
+      const auto time_ms = TimeKernel(args.num_runs, kernel, queue, device, global, local, true);
+
+      // Kernel run was not successful
+      if (time_ms == -1.0) {
+        continue;
+      }
+
+      // Compares the results
+      auto l2_error = 0.0;
+      for (const auto id : settings.outputs) {
+        device_buffers[id].Read(queue, buffer_sizes[id], result_buffers[id]);
+        for (auto index = size_t{0}; index<buffer_sizes[id]; ++index) {
+          const auto diff = SquaredDifference(result_buffers[id][index], reference_buffers[id][index]);
+          l2_error += diff;
+        }
+        l2_error /= static_cast<double>(buffer_sizes[id]);
+        if (std::isnan(l2_error) || l2_error > 1.0e-4) {
+          throw std::runtime_error("L2 error too large");
+        }
+      }
+      results.push_back(TuningResult{settings.kernel_name, time_ms, configuration});
+    }
+    catch (...) {
+    }
+  }
+
+  // Completed the tuning process
+  if (results.size() == 0) { return StatusCode::kUnexpectedError; }
+
+  // Computes the best results
+  auto comparison = [](const TuningResult& lhs, const TuningResult& rhs) { return lhs.score < rhs.score; };
+  const auto best_configuration = std::min_element(results.begin(), results.end(), comparison);
+  const auto best_time_ms = best_configuration->score;
+  if (best_time_ms == 0.0) { return StatusCode::kUnexpectedError; }
+
+  // Stores the best parameters
+  for (const auto config : best_configuration->config) {
+    parameters[config.first] = config.second;
+  }
+  return StatusCode::kSuccess;
+}
+
+// Compiles the above function
+template StatusCode TunerAPI<half>(Queue &queue, const Arguments<half> &args, const int V, const GetTunerDefaultsFunc GetTunerDefaults, const GetTunerSettingsFunc<half> GetTunerSettings, const TestValidArgumentsFunc<half> TestValidArguments, const SetConstraintsFunc SetConstraints, const SetArgumentsFunc<half> SetArguments, std::unordered_map<std::string,size_t>&);
+template StatusCode TunerAPI<float>(Queue &queue, const Arguments<float> &args, const int V, const GetTunerDefaultsFunc GetTunerDefaults, const GetTunerSettingsFunc<float> GetTunerSettings, const TestValidArgumentsFunc<float> TestValidArguments, const SetConstraintsFunc SetConstraints, const SetArgumentsFunc<float> SetArguments, std::unordered_map<std::string,size_t>&);
+template StatusCode TunerAPI<double>(Queue &queue, const Arguments<double> &args, const int V, const GetTunerDefaultsFunc GetTunerDefaults, const GetTunerSettingsFunc<double> GetTunerSettings, const TestValidArgumentsFunc<double> TestValidArguments, const SetConstraintsFunc SetConstraints, const SetArgumentsFunc<double> SetArguments, std::unordered_map<std::string,size_t>&);
+template StatusCode TunerAPI<float2>(Queue &queue, const Arguments<float2> &args, const int V, const GetTunerDefaultsFunc GetTunerDefaults, const GetTunerSettingsFunc<float2> GetTunerSettings, const TestValidArgumentsFunc<float2> TestValidArguments, const SetConstraintsFunc SetConstraints, const SetArgumentsFunc<float2> SetArguments, std::unordered_map<std::string,size_t>&);
+template StatusCode TunerAPI<double2>(Queue &queue, const Arguments<double2> &args, const int V, const GetTunerDefaultsFunc GetTunerDefaults, const GetTunerSettingsFunc<double2> GetTunerSettings, const TestValidArgumentsFunc<double2> TestValidArguments, const SetConstraintsFunc SetConstraints, const SetArgumentsFunc<double2> SetArguments, std::unordered_map<std::string,size_t>&);
+
+// =================================================================================================
+} // namespace clblast
diff --git a/src/utilities/timing.cpp b/src/utilities/timing.cpp
index af6a8ff2..1afb0d08 100644
--- a/src/utilities/timing.cpp
+++ b/src/utilities/timing.cpp
@@ -62,15 +62,16 @@ double RunKernelTimed(const size_t num_runs, Kernel &kernel, Queue &queue, const
 }
 
 double TimeKernel(const size_t num_runs, Kernel &kernel, Queue &queue, const Device &device,
-                  std::vector<size_t> global, const std::vector<size_t> &local) {
+                  std::vector<size_t> global, const std::vector<size_t> &local,
+                  const bool silent) {
   try {
     const auto time_ms = RunKernelTimed(num_runs, kernel, queue, device, global, local);
-    printf(" %9.2lf ms |", time_ms);
+    if (!silent) { printf(" %9.2lf ms |", time_ms); }
     return time_ms;
   }
   catch (...) {
     const auto status_code = DispatchExceptionCatchAll(true);
-    printf("  error %-5d |", static_cast<int>(status_code));
+    if (!silent) { printf("  error %-5d |", static_cast<int>(status_code)); }
     return -1.0; // invalid
   }
 }
diff --git a/src/utilities/timing.hpp b/src/utilities/timing.hpp
index c167cd5f..7761fd83 100644
--- a/src/utilities/timing.hpp
+++ b/src/utilities/timing.hpp
@@ -44,7 +44,8 @@ double RunKernelTimed(const size_t num_runs, Kernel &kernel, Queue &queue, const
                       std::vector<size_t> global, const std::vector<size_t> &local);
 
 double TimeKernel(const size_t num_runs, Kernel &kernel, Queue &queue, const Device &device,
-                  std::vector<size_t> global, const std::vector<size_t> &local);
+                  std::vector<size_t> global, const std::vector<size_t> &local,
+                  const bool silent = false);
 
 // =================================================================================================