Compare TUV-X Tridiagonal Solver Code with LAPACKE's implementation

.dockerignore

@@ -13,4 +13,5 @@
 !packaging/
 !test/
 !tool/
-!CMakeLists.txt
+!CMakeLists.txt
+!benchmark

.github/workflows/docker.yml

@@ -1,16 +1,21 @@
 name: Docker
 
-on: [push, pull_request]
+on:
+  push:
+    branches:
+      - main
+  pull_request:
+  workflow_dispatch:
 
 concurrency:
   group: ${{ github.workflow }}-${{ github.ref || github.run_id }}
   cancel-in-progress: true
 
 jobs:
   docker-build-and-test:
-    if: github.event_name != 'pull_request' || github.event.pull_request.head.repo.full_name != github.event.pull_request.base.repo.full_name
     name: Build and Test - ${{ matrix.dockerfile }}
     runs-on: ubuntu-latest
+    continue-on-error: true
     strategy:
       fail-fast: false
       matrix:

.github/workflows/ubuntu.yml

@@ -1,6 +1,11 @@
 name: Ubuntu
 
-on: [ push, pull_request ]
+on: 
+  push:
+    branches:
+      - main
+  pull_request:
+  workflow_dispatch:
 
 concurrency:
   group: ${{ github.workflow }}-${{ github.event.pull_request.number || github.ref }}
@@ -9,7 +14,7 @@ concurrency:
 jobs:
   gcc:
     runs-on: ubuntu-24.04
-    if: github.event_name != 'pull_request' || github.event.pull_request.head.repo.full_name != github.event.pull_request.base.repo.full_name
+    continue-on-error: true
     strategy:
       matrix:
         gcc_version: [12, 13, 14]
@@ -24,13 +29,13 @@ jobs:
     - name: Install dependencies
       run: |
         sudo apt-get update
-        sudo apt-get install -y libnetcdf-dev netcdf-bin libnetcdff-dev liblapack-dev
+        sudo apt-get install -y libnetcdf-dev netcdf-bin libnetcdff-dev liblapack-dev liblapacke-dev
         sudo apt-get install -y python3-numpy python3-scipy
     - name: Run Cmake 
       run: cmake -S . -B build 
     - name: Build
-      run: cmake --build build --parallel
+      run: cmake --build build 
     - name: Run tests
       run: |
         cd build
-        ctest --rerun-failed --output-on-failure . --verbose
+        ctest --rerun-failed --output-on-failure . --verbose

CMakeLists.txt

@@ -37,6 +37,7 @@ include(CMakeDependentOption)
 option(TUVX_ENABLE_MPI "Enable MPI parallel support" OFF)
 cmake_dependent_option(TUVX_ENABLE_OPENMP "Enable OpenMP support" OFF "TUVX_ENABLE_MPI" OFF)
 option(TUVX_ENABLE_TESTS "Build tests" ON)
+option(TUVX_ENABLE_BENCHMARK "Build benchmark examples" OFF)
 option(TUVX_ENABLE_COVERAGE "Enable code coverage output" OFF)
 option(TUVX_ENABLE_MEMCHECK "Enable memory checking in tests" OFF)
 option(TUVX_ENABLE_NC_CONFIG "Use nc-config to determine NetCDF libraries" OFF)
