From: Denis Davydov <davydden@gmail.com>
Date: Mon, 20 Nov 2017 09:11:56 +0000 (+0100)
Subject: find ScaLAPACK with cmake, document the new interface, add changelog
X-Git-Tag: v9.0.0-rc1~701^2~6
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=725d2633261d848f478041f19067869d2b530cc9;p=dealii.git

find ScaLAPACK with cmake, document the new interface, add changelog
---

diff --git a/cmake/configure/configure_scalapack.cmake b/cmake/configure/configure_scalapack.cmake
new file mode 100644
index 0000000000..a0b5099559
--- /dev/null
+++ b/cmake/configure/configure_scalapack.cmake
@@ -0,0 +1,32 @@
+## ---------------------------------------------------------------------
+##
+## Copyright (C) 2017 by the deal.II authors
+##
+## This file is part of the deal.II library.
+##
+## The deal.II library is free software; you can use it, redistribute
+## it, and/or modify it under the terms of the GNU Lesser General
+## Public License as published by the Free Software Foundation; either
+## version 2.1 of the License, or (at your option) any later version.
+## The full text of the license can be found in the file LICENSE at
+## the top level of the deal.II distribution.
+##
+## ---------------------------------------------------------------------
+
+#
+# Configuration for the SCALAPACK library:
+#
+
+SET(FEATURE_SCALAPACK_DEPENDS MPI LAPACK)
+
+
+MACRO(FEATURE_SCALAPACK_FIND_EXTERNAL var)
+  FIND_PACKAGE(SCALAPACK)
+
+  IF(SCALAPACK_FOUND)
+    SET(${var} TRUE)
+    CHECK_MPI_INTERFACE(SCALAPACK ${var})
+  ENDIF()
+ENDMACRO()
+
+CONFIGURE_FEATURE(SCALAPACK)
diff --git a/cmake/modules/FindSCALAPACK.cmake b/cmake/modules/FindSCALAPACK.cmake
new file mode 100644
index 0000000000..715bdbb689
--- /dev/null
+++ b/cmake/modules/FindSCALAPACK.cmake
@@ -0,0 +1,95 @@
+## ---------------------------------------------------------------------
+##
+## Copyright (C) 2017 by the deal.II authors
+##
+## This file is part of the deal.II library.
+##
+## The deal.II library is free software; you can use it, redistribute
+## it, and/or modify it under the terms of the GNU Lesser General
+## Public License as published by the Free Software Foundation; either
+## version 2.1 of the License, or (at your option) any later version.
+## The full text of the license can be found in the file LICENSE at
+## the top level of the deal.II distribution.
+##
+## ---------------------------------------------------------------------
+
+#
+# Try to find the SCALAPACK library
+#
+# This module exports
+#
+#   SCALAPACK_LIBRARIES
+#   SCALAPACK_LINKER_FLAGS
+#
+
+SET(SCALAPACK_DIR "" CACHE PATH "An optional hint to a SCALAPACK directory")
+SET(BLACS_DIR "" CACHE PATH "An optional hint to a BLACS directory")
+SET_IF_EMPTY(SCALAPACK_DIR "$ENV{SCALAPACK_DIR}")
+SET_IF_EMPTY(BLACS_DIR "$ENV{BLACS_DIR}")
+
+#
+# Search for scalapack:
+#
+
+DEAL_II_FIND_LIBRARY(SCALAPACK_LIBRARY NAMES scalapack scalapack-openmpi scalapack-pvm scalapack-mpi scalapack-mpich scalapack-mpich2 scalapack-lam
+  HINTS ${SCALAPACK_DIR}
+  PATH_SUFFIXES lib${LIB_SUFFIX} lib64 lib
+  )
+
+#
+# Well, depending on the version of scalapack and the distribution it might
+# be necessary to search for blacs, too.
+IF (SCALAPACK_LIBRARY)
+  MESSAGE(STATUS "Check if BLACS is embedded in ScaLAPACK library")
+
+  CLEAR_CMAKE_REQUIRED()
+  SET(CMAKE_REQUIRED_FLAGS ${MPI_CXX_COMPILE_FLAGS} ${MPI_CXX_LINK_FLAGS})
+  SET(CMAKE_REQUIRED_INCLUDES ${MPI_CXX_INCLUDE_PATH})
+  SET(CMAKE_REQUIRED_LIBRARIES ${MPI_LIBRARIES} ${SCALAPACK_LIBRARY} ${LAPACK_LIBRARIES})
+  CHECK_CXX_SOURCE_COMPILES(
+    "
+    #include <mpi.h>
+    extern \"C\"
+    {
+      int Csys2blacs_handle(MPI_Comm comm);
+    }
+    int main() {
+      const int res = Csys2blacs_handle(MPI_COMM_WORLD);
+    }
+    "
+    SCALAPACK_LIBRARY_HAS_BLACS
+  )
+  RESET_CMAKE_REQUIRED()
+
+  # If Blacs is not embedded, try to find it
+  IF (NOT SCALAPACK_LIBRARY_HAS_BLACS)
+    MESSAGE(STATUS "Try to find BLACS")
+    # FIXME: use instead:
+    #   FIND_PACKAGE(BLACS ${SCALAPACK_FIND_QUIETLY})
+    #   SET(SCALAPACK_LIBRARY ${SCALAPACK_LIBRARY} ${BLACS_LIBRARY})
+
+    # FIXME: this is potentially dangerous endeavour as we can easily pickup
+    # system-provided BLACS build with some MPI which is inconsistent with
+    # ScaLAPACK library found above.
+    FOREACH(_lib blacs blacsCinit blacsF77init)
+      STRING(TOUPPER "${_lib}" _lib_upper)
+      DEAL_II_FIND_LIBRARY(${_lib_upper}_LIBRARY
+        NAMES ${_lib} ${_lib}_MPI-LINUX-0 ${_lib}_MPI-DARWIN-0 ${_lib}-openmpi
+        HINTS ${BLACS_DIR} ${SCALAPACK_DIR} ${SCALAPACK_DIR}/../blacs/
+        PATH_SUFFIXES lib${LIB_SUFFIX} lib64 lib LIB
+      )
+    ENDFOREACH()
+  ENDIF()
+ENDIF()
+
+
+DEAL_II_PACKAGE_HANDLE(SCALAPACK
+  LIBRARIES
+    REQUIRED SCALAPACK_LIBRARY LAPACK_LIBRARIES
+    OPTIONAL BLACS_LIBRARY BLACSCINIT_LIBRARY BLACSF77INIT_LIBRARY MPI_Fortran_LIBRARIES
+  LINKER_FLAGS
+    OPTIONAL LAPACK_LINKER_FLAGS
+  CLEAR
+    SCALAPACK_LIBRARY
+    BLACS_LIBRARY BLACSCINIT_LIBRARY BLACSF77INIT_LIBRARY
+  )
diff --git a/doc/doxygen/options.dox.in b/doc/doxygen/options.dox.in
index 5aa46e195d..d5555be6b4 100644
--- a/doc/doxygen/options.dox.in
+++ b/doc/doxygen/options.dox.in
@@ -187,6 +187,7 @@ PREDEFINED             = DOXYGEN=1 \
                          DEAL_II_WITH_OPENCASCADE=1 \
                          DEAL_II_WITH_P4EST=1 \
                          DEAL_II_WITH_PETSC=1 \
