1 #pragma once 2 3 #include <petsc/private/matimpl.h> 4 #include <petscblaslapack.h> 5 6 typedef struct { 7 PetscBLASInt ictxt; /* process grid context */ 8 PetscBLASInt nprow, npcol; /* number of process rows and columns */ 9 PetscBLASInt myrow, mycol; /* coordinates of local process on the grid */ 10 PetscInt grid_refct; /* reference count */ 11 PetscBLASInt ictxrow, ictxcol; /* auxiliary 1d process grid contexts */ 12 } Mat_ScaLAPACK_Grid; 13 14 typedef struct { 15 Mat_ScaLAPACK_Grid *grid; /* process grid */ 16 PetscBLASInt desc[9]; /* ScaLAPACK descriptor */ 17 PetscBLASInt M, N; /* global dimensions, for rows and columns */ 18 PetscBLASInt locr, locc; /* dimensions of local array */ 19 PetscBLASInt mb, nb; /* block size, for rows and columns */ 20 PetscBLASInt rsrc, csrc; /* coordinates of process owning first row and column */ 21 PetscScalar *loc; /* pointer to local array */ 22 PetscBLASInt lld; /* local leading dimension */ 23 PetscBLASInt *pivots; /* pivots in LU factorization */ 24 PetscBool roworiented; /* if true, row-oriented input (default) */ 25 } Mat_ScaLAPACK; 26 27 PETSC_INTERN PetscErrorCode MatMatMultSymbolic_ScaLAPACK(Mat, Mat, PetscReal, Mat); 28 PETSC_INTERN PetscErrorCode MatMatMultNumeric_ScaLAPACK(Mat, Mat, Mat); 29 30 /* Macro to check nonzero info after ScaLAPACK call */ 31 #define PetscCheckScaLapackInfo(routine, info) \ 32 do { \ 33 PetscCheck(!info, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in ScaLAPACK subroutine %s: info=%d", routine, (int)info); \ 34 } while (0) 35 36 #define PETSC_PASTE4_(a, b, c, d) a##b##c##d 37 #define PETSC_PASTE4(a, b, c, d) PETSC_PASTE4_(a, b, c, d) 38 39 #if defined(PETSC_BLASLAPACK_CAPS) 40 #define PETSC_SCALAPACK_PREFIX_ P 41 #define PETSCBLASNOTYPE(x, X) PETSC_PASTE2(X, PETSC_BLASLAPACK_SUFFIX_) 42 #define PETSCSCALAPACK(x, X) PETSC_PASTE4(PETSC_SCALAPACK_PREFIX_, PETSC_BLASLAPACK_PREFIX_, X, PETSC_BLASLAPACK_SUFFIX_) 43 #else 44 #define PETSC_SCALAPACK_PREFIX_ p 45 #define PETSCBLASNOTYPE(x, X) PETSC_PASTE2(x, PETSC_BLASLAPACK_SUFFIX_) 46 #define PETSCSCALAPACK(x, X) PETSC_PASTE4(PETSC_SCALAPACK_PREFIX_, PETSC_BLASLAPACK_PREFIX_, x, PETSC_BLASLAPACK_SUFFIX_) 47 #endif 48 49 /* BLACS routines (C interface) */ 50 BLAS_EXTERN PetscBLASInt Csys2blacs_handle(MPI_Comm syscontext); 51 BLAS_EXTERN void Cblacs_pinfo(PetscBLASInt *mypnum, PetscBLASInt *nprocs); 52 BLAS_EXTERN void Cblacs_get(PetscBLASInt context, PetscBLASInt request, PetscBLASInt *value); 53 BLAS_EXTERN PetscBLASInt Cblacs_pnum(PetscBLASInt context, PetscBLASInt prow, PetscBLASInt pcol); 54 BLAS_EXTERN PetscBLASInt Cblacs_gridinit(PetscBLASInt *context, const char *order, PetscBLASInt np_row, PetscBLASInt np_col); 55 BLAS_EXTERN void Cblacs_gridinfo(PetscBLASInt context, PetscBLASInt *np_row, PetscBLASInt *np_col, PetscBLASInt *my_row, PetscBLASInt *my_col); 56 BLAS_EXTERN void Cblacs_gridexit(PetscBLASInt context); 57 BLAS_EXTERN void Cblacs_exit(PetscBLASInt error_code); 58 BLAS_EXTERN void Cdgebs2d(PetscBLASInt ctxt, const char *scope, const char *top, PetscBLASInt m, PetscBLASInt n, PetscScalar *A, PetscBLASInt lda); 59 BLAS_EXTERN void Cdgebr2d(PetscBLASInt ctxt, const char *scope, const char *top, PetscBLASInt m, PetscBLASInt n, PetscScalar *A, PetscBLASInt lda, PetscBLASInt rsrc, PetscBLASInt csrc); 60 BLAS_EXTERN void Cdgsum2d(PetscBLASInt ctxt, const char *scope, const char *top, PetscBLASInt m, PetscBLASInt n, PetscScalar *A, PetscBLASInt lda, PetscBLASInt rsrc, PetscBLASInt csrc); 61 62 /* PBLAS */ 63 #define PBLASgemv_ PETSCSCALAPACK(gemv, GEMV) 64 #define PBLASgemm_ PETSCSCALAPACK(gemm, GEMM) 65 #if defined(PETSC_USE_COMPLEX) 66 #define PBLAStran_ PETSCSCALAPACK(tranc, TRANC) 67 #else 68 #define PBLAStran_ PETSCSCALAPACK(tran, TRAN) 69 #endif 70 71 BLAS_EXTERN void PBLASgemv_(const char *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 72 BLAS_EXTERN void PBLASgemm_(const char *, const char *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 73 BLAS_EXTERN void PBLAStran_(PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 74 75 /* ScaLAPACK */ 76 #define SCALAPACKlange_ PETSCSCALAPACK(lange, LANGE) 77 #define SCALAPACKpotrf_ PETSCSCALAPACK(potrf, POTRF) 78 #define SCALAPACKpotrs_ PETSCSCALAPACK(potrs, POTRS) 79 #define SCALAPACKgetrf_ PETSCSCALAPACK(getrf, GETRF) 80 #define SCALAPACKgetrs_ PETSCSCALAPACK(getrs, GETRS) 81 82 BLAS_EXTERN PetscReal SCALAPACKlange_(const char *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *); 83 BLAS_EXTERN void SCALAPACKpotrf_(const char *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 84 BLAS_EXTERN void SCALAPACKpotrs_(const char *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 85 BLAS_EXTERN void SCALAPACKgetrf_(PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 86 BLAS_EXTERN void SCALAPACKgetrs_(const char *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 87 88 /* auxiliary routines */ 89 #define SCALAPACKnumroc_ PETSCBLASNOTYPE(numroc, NUMROC) 90 #define SCALAPACKdescinit_ PETSCBLASNOTYPE(descinit, DESCINIT) 91 #define SCALAPACKinfog2l_ PETSCBLASNOTYPE(infog2l, INFOG2L) 92 #define SCALAPACKgemr2d_ PETSCSCALAPACK(gemr2d, GEMR2D) 93 #define SCALAPACKmatadd_ PETSCSCALAPACK(matadd, MATADD) 94 #define SCALAPACKelset_ PETSCSCALAPACK(elset, ELSET) 95 #define SCALAPACKelget_ PETSCSCALAPACK(elget, ELGET) 96 97 BLAS_EXTERN PetscBLASInt SCALAPACKnumroc_(PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 98 BLAS_EXTERN void SCALAPACKdescinit_(PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 99 BLAS_EXTERN void SCALAPACKinfog2l_(PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 100 BLAS_EXTERN void SCALAPACKgemr2d_(PetscBLASInt *, PetscBLASInt *, const PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 101 BLAS_EXTERN void SCALAPACKmatadd_(PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 102 BLAS_EXTERN void SCALAPACKelset_(PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *, PetscScalar *); 103 BLAS_EXTERN void SCALAPACKelget_(const char *, const char *, PetscScalar *, PetscScalar *, PetscBLASInt *, PetscBLASInt *, PetscBLASInt *); 104 105 /* 106 Macros to test valid arguments 107 */ 108 #if !defined(PETSC_USE_DEBUG) 109 110 #define MatScaLAPACKCheckDistribution(a, arga, b, argb) \ 111 do { \ 112 (void)(a); \ 113 (void)(b); \ 114 } while (0) 115 116 #else 117 118 #define MatScaLAPACKCheckDistribution(a, arga, b, argb) \ 119 do { \ 120 Mat_ScaLAPACK *_aa = (Mat_ScaLAPACK *)(a)->data, *_bb = (Mat_ScaLAPACK *)(b)->data; \ 121 PetscCheck(_aa->mb == _bb->mb && _aa->nb == _bb->nb && _aa->rsrc == _bb->rsrc && _aa->csrc == _bb->csrc && _aa->grid->nprow == _bb->grid->nprow && _aa->grid->npcol == _bb->grid->npcol && _aa->grid->myrow == _bb->grid->myrow && \ 122 _aa->grid->mycol == _bb->grid->mycol, \ 123 PetscObjectComm((PetscObject)(a)), PETSC_ERR_ARG_INCOMP, "Arguments #%d and #%d have different ScaLAPACK distribution", arga, argb); \ 124 } while (0) 125 126 #endif 127