Actual source code: matstashspace.c

  2: #define PETSCMAT_DLL

 4:  #include include/private/matimpl.h

  6: /* Get new PetscMatStashSpace into the existing space */
  9: PetscErrorCode PetscMatStashSpaceGet(PetscInt bs2,PetscInt n,PetscMatStashSpace *space)
 10: {
 11:   PetscMatStashSpace a;
 12:   PetscErrorCode     ierr;
 13: 
 15:   if (!n) return(0);

 17:   PetscMalloc(sizeof(struct _MatStashSpace),&a);
 18:   PetscMalloc(n*(bs2*sizeof(MatScalar)+2*sizeof(PetscInt)),&(a->space_head));
 19:   a->val              = a->space_head;
 20:   a->idx              = (PetscInt*)(a->val + bs2*n);
 21:   a->idy              = (PetscInt*)(a->idx + n);
 22:   a->local_remaining  = n;
 23:   a->local_used       = 0;
 24:   a->total_space_size = 0;
 25:   a->next             = PETSC_NULL;

 27:   if (*space){
 28:     (*space)->next      = a;
 29:     a->total_space_size = (*space)->total_space_size;
 30:   }
 31:   a->total_space_size += n;
 32:   *space               = a;
 33:   return(0);
 34: }

 36: /* Copy the values in space into arrays val, idx and idy. Then destroy space */
 39: PetscErrorCode PetscMatStashSpaceContiguous(PetscInt bs2,PetscMatStashSpace *space,PetscScalar *val,PetscInt *idx,PetscInt *idy)
 40: {
 41:   PetscMatStashSpace a;
 42:   PetscErrorCode     ierr;

 45:   while ((*space) != PETSC_NULL){
 46:     a    = (*space)->next;
 47:     PetscMemcpy(val,(*space)->val,((*space)->local_used*bs2)*sizeof(MatScalar));
 48:     val += bs2*(*space)->local_used;
 49:     PetscMemcpy(idx,(*space)->idx,((*space)->local_used)*sizeof(PetscInt));
 50:     idx += (*space)->local_used;
 51:     PetscMemcpy(idy,(*space)->idy,((*space)->local_used)*sizeof(PetscInt));
 52:     idy += (*space)->local_used;
 53: 
 54:      PetscFree((*space)->space_head);
 55:      PetscFree(*space);
 56:     *space = a;
 57:   }
 58:   return(0);
 59: }

 63: PetscErrorCode PetscMatStashSpaceDestroy(PetscMatStashSpace space)
 64: {
 65:   PetscMatStashSpace a;
 66:   PetscErrorCode     ierr;

 69:   while (space != PETSC_NULL){
 70:     a     = space->next;
 71:     PetscFree(space->space_head);
 72:     PetscFree(space);
 73:     space = a;
 74:   }
 75:   return(0);
 76: }