+                         DEAL_II_WITH_SCALAPACK=1 \
                          DEAL_II_WITH_SLEPC=1 \
                          DEAL_II_WITH_SUNDIALS=1 \
                          DEAL_II_WITH_THREADS=1 \
diff --git a/doc/news/changes/major/20171120DenisDavydov-BenjaminBrands b/doc/news/changes/major/20171120DenisDavydov-BenjaminBrands
new file mode 100644
index 0000000000..e424c1da28
--- /dev/null
+++ b/doc/news/changes/major/20171120DenisDavydov-BenjaminBrands
@@ -0,0 +1,4 @@
+New: add ScaLAPACKMatrix and ProcessGrid wrappers for high-performance dense linear
+algebra routines for parallel distributed memory machines.
+<br>
+(Denis Davydov, Benjamin Brands, 2017/11/20)
diff --git a/doc/readme.html b/doc/readme.html
index d8218350a6..b66ae341e0 100644
--- a/doc/readme.html
+++ b/doc/readme.html
@@ -558,9 +558,9 @@
       <dt><a name="CUDA"></a><a href="https://developer.nvidia.com/cuda-zone">CUDA</a></dt>
       <dd>
   <p>
-    <a href="https://developer.nvidia.com/cuda-zone">CUDA</a> is a parallel 
-    computing platform and API model created by Nvidia. It allows software 
-    developers and software engineers to use CUDA-enabled GPU for general 
+    <a href="https://developer.nvidia.com/cuda-zone">CUDA</a> is a parallel
+    computing platform and API model created by Nvidia. It allows software
+    developers and software engineers to use CUDA-enabled GPU for general
     purpose processing.
 
     Details about compatibility and configuration can be found <a href=external-libs/cuda.html
