1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: ij.c,v 1.31 1999/05/04 20:31:42 balay Exp balay $"; 3 #endif 4 5 #include "src/mat/impls/aij/seq/aij.h" 6 7 #undef __FUNC__ 8 #define __FUNC__ "MatToSymmetricIJ_SeqAIJ" 9 /* 10 MatToSymmetricIJ_SeqAIJ - Convert a (generally nonsymmetric) sparse AIJ matrix 11 to IJ format (ignore the "A" part) Allocates the space needed. Uses only 12 the lower triangular part of the matrix. 13 14 Description: 15 Take the data in the row-oriented sparse storage and build the 16 IJ data for the Matrix. Return 0 on success, row + 1 on failure 17 at that row. Produces the ij for a symmetric matrix by only using 18 the lower triangular part of the matrix. 19 20 Input Parameters: 21 . Matrix - matrix to convert 22 . shiftin - the shift for the original matrix (0 or 1) 23 . shiftout - the shift required for the ordering routine (0 or 1) 24 25 Output Parameters: 26 . ia - ia part of IJ representation (row information) 27 . ja - ja part (column indices) 28 29 Notes: 30 Both ia and ja may be freed with PetscFree(); 31 This routine is provided for ordering routines that require a 32 symmetric structure. It is required since those routines call 33 SparsePak routines that expect a symmetric matrix. 34 */ 35 int MatToSymmetricIJ_SeqAIJ(int m,int *ai,int *aj,int shiftin, int shiftout, 36 int **iia, int **jja ) 37 { 38 int *work,*ia,*ja,*j,i, nz, row, col,ierr; 39 40 PetscFunctionBegin; 41 /* allocate space for row pointers */ 42 *iia = ia = (int *) PetscMalloc( (m+1)*sizeof(int) );CHKPTRQ(ia); 43 ierr = PetscMemzero(ia,(m+1)*sizeof(int));CHKERRQ(ierr); 44 work = (int *) PetscMalloc( (m+1)*sizeof(int) );CHKPTRQ(work); 45 46 /* determine the number of columns in each row */ 47 ia[0] = shiftout; 48 for (row = 0; row < m; row++) { 49 nz = ai[row+1] - ai[row]; 50 j = aj + ai[row] + shiftin; 51 while (nz--) { 52 col = *j++ + shiftin; 53 if (col > row) { break;} 54 if (col != row) ia[row+1]++; 55 ia[col+1]++; 56 } 57 } 58 59 /* shiftin ia[i] to point to next row */ 60 for ( i=1; i<m+1; i++ ) { 61 row = ia[i-1]; 62 ia[i] += row; 63 work[i-1] = row - shiftout; 64 } 65 66 /* allocate space for column pointers */ 67 nz = ia[m] + (!shiftin); 68 *jja = ja = (int *) PetscMalloc( nz*sizeof(int) );CHKPTRQ(ja); 69 70 /* loop over lower triangular part putting into ja */ 71 for (row = 0; row < m; row++) { 72 nz = ai[row+1] - ai[row]; 73 j = aj + ai[row] + shiftin; 74 while (nz--) { 75 col = *j++ + shiftin; 76 if (col > row) { break;} 77 if (col != row) {ja[work[col]++] = row + shiftout; } 78 ja[work[row]++] = col + shiftout; 79 } 80 } 81 ierr = PetscFree(work);CHKERRQ(ierr); 82 PetscFunctionReturn(0); 83 } 84 85 86 87