@@ -105,6 +106,9 @@ endif()
 target_link_libraries(tuv-x 
   PUBLIC 
     musica::tuvx
+    yaml-cpp::yaml-cpp
+    ${BLAS_LIBRARIES}
+    ${LAPACK_LIBRARIES}
 )
 
 target_include_directories(tuv-x
@@ -151,4 +155,11 @@ if(PROJECT_IS_TOP_LEVEL)
   add_subdirectory(packaging)
 endif()
 
-################################################################################
+################################################################################
+# benchmarking
+
+if(TUVX_ENABLE_BENCHMARK)
+  add_subdirectory(benchmark)
+endif()
+
+################################################################################

benchmark/CMakeLists.txt

@@ -0,0 +1,12 @@
+add_executable(benchmark_tridiagonal_solver benchmark_tridiagonal_solver.cpp)
+target_include_directories(
+  benchmark_tridiagonal_solver PUBLIC ${OpenBLAS_INCLUDE_DIRS}
+                                      ${LAPACK_INCLUDE_DIRS})
+
+target_link_libraries(benchmark_tridiagonal_solver
+  PUBLIC 
+    LAPACK::LAPACK 
+    ${LAPACKE_LIBRARIES}
+    benchmark::benchmark 
+    musica::tuvx
+)

benchmark/benchmark_tridiagonal_solver.cpp

@@ -0,0 +1,126 @@
+
+#include <tuvx/linear_algebra/linear_algebra.hpp>
+
+#include <benchmark/benchmark.h>
+
+#include <random>
+#include <vector>
+
+#ifdef TUVX_COMPILE_WITH_INTEL
+  #include <mkl_lapacke.h>
+#elif TUVX_COMPILE_WITH_GCC
+  #include <lapacke.h>
+#endif
+
+const std::size_t SYSTEM_SIZE = 1e6;
+
+const bool MAKE_DIAGONALLY_DOMINANT = true;
+
+const unsigned RANDOM_NUMBER_SEED = 1;
+
+/// @brief This function benchmarks the lapacke tridiagonal matrix solver for single precision
+/// @param state Benchmarking argument
+static void BM_LAPACKE_SINGLE_PRECISISON(benchmark::State& state)
+{
+  std::vector<float> x(SYSTEM_SIZE);
+  std::vector<float> b(SYSTEM_SIZE);
+  tuvx::TridiagonalMatrix<float> A(SYSTEM_SIZE);
+  std::mt19937 random_device(RANDOM_NUMBER_SEED);
+  for (auto _ : state)
+  {
+    state.PauseTiming();
+    tuvx::FillRandom<float>(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT);
+    tuvx::FillRandom<float>(x, RANDOM_NUMBER_SEED);
+    b = tuvx::Dot<float>(A, x);
+    state.ResumeTiming();
+
+    LAPACKE_sgtsv(
+        LAPACK_ROW_MAJOR,
+        SYSTEM_SIZE,
+        1,
+        A.lower_diagonal_.data(),
+        A.main_diagonal_.data(),
+        A.upper_diagonal_.data(),
+        b.data(),
+        1);
+  }
+}
+
+/// @brief This function benchmarks the lapacke tridiagonal matrix solver for double precision
+/// @param state Benchmarking argument
+static void BM_LAPACKE_DOUBLE_PRECISISON(benchmark::State& state)
+{
+  std::vector<double> x(SYSTEM_SIZE);
+  std::vector<double> b(SYSTEM_SIZE);
+  tuvx::TridiagonalMatrix<double> A(SYSTEM_SIZE);
+
+  std::mt19937 random_device(RANDOM_NUMBER_SEED);
+  // Perform setup here
+  for (auto _ : state)
+  {
+    state.PauseTiming();
+    tuvx::FillRandom<double>(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT);
+    tuvx::FillRandom<double>(x, RANDOM_NUMBER_SEED);
+    b = tuvx::Dot<double>(A, x);
+    state.ResumeTiming();
+
+    LAPACKE_dgtsv(
+        LAPACK_ROW_MAJOR,
+        SYSTEM_SIZE,
+        1,
+        A.lower_diagonal_.data(),
+        A.main_diagonal_.data(),
+        A.upper_diagonal_.data(),
+        b.data(),
+        1);
+  }
+}
+
+/// @brief This function benchmarks the tuvx tridiagonal matrix solver for single precision
+/// @param state Benchmarking argument
+static void BM_TUVX_DOUBLE_PRECISISON(benchmark::State& state)
+{
+  std::vector<double> x(SYSTEM_SIZE);
+  std::vector<double> b(SYSTEM_SIZE);
+  tuvx::TridiagonalMatrix<double> A(SYSTEM_SIZE);
+  std::vector<double> x_approx(SYSTEM_SIZE);
+
+  for (auto _ : state)
+  {
+    state.PauseTiming();
+    tuvx::FillRandom<double>(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT);
+    tuvx::FillRandom<double>(x, RANDOM_NUMBER_SEED);
+    state.ResumeTiming();
+    tuvx::Solve<double>(A, b);
+  }
+}
+
+/// @brief This function benchmarks the tuvx tridiagonal matrix solver for double precision
+/// @param state Benchmarking argument
+static void BM_TUVX_SINGLE_PRECISISON(benchmark::State& state)
+{
+  std::vector<float> x(SYSTEM_SIZE);
+  std::vector<float> b(SYSTEM_SIZE);
+  tuvx::TridiagonalMatrix<float> A(SYSTEM_SIZE);
+  std::vector<float> x_approx(SYSTEM_SIZE);
+
+  // Perform setup here
+  for (auto _ : state)
+  {
+    state.PauseTiming();
+    tuvx::FillRandom<float>(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT);
+    tuvx::FillRandom<float>(x, RANDOM_NUMBER_SEED);
+    b = tuvx::Dot<float>(A, x);
+    state.ResumeTiming();
+    tuvx::Solve<float>(A, b);
+  }
+}
+
+/// @brief Register the functions defined above as a benchmark
+BENCHMARK(BM_LAPACKE_DOUBLE_PRECISISON);
+BENCHMARK(BM_LAPACKE_SINGLE_PRECISISON);
+BENCHMARK(BM_TUVX_DOUBLE_PRECISISON);
+BENCHMARK(BM_TUVX_SINGLE_PRECISISON);
+
+/// @brief Run all benchmarks
+BENCHMARK_MAIN();