@@ -715,6 +715,23 @@
       </dd>
 
 
+      <dt>
+	<a name="scalapack"></a>
+	<a href="http://www.netlib.org/scalapack/">scalapack</a>
+      </dt>
+      <dd>
+	<p>
+	  <a href="http://www.netlib.org/scalapack/">scalapack</a>
+	  is a library of high-performance linear algebra routines for parallel
+    distributed memory machines. ScaLAPACK solves dense and banded linear
+    systems, least squares problems, eigenvalue problems, and singular value
+    problems. In order to enable this feature, add
+    <code>-DSCALAPACK_DIR=/path/to/scalapack</code> to the argument list for
+    <code>cmake</code>.
+	</p>
+      </dd>
+
+
       <dt><a name="slepc"></a><a href="http://www.grycap.upv.es/slepc/" target="_top">SLEPc</a></dt>
       <dd>
 	<p>
diff --git a/doc/users/cmake.html b/doc/users/cmake.html
index 441f1024d8..ec370a953f 100644
--- a/doc/users/cmake.html
+++ b/doc/users/cmake.html
@@ -473,6 +473,7 @@ DEAL_II_WITH_NETCDF
 DEAL_II_WITH_OPENCASCADE
 DEAL_II_WITH_P4EST
 DEAL_II_WITH_PETSC
+DEAL_II_WITH_SCALAPACK
 DEAL_II_WITH_SLEPC
 DEAL_II_WITH_SUNDIALS
 DEAL_II_WITH_THREADS
@@ -626,6 +627,7 @@ METIS_DIR,
 MUPARSER_DIR,
 P4EST_DIR (and SC_DIR),
 PETSC_DIR and PETSC_ARCH (forming ${PETSC_DIR}/${PETSC_ARCH}),
+SCALAPACK_DIR (and optionally BLACS_DIR)
 SLEPC_DIR (forming ${SLEPC_DIR}/${PETSC_ARCH}),
 TBB_DIR,
 TRILINOS_DIR,
diff --git a/doc/users/cmakelists.html b/doc/users/cmakelists.html
index cd42b5040b..c77b1d99a3 100644
--- a/doc/users/cmakelists.html
+++ b/doc/users/cmakelists.html
@@ -860,6 +860,7 @@ DEAL_II_WITH_NETCDF
 DEAL_II_WITH_OPENCASCADE
 DEAL_II_WITH_P4EST
 DEAL_II_WITH_PETSC
+DEAL_II_WITH_SCALAPACK
 DEAL_II_WITH_SLEPC
 DEAL_II_WITH_SUNDIALS
 DEAL_II_WITH_THREADS
diff --git a/include/deal.II/base/config.h.in b/include/deal.II/base/config.h.in
index 95da58f00e..c56784bb77 100644
--- a/include/deal.II/base/config.h.in
+++ b/include/deal.II/base/config.h.in
@@ -51,6 +51,7 @@
 #cmakedefine DEAL_II_WITH_OPENCASCADE
 #cmakedefine DEAL_II_WITH_P4EST
 #cmakedefine DEAL_II_WITH_PETSC
