Actual source code: aijbaij.c
2: #include src/mat/impls/baij/seq/baij.h
7: PetscErrorCode MatConvert_SeqBAIJ_SeqAIJ(Mat A,const MatType newtype,Mat *newmat)
8: {
9: Mat B;
10: Mat_SeqBAIJ *a = (Mat_SeqBAIJ*)A->data;
12: PetscInt bs = A->bs,*ai = a->i,*aj = a->j,n = A->M/bs,i,j,k;
13: PetscInt *rowlengths,*rows,*cols,maxlen = 0,ncols;
14: PetscScalar *aa = a->a;
17: PetscMalloc(n*bs*sizeof(int),&rowlengths);
18: for (i=0; i<n; i++) {
19: maxlen = PetscMax(maxlen,(ai[i+1] - ai[i]));
20: for (j=0; j<bs; j++) {
21: rowlengths[i*bs+j] = bs*(ai[i+1] - ai[i]);
22: }
23: }
24: MatCreate(A->comm,A->m,A->n,A->m,A->n,&B);
25: MatSetType(B,newtype);
26: MatSeqAIJSetPreallocation(B,0,rowlengths);
27: MatSetOption(B,MAT_COLUMN_ORIENTED);
28: MatSetOption(B,MAT_ROWS_SORTED);
29: MatSetOption(B,MAT_COLUMNS_SORTED);
30: PetscFree(rowlengths);
32: PetscMalloc(bs*sizeof(int),&rows);
33: PetscMalloc(bs*maxlen*sizeof(int),&cols);
34: for (i=0; i<n; i++) {
35: for (j=0; j<bs; j++) {
36: rows[j] = i*bs+j;
37: }
38: ncols = ai[i+1] - ai[i];
39: for (k=0; k<ncols; k++) {
40: for (j=0; j<bs; j++) {
41: cols[k*bs+j] = bs*(*aj) + j;
42: }
43: aj++;
44: }
45: MatSetValues(B,bs,rows,bs*ncols,cols,aa,INSERT_VALUES);
46: aa += ncols*bs*bs;
47: }
48: PetscFree(cols);
49: PetscFree(rows);
50: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
51: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
52:
53: B->bs = A->bs;
55: /* Fake support for "inplace" convert. */
56: if (*newmat == A) {
57: MatDestroy(A);
58: }
59: *newmat = B;
60: return(0);
61: }
64: #include src/mat/impls/aij/seq/aij.h
69: PetscErrorCode MatConvert_SeqAIJ_SeqBAIJ(Mat A,const MatType newtype,Mat *newmat)
70: {
71: Mat B;
72: Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data;
73: Mat_SeqBAIJ *b;
75: PetscInt *ai=a->i,m=A->M,n=A->N,i,*rowlengths;
78: if (n != m) SETERRQ(PETSC_ERR_ARG_WRONG,"Matrix must be square");
80: PetscMalloc(m*sizeof(int),&rowlengths);
81: for (i=0; i<m; i++) {
82: rowlengths[i] = ai[i+1] - ai[i];
83: }
84: MatCreate(A->comm,m,n,m,n,&B);
85: MatSetType(B,newtype);
86: MatSeqBAIJSetPreallocation(B,1,0,rowlengths);
87: PetscFree(rowlengths);
89: MatSetOption(B,MAT_ROW_ORIENTED);
90: MatSetOption(B,MAT_ROWS_SORTED);
91: MatSetOption(B,MAT_COLUMNS_SORTED);
92:
93: b = (Mat_SeqBAIJ*)(B->data);
95: PetscMemcpy(b->i,a->i,(m+1)*sizeof(int));
96: PetscMemcpy(b->ilen,a->ilen,m*sizeof(int));
97: PetscMemcpy(b->j,a->j,a->nz*sizeof(int));
98: PetscMemcpy(b->a,a->a,a->nz*sizeof(MatScalar));
99:
100: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
101: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
103: /* Fake support for "inplace" convert. */
104: if (*newmat == A) {
105: MatDestroy(A);
106: }
107: *newmat = B;
108: return(0);
109: }