Actual source code: matelemimpl.h
petsc-3.4.5 2014-06-29
1: #if !defined(_matelemimpl_h)
2: #define _matelemimpl_h
4: #include <elemental.hpp>
5: #include <petsc-private/matimpl.h>
7: #if defined(PETSC_USE_COMPLEX)
8: typedef elem::Complex<PetscReal> PetscElemScalar;
9: #else
10: typedef PetscScalar PetscElemScalar;
11: #endif
13: typedef struct {
14: PetscInt commsize;
15: PetscInt m[2]; /* Number of entries in a local block of the row (column) space */
16: PetscInt mr[2]; /* First incomplete/ragged rank of (row) column space.
17: We expose a blocked ordering to the user because that is what all other PETSc infrastructure uses.
18: With the blocked ordering when the number of processes do not evenly divide the vector size,
19: we still need to be able to convert from PETSc/blocked ordering to VC/VR ordering. */
20: elem::Grid *grid;
21: elem::DistMatrix<PetscElemScalar> *emat;
22: elem::Matrix<PetscElemScalar> *esubmat; /* Used for adding off-proc matrix entries */
23: elem::AxpyInterface<PetscElemScalar> *interface;
24: elem::DistMatrix<PetscInt,elem::VC,elem::STAR> *pivot; /* pivot vector representing the pivot matrix P in PA = LU */
25: } Mat_Elemental;
27: typedef struct {
28: elem::Grid *grid;
29: PetscInt grid_refct;
30: } Mat_Elemental_Grid;
32: PETSC_STATIC_INLINE void P2RO(Mat A,PetscInt rc,PetscInt p,PetscInt *rank,PetscInt *offset)
33: {
34: Mat_Elemental *a = (Mat_Elemental*)A->data;
35: PetscInt critical = a->m[rc]*a->mr[rc];
36: if (p < critical) {
37: *rank = p / a->m[rc];
38: *offset = p % a->m[rc];
39: } else {
40: *rank = a->mr[rc] + (p - critical) / (a->m[rc] - 1);
41: *offset = (p - critical) % (a->m[rc] - 1);
42: }
43: }
44: PETSC_STATIC_INLINE void RO2P(Mat A,PetscInt rc,PetscInt rank,PetscInt offset,PetscInt *p)
45: {
46: Mat_Elemental *a = (Mat_Elemental*)A->data;
47: if (rank < a->mr[rc]) {
48: *p = rank*a->m[rc] + offset;
49: } else {
50: *p = a->mr[rc]*a->m[rc] + (rank - a->mr[rc])*(a->m[rc]-1) + offset;
51: }
52: }
54: PETSC_STATIC_INLINE void E2RO(Mat A,PetscInt rc,PetscInt p,PetscInt *rank,PetscInt *offset)
55: {
56: Mat_Elemental *a = (Mat_Elemental*)A->data;
57: *rank = p % a->commsize;
58: *offset = p / a->commsize;
59: }
60: PETSC_STATIC_INLINE void RO2E(Mat A,PetscInt rc,PetscInt rank,PetscInt offset,PetscInt *e)
61: {
62: Mat_Elemental *a = (Mat_Elemental*)A->data;
63: *e = offset * a->commsize + rank;
64: }
66: #endif