Actual source code: matpreallocator.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  1: #include <petsc/private/matimpl.h>      /*I "petscmat.h" I*/
  2: #include <../src/sys/utils/hash.h>

  4: typedef struct {
  5:   PetscHashJK ht;
  6:   PetscInt   *dnz, *onz;
  7: } Mat_Preallocator;

 11: PetscErrorCode MatDestroy_Preallocator(Mat A)
 12: {
 13:   Mat_Preallocator *p = (Mat_Preallocator *) A->data;
 14:   PetscErrorCode    ierr;

 17:   MatStashDestroy_Private(&A->stash);
 18:   PetscHashJKDestroy(&p->ht);
 19:   PetscFree2(p->dnz, p->onz);
 20:   PetscFree(A->data);
 21:   PetscObjectChangeTypeName((PetscObject) A, 0);
 22:   PetscObjectComposeFunction((PetscObject) A, "MatPreallocatorPreallocate_C", NULL);
 23:   return(0);
 24: }

 28: PetscErrorCode MatSetUp_Preallocator(Mat A)
 29: {
 30:   Mat_Preallocator *p = (Mat_Preallocator *) A->data;
 31:   PetscInt          m, bs;
 32:   PetscErrorCode    ierr;

 35:   PetscLayoutSetUp(A->rmap);
 36:   PetscLayoutSetUp(A->cmap);
 37:   MatGetLocalSize(A, &m, NULL);
 38:   PetscHashJKCreate(&p->ht);
 39:   MatGetBlockSize(A, &bs);
 40:   MatStashCreate_Private(PetscObjectComm((PetscObject) A), bs, &A->stash);
 41:   PetscCalloc2(m, &p->dnz, m, &p->onz);
 42:   return(0);
 43: }

 47: PetscErrorCode MatSetValues_Preallocator(Mat A, PetscInt m, const PetscInt *rows, PetscInt n, const PetscInt *cols, const PetscScalar *values, InsertMode addv)
 48: {
 49:   Mat_Preallocator *p = (Mat_Preallocator *) A->data;
 50:   PetscInt          rStart, rEnd, r, cStart, cEnd, c;
 51:   PetscErrorCode    ierr;

 54:   /* TODO: Handle blocksize */
 55:   MatGetOwnershipRange(A, &rStart, &rEnd);
 56:   MatGetOwnershipRangeColumn(A, &cStart, &cEnd);
 57:   for (r = 0; r < m; ++r) {
 58:     PetscHashJKKey  key;
 59:     PetscHashJKIter missing, iter;

 61:     key.j = rows[r];
 62:     if (key.j < 0) continue;
 63:     if ((key.j < rStart) || (key.j >= rEnd)) {
 64:       MatStashValuesRow_Private(&A->stash, key.j, n, cols, values, PETSC_FALSE);
 65:     } else {
 66:       for (c = 0; c < n; ++c) {
 67:         key.k = cols[c];
 68:         if (key.k < 0) continue;
 69:         PetscHashJKPut(p->ht, key, &missing, &iter);
 70:         if (missing) {
 71:           PetscHashJKSet(p->ht, iter, 1);
 72:           if ((key.k >= cStart) && (key.k < cEnd)) ++p->dnz[key.j-rStart];
 73:           else                                     ++p->onz[key.j-rStart];
 74:         }
 75:       }
 76:     }
 77:   }
 78:   return(0);
 79: }

 83: PetscErrorCode MatAssemblyBegin_Preallocator(Mat A, MatAssemblyType type)
 84: {
 85:   PetscInt       nstash, reallocs;

 89:   PetscLayoutSetUp(A->rmap);
 90:   MatStashScatterBegin_Private(A, &A->stash, A->rmap->range);
 91:   MatStashGetInfo_Private(&A->stash, &nstash, &reallocs);
 92:   PetscInfo2(A, "Stash has %D entries, uses %D mallocs.\n", nstash, reallocs);
 93:   return(0);
 94: }

 98: PetscErrorCode MatAssemblyEnd_Preallocator(Mat A, MatAssemblyType type)
 99: {
100:   PetscScalar   *val;
101:   PetscInt      *row, *col;
102:   PetscInt       i, j, rstart, ncols, flg;
103:   PetscMPIInt    n;

107:   while (1) {
108:     MatStashScatterGetMesg_Private(&A->stash, &n, &row, &col, &val, &flg);
109:     if (!flg) break;

111:     for (i = 0; i < n; ) {
112:       /* Now identify the consecutive vals belonging to the same row */
113:       for (j = i, rstart = row[j]; j < n; j++) {
114:         if (row[j] != rstart) break;
115:       }
116:       if (j < n) ncols = j-i;
117:       else       ncols = n-i;
118:       /* Now assemble all these values with a single function call */
119:       MatSetValues_Preallocator(A, 1, row+i, ncols, col+i, val+i, INSERT_VALUES);
120:       i = j;
121:     }
122:   }
123:   MatStashScatterEnd_Private(&A->stash);
124:   return(0);
125: }