+#cmakedefine DEAL_II_WITH_SCALAPACK
 #cmakedefine DEAL_II_WITH_SLEPC
 #cmakedefine DEAL_II_WITH_SUNDIALS
 #cmakedefine DEAL_II_WITH_THREADS
diff --git a/include/deal.II/lac/scalapack.h b/include/deal.II/lac/scalapack.h
index 48fb3d6c22..d216e74932 100644
--- a/include/deal.II/lac/scalapack.h
+++ b/include/deal.II/lac/scalapack.h
@@ -18,7 +18,7 @@
 
 #include <deal.II/base/config.h>
 
-// FIXME: #ifdef DEAL_II_WITH_SCALAPACK (<--- Lapack+MPI+Scalapack)
+#ifdef DEAL_II_WITH_SCALAPACK
 
 #include <deal.II/base/exceptions.h>
 #include <deal.II/lac/full_matrix.h>
@@ -492,4 +492,6 @@ ScaLAPACKMatrix<NumberType>::local_el(const int loc_row, const int loc_column)
 
 DEAL_II_NAMESPACE_CLOSE
 
+#endif // DEAL_II_WITH_SCALAPACK
+
 #endif
diff --git a/source/lac/scalapack.cc b/source/lac/scalapack.cc
index a2283156de..e4d26217e7 100644
--- a/source/lac/scalapack.cc
+++ b/source/lac/scalapack.cc
@@ -15,6 +15,9 @@
 
 
 #include <deal.II/lac/scalapack.h>
+
+#ifdef DEAL_II_WITH_SCALAPACK
+
 #include <deal.II/base/mpi.h>
 
 #include <deal.II/base/conditional_ostream.h>
@@ -937,3 +940,5 @@ int ScaLAPACKMatrix<NumberType>::get_process_grid_columns() const
 template class ScaLAPACKMatrix<double>;
 
 DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_SCALAPACK
diff --git a/tests/run_testsuite.cmake b/tests/run_testsuite.cmake
index 6e53bb50a6..6c1b827fd1 100644
--- a/tests/run_testsuite.cmake
+++ b/tests/run_testsuite.cmake
@@ -270,7 +270,7 @@ FOREACH(_var ${_variables})
   IF( _var MATCHES "^(TEST|DEAL_II|ALLOW|WITH|FORCE|COMPONENT)_" OR
       _var MATCHES "^(DOCUMENTATION|EXAMPLES)" OR
       _var MATCHES "^(ADOLC|ARPACK|BOOST|OPENCASCADE|MUPARSER|HDF5|METIS|MPI)_" OR
-      _var MATCHES "^(NETCDF|P4EST|PETSC|SLEPC|THREADS|TBB|TRILINOS)_" OR
+      _var MATCHES "^(NETCDF|P4EST|PETSC|SCALAPACK|SLEPC|THREADS|TBB|TRILINOS)_" OR
       _var MATCHES "^(UMFPACK|ZLIB|LAPACK|MUPARSER|CUDA)_" OR
       _var MATCHES "^(CMAKE|DEAL_II)_(C|CXX|Fortran|BUILD)_(COMPILER|FLAGS)" OR
       _var MATCHES "^CMAKE_BUILD_TYPE$" OR
@@ -483,7 +483,7 @@ IF("${_res}" STREQUAL "0")
       --build . --target setup_tests
       -- ${MAKEOPTS}
       WORKING_DIRECTORY ${CTEST_BINARY_DIRECTORY}
-      OUTPUT_QUIET 
+      OUTPUT_QUIET
       RESULT_VARIABLE _res
       )
 
diff --git a/tests/scalapack/CMakeLists.txt b/tests/scalapack/CMakeLists.txt
index d33eafa20b..b119a5be8f 100644
--- a/tests/scalapack/CMakeLists.txt
+++ b/tests/scalapack/CMakeLists.txt
@@ -1,4 +1,6 @@
 CMAKE_MINIMUM_REQUIRED(VERSION 2.8.12)
 INCLUDE(../setup_testsubproject.cmake)
 PROJECT(testsuite CXX)
-DEAL_II_PICKUP_TESTS()
+IF(DEAL_II_WITH_SCALAPACK)
+  DEAL_II_PICKUP_TESTS()
+ENDIF()