Actual source code: petscelemental.h

  1: #ifndef _PETSCELEMENTAL_H
  2: #define _PETSCELEMENTAL_H

  4: #include <petsc/private/matimpl.h>
  5: #include <petscmatelemental.h>

  7: typedef struct {
  8:   PetscInt                         commsize;
  9:   PetscInt                         m[2];  /* Number of entries in a local block of the row (column) space */
 10:   PetscInt                         mr[2]; /* First incomplete/ragged rank of (row) column space.
 11:                           We expose a blocked ordering to the user because that is what all other PETSc infrastructure uses.
 12:                           With the blocked ordering when the number of processes do not evenly divide the vector size,
 13:                           we still need to be able to convert from PETSc/blocked ordering to VC/VR ordering. */
 14:   El::Grid                        *grid;
 15:   El::DistMatrix<PetscElemScalar> *emat;
 16:   PetscInt                         pivoting; /* 0: no pivoting; 1: partial pivoting; 2: full pivoting */
 17:   El::DistPermutation             *P, *Q;
 18:   PetscBool                        roworiented; /* if true, row oriented input (default) */
 19: } Mat_Elemental;

 21: typedef struct {
 22:   El::Grid *grid;
 23:   PetscInt  grid_refct;
 24: } Mat_Elemental_Grid;

 26: PETSC_INTERN PetscErrorCode MatMatMultSymbolic_Elemental(Mat, Mat, PetscReal, Mat);
 27: PETSC_INTERN PetscErrorCode MatMatMultNumeric_Elemental(Mat, Mat, Mat);

 29: /*
 30:  P2RO, RO2P, E2RO and RO2E convert indices between PETSc <-> (Rank,Offset) <-> Elemental
 31:  (PETSc parallel vectors can be used with Elemental matries without changes)
 32: */
 33: static inline void P2RO(Mat A, PetscInt rc, PetscInt p, PetscInt *rank, PetscInt *offset)
 34: {
 35:   Mat_Elemental *a        = (Mat_Elemental *)A->data;
 36:   PetscInt       critical = a->m[rc] * a->mr[rc];
 37:   if (p < critical) {
 38:     *rank   = p / a->m[rc];
 39:     *offset = p % a->m[rc];
 40:   } else {
 41:     *rank   = a->mr[rc] + (p - critical) / (a->m[rc] - 1);
 42:     *offset = (p - critical) % (a->m[rc] - 1);
 43:   }
 44: }
 45: static inline void RO2P(Mat A, PetscInt rc, PetscInt rank, PetscInt offset, PetscInt *p)
 46: {
 47:   Mat_Elemental *a = (Mat_Elemental *)A->data;
 48:   if (rank < a->mr[rc]) {
 49:     *p = rank * a->m[rc] + offset;
 50:   } else {
 51:     *p = a->mr[rc] * a->m[rc] + (rank - a->mr[rc]) * (a->m[rc] - 1) + offset;
 52:   }
 53: }

 55: static inline void E2RO(Mat A, PetscInt, PetscInt p, PetscInt *rank, PetscInt *offset)
 56: {
 57:   Mat_Elemental *a = (Mat_Elemental *)A->data;
 58:   *rank            = p % a->commsize;
 59:   *offset          = p / a->commsize;
 60: }
 61: static inline void RO2E(Mat A, PetscInt, PetscInt rank, PetscInt offset, PetscInt *e)
 62: {
 63:   Mat_Elemental *a = (Mat_Elemental *)A->data;
 64:   *e               = offset * a->commsize + rank;
 65: }

 67: #endif