diff --git a/.dockerignore b/.dockerignore index 69a3598b..5641c612 100644 --- a/.dockerignore +++ b/.dockerignore @@ -13,4 +13,5 @@ !packaging/ !test/ !tool/ -!CMakeLists.txt \ No newline at end of file +!CMakeLists.txt +!benchmark \ No newline at end of file diff --git a/.github/workflows/docker.yml b/.github/workflows/docker.yml index 15d7a717..de4873c7 100644 --- a/.github/workflows/docker.yml +++ b/.github/workflows/docker.yml @@ -1,6 +1,11 @@ name: Docker -on: [push, pull_request] +on: + push: + branches: + - main + pull_request: + workflow_dispatch: concurrency: group: ${{ github.workflow }}-${{ github.ref || github.run_id }} @@ -8,9 +13,9 @@ concurrency: 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: diff --git a/.github/workflows/ubuntu.yml b/.github/workflows/ubuntu.yml index f01589f7..8f7db69d 100644 --- a/.github/workflows/ubuntu.yml +++ b/.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 \ No newline at end of file + ctest --rerun-failed --output-on-failure . --verbose diff --git a/CMakeLists.txt b/CMakeLists.txt index e1af50cb..9fa1ab4d 100644 --- a/CMakeLists.txt +++ b/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() -################################################################################ \ No newline at end of file +################################################################################ +# benchmarking + +if(TUVX_ENABLE_BENCHMARK) + add_subdirectory(benchmark) +endif() + +################################################################################ diff --git a/benchmark/CMakeLists.txt b/benchmark/CMakeLists.txt new file mode 100644 index 00000000..97c1247a --- /dev/null +++ b/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 +) diff --git a/benchmark/benchmark_tridiagonal_solver.cpp b/benchmark/benchmark_tridiagonal_solver.cpp new file mode 100644 index 00000000..c0496262 --- /dev/null +++ b/benchmark/benchmark_tridiagonal_solver.cpp @@ -0,0 +1,126 @@ + +#include + +#include + +#include +#include + +#ifdef TUVX_COMPILE_WITH_INTEL + #include +#elif TUVX_COMPILE_WITH_GCC + #include +#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 x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + std::mt19937 random_device(RANDOM_NUMBER_SEED); + for (auto _ : state) + { + state.PauseTiming(); + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(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 x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + + std::mt19937 random_device(RANDOM_NUMBER_SEED); + // Perform setup here + for (auto _ : state) + { + state.PauseTiming(); + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(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 x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + std::vector x_approx(SYSTEM_SIZE); + + for (auto _ : state) + { + state.PauseTiming(); + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + state.ResumeTiming(); + tuvx::Solve(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 x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + std::vector x_approx(SYSTEM_SIZE); + + // Perform setup here + for (auto _ : state) + { + state.PauseTiming(); + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(A, x); + state.ResumeTiming(); + tuvx::Solve(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(); diff --git a/cmake/FindLAPACKE.cmake b/cmake/FindLAPACKE.cmake new file mode 100644 index 00000000..a53891eb --- /dev/null +++ b/cmake/FindLAPACKE.cmake @@ -0,0 +1,421 @@ +### +# +# @copyright (c) 2009-2014 The University of Tennessee and The University +# of Tennessee Research Foundation. +# All rights reserved. +# @copyright (c) 2012-2016 Inria. All rights reserved. +# @copyright (c) 2012-2014 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, Univ. Bordeaux. All rights reserved. +# +### +# +# - Find LAPACKE include dirs and libraries +# Use this module by invoking find_package with the form: +# find_package(LAPACKE +# [REQUIRED] # Fail with error if lapacke is not found +# [COMPONENTS ...] # dependencies +# ) +# +# LAPACKE depends on the following libraries: +# - LAPACK +# +# This module finds headers and lapacke library. +# Results are reported in variables: +# LAPACKE_FOUND - True if headers and requested libraries were found +# LAPACKE_LINKER_FLAGS - list of required linker flags (excluding -l and -L) +# LAPACKE_INCLUDE_DIRS - lapacke include directories +# LAPACKE_LIBRARY_DIRS - Link directories for lapacke libraries +# LAPACKE_LIBRARIES - lapacke component libraries to be linked +# LAPACKE_INCLUDE_DIRS_DEP - lapacke + dependencies include directories +# LAPACKE_LIBRARY_DIRS_DEP - lapacke + dependencies link directories +# LAPACKE_LIBRARIES_DEP - lapacke libraries + dependencies +# +# The user can give specific paths where to find the libraries adding cmake +# options at configure (ex: cmake path/to/project -DLAPACKE_DIR=path/to/lapacke): +# LAPACKE_DIR - Where to find the base directory of lapacke +# LAPACKE_INCDIR - Where to find the header files +# LAPACKE_LIBDIR - Where to find the library files +# The module can also look for the following environment variables if paths +# are not given as cmake variable: LAPACKE_DIR, LAPACKE_INCDIR, LAPACKE_LIBDIR +# +# If the static version of the LAPACKE libraries is required, please add the +# following in your CMakeLists.txt before calling find_package(LAPACKE): +# set(LAPACKE_STATIC TRUE) +# +# LAPACKE could be directly embedded in LAPACK library (ex: Intel MKL) so that +# we test a lapacke function with the lapack libraries found and set LAPACKE +# variables to LAPACK ones if test is successful. To skip this feature and +# look for a stand alone lapacke, please add the following in your +# CMakeLists.txt before to call find_package(LAPACKE): +# set(LAPACKE_STANDALONE TRUE) + +#============================================================================= +# Copyright 2012-2013 Inria +# Copyright 2012-2013 Emmanuel Agullo +# Copyright 2012-2013 Mathieu Faverge +# Copyright 2012 Cedric Castagnede +# Copyright 2013-2016 Florent Pruvost +# +# This file is part of a computer program whose purpose is to process +# Matrices Over Runtime Systems @ Exascale (MORSE). More information +# can be found on the following website: http://www.inria.fr/en/teams/morse. +# +# This software is governed by the CeCILL-C license under French law and +# abiding by the rules of distribution of free software. You can use, +# modify and/ or redistribute the software under the terms of the CeCILL-C +# license as circulated by CEA, CNRS and INRIA at the following URL +# "http://www.cecill.info". +# +# As a counterpart to the access to the source code and rights to copy, +# modify and redistribute granted by the license, users are provided only +# with a limited warranty and the software's author, the holder of the +# economic rights, and the successive licensors have only limited +# liability. +# +# In this respect, the user's attention is drawn to the risks associated +# with loading, using, modifying and/or developing or reproducing the +# software by the user in light of its specific status of free software, +# that may mean that it is complicated to manipulate, and that also +# therefore means that it is reserved for developers and experienced +# professionals having in-depth computer knowledge. Users are therefore +# encouraged to load and test the software's suitability as regards their +# requirements in conditions enabling the security of their systems and/or +# data to be ensured and, more generally, to use and operate it in the +# same conditions as regards security. +# +# The fact that you are presently reading this means that you have had +# knowledge of the CeCILL-C license and that you accept its terms. +# +# This software is distributed WITHOUT ANY WARRANTY; without even the +# implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. +# See the License for more information. +#============================================================================= + +if (NOT LAPACKE_FOUND) + set(LAPACKE_DIR "" CACHE PATH "Installation directory of LAPACKE library") + if (NOT LAPACKE_FIND_QUIETLY) + message(STATUS "A cache variable, namely LAPACKE_DIR, has been set to specify the install directory of LAPACKE") + endif() +endif() + +# LAPACKE depends on LAPACK anyway, try to find it +if (NOT LAPACK_FOUND) + if(LAPACKE_FIND_REQUIRED) + find_package(LAPACKEXT REQUIRED) + else() + find_package(LAPACKEXT) + endif() +endif() + +# LAPACKE depends on LAPACK +if (LAPACK_FOUND) + + if (NOT LAPACKE_STANDALONE) + # check if a lapacke function exists in the LAPACK lib + include(CheckFunctionExists) + set(CMAKE_REQUIRED_LIBRARIES "${LAPACK_LINKER_FLAGS};${LAPACK_LIBRARIES}") + unset(LAPACKE_WORKS CACHE) + check_function_exists(LAPACKE_dgeqrf LAPACKE_WORKS) + mark_as_advanced(LAPACKE_WORKS) + set(CMAKE_REQUIRED_LIBRARIES) + + if(LAPACKE_WORKS) + if(NOT LAPACKE_FIND_QUIETLY) + message(STATUS "Looking for lapacke: test with lapack succeeds") + endif() + # test succeeds: LAPACKE is in LAPACK + set(LAPACKE_LIBRARIES "${LAPACK_LIBRARIES}") + set(LAPACKE_LIBRARIES_DEP "${LAPACK_LIBRARIES}") + if (LAPACK_LIBRARY_DIRS) + set(LAPACKE_LIBRARY_DIRS "${LAPACK_LIBRARY_DIRS}") + endif() + if(LAPACK_INCLUDE_DIRS) + set(LAPACKE_INCLUDE_DIRS "${LAPACK_INCLUDE_DIRS}") + set(LAPACKE_INCLUDE_DIRS_DEP "${LAPACK_INCLUDE_DIRS}") + endif() + if (LAPACK_LINKER_FLAGS) + set(LAPACKE_LINKER_FLAGS "${LAPACK_LINKER_FLAGS}") + endif() + endif() + endif (NOT LAPACKE_STANDALONE) + + if (LAPACKE_STANDALONE OR NOT LAPACKE_WORKS) + + if(NOT LAPACKE_WORKS AND NOT LAPACKE_FIND_QUIETLY) + message(STATUS "Looking for lapacke : test with lapack fails") + endif() + # test fails: try to find LAPACKE lib exterior to LAPACK + + # Try to find LAPACKE lib + ####################### + + # Looking for include + # ------------------- + + # Add system include paths to search include + # ------------------------------------------ + unset(_inc_env) + set(ENV_LAPACKE_DIR "$ENV{LAPACKE_DIR}") + set(ENV_LAPACKE_INCDIR "$ENV{LAPACKE_INCDIR}") + if(ENV_LAPACKE_INCDIR) + list(APPEND _inc_env "${ENV_LAPACKE_INCDIR}") + elseif(ENV_LAPACKE_DIR) + list(APPEND _inc_env "${ENV_LAPACKE_DIR}") + list(APPEND _inc_env "${ENV_LAPACKE_DIR}/include") + list(APPEND _inc_env "${ENV_LAPACKE_DIR}/include/lapacke") + else() + if(WIN32) + string(REPLACE ":" ";" _inc_env "$ENV{INCLUDE}") + else() + string(REPLACE ":" ";" _path_env "$ENV{INCLUDE}") + list(APPEND _inc_env "${_path_env}") + string(REPLACE ":" ";" _path_env "$ENV{C_INCLUDE_PATH}") + list(APPEND _inc_env "${_path_env}") + string(REPLACE ":" ";" _path_env "$ENV{CPATH}") + list(APPEND _inc_env "${_path_env}") + string(REPLACE ":" ";" _path_env "$ENV{INCLUDE_PATH}") + list(APPEND _inc_env "${_path_env}") + endif() + endif() + list(APPEND _inc_env "${CMAKE_PLATFORM_IMPLICIT_INCLUDE_DIRECTORIES}") + list(APPEND _inc_env "${CMAKE_C_IMPLICIT_INCLUDE_DIRECTORIES}") + list(REMOVE_DUPLICATES _inc_env) + + + # Try to find the lapacke header in the given paths + # ------------------------------------------------- + # call cmake macro to find the header path + if(LAPACKE_INCDIR) + set(LAPACKE_lapacke.h_DIRS "LAPACKE_lapacke.h_DIRS-NOTFOUND") + find_path(LAPACKE_lapacke.h_DIRS + NAMES lapacke.h + HINTS ${LAPACKE_INCDIR}) + else() + if(LAPACKE_DIR) + set(LAPACKE_lapacke.h_DIRS "LAPACKE_lapacke.h_DIRS-NOTFOUND") + find_path(LAPACKE_lapacke.h_DIRS + NAMES lapacke.h + HINTS ${LAPACKE_DIR} + PATH_SUFFIXES "include" "include/lapacke") + else() + set(LAPACKE_lapacke.h_DIRS "LAPACKE_lapacke.h_DIRS-NOTFOUND") + find_path(LAPACKE_lapacke.h_DIRS + NAMES lapacke.h + HINTS ${_inc_env}) + endif() + endif() + mark_as_advanced(LAPACKE_lapacke.h_DIRS) + + # If found, add path to cmake variable + # ------------------------------------ + if (LAPACKE_lapacke.h_DIRS) + set(LAPACKE_INCLUDE_DIRS "${LAPACKE_lapacke.h_DIRS}") + else () + set(LAPACKE_INCLUDE_DIRS "LAPACKE_INCLUDE_DIRS-NOTFOUND") + if(NOT LAPACKE_FIND_QUIETLY) + message(STATUS "Looking for lapacke -- lapacke.h not found") + endif() + endif() + + + # Looking for lib + # --------------- + + # Add system library paths to search lib + # -------------------------------------- + unset(_lib_env) + set(ENV_LAPACKE_LIBDIR "$ENV{LAPACKE_LIBDIR}") + if(ENV_LAPACKE_LIBDIR) + list(APPEND _lib_env "${ENV_LAPACKE_LIBDIR}") + elseif(ENV_LAPACKE_DIR) + list(APPEND _lib_env "${ENV_LAPACKE_DIR}") + list(APPEND _lib_env "${ENV_LAPACKE_DIR}/lib") + else() + if(WIN32) + string(REPLACE ":" ";" _lib_env "$ENV{LIB}") + else() + if(APPLE) + string(REPLACE ":" ";" _lib_env "$ENV{DYLD_LIBRARY_PATH}") + else() + string(REPLACE ":" ";" _lib_env "$ENV{LD_LIBRARY_PATH}") + endif() + list(APPEND _lib_env "${CMAKE_PLATFORM_IMPLICIT_LINK_DIRECTORIES}") + list(APPEND _lib_env "${CMAKE_C_IMPLICIT_LINK_DIRECTORIES}") + endif() + endif() + list(REMOVE_DUPLICATES _lib_env) + + # Try to find the lapacke lib in the given paths + # ---------------------------------------------- + + # name of the lapacke library + set(LAPACKE_lapacke_NAMES "lapacke") + if(LAPACKE_STATIC) + if(WIN32) + set(LAPACKE_lapacke_NAMES "liblapacke.lib") + endif() + + if(UNIX) + set(LAPACKE_lapacke_NAMES "liblapacke.a") + endif() + endif() + + # call cmake macro to find the lib path + if(LAPACKE_LIBDIR) + set(LAPACKE_lapacke_LIBRARY "LAPACKE_lapacke_LIBRARY-NOTFOUND") + find_library(LAPACKE_lapacke_LIBRARY + NAMES ${LAPACKE_lapacke_NAMES} + HINTS ${LAPACKE_LIBDIR}) + else() + if(LAPACKE_DIR) + set(LAPACKE_lapacke_LIBRARY "LAPACKE_lapacke_LIBRARY-NOTFOUND") + find_library(LAPACKE_lapacke_LIBRARY + NAMES ${LAPACKE_lapacke_NAMES} + HINTS ${LAPACKE_DIR} + PATH_SUFFIXES lib lib32 lib64) + else() + set(LAPACKE_lapacke_LIBRARY "LAPACKE_lapacke_LIBRARY-NOTFOUND") + find_library(LAPACKE_lapacke_LIBRARY + NAMES ${LAPACKE_lapacke_NAMES} + HINTS ${_lib_env}) + endif() + endif() + mark_as_advanced(LAPACKE_lapacke_LIBRARY) + + # If found, add path to cmake variable + # ------------------------------------ + if (LAPACKE_lapacke_LIBRARY) + get_filename_component(lapacke_lib_path "${LAPACKE_lapacke_LIBRARY}" PATH) + # set cmake variables + set(LAPACKE_LIBRARIES "${LAPACKE_lapacke_LIBRARY}") + set(LAPACKE_LIBRARY_DIRS "${lapacke_lib_path}") + else () + set(LAPACKE_LIBRARIES "LAPACKE_LIBRARIES-NOTFOUND") + set(LAPACKE_LIBRARY_DIRS "LAPACKE_LIBRARY_DIRS-NOTFOUND") + if (NOT LAPACKE_FIND_QUIETLY) + message(STATUS "Looking for lapacke -- lib lapacke not found") + endif() + endif () + + # check a function to validate the find + if(LAPACKE_LIBRARIES) + + set(REQUIRED_LDFLAGS) + set(REQUIRED_INCDIRS) + set(REQUIRED_LIBDIRS) + set(REQUIRED_LIBS) + + # LAPACKE + if (LAPACKE_INCLUDE_DIRS) + set(REQUIRED_INCDIRS "${LAPACKE_INCLUDE_DIRS}") + endif() + if (LAPACKE_LIBRARY_DIRS) + set(REQUIRED_LIBDIRS "${LAPACKE_LIBRARY_DIRS}") + endif() + set(REQUIRED_LIBS "${LAPACKE_LIBRARIES}") + # LAPACK + if (LAPACK_INCLUDE_DIRS) + list(APPEND REQUIRED_INCDIRS "${LAPACK_INCLUDE_DIRS}") + endif() + if (LAPACK_LIBRARY_DIRS) + list(APPEND REQUIRED_LIBDIRS "${LAPACK_LIBRARY_DIRS}") + endif() + list(APPEND REQUIRED_LIBS "${LAPACK_LIBRARIES}") + if (LAPACK_LINKER_FLAGS) + list(APPEND REQUIRED_LDFLAGS "${LAPACK_LINKER_FLAGS}") + endif() + # Fortran + if (CMAKE_C_COMPILER_ID MATCHES "GNU") + find_library( + FORTRAN_gfortran_LIBRARY + NAMES gfortran + HINTS ${_lib_env} + ) + mark_as_advanced(FORTRAN_gfortran_LIBRARY) + if (FORTRAN_gfortran_LIBRARY) + list(APPEND REQUIRED_LIBS "${FORTRAN_gfortran_LIBRARY}") + endif() + elseif (CMAKE_C_COMPILER_ID MATCHES "Intel") + find_library( + FORTRAN_ifcore_LIBRARY + NAMES ifcore + HINTS ${_lib_env} + ) + mark_as_advanced(FORTRAN_ifcore_LIBRARY) + if (FORTRAN_ifcore_LIBRARY) + list(APPEND REQUIRED_LIBS "${FORTRAN_ifcore_LIBRARY}") + endif() + endif() + # m + find_library(M_LIBRARY NAMES m HINTS ${_lib_env}) + mark_as_advanced(M_LIBRARY) + if(M_LIBRARY) + list(APPEND REQUIRED_LIBS "-lm") + endif() + # set required libraries for link + set(CMAKE_REQUIRED_INCLUDES "${REQUIRED_INCDIRS}") + set(CMAKE_REQUIRED_LIBRARIES) + list(APPEND CMAKE_REQUIRED_LIBRARIES "${REQUIRED_LDFLAGS}") + foreach(lib_dir ${REQUIRED_LIBDIRS}) + list(APPEND CMAKE_REQUIRED_LIBRARIES "-L${lib_dir}") + endforeach() + list(APPEND CMAKE_REQUIRED_LIBRARIES "${REQUIRED_LIBS}") + string(REGEX REPLACE "^ -" "-" CMAKE_REQUIRED_LIBRARIES "${CMAKE_REQUIRED_LIBRARIES}") + + # test link + unset(LAPACKE_WORKS CACHE) + include(CheckFunctionExists) + check_function_exists(LAPACKE_dgeqrf LAPACKE_WORKS) + mark_as_advanced(LAPACKE_WORKS) + + if(LAPACKE_WORKS) + # save link with dependencies + set(LAPACKE_LIBRARIES_DEP "${REQUIRED_LIBS}") + set(LAPACKE_LIBRARY_DIRS_DEP "${REQUIRED_LIBDIRS}") + set(LAPACKE_INCLUDE_DIRS_DEP "${REQUIRED_INCDIRS}") + set(LAPACKE_LINKER_FLAGS "${REQUIRED_LDFLAGS}") + list(REMOVE_DUPLICATES LAPACKE_LIBRARY_DIRS_DEP) + list(REMOVE_DUPLICATES LAPACKE_INCLUDE_DIRS_DEP) + list(REMOVE_DUPLICATES LAPACKE_LINKER_FLAGS) + else() + if(NOT LAPACKE_FIND_QUIETLY) + message(STATUS "Looking for lapacke: test of LAPACKE_dgeqrf with lapacke and lapack libraries fails") + message(STATUS "CMAKE_REQUIRED_LIBRARIES: ${CMAKE_REQUIRED_LIBRARIES}") + message(STATUS "CMAKE_REQUIRED_INCLUDES: ${CMAKE_REQUIRED_INCLUDES}") + message(STATUS "Check in CMakeFiles/CMakeError.log to figure out why it fails") + endif() + endif() + set(CMAKE_REQUIRED_INCLUDES) + set(CMAKE_REQUIRED_FLAGS) + set(CMAKE_REQUIRED_LIBRARIES) + endif(LAPACKE_LIBRARIES) + + endif (LAPACKE_STANDALONE OR NOT LAPACKE_WORKS) + +else(LAPACK_FOUND) + + if (NOT LAPACKE_FIND_QUIETLY) + message(STATUS "LAPACKE requires LAPACK but LAPACK has not been found." + "Please look for LAPACK first.") + endif() + +endif(LAPACK_FOUND) + +if (LAPACKE_LIBRARIES) + list(GET LAPACKE_LIBRARIES 0 first_lib) + get_filename_component(first_lib_path "${first_lib}" PATH) + if (${first_lib_path} MATCHES "(/lib(32|64)?$)|(/lib/intel64$|/lib/ia32$)") + string(REGEX REPLACE "(/lib(32|64)?$)|(/lib/intel64$|/lib/ia32$)" "" not_cached_dir "${first_lib_path}") + set(LAPACKE_DIR_FOUND "${not_cached_dir}" CACHE PATH "Installation directory of LAPACKE library" FORCE) + else() + set(LAPACKE_DIR_FOUND "${first_lib_path}" CACHE PATH "Installation directory of LAPACKE library" FORCE) + endif() +endif() +mark_as_advanced(LAPACKE_DIR) +mark_as_advanced(LAPACKE_DIR_FOUND) + +# check that LAPACKE has been found +# --------------------------------- +include(FindPackageHandleStandardArgs) +find_package_handle_standard_args(LAPACKE DEFAULT_MSG + LAPACKE_LIBRARIES + LAPACKE_WORKS) \ No newline at end of file diff --git a/cmake/dependencies.cmake b/cmake/dependencies.cmake index 29371bed..d451d662 100644 --- a/cmake/dependencies.cmake +++ b/cmake/dependencies.cmake @@ -6,6 +6,7 @@ include(FetchContent) find_package(BLAS) find_package(LAPACK) +find_package(LAPACKE) # ############################################################################## # Memory check @@ -68,7 +69,20 @@ if(TUVX_BUILD_DOCS) endif() # ############################################################################## -# google test +# google test and benchmark + +if(TUVX_ENABLE_BENCHMARK) + FetchContent_Declare( + googlebenchmark + GIT_REPOSITORY https://github.com/google/benchmark.git + GIT_TAG v1.8.3) + + set(BENCHMARK_DOWNLOAD_DEPENDENCIES ON) + set(BENCHMARK_ENABLE_GTEST_TESTS OFF) + set(BENCHMARK_ENABLE_ASSEMBLY_TESTS OFF) + set(BENCHMARK_ENABLE_TESTING OFF) + FetchContent_MakeAvailable(googlebenchmark) +endif() if(TUVX_ENABLE_TESTS) FetchContent_Declare( @@ -89,3 +103,5 @@ if(TUVX_ENABLE_TESTS) set_target_properties(gtest PROPERTIES CXX_CLANG_TIDY "") set_target_properties(gtest_main PROPERTIES CXX_CLANG_TIDY "") endif() + +# ############################################################################## diff --git a/cmake/test_util.cmake b/cmake/test_util.cmake index d77af674..7cdade62 100644 --- a/cmake/test_util.cmake +++ b/cmake/test_util.cmake @@ -50,8 +50,9 @@ function(create_standard_cxx_test) add_executable(test_${TEST_NAME} ${TEST_SOURCES}) - target_link_libraries(test_${TEST_NAME} PUBLIC musica::tuvx GTest::gtest_main) - + target_include_directories(test_${TEST_NAME} PUBLIC ${LAPACK_INCLUDE_DIRS}) + target_link_libraries(test_${TEST_NAME} PUBLIC LAPACK::LAPACK ${LAPACKE_LIBRARIES} musica::tuvx GTest::gtest_main) + # link additional libraries foreach(library ${TEST_LIBRARIES}) target_link_libraries(test_${TEST_NAME} PUBLIC ${library}) diff --git a/docker/Dockerfile b/docker/Dockerfile index 85874c0a..fac338b8 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -2,19 +2,19 @@ FROM fedora:37 RUN dnf -y update \ && dnf -y install \ - gcc-fortran \ - gcc-c++ \ + cmake \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ - valgrind \ + make \ + netcdf-fortran-devel \ python3 \ python3-pip \ - lapack-devel \ + valgrind \ yaml-cpp-devel \ && dnf clean all @@ -23,13 +23,14 @@ RUN pip3 install numpy scipy ENV LD_LIBRARY_PATH=/usr/local/lib64 # build the tuv-x tool -COPY . /tuv-x/ -RUN mkdir /build \ - && cd /build \ +COPY . /tuv-x +RUN cd /tuv-x \ + && mkdir build \ + && cd build \ && cmake -D CMAKE_BUILD_TYPE=release \ -D TUVX_INSTALL_INCLUDE_DIR=/usr/local/include \ -D TUVX_INSTALL_MOD_DIR=/usr/local/include \ - /tuv-x \ + .. \ && make install -j 8 -WORKDIR /build +WORKDIR /tuv-x/build diff --git a/docker/Dockerfile.coverage b/docker/Dockerfile.coverage index e06e0a9e..8f5363cb 100644 --- a/docker/Dockerfile.coverage +++ b/docker/Dockerfile.coverage @@ -2,19 +2,19 @@ FROM fedora:37 RUN dnf -y update \ && dnf -y install \ - gcc-fortran \ - gcc-c++ \ + cmake \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ - valgrind \ + make \ + netcdf-fortran-devel \ python3 \ python3-pip \ - lapack-devel \ + valgrind \ yaml-cpp-devel \ && dnf clean all diff --git a/docker/Dockerfile.docs b/docker/Dockerfile.docs index 2c82b0a6..67b84af1 100644 --- a/docker/Dockerfile.docs +++ b/docker/Dockerfile.docs @@ -2,20 +2,20 @@ FROM fedora:37 RUN dnf -y update \ && dnf -y install \ + cmake \ doxygen \ - gcc-fortran \ - gcc-c++ \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ - valgrind \ + make \ + netcdf-fortran-devel \ python3 \ python3-pip \ - lapack-devel \ + valgrind \ yaml-cpp-devel \ && dnf clean all @@ -31,9 +31,11 @@ RUN echo "The suffix is '$SWITCHER_SUFFIX'" RUN mkdir /build \ && cd /build \ - && cmake -D TUVX_ENABLE_TESTS=OFF \ - -D TUVX_BUILD_DOCS=ON \ - /tuv-x \ + && cmake \ + -D TUVX_ENABLE_TESTS=OFF \ + -D TUVX_ENABLE_BENCHMARK=OFF \ + -D TUVX_BUILD_DOCS=ON \ + /tuv-x \ && make docs WORKDIR /build diff --git a/docker/Dockerfile.memcheck b/docker/Dockerfile.memcheck index 3c640d0e..8a2127ca 100644 --- a/docker/Dockerfile.memcheck +++ b/docker/Dockerfile.memcheck @@ -2,19 +2,19 @@ FROM fedora:37 RUN dnf -y update \ && dnf -y install \ - gcc-fortran \ - gcc-c++ \ + cmake \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ - valgrind \ + make \ + netcdf-fortran-devel \ python3 \ python3-pip \ - lapack-devel \ + valgrind \ yaml-cpp-devel \ && dnf clean all diff --git a/docker/Dockerfile.mpi b/docker/Dockerfile.mpi index a891c312..7f3a382a 100644 --- a/docker/Dockerfile.mpi +++ b/docker/Dockerfile.mpi @@ -10,20 +10,20 @@ USER test_user WORKDIR /home/test_user RUN sudo dnf -y install \ - openmpi-devel \ - gcc-fortran \ - gcc-c++ \ + cmake \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ + make \ + netcdf-fortran-devel \ + openmpi-devel \ python3 \ python3-pip \ valgrind-openmpi \ - lapack-devel \ yaml-cpp-devel \ && sudo dnf clean all diff --git a/docker/Dockerfile.mpi.memcheck b/docker/Dockerfile.mpi.memcheck index 7e04de53..88e0c44b 100644 --- a/docker/Dockerfile.mpi.memcheck +++ b/docker/Dockerfile.mpi.memcheck @@ -10,20 +10,20 @@ USER test_user WORKDIR /home/test_user RUN sudo dnf -y install \ - openmpi-devel \ - gcc-fortran \ - gcc-c++ \ + cmake \ gcc \ + gcc-c++ \ + gcc-fortran \ gdb \ git \ - netcdf-fortran-devel \ - cmake \ - make \ + lapack-devel \ lcov \ + make \ + netcdf-fortran-devel \ + openmpi-devel \ python3 \ python3-pip \ valgrind-openmpi \ - lapack-devel \ yaml-cpp-devel \ && sudo dnf clean all diff --git a/include/tuvx/linear_algebra/linear_algebra.hpp b/include/tuvx/linear_algebra/linear_algebra.hpp index 33bfd888..29c7f5b3 100644 --- a/include/tuvx/linear_algebra/linear_algebra.hpp +++ b/include/tuvx/linear_algebra/linear_algebra.hpp @@ -9,10 +9,7 @@ namespace tuvx { - /// @struct Tridiagonal Matrix Data Structure. - /// @brief Given a matrix of shape (n, n), the - /// data structure holds three std::vectors - /// storing the diagonal values. + /// @brief Tridiagonal Matrix Data Structure. template struct TridiagonalMatrix { @@ -20,59 +17,56 @@ namespace tuvx std::vector upper_diagonal_; // upper diagonal std::vector lower_diagonal_; // lower diagonal std::vector main_diagonal_; // main diagonal - /// @fn tridiagonal matrix structure constructor - /// @brief initializes the three internal vectors based on + /// @brief Initializes the three internal vectors based on /// size. Main diagonal will be of size $n$, while the upper/lower - /// diagonals will be of size $n-1# + /// diagonals will be of size $n-1$ + /// @param size Size of the matrix to be initialized TridiagonalMatrix(std::size_t size); }; - /// @fn Random vector function - /// @brief Fill a std::vector with random values + /// @brief Fills a std::vector with random values /// @param x Vector to allocate and fill + /// @param seed Seed for random number generation template - void FillRandom(std::vector &x); + void FillRandom(std::vector &x, const unsigned &seed); - /// @fn Initialize random matrix function /// @brief Fills a matrix with uniformly distributed random values. - /// @param A tridiagonal matrix to allocate and fill + /// @param A Tridiagonal matrix to allocate and fill + /// @param seed Seed for random number generation + /// @param make_diagonally_dominant Make the tridiagonal matrix diagonally dominant (diagonal value is greater than sum of + /// the corrosponding row values) template - void FillRandom(TridiagonalMatrix &A); + void FillRandom(TridiagonalMatrix &A, const unsigned &seed, const bool &make_diagonally_dominant = false); - /// @fn print vector function - /// @brief displays the data stored inside a std::vector + /// @brief Displays the data stored inside a std::vector /// @param x Vector to print template void Print(const std::vector &x); - /// @fn print matrix function - /// @brief displays the data stored inside a tridiagonal matrix. + /// @brief Displays the data stored inside a tridiagonal matrix. /// For now, this function simply calls printvec() on the three diagonals. - /// @param A tridiagonal matrix to print + /// @param A Tridiagonal matrix to print template void Print(const TridiagonalMatrix &A); - /// @fn Tridiagonal Linear System Solver - /// @brief Thomas' algorithm for solving tridiagonal linear system (A x = b) - /// @param A tridiagonal coeffcient matrix - /// @param b right hand side vector of the tridiagonal system. - /// @returns x solution that solves A x = b. + /// @brief Thomas' algorithm for solving tridiagonal linear system (A x = b). + /// store solution in b. + /// @param A Tridiagonal coeffcient matrix + /// @param b Right hand side vector of the tridiagonal system. template - std::vector Solve(TridiagonalMatrix A, std::vector b); + void Solve(TridiagonalMatrix &A, std::vector &b); - /// @fn Tridiagonal Matrix - vector dot product /// @brief Specialized dot product function for tridiagonal matrices. - /// @param A tridiagonal matrix - /// @param x vector to multiply the matrix with - /// @returns dot product between A and x + /// @param A Tridiagonal matrix + /// @param x Vector to multiply the matrix with + /// @returns Dot product between A and x template std::vector Dot(const TridiagonalMatrix &A, const std::vector &b); - /// @fn residual vector function - /// @brief Computes the LP error norm between two vectors. Used for computing + /// @brief Computes the relative error between two vectors. Used for computing /// approximation errors. - /// @param x true solution - /// @param x_approx approximated solution + /// @param x True solution + /// @param x_approx Approximated solution template T ComputeError(const std::vector &x, const std::vector &x_approx); diff --git a/include/tuvx/linear_algebra/linear_algebra.inl b/include/tuvx/linear_algebra/linear_algebra.inl index 0495bb58..a04caa09 100644 --- a/include/tuvx/linear_algebra/linear_algebra.inl +++ b/include/tuvx/linear_algebra/linear_algebra.inl @@ -1,99 +1,120 @@ -namespace tuvx { +#include +namespace tuvx +{ -template -inline TridiagonalMatrix::TridiagonalMatrix(std::size_t size) { - this->size_ = size; - this->main_diagonal_ = std::vector(size); - this->upper_diagonal_ = std::vector(size - 1); - this->lower_diagonal_ = std::vector(size - 1); -} + template + inline TridiagonalMatrix::TridiagonalMatrix(std::size_t size) + { + this->size_ = size; + this->main_diagonal_ = std::vector(size); + this->upper_diagonal_ = std::vector(size - 1); + this->lower_diagonal_ = std::vector(size - 1); + } -template -inline std::vector Dot(const TridiagonalMatrix &A, - const std::vector &x) { - std::size_t size = x.size(); - std::vector v(size); - v[0] = A.main_diagonal_[0] * x[0] + A.upper_diagonal_[0] * x[1]; + template + inline std::vector Dot(const TridiagonalMatrix &A, const std::vector &x) + { + std::size_t size = x.size(); + std::vector v(size); + v[0] = A.main_diagonal_[0] * x[0] + A.upper_diagonal_[0] * x[1]; - std::size_t i = 0; - for (i = 1; i < size - 1; i++) { - v[i] = A.main_diagonal_[i] * x[i] + A.upper_diagonal_[i] * x[i + 1] + - A.lower_diagonal_[i - 1] * x[i - 1]; + std::size_t i = 0; + for (i = 1; i < size - 1; i++) + { + v[i] = A.main_diagonal_[i] * x[i] + A.upper_diagonal_[i] * x[i + 1] + A.lower_diagonal_[i - 1] * x[i - 1]; + } + v[i] = A.main_diagonal_[i] * x[i] + A.lower_diagonal_[i - 1] * x[i - 1]; + return v; } - v[i] = A.main_diagonal_[i] * x[i] + A.lower_diagonal_[i - 1] * x[i - 1]; - return v; -} -template -inline std::vector Solve(TridiagonalMatrix A, std::vector b) { - T temp; - std::size_t N = b.size(); - std::vector x(N); - // forward pass - for (std::size_t i = 1; i < N; i++) { - temp = A.lower_diagonal_[i - 1] / A.main_diagonal_[i - 1]; - A.main_diagonal_[i] -= temp * A.upper_diagonal_[i - 1]; - b[i] -= temp * b[i - 1]; - } - // back substitution - x[N - 1] = b[N - 1] / A.main_diagonal_[N - 1]; - for (std::size_t i = N - 2;; i--) { - x[i] = (b[i] - A.upper_diagonal_[i] * x[i + 1]) / A.main_diagonal_[i]; - if (i == 0) { - break; + template + inline void Solve(TridiagonalMatrix &A, std::vector &b) + { + T temp; + std::size_t N = b.size(); + // forward pass + for (std::size_t i = 1; i < N; i++) + { + temp = A.lower_diagonal_[i - 1] / A.main_diagonal_[i - 1]; + A.main_diagonal_[i] -= temp * A.upper_diagonal_[i - 1]; + b[i] -= temp * b[i - 1]; + } + // back substitution + b[N - 1] = b[N - 1] / A.main_diagonal_[N - 1]; + for (std::size_t i = N - 2;; i--) + { + b[i] = (b[i] - A.upper_diagonal_[i] * b[i + 1]) / A.main_diagonal_[i]; + if (i == 0) + { + break; + } } } - return x; -} -template inline void FillRandom(std::vector &x) { - std::random_device dev; - std::mt19937 rng(dev()); - std::uniform_int_distribution dist6(1, 6); - for (std::size_t i = 0; i < x.size(); i++) { - x[i] = (T)dist6(rng) + 1; + template + inline void FillRandom(std::vector &x, const unsigned &seed) + { + // sample from normal distribution + + std::mt19937 random_device(seed); + std::normal_distribution distribution(5.0, 1.0); + for (std::size_t i = 0; i < x.size(); i++) + { + x[i] = (T)distribution(random_device); + } } -} -template inline void FillRandom(TridiagonalMatrix &A) { - FillRandom(A.main_diagonal_); - FillRandom(A.lower_diagonal_); - FillRandom(A.upper_diagonal_); - // make diagonally dominant (diagonal value greater than sum of its row) - std::size_t i = 0; - A.main_diagonal_[i] += A.upper_diagonal_[i] + (T)100; - for (i = 1; i < A.size_ - 1; i++) { - A.main_diagonal_[i] += - A.lower_diagonal_[i - 1] + A.upper_diagonal_[i] + (T)100; + template + inline void FillRandom(TridiagonalMatrix &A, const unsigned &seed, const bool &make_diagonally_dominant) + { + FillRandom(A.main_diagonal_, seed); + FillRandom(A.lower_diagonal_, seed); + FillRandom(A.upper_diagonal_, seed); + + if (make_diagonally_dominant) + { + // make diagonally dominant (diagonal value greater than sum of its row) + std::size_t i = 0; + A.main_diagonal_[i] += A.upper_diagonal_[i]; + for (i = 1; i < A.size_ - 1; i++) + { + A.main_diagonal_[i] += A.lower_diagonal_[i - 1] + A.upper_diagonal_[i]; + } + A.main_diagonal_[i] += A.lower_diagonal_[i - 1]; + } } - A.main_diagonal_[i] += A.lower_diagonal_[i - 1] + (T)100; -} -template inline void Print(const std::vector &x) { - std::cout << std::endl; - for (std::size_t i = 0; i < x.size(); i++) { - std::cout << x.at(i) << std::endl; + template + inline void Print(const std::vector &x) + { + std::cout << std::endl; + for (std::size_t i = 0; i < x.size(); i++) + { + std::cout << x.at(i) << std::endl; + } + std::cout << std::endl; } - std::cout << std::endl; -} -template -T ComputeError(const std::vector &x, const std::vector &x_approx) { - T error = 0; - for (std::size_t i = 0; i < x.size(); i++) { - error += - (abs(x[i] - x_approx[i]) / std::max(x[i], x_approx[i])) / (T)x.size(); + template + T ComputeError(const std::vector &x, const std::vector &x_approx) + { + T error = 0; + for (std::size_t i = 0; i < x.size(); i++) + { + error += (std::abs(x[i] - x_approx[i]) / std::max(x[i], x_approx[i])) / (T)x.size(); + } + return error; } - return error; -} -template inline void Print(const TridiagonalMatrix &x) { - std::cout << "----" << std::endl; - Print(x.upper_diagonal_); - Print(x.main_diagonal_); - Print(x.lower_diagonal_); - std::cout << "----" << std::endl; -} + template + inline void Print(const TridiagonalMatrix &x) + { + std::cout << "----" << std::endl; + Print(x.upper_diagonal_); + Print(x.main_diagonal_); + Print(x.lower_diagonal_); + std::cout << "----" << std::endl; + } -} // namespace tuvx +} // namespace tuvx diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index b11530c1..6dd3a25b 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -6,108 +6,100 @@ add_library(tuvx_object OBJECT) target_compile_features(tuvx_object PUBLIC cxx_std_11) -set_target_properties(tuvx_object PROPERTIES - Fortran_MODULE_DIRECTORY ${TUVX_MOD_DIR} -) +set_target_properties(tuvx_object PROPERTIES Fortran_MODULE_DIRECTORY + ${TUVX_MOD_DIR}) -target_include_directories(tuvx_object - PUBLIC - $ - $ - $ - $ -) +target_include_directories( + tuvx_object + PUBLIC $ + $ + $ + $) message(STATUS "blas libraries: ${BLAS_LIBRARIES}") message(STATUS "lapack libraries: ${LAPACK_LIBRARIES}") -target_link_libraries(tuvx_object - PRIVATE - PkgConfig::netcdff - PkgConfig::netcdfc - yaml-cpp::yaml-cpp -) +target_link_libraries(tuvx_object PRIVATE PkgConfig::netcdff + PkgConfig::netcdfc + yaml-cpp::yaml-cpp) + +if("${CMAKE_CXX_COMPILER_ID}" STREQUAL "IntelLLVM") + target_compile_definitions(tuvx_object PUBLIC TUVX_COMPILE_WITH_INTEL) +elseif("${CMAKE_CXX_COMPILER_ID}" STREQUAL "GNU") + target_compile_definitions(tuvx_object PUBLIC TUVX_COMPILE_WITH_GCC) +endif() if(BLAS_LIBRARIES) - target_link_libraries(tuvx_object - PRIVATE - ${BLAS_LIBRARIES} - ) + target_link_libraries(tuvx_object PRIVATE ${BLAS_LIBRARIES}) else() message(FATAL_ERROR "BLAS Libraries not found.") endif() if(LAPACK_LIBRARIES) - target_link_libraries(tuvx_object - PRIVATE - ${LAPACK_LIBRARIES} - ) + target_link_libraries(tuvx_object PRIVATE ${LAPACK_LIBRARIES} LAPACK::LAPACK) else() message(FATAL_ERROR "LAPACK libraries not found.") endif() - # tuvx library add_library(tuvx $) add_library(musica::tuvx ALIAS tuvx) if(NOT BUILD_SHARED_LIBS) - set_target_properties(tuvx_object - PROPERTIES - POSITION_INDEPENDENT_CODE ON - ) + set_target_properties(tuvx_object PROPERTIES POSITION_INDEPENDENT_CODE ON) endif() -set_target_properties(tuvx PROPERTIES - ARCHIVE_OUTPUT_DIRECTORY ${TUVX_LIB_DIR} - VERSION ${PROJECT_VERSION} - SOVERSION ${PROJECT_VERSION_MAJOR} -) +set_target_properties( + tuvx + PROPERTIES ARCHIVE_OUTPUT_DIRECTORY ${TUVX_LIB_DIR} + VERSION ${PROJECT_VERSION} + SOVERSION ${PROJECT_VERSION_MAJOR}) -target_link_libraries(tuvx - PRIVATE - tuvx_object -) +target_link_libraries(tuvx PUBLIC tuvx_object) -message(STATUS "TUV-x build include directories: ${TUVX_MOD_DIR};${CMAKE_CURRENT_SOURCE_DIR}/../include") -message(STATUS "TUV-x install include directories: ${TUVX_INSTALL_MOD_DIR};${TUVX_INSTALL_INCLUDE_DIR}") -target_include_directories(tuvx - PUBLIC - $ - $ - $ - $ +message( + STATUS + "TUV-x build include directories: ${TUVX_MOD_DIR};${CMAKE_CURRENT_SOURCE_DIR}/../include" ) - -target_sources(tuvx_object - PRIVATE - constants.F90 - core.F90 - cross_section.F90 - cross_section_factory.F90 - cross_section_warehouse.F90 - diagnostic_util.F90 - dose_rates.F90 - grid.F90 - grid_factory.F90 - grid_warehouse.F90 - heating_rates.F90 - interpolate.F90 - la_sr_bands.F90 - linear_algebra.F90 - netcdf.F90 - output.F90 - photolysis_rates.F90 - profile.F90 - profile_factory.F90 - profile_warehouse.F90 - quantum_yield.F90 - quantum_yield_factory.F90 - spectral_weight.F90 - spectral_weight_factory.F90 - spherical_geometry.F90 - util.F90 +message( + STATUS + "TUV-x install include directories: ${TUVX_INSTALL_MOD_DIR};${TUVX_INSTALL_INCLUDE_DIR}" ) +target_include_directories( + tuvx + PUBLIC $ + $ + $ + $) + +target_sources( + tuvx_object + PRIVATE constants.F90 + core.F90 + cross_section.F90 + cross_section_factory.F90 + cross_section_warehouse.F90 + diagnostic_util.F90 + dose_rates.F90 + grid.F90 + grid_factory.F90 + grid_warehouse.F90 + heating_rates.F90 + interpolate.F90 + la_sr_bands.F90 + linear_algebra.F90 + netcdf.F90 + output.F90 + photolysis_rates.F90 + profile.F90 + profile_factory.F90 + profile_warehouse.F90 + quantum_yield.F90 + quantum_yield_factory.F90 + spectral_weight.F90 + spectral_weight_factory.F90 + spherical_geometry.F90 + util.F90) add_subdirectory(cross_sections) add_subdirectory(grids) diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index a61c517e..52649f39 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -2,7 +2,7 @@ # Test utilities add_library(tuvx_test_utils unit/test_utils.F90) -target_link_libraries(tuvx_test_utils musica::tuvx) +target_link_libraries(tuvx_test_utils PUBLIC LAPACK::LAPACK musica::tuvx) set_target_properties(tuvx_test_utils PROPERTIES Fortran_MODULE_DIRECTORY ${PROJECT_BINARY_DIR}/test_include diff --git a/test/benchmark/linalg/benchmark_tridiag.cpp b/test/benchmark/linalg/benchmark_tridiag.cpp deleted file mode 100644 index e69de29b..00000000 diff --git a/test/unit/linear_algebra/test_error_function.cpp b/test/unit/linear_algebra/test_error_function.cpp index 9c40605b..3c23e66c 100644 --- a/test/unit/linear_algebra/test_error_function.cpp +++ b/test/unit/linear_algebra/test_error_function.cpp @@ -8,29 +8,23 @@ #include #include -const std::size_t size = 10; // size of the system to test -const double tol_dp = 1e-14; // tolorance for double -const float tol_sp = 1e-5; // tolorance for single point floating precision +const std::size_t SIZE = 10; // size of the system to test +const double TOL_DP = 1e-14; // tolorance for double +const float TOL_SP = 1e-5; // tolorance for single point floating precision -using namespace tuvx; +const unsigned RANDOM_NUMBER_SEED = 1; -typedef TridiagonalMatrix trid_matd; -typedef std::vector vecd; - -typedef TridiagonalMatrix trid_matf; -typedef std::vector vecf; - -/// @test Error function test /// @brief Test the correctness of the error function used for -/// testing the Linear approximation solvers +/// testing the Linear approximation solvers with double precision data types. +/// @param ErrorFunctionTest Name of the test suite TEST(ErrorFunctionTest, DoublePrecision) { // same vector should return 0 error - vecd x(size); - vecd x1(size); - FillRandom(x); + std::vector x(SIZE); + std::vector x1(SIZE); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); x1 = x; - double error = ComputeError(x, x1); + double error = tuvx::ComputeError(x, x1); EXPECT_EQ(error, 0.0); // L1 norm between x and (x[1]+0.1) should be 0.1/size; @@ -38,21 +32,20 @@ TEST(ErrorFunctionTest, DoublePrecision) x1 = x; x[0] += 0.1; - error = ComputeError(x, x1); - EXPECT_LE(error - 0.1 / size, tol_dp); + error = tuvx::ComputeError(x, x1); + EXPECT_LE(error - 0.1 / SIZE, TOL_SP); } -/// @test Error function test /// @brief Test the correctness of the error function used for -/// testing the Linear approximation solvers +/// testing the Linear approximation solvers with single precision data types. TEST(ErrorFunctionTest, SinglePrecision) { // same vector should return 0 error - vecf x(size); - vecf x1(size); - FillRandom(x); + std::vector x(SIZE); + std::vector x1(SIZE); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); x1 = x; - float error = ComputeError(x, x1); + float error = tuvx::ComputeError(x, x1); EXPECT_EQ(error, 0.0f); // L1 norm between x and (x[1]+0.1) should be 0.1/size; @@ -60,6 +53,6 @@ TEST(ErrorFunctionTest, SinglePrecision) x1 = x; x[0] += 0.1f; - error = ComputeError(x, x1); - EXPECT_LE(error - 0.1f / size, tol_sp); + error = tuvx::ComputeError(x, x1); + EXPECT_LE(error - 0.1f / SIZE, TOL_SP); } diff --git a/test/unit/linear_algebra/test_tridiagonal_solver.cpp b/test/unit/linear_algebra/test_tridiagonal_solver.cpp index 8765e40d..70035cb1 100644 --- a/test/unit/linear_algebra/test_tridiagonal_solver.cpp +++ b/test/unit/linear_algebra/test_tridiagonal_solver.cpp @@ -5,67 +5,141 @@ #include #include #include -#include #include #include -using namespace tuvx; +#ifdef TUVX_COMPILE_WITH_INTEL + #include +#elif TUVX_COMPILE_WITH_GCC + #include +#endif -typedef TridiagonalMatrix trid_matd; -typedef std::vector vecd; +const double TOL_DP = std::numeric_limits::epsilon(); +const float TOL_SP = std::numeric_limits::epsilon(); -typedef TridiagonalMatrix trid_matf; -typedef std::vector vecf; +const std::size_t NUMBER_OF_RUNS = 20; +const std::size_t SYSTEM_SIZE = 10; +const bool MAKE_DIAGONALLY_DOMINANT = true; -const double tol_dp = std::numeric_limits::epsilon(); -const float tol_sp = std::numeric_limits::epsilon(); - -/// @test Tridiagonal Solver Test for single Precision Floats. -/// @brief Generate random tridiagonal matrix $A$ and vector $x$, +const unsigned RANDOM_NUMBER_SEED = 1; +/// @brief Tridiagonal Solver Test for single Precision Floats. +/// Generate random tridiagonal matrix $A$ and vector $x$, /// compute $b=A \cdot x$, and check if solution is reconstructed /// accurately using L2 norm (single precision). Check for different /// sizes to check consistency TEST(TridiagSolveTest, SinglePrecision) { - std::size_t sizes[5] = { 5, 10, 1000, 100000, 10000000 }; float error = 0; - for (std::size_t i = 0; i < 5; i++) + + for (std::size_t j = 0; j < NUMBER_OF_RUNS; j++) { - vecf x(sizes[i]); - vecf b(sizes[i]); - trid_matf A(sizes[i]); + std::vector x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); - FillRandom(A); - FillRandom(x); - b = Dot(A, x); + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(A, x); + tuvx::Solve(A, b); - vecf x_approx = Solve(A, b); - error = ComputeError(x, x_approx); - EXPECT_LE(error, tol_sp); + error += tuvx::ComputeError(x, b); } + + error /= NUMBER_OF_RUNS; + + EXPECT_LE(error, TOL_SP); } -/// @test Tridiagonal Solver Test for Double Precision Floats. /// @brief Generate random tridiagonal matrix $A$ and vector $x$, +/// Tridiagonal Solver Test for Double Precision Floats. /// compute $b=A \cdot x$, and check if solution is reconstructed /// accurately using L2 norm (double precision). Check for different /// sizes to check consistency TEST(TridiagSolveTest, DoublePrecision) { - std::size_t sizes[5] = { 5, 10, 1000, 100000, 10000000 }; double error = 0; - for (std::size_t i = 0; i < 5; i++) + for (std::size_t j = 0; j < NUMBER_OF_RUNS; j++) + { + std::vector x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(A, x); + tuvx::Solve(A, b); + + error += tuvx::ComputeError(x, b); + } + error /= NUMBER_OF_RUNS; + EXPECT_LE(error, TOL_DP); +} + +/// @brief LAPACKE Tridiagonal Solver Test for single Precision Floats. +/// Generate random tridiagonal matrix $A$ and vector $x$, +/// compute $b=A \cdot x$, and check if solution is reconstructed +/// accurately using L2 norm (single precision). Check for different +/// sizes to check consistency +TEST(LapackeTest, SinglePrecision) +{ + float error = 0; + + for (std::size_t j = 0; j < NUMBER_OF_RUNS; j++) + { + std::vector x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(A, x); + + LAPACKE_sgtsv( + LAPACK_ROW_MAJOR, + SYSTEM_SIZE, + 1, + A.lower_diagonal_.data(), + A.main_diagonal_.data(), + A.upper_diagonal_.data(), + b.data(), + 1); + + error += tuvx::ComputeError(x, b); + } + error /= NUMBER_OF_RUNS; + EXPECT_LE(error, TOL_SP); +} + +/// @brief LAPACKE Tridiagonal Solver Test for double Precision Floats. +/// Generate random tridiagonal matrix $A$ and vector $x$, +/// compute $b=A \cdot x$, and check if solution is reconstructed +/// accurately using L2 norm (double precision). Check for different +/// sizes to check consistency +TEST(LapackeTest, DoublePrecision) +{ + double error = 0; + for (std::size_t j = 0; j < NUMBER_OF_RUNS; j++) { - vecd x(sizes[i]); - vecd b(sizes[i]); - trid_matd A(sizes[i]); + std::vector x(SYSTEM_SIZE); + std::vector b(SYSTEM_SIZE); + tuvx::TridiagonalMatrix A(SYSTEM_SIZE); + + tuvx::FillRandom(A, RANDOM_NUMBER_SEED, MAKE_DIAGONALLY_DOMINANT); + tuvx::FillRandom(x, RANDOM_NUMBER_SEED); + b = tuvx::Dot(A, x); - FillRandom(A); - FillRandom(x); - b = Dot(A, x); + LAPACKE_dgtsv( + LAPACK_ROW_MAJOR, + SYSTEM_SIZE, + 1, + A.lower_diagonal_.data(), + A.main_diagonal_.data(), + A.upper_diagonal_.data(), + b.data(), + 1); - vecd x_approx = Solve(A, b); - error = ComputeError(x, x_approx); - EXPECT_LE(error, tol_sp); + error += tuvx::ComputeError(x, b); } + error /= NUMBER_OF_RUNS; + EXPECT_LE(error, TOL_DP); }