1 #define PETSCMAT_DLL 2 3 #include "src/mat/impls/baij/seq/baij.h" 4 5 EXTERN_C_BEGIN 6 #undef __FUNCT__ 7 #define __FUNCT__ "MatConvert_SeqBAIJ_SeqAIJ" 8 PetscErrorCode PETSCMAT_DLLEXPORT MatConvert_SeqBAIJ_SeqAIJ(Mat A, MatType newtype,MatReuse reuse,Mat *newmat) 9 { 10 Mat B; 11 Mat_SeqBAIJ *a = (Mat_SeqBAIJ*)A->data; 12 PetscErrorCode ierr; 13 PetscInt bs = A->bs,*ai = a->i,*aj = a->j,n = A->M/bs,i,j,k; 14 PetscInt *rowlengths,*rows,*cols,maxlen = 0,ncols; 15 PetscScalar *aa = a->a; 16 17 PetscFunctionBegin; 18 ierr = PetscMalloc(n*bs*sizeof(PetscInt),&rowlengths);CHKERRQ(ierr); 19 for (i=0; i<n; i++) { 20 maxlen = PetscMax(maxlen,(ai[i+1] - ai[i])); 21 for (j=0; j<bs; j++) { 22 rowlengths[i*bs+j] = bs*(ai[i+1] - ai[i]); 23 } 24 } 25 ierr = MatCreate(A->comm,&B);CHKERRQ(ierr); 26 ierr = MatSetSizes(B,A->m,A->n,A->m,A->n);CHKERRQ(ierr); 27 ierr = MatSetType(B,newtype);CHKERRQ(ierr); 28 ierr = MatSeqAIJSetPreallocation(B,0,rowlengths);CHKERRQ(ierr); 29 ierr = MatSetOption(B,MAT_COLUMN_ORIENTED);CHKERRQ(ierr); 30 ierr = MatSetOption(B,MAT_ROWS_SORTED);CHKERRQ(ierr); 31 ierr = MatSetOption(B,MAT_COLUMNS_SORTED);CHKERRQ(ierr); 32 ierr = PetscFree(rowlengths);CHKERRQ(ierr); 33 34 ierr = PetscMalloc(bs*sizeof(PetscInt),&rows);CHKERRQ(ierr); 35 ierr = PetscMalloc(bs*maxlen*sizeof(PetscInt),&cols);CHKERRQ(ierr); 36 for (i=0; i<n; i++) { 37 for (j=0; j<bs; j++) { 38 rows[j] = i*bs+j; 39 } 40 ncols = ai[i+1] - ai[i]; 41 for (k=0; k<ncols; k++) { 42 for (j=0; j<bs; j++) { 43 cols[k*bs+j] = bs*(*aj) + j; 44 } 45 aj++; 46 } 47 ierr = MatSetValues(B,bs,rows,bs*ncols,cols,aa,INSERT_VALUES);CHKERRQ(ierr); 48 aa += ncols*bs*bs; 49 } 50 ierr = PetscFree(cols);CHKERRQ(ierr); 51 ierr = PetscFree(rows);CHKERRQ(ierr); 52 ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 53 ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 54 55 B->bs = A->bs; 56 57 if (reuse == MAT_REUSE_MATRIX) { 58 ierr = MatHeaderReplace(A,B);CHKERRQ(ierr); 59 } else { 60 *newmat = B; 61 } 62 PetscFunctionReturn(0); 63 } 64 EXTERN_C_END 65 66 #include "src/mat/impls/aij/seq/aij.h" 67 68 EXTERN_C_BEGIN 69 #undef __FUNCT__ 70 #define __FUNCT__ "MatConvert_SeqAIJ_SeqBAIJ" 71 PetscErrorCode PETSCMAT_DLLEXPORT MatConvert_SeqAIJ_SeqBAIJ(Mat A, MatType newtype,MatReuse reuse,Mat *newmat) 72 { 73 Mat B; 74 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 75 Mat_SeqBAIJ *b; 76 PetscErrorCode ierr; 77 PetscInt *ai=a->i,m=A->M,n=A->N,i,*rowlengths; 78 79 PetscFunctionBegin; 80 if (n != m) SETERRQ(PETSC_ERR_ARG_WRONG,"Matrix must be square"); 81 82 ierr = PetscMalloc(m*sizeof(PetscInt),&rowlengths);CHKERRQ(ierr); 83 for (i=0; i<m; i++) { 84 rowlengths[i] = ai[i+1] - ai[i]; 85 } 86 ierr = MatCreate(A->comm,&B);CHKERRQ(ierr); 87 ierr = MatSetSizes(B,m,n,m,n);CHKERRQ(ierr); 88 ierr = MatSetType(B,newtype);CHKERRQ(ierr); 89 ierr = MatSeqBAIJSetPreallocation_SeqBAIJ(B,1,0,rowlengths);CHKERRQ(ierr); 90 ierr = PetscFree(rowlengths);CHKERRQ(ierr); 91 92 ierr = MatSetOption(B,MAT_ROW_ORIENTED);CHKERRQ(ierr); 93 ierr = MatSetOption(B,MAT_ROWS_SORTED);CHKERRQ(ierr); 94 ierr = MatSetOption(B,MAT_COLUMNS_SORTED);CHKERRQ(ierr); 95 96 b = (Mat_SeqBAIJ*)(B->data); 97 98 ierr = PetscMemcpy(b->i,a->i,(m+1)*sizeof(PetscInt));CHKERRQ(ierr); 99 ierr = PetscMemcpy(b->ilen,a->ilen,m*sizeof(PetscInt));CHKERRQ(ierr); 100 ierr = PetscMemcpy(b->j,a->j,a->nz*sizeof(PetscInt));CHKERRQ(ierr); 101 ierr = PetscMemcpy(b->a,a->a,a->nz*sizeof(MatScalar));CHKERRQ(ierr); 102 103 ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 104 ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 105 106 if (reuse == MAT_REUSE_MATRIX) { 107 ierr = MatHeaderReplace(A,B);CHKERRQ(ierr); 108 } else { 109 *newmat = B; 110 } 111 PetscFunctionReturn(0); 112 } 113 EXTERN_C_END 114