Actual source code: aijbaij.c
1: /*$Id: aijbaij.c,v 1.9 2001/08/07 03:02:55 balay Exp $*/
3: #include src/mat/impls/baij/seq/baij.h
5: EXTERN_C_BEGIN
8: int MatConvert_SeqBAIJ_SeqAIJ(Mat A,const MatType newtype,Mat *newmat) {
9: Mat B;
10: Mat_SeqBAIJ *a = (Mat_SeqBAIJ*)A->data;
11: int ierr,bs = a->bs,*ai = a->i,*aj = a->j,n = A->M/bs,i,j,k;
12: int *rowlengths,*rows,*cols,maxlen = 0,ncols;
13: PetscScalar *aa = a->a;
16: PetscMalloc(n*bs*sizeof(int),&rowlengths);
17: for (i=0; i<n; i++) {
18: maxlen = PetscMax(maxlen,(ai[i+1] - ai[i]));
19: for (j=0; j<bs; j++) {
20: rowlengths[i*bs+j] = bs*(ai[i+1] - ai[i]);
21: }
22: }
23: MatCreate(A->comm,A->m,A->n,A->m,A->n,&B);
24: MatSetType(B,newtype);
25: MatSeqAIJSetPreallocation(B,0,rowlengths);
26: MatSetOption(B,MAT_COLUMN_ORIENTED);
27: MatSetOption(B,MAT_ROWS_SORTED);
28: MatSetOption(B,MAT_COLUMNS_SORTED);
29: PetscFree(rowlengths);
31: PetscMalloc(bs*sizeof(int),&rows);
32: PetscMalloc(bs*maxlen*sizeof(int),&cols);
33: for (i=0; i<n; i++) {
34: for (j=0; j<bs; j++) {
35: rows[j] = i*bs+j;
36: }
37: ncols = ai[i+1] - ai[i];
38: for (k=0; k<ncols; k++) {
39: for (j=0; j<bs; j++) {
40: cols[k*bs+j] = bs*(*aj) + j;
41: }
42: aj++;
43: }
44: MatSetValues(B,bs,rows,bs*ncols,cols,aa,INSERT_VALUES);
45: aa += ncols*bs*bs;
46: }
47: PetscFree(cols);
48: PetscFree(rows);
49: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
50: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
52: /* Fake support for "inplace" convert. */
53: if (*newmat == A) {
54: MatDestroy(A);
55: }
56: *newmat = B;
57: return(0);
58: }
59: EXTERN_C_END
61: #include src/mat/impls/aij/seq/aij.h
63: EXTERN_C_BEGIN
66: int MatConvert_SeqAIJ_SeqBAIJ(Mat A,const MatType newtype,Mat *newmat) {
67: Mat B;
68: Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data;
69: Mat_SeqBAIJ *b;
70: int ierr,*ai=a->i,m=A->M,n=A->N,i,*rowlengths;
73: if (n != m) SETERRQ(PETSC_ERR_ARG_WRONG,"Matrix must be square");
75: PetscMalloc(m*sizeof(int),&rowlengths);
76: for (i=0; i<m; i++) {
77: rowlengths[i] = ai[i+1] - ai[i];
78: }
79: MatCreate(A->comm,m,n,m,n,&B);
80: MatSetType(B,newtype);
81: MatSeqBAIJSetPreallocation(B,1,0,rowlengths);
82: PetscFree(rowlengths);
84: MatSetOption(B,MAT_ROW_ORIENTED);
85: MatSetOption(B,MAT_ROWS_SORTED);
86: MatSetOption(B,MAT_COLUMNS_SORTED);
87:
88: b = (Mat_SeqBAIJ*)(B->data);
90: PetscMemcpy(b->i,a->i,(m+1)*sizeof(int));
91: PetscMemcpy(b->ilen,a->ilen,m*sizeof(int));
92: PetscMemcpy(b->j,a->j,a->nz*sizeof(int));
93: PetscMemcpy(b->a,a->a,a->nz*sizeof(MatScalar));
94:
95: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
96: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
98: /* Fake support for "inplace" convert. */
99: if (*newmat == A) {
100: MatDestroy(A);
101: }
102: *newmat = B;
103: return(0);
104: }
105: EXTERN_C_END