Actual source code: matelemimpl.h
petsc-3.13.6 2020-09-29
1: #if !defined(_matelemimpl_h)
2: #define _matelemimpl_h
4: #include <El.hpp>
5: #include <petsc/private/matimpl.h>
6: #include <petscmatelemental.h>
8: typedef struct {
9: PetscInt commsize;
10: PetscInt m[2]; /* Number of entries in a local block of the row (column) space */
11: PetscInt mr[2]; /* First incomplete/ragged rank of (row) column space.
12: We expose a blocked ordering to the user because that is what all other PETSc infrastructure uses.
13: With the blocked ordering when the number of processes do not evenly divide the vector size,
14: we still need to be able to convert from PETSc/blocked ordering to VC/VR ordering. */
15: El::Grid *grid;
16: El::DistMatrix<PetscElemScalar> *emat;
17: PetscInt pivoting; /* 0: no pivoting; 1: partial pivoting; 2: full pivoting */
18: El::DistPermutation *P,*Q;
19: PetscBool roworiented; /* if true, row oriented input (default) */
20: } Mat_Elemental;
22: typedef struct {
23: El::Grid *grid;
24: PetscInt grid_refct;
25: } Mat_Elemental_Grid;
27: PETSC_INTERN PetscErrorCode MatMatMultSymbolic_Elemental(Mat,Mat,PetscReal,Mat);
28: PETSC_INTERN PetscErrorCode MatMatMultNumeric_Elemental(Mat,Mat,Mat);
30: /*
31: P2RO, RO2P, E2RO and RO2E convert indices between PETSc <-> (Rank,Offset) <-> Elemental
32: (PETSc parallel vectors can be used with Elemental matries without changes)
33: */
34: PETSC_STATIC_INLINE void P2RO(Mat A,PetscInt rc,PetscInt p,PetscInt *rank,PetscInt *offset)
35: {
36: Mat_Elemental *a = (Mat_Elemental*)A->data;
37: PetscInt critical = a->m[rc]*a->mr[rc];
38: if (p < critical) {
39: *rank = p / a->m[rc];
40: *offset = p % a->m[rc];
41: } else {
42: *rank = a->mr[rc] + (p - critical) / (a->m[rc] - 1);
43: *offset = (p - critical) % (a->m[rc] - 1);
44: }
45: }
46: PETSC_STATIC_INLINE void RO2P(Mat A,PetscInt rc,PetscInt rank,PetscInt offset,PetscInt *p)
47: {
48: Mat_Elemental *a = (Mat_Elemental*)A->data;
49: if (rank < a->mr[rc]) {
50: *p = rank*a->m[rc] + offset;
51: } else {
52: *p = a->mr[rc]*a->m[rc] + (rank - a->mr[rc])*(a->m[rc]-1) + offset;
53: }
54: }
56: PETSC_STATIC_INLINE void E2RO(Mat A,PetscInt rc,PetscInt p,PetscInt *rank,PetscInt *offset)
57: {
58: Mat_Elemental *a = (Mat_Elemental*)A->data;
59: *rank = p % a->commsize;
60: *offset = p / a->commsize;
61: }
62: PETSC_STATIC_INLINE void RO2E(Mat A,PetscInt rc,PetscInt rank,PetscInt offset,PetscInt *e)
63: {
64: Mat_Elemental *a = (Mat_Elemental*)A->data;
65: *e = offset * a->commsize + rank;
66: }
68: #endif