1*75cae7c1SHong Zhang 2*75cae7c1SHong Zhang #define PETSCMAT_DLL 3*75cae7c1SHong Zhang 4*75cae7c1SHong Zhang #include "src/mat/matimpl.h" 5*75cae7c1SHong Zhang #include "src/mat/utils/matstashspace.h" 6*75cae7c1SHong Zhang 7*75cae7c1SHong Zhang /* Get new PetscMatStashSpace into the existing space */ 8*75cae7c1SHong Zhang #undef __FUNCT__ 9*75cae7c1SHong Zhang #define __FUNCT__ "PetscMatStashSpaceGet" 10*75cae7c1SHong Zhang PetscErrorCode PetscMatStashSpaceGet(PetscInt bs2,PetscInt n,PetscMatStashSpace *space) 11*75cae7c1SHong Zhang { 12*75cae7c1SHong Zhang PetscMatStashSpace a; 13*75cae7c1SHong Zhang PetscErrorCode ierr; 14*75cae7c1SHong Zhang 15*75cae7c1SHong Zhang PetscFunctionBegin; 16*75cae7c1SHong Zhang if (!n) PetscFunctionReturn(0); 17*75cae7c1SHong Zhang 18*75cae7c1SHong Zhang ierr = PetscMalloc(sizeof(struct _MatStashSpace),&a);CHKERRQ(ierr); 19*75cae7c1SHong Zhang ierr = PetscMalloc(n*(bs2*sizeof(MatScalar)+2*sizeof(PetscInt)),&(a->space_head));CHKERRQ(ierr); 20*75cae7c1SHong Zhang a->val = a->space_head; 21*75cae7c1SHong Zhang a->idx = (PetscInt*)(a->val + bs2*n); 22*75cae7c1SHong Zhang a->idy = (PetscInt*)(a->idx + n); 23*75cae7c1SHong Zhang a->local_remaining = n; 24*75cae7c1SHong Zhang a->local_used = 0; 25*75cae7c1SHong Zhang a->total_space_size = 0; 26*75cae7c1SHong Zhang a->next = PETSC_NULL; 27*75cae7c1SHong Zhang 28*75cae7c1SHong Zhang if (*space){ 29*75cae7c1SHong Zhang (*space)->next = a; 30*75cae7c1SHong Zhang a->total_space_size = (*space)->total_space_size; 31*75cae7c1SHong Zhang } 32*75cae7c1SHong Zhang a->total_space_size += n; 33*75cae7c1SHong Zhang *space = a; 34*75cae7c1SHong Zhang /* printf(" 2. total_space_size : %d\n",a->total_space_size); */ 35*75cae7c1SHong Zhang /* 36*75cae7c1SHong Zhang PetscMPIInt rank; 37*75cae7c1SHong Zhang ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 38*75cae7c1SHong Zhang printf("[%d] MatStashSpaceCreate space %p, next %p, total_space_size: %d, space->val %p\n",rank,*space,(*space)->next,(*space)->total_space_size,(*space)->val); 39*75cae7c1SHong Zhang */ 40*75cae7c1SHong Zhang PetscFunctionReturn(0); 41*75cae7c1SHong Zhang } 42*75cae7c1SHong Zhang 43*75cae7c1SHong Zhang /* Copy the values in space into arrays val, idx and idy. Then destroy space */ 44*75cae7c1SHong Zhang #undef __FUNCT__ 45*75cae7c1SHong Zhang #define __FUNCT__ "PetscMatStashSpaceContiguous" 46*75cae7c1SHong Zhang PetscErrorCode PetscMatStashSpaceContiguous(PetscInt bs2,PetscMatStashSpace *space,PetscScalar *val,PetscInt *idx,PetscInt *idy) 47*75cae7c1SHong Zhang { 48*75cae7c1SHong Zhang PetscMatStashSpace a; 49*75cae7c1SHong Zhang PetscErrorCode ierr; 50*75cae7c1SHong Zhang 51*75cae7c1SHong Zhang PetscFunctionBegin; 52*75cae7c1SHong Zhang while ((*space) != PETSC_NULL){ 53*75cae7c1SHong Zhang a = (*space)->next; 54*75cae7c1SHong Zhang ierr = PetscMemcpy(val,(*space)->val,((*space)->local_used*bs2)*sizeof(MatScalar));CHKERRQ(ierr); 55*75cae7c1SHong Zhang val += bs2*(*space)->local_used; 56*75cae7c1SHong Zhang ierr = PetscMemcpy(idx,(*space)->idx,((*space)->local_used)*sizeof(PetscInt));CHKERRQ(ierr); 57*75cae7c1SHong Zhang idx += (*space)->local_used; 58*75cae7c1SHong Zhang ierr = PetscMemcpy(idy,(*space)->idy,((*space)->local_used)*sizeof(PetscInt));CHKERRQ(ierr); 59*75cae7c1SHong Zhang idy += (*space)->local_used; 60*75cae7c1SHong Zhang 61*75cae7c1SHong Zhang ierr = PetscFree((*space)->space_head);CHKERRQ(ierr); 62*75cae7c1SHong Zhang ierr = PetscFree(*space);CHKERRQ(ierr); 63*75cae7c1SHong Zhang *space = a; 64*75cae7c1SHong Zhang } 65*75cae7c1SHong Zhang PetscFunctionReturn(0); 66*75cae7c1SHong Zhang } 67*75cae7c1SHong Zhang 68*75cae7c1SHong Zhang #undef __FUNCT__ 69*75cae7c1SHong Zhang #define __FUNCT__ "PetscMatStashSpaceDestroy" 70*75cae7c1SHong Zhang PetscErrorCode PetscMatStashSpaceDestroy(PetscMatStashSpace space) 71*75cae7c1SHong Zhang { 72*75cae7c1SHong Zhang PetscMatStashSpace a; 73*75cae7c1SHong Zhang PetscErrorCode ierr; 74*75cae7c1SHong Zhang 75*75cae7c1SHong Zhang PetscFunctionBegin; 76*75cae7c1SHong Zhang while (space != PETSC_NULL){ 77*75cae7c1SHong Zhang a = space->next; 78*75cae7c1SHong Zhang ierr = PetscFree(space->space_head);CHKERRQ(ierr); 79*75cae7c1SHong Zhang ierr = PetscFree(space);CHKERRQ(ierr); 80*75cae7c1SHong Zhang space = a; 81*75cae7c1SHong Zhang } 82*75cae7c1SHong Zhang PetscFunctionReturn(0); 83*75cae7c1SHong Zhang } 84