]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Mark up the scalapack wrapper file for module scripts. 18494/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 23 May 2025 16:15:13 +0000 (10:15 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 23 May 2025 16:16:51 +0000 (10:16 -0600)
include/deal.II/lac/scalapack.templates.h

index 2533f5dc9bb1a7e1b70bb9e72410943cd6cc6bad..7dc9ce3c5def744d93b5980d660143efdcce0545 100644 (file)
 
 #include <deal.II/base/config.h>
 
+// This file does not actually import anything into namespace dealii,
+// but to avoid it being completely empty to some of our scripts, we
+// need to make sure it opens and closes the namespace at least once.
+DEAL_II_NAMESPACE_OPEN
+DEAL_II_NAMESPACE_CLOSE // Do not convert for module purposes
+
+
 #ifdef DEAL_II_WITH_SCALAPACK
 
 #  include <deal.II/base/mpi.h>
 #    include <cfenv>
 #  endif
 
-// useful examples:
-// https://stackoverflow.com/questions/14147705/cholesky-decomposition-scalapack-error/14203864
-// http://icl.cs.utk.edu/lapack-forum/viewtopic.php?t=139   // second post by
-// Julien Langou
-// https://andyspiros.wordpress.com/2011/07/08/an-example-of-blacs-with-c/
-// http://qboxcode.org/trac/browser/qb/tags/rel1_63_4/src/Matrix.C
-// https://gitlab.phys.ethz.ch/lwossnig/lecture/blob/a534f562dfb2ad5c564abe5c2356d5d956fb7218/examples/mpi/scalapack.cpp
-// https://github.com/elemental/Elemental/blob/master/src/core/imports/scalapack.cpp
-// https://scicomp.stackexchange.com/questions/7766/performance-optimization-or-tuning-possible-for-scalapack-gemm
-//
-// info:
-// http://www.netlib.org/scalapack/slug/index.html       // User guide
-// http://www.netlib.org/scalapack/slug/node135.html // How to Measure Errors
 
-extern "C"
+  // useful examples:
+  // https://stackoverflow.com/questions/14147705/cholesky-decomposition-scalapack-error/14203864
+  // http://icl.cs.utk.edu/lapack-forum/viewtopic.php?t=139   // second post by
+  // Julien Langou
+  // https://andyspiros.wordpress.com/2011/07/08/an-example-of-blacs-with-c/
+  // http://qboxcode.org/trac/browser/qb/tags/rel1_63_4/src/Matrix.C
+  // https://gitlab.phys.ethz.ch/lwossnig/lecture/blob/a534f562dfb2ad5c564abe5c2356d5d956fb7218/examples/mpi/scalapack.cpp
+  // https://github.com/elemental/Elemental/blob/master/src/core/imports/scalapack.cpp
+  // https://scicomp.stackexchange.com/questions/7766/performance-optimization-or-tuning-possible-for-scalapack-gemm
+  //
+  // info:
+  // http://www.netlib.org/scalapack/slug/index.html       // User guide
+  // http://www.netlib.org/scalapack/slug/node135.html // How to Measure Errors
+
+  extern "C"
 {
   /* Basic Linear Algebra Communication Subprograms (BLACS) declarations */
   // https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dinitb.htm#dinitb
@@ -51,8 +59,7 @@ extern "C"
    *
    * https://www.ibm.com/support/knowledgecenter/en/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dbpnf.htm
    */
-  void
-  Cblacs_pinfo(int *rank, int *nprocs);
+  void Cblacs_pinfo(int *rank, int *nprocs);
 
   /**
    * Return internal BLACS value in @p val based on the input @p what and @p icontxt.
@@ -61,8 +68,7 @@ extern "C"
    *
    * https://www.ibm.com/support/knowledgecenter/en/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dbget.htm
    */
-  void
-  Cblacs_get(int icontxt, int what, int *val);
+  void Cblacs_get(int icontxt, int what, int *val);
 
   /**
    * Map the processes sequentially in row-major or column-major order
@@ -74,59 +80,53 @@ extern "C"
    *
    * https://www.ibm.com/support/knowledgecenter/en/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dbint.htm
    */
-  void
-  Cblacs_gridinit(int        *context,
-                  const char *order,
-                  int         grid_height,
-                  int         grid_width);
+  void Cblacs_gridinit(int        *context,
+                       const char *order,
+                       int         grid_height,
+                       int         grid_width);
 
   /**
    * Return the process row and column index.
    *
    * https://www.ibm.com/support/knowledgecenter/en/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dbinfo.htm
    */
-  void
-  Cblacs_gridinfo(int  context,
-                  int *grid_height,
-                  int *grid_width,
-                  int *grid_row,
-                  int *grid_col);
+  void Cblacs_gridinfo(int  context,
+                       int *grid_height,
+                       int *grid_width,
+                       int *grid_row,
+                       int *grid_col);
 
   /**
    * Given the system process number, return the row and column coordinates in
    * the BLACS' process grid.
    */
-  void
-  Cblacs_pcoord(int ictxt, int pnum, int *prow, int *pcol);
+  void Cblacs_pcoord(int ictxt, int pnum, int *prow, int *pcol);
 
   /**
    * Release a BLACS context.
    */
-  void
-  Cblacs_gridexit(int context);
+  void Cblacs_gridexit(int context);
 
   /**
    * This routines holds up execution of all processes within the indicated
    * scope until they have all called the routine.
    */
-  void
-  Cblacs_barrier(int, const char *);
+  void Cblacs_barrier(int, const char *);
 
   /**
    * Free all BLACS contexts and releases all allocated memory.
    */
-  void
-  Cblacs_exit(int error_code);
+  void Cblacs_exit(int error_code);
 
   /**
    * Receives a message from a process @prsrc, @p csrc into a general rectangular matrix.
    *
    * https://software.intel.com/en-us/mkl-developer-reference-c-gerv2d
    */
-  void
-  Cdgerv2d(int context, int M, int N, double *A, int lda, int rsrc, int csrc);
-  void
-  Csgerv2d(int context, int M, int N, float *A, int lda, int rsrc, int csrc);
+  void Cdgerv2d(
+    int context, int M, int N, double *A, int lda, int rsrc, int csrc);
+  void Csgerv2d(
+    int context, int M, int N, float *A, int lda, int rsrc, int csrc);
 
   /**
    * Sends the general rectangular matrix A to the destination
@@ -134,16 +134,15 @@ extern "C"
    *
    * https://software.intel.com/en-us/mkl-developer-reference-c-2018-beta-gesd2d
    */
-  void
-  Cdgesd2d(int context, int M, int N, double *A, int lda, int rdest, int cdest);
-  void
-  Csgesd2d(int context, int M, int N, float *A, int lda, int rdest, int cdest);
+  void Cdgesd2d(
+    int context, int M, int N, double *A, int lda, int rdest, int cdest);
+  void Csgesd2d(
+    int context, int M, int N, float *A, int lda, int rdest, int cdest);
 
   /**
    * Get BLACS context from MPI @p comm.
    */
-  int
-  Csys2blacs_handle(MPI_Comm comm);
+  int Csys2blacs_handle(MPI_Comm comm);
 
   /**
    * Compute how many rows and columns each process owns (NUMber of Rows Or
@@ -151,12 +150,11 @@ extern "C"
    *
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_dnumy.htm
    */
-  int
-  numroc_(const int *n,
-          const int *nb,
-          const int *iproc,
-          const int *isproc,
-          const int *nprocs);
+  int numroc_(const int *n,
+              const int *nb,
+              const int *iproc,
+              const int *isproc,
+              const int *nprocs);
 
   /**
    * Compute the Cholesky factorization of an N-by-N real
@@ -166,22 +164,20 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/d5/d9e/pdpotrf_8f_source.html
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lpotrf.htm
    */
-  void
-  pdpotrf_(const char *UPLO,
-           const int  *N,
-           double     *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
-  void
-  pspotrf_(const char *UPLO,
-           const int  *N,
-           float      *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
+  void pdpotrf_(const char *UPLO,
+                const int  *N,
+                double     *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
+  void pspotrf_(const char *UPLO,
+                const int  *N,
+                float      *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
 
   /**
    * Computes an LU factorization of a general distributed matrix sub( A )
@@ -190,24 +186,22 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/df/dfe/pdgetrf_8f_source.html
    * https://www.ibm.com/support/knowledgecenter/en/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lgetrf.htm
    */
-  void
-  pdgetrf_(const int *m,
-           const int *n,
-           double    *A,
-           const int *IA,
-           const int *JA,
-           const int *DESCA,
-           int       *ipiv,
-           int       *INFO);
-  void
-  psgetrf_(const int *m,
-           const int *n,
-           float     *A,
-           const int *IA,
-           const int *JA,
-           const int *DESCA,
-           int       *ipiv,
-           int       *INFO);
+  void pdgetrf_(const int *m,
+                const int *n,
+                double    *A,
+                const int *IA,
+                const int *JA,
+                const int *DESCA,
+                int       *ipiv,
+                int       *INFO);
+  void psgetrf_(const int *m,
+                const int *n,
+                float     *A,
+                const int *IA,
+                const int *JA,
+                const int *DESCA,
+                int       *ipiv,
+                int       *INFO);
 
   /**
    * Compute the inverse of a real symmetric positive definite
@@ -219,22 +213,20 @@ extern "C"
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lpotri.htm
    * https://software.intel.com/en-us/mkl-developer-reference-c-p-potri
    */
-  void
-  pdpotri_(const char *UPLO,
-           const int  *N,
-           double     *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
-  void
-  pspotri_(const char *UPLO,
-           const int  *N,
-           float      *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
+  void pdpotri_(const char *UPLO,
+                const int  *N,
+                double     *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
+  void pspotri_(const char *UPLO,
+                const int  *N,
+                float      *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
 
   /**
    * PDGETRI computes the inverse of a distributed matrix using the LU
@@ -245,30 +237,28 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/d3/df3/pdgetri_8f_source.html
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lgetri.htm
    */
-  void
-  pdgetri_(const int *N,
-           double    *A,
-           const int *IA,
-           const int *JA,
-           const int *DESCA,
-           const int *ipiv,
-           double    *work,
-           int       *lwork,
-           int       *iwork,
-           int       *liwork,
-           int       *info);
-  void
-  psgetri_(const int *N,
-           float     *A,
-           const int *IA,
-           const int *JA,
-           const int *DESCA,
-           const int *ipiv,
-           float     *work,
-           int       *lwork,
-           int       *iwork,
-           int       *liwork,
-           int       *info);
+  void pdgetri_(const int *N,
+                double    *A,
+                const int *IA,
+                const int *JA,
+                const int *DESCA,
+                const int *ipiv,
+                double    *work,
+                int       *lwork,
+                int       *iwork,
+                int       *liwork,
+                int       *info);
+  void psgetri_(const int *N,
+                float     *A,
+                const int *IA,
+                const int *JA,
+                const int *DESCA,
+                const int *ipiv,
+                float     *work,
+                int       *lwork,
+                int       *iwork,
+                int       *liwork,
+                int       *info);
 
 
   /**
@@ -279,24 +269,22 @@ extern "C"
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lpdtri.htm
    * https://software.intel.com/en-us/mkl-developer-reference-c-p-trtri
    */
-  void
-  pdtrtri_(const char *UPLO,
-           const char *DIAG,
-           const int  *N,
-           double     *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
-  void
-  pstrtri_(const char *UPLO,
-           const char *DIAG,
-           const int  *N,
-           float      *A,
-           const int  *IA,
-           const int  *JA,
-           const int  *DESCA,
-           int        *INFO);
+  void pdtrtri_(const char *UPLO,
+                const char *DIAG,
+                const int  *N,
+                double     *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
+  void pstrtri_(const char *UPLO,
+                const char *DIAG,
+                const int  *N,
+                float      *A,
+                const int  *IA,
+                const int  *JA,
+                const int  *DESCA,
+                int        *INFO);
 
   /**
    * Estimate the reciprocal of the condition number (in the
@@ -307,34 +295,32 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/d4/df7/pdpocon_8f.html
    * https://software.intel.com/en-us/mkl-developer-reference-fortran-pocon
    */
-  void
-  pdpocon_(const char   *uplo,
-           const int    *N,
-           const double *A,
-           const int    *IA,
-           const int    *JA,
-           const int    *DESCA,
-           const double *ANORM,
-           double       *RCOND,
-           double       *WORK,
-           const int    *LWORK,
-           int          *IWORK,
-           const int    *LIWORK,
-           int          *INFO);
-  void
-  pspocon_(const char  *uplo,
-           const int   *N,
-           const float *A,
-           const int   *IA,
-           const int   *JA,
-           const int   *DESCA,
-           const float *ANORM,
-           float       *RCOND,
-           float       *WORK,
-           const int   *LWORK,
-           int         *IWORK,
-           const int   *LIWORK,
-           int         *INFO);
+  void pdpocon_(const char   *uplo,
+                const int    *N,
+                const double *A,
+                const int    *IA,
+                const int    *JA,
+                const int    *DESCA,
+                const double *ANORM,
+                double       *RCOND,
+                double       *WORK,
+                const int    *LWORK,
+                int          *IWORK,
+                const int    *LIWORK,
+                int          *INFO);
+  void pspocon_(const char  *uplo,
+                const int   *N,
+                const float *A,
+                const int   *IA,
+                const int   *JA,
+                const int   *DESCA,
+                const float *ANORM,
+                float       *RCOND,
+                float       *WORK,
+                const int   *LWORK,
+                int         *IWORK,
+                const int   *LIWORK,
+                int         *INFO);
 
   /**
    * Norm of a real symmetric matrix
@@ -342,24 +328,22 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/dd/d12/pdlansy_8f_source.html
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_pdlansy.htm#pdlansy
    */
-  double
-  pdlansy_(const char   *norm,
-           const char   *uplo,
-           const int    *N,
-           const double *A,
-           const int    *IA,
-           const int    *JA,
-           const int    *DESCA,
-           double       *work);
-  float
-  pslansy_(const char  *norm,
-           const char  *uplo,
-           const int   *N,
-           const float *A,
-           const int   *IA,
-           const int   *JA,
-           const int   *DESCA,
-           float       *work);
+  double pdlansy_(const char   *norm,
+                  const char   *uplo,
+                  const int    *N,
+                  const double *A,
+                  const int    *IA,
+                  const int    *JA,
+                  const int    *DESCA,
+                  double       *work);
+  float  pslansy_(const char  *norm,
+                  const char  *uplo,
+                  const int   *N,
+                  const float *A,
+                  const int   *IA,
+                  const int   *JA,
+                  const int   *DESCA,
+                  float       *work);
 
   /**
    * Compute the Least Common Multiple (LCM) of two positive integers @p M and @p N.
@@ -368,31 +352,28 @@ extern "C"
    *
    * http://www.netlib.org/scalapack/explore-html/d0/d9b/ilcm_8f_source.html
    */
-  int
-  ilcm_(const int *M, const int *N);
+  int ilcm_(const int *M, const int *N);
 
   /**
    * Return the ceiling of the division of two integers.
    *
    * http://www.netlib.org/scalapack/explore-html/df/d07/iceil_8f_source.html
    */
-  int
-  iceil_(const int *i1, const int *i2);
+  int iceil_(const int *i1, const int *i2);
 
   /**
    * Initialize the descriptor vector with the 8 input arguments
    */
-  void
-  descinit_(int       *desc,
-            const int *m,
-            const int *n,
-            const int *mb,
-            const int *nb,
-            const int *irsrc,
-            const int *icsrc,
-            const int *ictxt,
-            const int *lld,
-            int       *info);
+  void descinit_(int       *desc,
+                 const int *m,
+                 const int *n,
+                 const int *mb,
+                 const int *nb,
+                 const int *irsrc,
+                 const int *icsrc,
+                 const int *ictxt,
+                 const int *lld,
+                 int       *info);
 
   /**
    * Compute the global index of a distributed matrix entry
@@ -409,42 +390,39 @@ extern "C"
    * @param nprocs The total number processes over which the distributed matrix
    * is distributed
    */
-  int
-  indxl2g_(const int *indxloc,
-           const int *nb,
-           const int *iproc,
-           const int *isrcproc,
-           const int *nprocs);
+  int indxl2g_(const int *indxloc,
+               const int *nb,
+               const int *iproc,
+               const int *isrcproc,
+               const int *nprocs);
 
   /**
    * Compute the solution to a real system of linear equations
    */
-  void
-  pdgesv_(const int *n,
-          const int *nrhs,
-          double    *A,
-          const int *ia,
-          const int *ja,
-          const int *desca,
-          int       *ipiv,
-          double    *B,
-          const int *ib,
-          const int *jb,
-          const int *descb,
-          int       *info);
-  void
-  psgesv_(const int *n,
-          const int *nrhs,
-          float     *A,
-          const int *ia,
-          const int *ja,
-          const int *desca,
-          int       *ipiv,
-          float     *B,
-          const int *ib,
-          const int *jb,
-          const int *descb,
-          int       *info);
+  void pdgesv_(const int *n,
+               const int *nrhs,
+               double    *A,
+               const int *ia,
+               const int *ja,
+               const int *desca,
+               int       *ipiv,
+               double    *B,
+               const int *ib,
+               const int *jb,
+               const int *descb,
+               int       *info);
+  void psgesv_(const int *n,
+               const int *nrhs,
+               float     *A,
+               const int *ia,
+               const int *ja,
+               const int *desca,
+               int       *ipiv,
+               float     *B,
+               const int *ib,
+               const int *jb,
+               const int *descb,
+               int       *info);
 
   /**
    * Perform one of the matrix-matrix operations:
@@ -458,80 +436,75 @@ extern "C"
    * $\mathrm{sub}(C)$ denotes C(IC:IC+M-1,JC:JC+N-1),  and, $op(X)$ is one of
    * $op(X) = X$ or $op(X) = X^T$.
    */
-  void
-  pdgemm_(const char   *transa,
-          const char   *transb,
-          const int    *m,
-          const int    *n,
-          const int    *k,
-          const double *alpha,
-          const double *A,
-          const int    *IA,
-          const int    *JA,
-          const int    *DESCA,
-          const double *B,
-          const int    *IB,
-          const int    *JB,
-          const int    *DESCB,
-          const double *beta,
-          double       *C,
-          const int    *IC,
-          const int    *JC,
-          const int    *DESCC);
-  void
-  psgemm_(const char  *transa,
-          const char  *transb,
-          const int   *m,
-          const int   *n,
-          const int   *k,
-          const float *alpha,
-          const float *A,
-          const int   *IA,
-          const int   *JA,
-          const int   *DESCA,
-          const float *B,
-          const int   *IB,
-          const int   *JB,
-          const int   *DESCB,
-          const float *beta,
-          float       *C,
-          const int   *IC,
-          const int   *JC,
-          const int   *DESCC);
+  void pdgemm_(const char   *transa,
+               const char   *transb,
+               const int    *m,
+               const int    *n,
+               const int    *k,
+               const double *alpha,
+               const double *A,
+               const int    *IA,
+               const int    *JA,
+               const int    *DESCA,
+               const double *B,
+               const int    *IB,
+               const int    *JB,
+               const int    *DESCB,
+               const double *beta,
+               double       *C,
+               const int    *IC,
+               const int    *JC,
+               const int    *DESCC);
+  void psgemm_(const char  *transa,
+               const char  *transb,
+               const int   *m,
+               const int   *n,
+               const int   *k,
+               const float *alpha,
+               const float *A,
+               const int   *IA,
+               const int   *JA,
+               const int   *DESCA,
+               const float *B,
+               const int   *IB,
+               const int   *JB,
+               const int   *DESCB,
+               const float *beta,
+               float       *C,
+               const int   *IC,
+               const int   *JC,
+               const int   *DESCC);
 
   /**
    * Return the value of the one norm, or the Frobenius norm, or the infinity
    * norm, or the element of largest absolute value of a distributed matrix
    */
-  double
-  pdlange_(const char   *norm,
-           const int    *m,
-           const int    *n,
-           const double *A,
-           const int    *ia,
-           const int    *ja,
-           const int    *desca,
-           double       *work);
-  float
-  pslange_(const char  *norm,
-           const int   *m,
-           const int   *n,
-           const float *A,
-           const int   *ia,
-           const int   *ja,
-           const int   *desca,
-           float       *work);
+  double pdlange_(const char   *norm,
+                  const int    *m,
+                  const int    *n,
+                  const double *A,
+                  const int    *ia,
+                  const int    *ja,
+                  const int    *desca,
+                  double       *work);
+  float  pslange_(const char  *norm,
+                  const int   *m,
+                  const int   *n,
+                  const float *A,
+                  const int   *ia,
+                  const int   *ja,
+                  const int   *desca,
+                  float       *work);
 
   /**
    * Compute the process coordinate which possesses the entry of a
    * distributed matrix specified by a global index
    */
-  int
-  indxg2p_(const int *glob,
-           const int *nb,
-           const int *iproc,
-           const int *isproc,
-           const int *nprocs);
+  int indxg2p_(const int *glob,
+               const int *nb,
+               const int *iproc,
+               const int *isproc,
+               const int *nprocs);
 
   /**
    * Compute all eigenvalues and, optionally, eigenvectors of a real symmetric
@@ -544,38 +517,36 @@ extern "C"
    * http://www.netlib.org/scalapack/explore-html/d0/d1a/pdsyev_8f.html
    * https://www.ibm.com/support/knowledgecenter/SSNR5K_4.2.0/com.ibm.cluster.pessl.v4r2.pssl100.doc/am6gr_lsyev.htm#lsyev
    */
-  void
-  pdsyev_(const char *jobz,
-          const char *uplo,
-          const int  *m,
-          double     *A,
-          const int  *ia,
-          const int  *ja,
-          int        *desca,
-          double     *w,
-          double     *z,
-          const int  *iz,
-          const int  *jz,
-          int        *descz,
-          double     *work,
-          const int  *lwork,
-          int        *info);
-  void
-  pssyev_(const char *jobz,
-          const char *uplo,
-          const int  *m,
-          float      *A,
-          const int  *ia,
-          const int  *ja,
-          int        *desca,
-          float      *w,
-          float      *z,
-          const int  *iz,
-          const int  *jz,
-          int        *descz,
-          float      *work,
-          const int  *lwork,
-          int        *info);
+  void pdsyev_(const char *jobz,
+               const char *uplo,
+               const int  *m,
+               double     *A,
+               const int  *ia,
+               const int  *ja,
+               int        *desca,
+               double     *w,
+               double     *z,
+               const int  *iz,
+               const int  *jz,
+               int        *descz,
+               double     *work,
+               const int  *lwork,
+               int        *info);
+  void pssyev_(const char *jobz,
+               const char *uplo,
+               const int  *m,
+               float      *A,
+               const int  *ia,
+               const int  *ja,
+               int        *desca,
+               float      *w,
+               float      *z,
+               const int  *iz,
+               const int  *jz,
+               int        *descz,
+               float      *work,
+               const int  *lwork,
+               int        *info);
 
   /**
    * Copy all or a part of a distributed matrix A to another distributed matrix
@@ -584,30 +555,28 @@ extern "C"
    * denotes $A(ia:ia+m-1, ja:ja+n-1)$ and $\mathrm{sub}(B)$ denotes
    * $B(ib:ib+m-1, jb:jb+n-1)$.
    */
-  void
-  pdlacpy_(const char   *uplo,
-           const int    *m,
-           const int    *n,
-           const double *A,
-           const int    *ia,
-           const int    *ja,
-           const int    *desca,
-           double       *B,
-           const int    *ib,
-           const int    *jb,
-           const int    *descb);
-  void
-  pslacpy_(const char  *uplo,
-           const int   *m,
-           const int   *n,
-           const float *A,
-           const int   *ia,
-           const int   *ja,
-           const int   *desca,
-           float       *B,
-           const int   *ib,
-           const int   *jb,
-           const int   *descb);
+  void pdlacpy_(const char   *uplo,
+                const int    *m,
+                const int    *n,
+                const double *A,
+                const int    *ia,
+                const int    *ja,
+                const int    *desca,
+                double       *B,
+                const int    *ib,
+                const int    *jb,
+                const int    *descb);
+  void pslacpy_(const char  *uplo,
+                const int   *m,
+                const int   *n,
+                const float *A,
+                const int   *ia,
+                const int   *ja,
+                const int   *desca,
+                float       *B,
+                const int   *ib,
+                const int   *jb,
+                const int   *descb);
 
   /**
    * Copies the content of a general rectangular distributed matrix @p A to another distributed matrix @p B
@@ -617,38 +586,34 @@ extern "C"
    * @p ictxt is a context which is at least a union of all processes in context
    * A and B
    */
-  void
-  pdgemr2d_(const int    *m,
-            const int    *n,
-            const double *A,
-            const int    *ia,
-            const int    *ja,
-            const int    *desca,
-            double       *B,
-            const int    *ib,
-            const int    *jb,
-            const int    *descb,
-            const int    *ictxt);
-  void
-  psgemr2d_(const int   *m,
-            const int   *n,
-            const float *A,
-            const int   *ia,
-            const int   *ja,
-            const int   *desca,
-            float       *B,
-            const int   *ib,
-            const int   *jb,
-            const int   *descb,
-            const int   *ictxt);
+  void pdgemr2d_(const int    *m,
+                 const int    *n,
+                 const double *A,
+                 const int    *ia,
+                 const int    *ja,
+                 const int    *desca,
+                 double       *B,
+                 const int    *ib,
+                 const int    *jb,
+                 const int    *descb,
+                 const int    *ictxt);
+  void psgemr2d_(const int   *m,
+                 const int   *n,
+                 const float *A,
+                 const int   *ia,
+                 const int   *ja,
+                 const int   *desca,
+                 float       *B,
+                 const int   *ib,
+                 const int   *jb,
+                 const int   *descb,
+                 const int   *ictxt);
 
   /**
    * helper routines determining machine precision
    */
-  double
-  pdlamch_(const int *ictxt, const char *cmach);
-  float
-  pslamch_(const int *ictxt, const char *cmach);
+  double pdlamch_(const int *ictxt, const char *cmach);
+  float  pslamch_(const int *ictxt, const char *cmach);
 
 
   /**
@@ -657,152 +622,146 @@ extern "C"
    * specifying a range of values or a range of indices for the desired
    * eigenvalues.
    */
-  void
-  pdsyevx_(const char   *jobz,
-           const char   *range,
-           const char   *uplo,
-           const int    *n,
-           double       *A,
-           const int    *ia,
-           const int    *ja,
-           const int    *desca,
-           const double *VL,
-           const double *VU,
-           const int    *il,
-           const int    *iu,
-           const double *abstol,
-           const int    *m,
-           const int    *nz,
-           double       *w,
-           double       *orfac,
-           double       *Z,
-           const int    *iz,
-           const int    *jz,
-           const int    *descz,
-           double       *work,
-           int          *lwork,
-           int          *iwork,
-           int          *liwork,
-           int          *ifail,
-           int          *iclustr,
-           double       *gap,
-           int          *info);
-  void
-  pssyevx_(const char  *jobz,
-           const char  *range,
-           const char  *uplo,
-           const int   *n,
-           float       *A,
-           const int   *ia,
-           const int   *ja,
-           const int   *desca,
-           const float *VL,
-           const float *VU,
-           const int   *il,
-           const int   *iu,
-           const float *abstol,
-           const int   *m,
-           const int   *nz,
-           float       *w,
-           float       *orfac,
-           float       *Z,
-           const int   *iz,
-           const int   *jz,
-           const int   *descz,
-           float       *work,
-           int         *lwork,
-           int         *iwork,
-           int         *liwork,
-           int         *ifail,
-           int         *iclustr,
-           float       *gap,
-           int         *info);
+  void pdsyevx_(const char   *jobz,
+                const char   *range,
+                const char   *uplo,
+                const int    *n,
+                double       *A,
+                const int    *ia,
+                const int    *ja,
+                const int    *desca,
+                const double *VL,
+                const double *VU,
+                const int    *il,
+                const int    *iu,
+                const double *abstol,
+                const int    *m,
+                const int    *nz,
+                double       *w,
+                double       *orfac,
+                double       *Z,
+                const int    *iz,
+                const int    *jz,
+                const int    *descz,
+                double       *work,
+                int          *lwork,
+                int          *iwork,
+                int          *liwork,
+                int          *ifail,
+                int          *iclustr,
+                double       *gap,
+                int          *info);
+  void pssyevx_(const char  *jobz,
+                const char  *range,
+                const char  *uplo,
+                const int   *n,
+                float       *A,
+                const int   *ia,
+                const int   *ja,
+                const int   *desca,
+                const float *VL,
+                const float *VU,
+                const int   *il,
+                const int   *iu,
+                const float *abstol,
+                const int   *m,
+                const int   *nz,
+                float       *w,
+                float       *orfac,
+                float       *Z,
+                const int   *iz,
+                const int   *jz,
+                const int   *descz,
+                float       *work,
+                int         *lwork,
+                int         *iwork,
+                int         *liwork,
+                int         *ifail,
+                int         *iclustr,
+                float       *gap,
+                int         *info);
 
   /*
    * PDGESVD computes the singular value decomposition (SVD) of an
    * M-by-N matrix A, optionally computing the left and/or right
    * singular vectors
    */
-  void
-  pdgesvd_(const char *jobu,
-           const char *jobvt,
-           const int  *m,
-           const int  *n,
-           double     *A,
-           const int  *ia,
-           const int  *ja,
-           const int  *desca,
-           double     *S,
-           double     *U,
-           const int  *iu,
-           const int  *ju,
-           const int  *descu,
-           double     *VT,
-           const int  *ivt,
-           const int  *jvt,
-           const int  *descvt,
-           double     *work,
-           int        *lwork,
-           int        *info);
-  void
-  psgesvd_(const char *jobu,
-           const char *jobvt,
-           const int  *m,
-           const int  *n,
-           float      *A,
-           const int  *ia,
-           const int  *ja,
-           const int  *desca,
-           float      *S,
-           float      *U,
-           const int  *iu,
-           const int  *ju,
-           const int  *descu,
-           float      *VT,
-           const int  *ivt,
-           const int  *jvt,
-           const int  *descvt,
-           float      *work,
-           int        *lwork,
-           int        *info);
+  void pdgesvd_(const char *jobu,
+                const char *jobvt,
+                const int  *m,
+                const int  *n,
+                double     *A,
+                const int  *ia,
+                const int  *ja,
+                const int  *desca,
+                double     *S,
+                double     *U,
+                const int  *iu,
+                const int  *ju,
+                const int  *descu,
+                double     *VT,
+                const int  *ivt,
+                const int  *jvt,
+                const int  *descvt,
+                double     *work,
+                int        *lwork,
+                int        *info);
+  void psgesvd_(const char *jobu,
+                const char *jobvt,
+                const int  *m,
+                const int  *n,
+                float      *A,
+                const int  *ia,
+                const int  *ja,
+                const int  *desca,
+                float      *S,
+                float      *U,
+                const int  *iu,
+                const int  *ju,
+                const int  *descu,
+                float      *VT,
+                const int  *ivt,
+                const int  *jvt,
+                const int  *descvt,
+                float      *work,
+                int        *lwork,
+                int        *info);
 
   /*
    * P_GELS solves overdetermined or underdetermined real linear
    * systems involving an M-by-N matrix A, or its transpose,
    * using a QR or LQ factorization of A.  It is assumed that A has full rank.
    */
-  void
-  pdgels_(const char *trans,
-          const int  *m,
-          const int  *n,
-          const int  *nrhs,
-          double     *A,
-          const int  *ia,
-          const int  *ja,
-          const int  *desca,
-          double     *B,
-          const int  *ib,
-          const int  *jb,
-          const int  *descb,
-          double     *work,
-          int        *lwork,
-          int        *info);
-  void
-  psgels_(const char *trans,
-          const int  *m,
-          const int  *n,
-          const int  *nrhs,
-          float      *A,
-          const int  *ia,
-          const int  *ja,
-          const int  *desca,
-          float      *B,
-          const int  *ib,
-          const int  *jb,
-          const int  *descb,
-          float      *work,
-          int        *lwork,
-          int        *info);
+  void pdgels_(const char *trans,
+               const int  *m,
+               const int  *n,
+               const int  *nrhs,
+               double     *A,
+               const int  *ia,
+               const int  *ja,
+               const int  *desca,
+               double     *B,
+               const int  *ib,
+               const int  *jb,
+               const int  *descb,
+               double     *work,
+               int        *lwork,
+               int        *info);
+  void psgels_(const char *trans,
+               const int  *m,
+               const int  *n,
+               const int  *nrhs,
+               float      *A,
+               const int  *ia,
+               const int  *ja,
+               const int  *desca,
+               float      *B,
+               const int  *ib,
+               const int  *jb,
+               const int  *descb,
+               float      *work,
+               int        *lwork,
+               int        *info);
 
   /*
    * Perform matrix sum:
@@ -811,65 +770,61 @@ extern "C"
    * @f
    * where $op(A)$ denotes either $op(A) = A$ or $op(A)=A^T$.
    */
-  void
-  pdgeadd_(const char   *transa,
-           const int    *m,
-           const int    *n,
-           const double *alpha,
-           const double *A,
-           const int    *IA,
-           const int    *JA,
-           const int    *DESCA,
-           const double *beta,
-           double       *C,
-           const int    *IC,
-           const int    *JC,
-           const int    *DESCC);
-  void
-  psgeadd_(const char  *transa,
-           const int   *m,
-           const int   *n,
-           const float *alpha,
-           const float *A,
-           const int   *IA,
-           const int   *JA,
-           const int   *DESCA,
-           const float *beta,
-           float       *C,
-           const int   *IC,
-           const int   *JC,
-           const int   *DESCC);
+  void pdgeadd_(const char   *transa,
+                const int    *m,
+                const int    *n,
+                const double *alpha,
+                const double *A,
+                const int    *IA,
+                const int    *JA,
+                const int    *DESCA,
+                const double *beta,
+                double       *C,
+                const int    *IC,
+                const int    *JC,
+                const int    *DESCC);
+  void psgeadd_(const char  *transa,
+                const int   *m,
+                const int   *n,
+                const float *alpha,
+                const float *A,
+                const int   *IA,
+                const int   *JA,
+                const int   *DESCA,
+                const float *beta,
+                float       *C,
+                const int   *IC,
+                const int   *JC,
+                const int   *DESCC);
 
   /**
    * Routine to transpose a matrix:
    * C = beta C + alpha A^T
    */
-  void
-  pdtran_(const int    *m,
-          const int    *n,
-          const double *alpha,
-          const double *A,
-          const int    *IA,
-          const int    *JA,
-          const int    *DESCA,
-          const double *beta,
-          double       *C,
-          const int    *IC,
-          const int    *JC,
-          const int    *DESCC);
-  void
-  pstran_(const int   *m,
-          const int   *n,
-          const float *alpha,
-          const float *A,
-          const int   *IA,
-          const int   *JA,
-          const int   *DESCA,
-          const float *beta,
-          float       *C,
-          const int   *IC,
-          const int   *JC,
-          const int   *DESCC);
+  void pdtran_(const int    *m,
+               const int    *n,
+               const double *alpha,
+               const double *A,
+               const int    *IA,
+               const int    *JA,
+               const int    *DESCA,
+               const double *beta,
+               double       *C,
+               const int    *IC,
+               const int    *JC,
+               const int    *DESCC);
+  void pstran_(const int   *m,
+               const int   *n,
+               const float *alpha,
+               const float *A,
+               const int   *IA,
+               const int   *JA,
+               const int   *DESCA,
+               const float *beta,
+               float       *C,
+               const int   *IC,
+               const int   *JC,
+               const int   *DESCC);
 
   /**
    * psyevr computes selected eigenvalues and, optionally, eigenvectors
@@ -877,56 +832,54 @@ extern "C"
    * algorithm. Eigenvalues/vectors can be selected by specifying a range of
    * values or a range of indices for the desired eigenvalues.
    */
-  void
-  pdsyevr_(const char   *jobz,
-           const char   *range,
-           const char   *uplo,
-           const int    *n,
-           double       *A,
-           const int    *IA,
-           const int    *JA,
-           const int    *DESCA,
-           const double *VL,
-           const double *VU,
-           const int    *IL,
-           const int    *IU,
-           int          *m,
-           int          *nz,
-           double       *w,
-           double       *Z,
-           const int    *IZ,
-           const int    *JZ,
-           const int    *DESCZ,
-           double       *work,
-           int          *lwork,
-           int          *iwork,
-           int          *liwork,
-           int          *info);
-  void
-  pssyevr_(const char  *jobz,
-           const char  *range,
-           const char  *uplo,
-           const int   *n,
-           float       *A,
-           const int   *IA,
-           const int   *JA,
-           const int   *DESCA,
-           const float *VL,
-           const float *VU,
-           const int   *IL,
-           const int   *IU,
-           int         *m,
-           int         *nz,
-           float       *w,
-           float       *Z,
-           const int   *IZ,
-           const int   *JZ,
-           const int   *DESCZ,
-           float       *work,
-           int         *lwork,
-           int         *iwork,
-           int         *liwork,
-           int         *info);
+  void pdsyevr_(const char   *jobz,
+                const char   *range,
+                const char   *uplo,
+                const int    *n,
+                double       *A,
+                const int    *IA,
+                const int    *JA,
+                const int    *DESCA,
+                const double *VL,
+                const double *VU,
+                const int    *IL,
+                const int    *IU,
+                int          *m,
+                int          *nz,
+                double       *w,
+                double       *Z,
+                const int    *IZ,
+                const int    *JZ,
+                const int    *DESCZ,
+                double       *work,
+                int          *lwork,
+                int          *iwork,
+                int          *liwork,
+                int          *info);
+  void pssyevr_(const char  *jobz,
+                const char  *range,
+                const char  *uplo,
+                const int   *n,
+                float       *A,
+                const int   *IA,
+                const int   *JA,
+                const int   *DESCA,
+                const float *VL,
+                const float *VU,
+                const int   *IL,
+                const int   *IU,
+                int         *m,
+                int         *nz,
+                float       *w,
+                float       *Z,
+                const int   *IZ,
+                const int   *JZ,
+                const int   *DESCZ,
+                float       *work,
+                int         *lwork,
+                int         *iwork,
+                int         *liwork,
+                int         *info);
 }
 
 
@@ -2294,4 +2247,11 @@ psyevr(const char  *jobz,
 
 #endif // DEAL_II_WITH_SCALAPACK
 
+// This file does not actually import anything into namespace dealii,
+// but to avoid it being completely empty to some of our scripts, we
+// need to make sure it opens and closes the namespace at least once.
+DEAL_II_NAMESPACE_OPEN // Do not convert for module purposes
+  DEAL_II_NAMESPACE_CLOSE
+
+
 #endif // dealii_scalapack_templates_h

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.