Actual source code: matstashspace.c

  2: #define PETSCMAT_DLL

 4:  #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:   PetscMalloc3(n*bs2,PetscScalar,&(a->space_head),n,PetscInt,&a->idx,n,PetscInt,&a->idy);
 19:   a->val              = a->space_head;
 20:   a->local_remaining  = n;
 21:   a->local_used       = 0;
 22:   a->total_space_size = 0;
 23:   a->next             = PETSC_NULL;

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

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

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

 61: PetscErrorCode PetscMatStashSpaceDestroy(PetscMatStashSpace space)
 62: {
 63:   PetscMatStashSpace a;
 64:   PetscErrorCode     ierr;

 67:   while (space != PETSC_NULL){
 68:     a     = space->next;
 69:     PetscFree3(space->space_head,space->idx,space->idy);
 70:     PetscFree(space);
 71:     space = a;
 72:   }
 73:   return(0);
 74: }