129: PetscErrorCode MatView_Preallocator(Mat A, PetscViewer viewer)
130: {
132:   return(0);
133: }

137: PetscErrorCode MatSetOption_Preallocator(Mat A, MatOption op, PetscBool flg)
138: {
140:   return(0);
141: }

145: PetscErrorCode MatPreallocatorPreallocate_Preallocator(Mat mat, PetscBool fill, Mat A)
146: {
147:   Mat_Preallocator *p = (Mat_Preallocator *) mat->data;
148:   PetscInt         *udnz = NULL, *uonz = NULL;
149:   PetscInt          bs;
150:   PetscErrorCode    ierr;

153:   MatGetBlockSize(mat, &bs);
154:   MatXAIJSetPreallocation(A, bs, p->dnz, p->onz, udnz, uonz);
155:   MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
156:   return(0);
157: }

161: /*@
162:   MatPreallocatorPreallocate - Preallocates the input matrix, optionally filling it with zeros

164:   Input Parameter:
165: + mat  - the preallocator
166: - fill - fill the matrix with zeros

168:   Output Parameter:
169: . A    - the matrix

171:   Level: advanced

173: .seealso: MATPREALLOCATOR
174: @*/
175: PetscErrorCode MatPreallocatorPreallocate(Mat mat, PetscBool fill, Mat A)
176: {

182:   PetscUseMethod(mat, "MatPreallocatorPreallocate_C", (Mat,PetscBool,Mat),(mat,fill,A));
183:   return(0);
184: }

186: /*MC
187:    MATPREALLOCATOR - MATPREALLOCATOR = "preallocator" - A matrix type to be used for computing a matrix preallocation.

189:    Operations Provided:
190: .  MatSetValues()

192:    Options Database Keys:
193: . -mat_type preallocator - sets the matrix type to "preallocator" during a call to MatSetFromOptions()

195:   Level: advanced

197: .seealso: Mat

199: M*/

203: PETSC_EXTERN PetscErrorCode MatCreate_Preallocator(Mat A)
204: {
205:   Mat_Preallocator *p;
206:   PetscErrorCode    ierr;

209:   PetscNewLog(A, &p);
210:   A->data = (void *) p;

212:   p->ht  = NULL;
213:   p->dnz = NULL;
214:   p->onz = NULL;

216:   /* matrix ops */
217:   PetscMemzero(A->ops, sizeof(struct _MatOps));
218:   A->ops->destroy                 = MatDestroy_Preallocator;
219:   A->ops->setup                   = MatSetUp_Preallocator;
220:   A->ops->setvalues               = MatSetValues_Preallocator;
221:   A->ops->assemblybegin           = MatAssemblyBegin_Preallocator;
222:   A->ops->assemblyend             = MatAssemblyEnd_Preallocator;
223:   A->ops->view                    = MatView_Preallocator;
224:   A->ops->setoption               = MatSetOption_Preallocator;

226:   /* special MATPREALLOCATOR functions */
227:   PetscObjectComposeFunction((PetscObject) A, "MatPreallocatorPreallocate_C", MatPreallocatorPreallocate_Preallocator);
228:   PetscObjectChangeTypeName((PetscObject) A, MATPREALLOCATOR);
229:   return(0);
230: }