1 /*$Id: aijbaij.c,v 1.1 2000/06/06 19:00:44 bsmith Exp bsmith $*/ 2 3 #include "src/mat/impls/baij/seq/baij.h" 4 5 #undef __FUNC__ 6 #define __FUNC__ /*<a name="MatConvert_SeqBAIJ_SeqAIJ"></a>*/"MatConvert_SeqBAI_SeqAIJ" 7 int MatConvert_SeqBAIJ_SeqAIJ(Mat A,MatType newtype,Mat *B) 8 { 9 Mat_SeqBAIJ *a = (Mat_SeqBAIJ*)A->data; 10 int ierr,bs = a->bs,*ai = a->i,*aj = a->j,n = A->M/bs,i,j,k; 11 int *rowlengths,*rows,*cols,maxlen = 0,ncols; 12 Scalar *aa = a->a; 13 14 PetscFunctionBegin; 15 rowlengths = (int*)PetscMalloc(n*bs*sizeof(int));CHKPTRQ(rowlengths); 16 for (i=0; i<n; i++) { 17 maxlen = PetscMax(maxlen,(ai[i+1] - ai[i])); 18 for (j=0; j<bs; j++) { 19 rowlengths[i*bs+j] = bs*(ai[i+1] - ai[i]); 20 } 21 } 22 ierr = MatCreateSeqAIJ(PETSC_COMM_SELF,A->m,A->n,0,rowlengths,B);CHKERRQ(ierr); 23 ierr = MatSetOption(*B,MAT_COLUMN_ORIENTED);CHKERRQ(ierr); 24 ierr = MatSetOption(*B,MAT_ROWS_SORTED);CHKERRQ(ierr); 25 ierr = MatSetOption(*B,MAT_COLUMNS_SORTED);CHKERRQ(ierr); 26 ierr = PetscFree(rowlengths);CHKERRQ(ierr); 27 28 rows = (int*)PetscMalloc(bs*sizeof(int));CHKPTRQ(rows); 29 cols = (int*)PetscMalloc(bs*maxlen*sizeof(int));CHKPTRQ(cols); 30 for (i=0; i<n; i++) { 31 for (j=0; j<bs; j++) { 32 rows[j] = i*bs+j; 33 } 34 ncols = ai[i+1] - ai[i]; 35 for (k=0; k<ncols; k++) { 36 for (j=0; j<bs; j++) { 37 cols[k*bs+j] = bs*(*aj) + j; 38 } 39 aj++; 40 } 41 ierr = MatSetValues(*B,bs,rows,bs*ncols,cols,aa,INSERT_VALUES); 42 aa += ncols*bs*bs; 43 } 44 ierr = PetscFree(cols);CHKERRQ(ierr); 45 ierr = PetscFree(rows);CHKERRQ(ierr); 46 ierr = MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 47 ierr = MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 48 PetscFunctionReturn(0); 49 } 50 51 52