xref: /petsc/src/mat/impls/aij/seq/ij.c (revision 7a4fe282d1b349e95b3be72d69d8dd3d3bcd7bc6)
1 #define PETSCMAT_DLL
2 
3 #include "../src/mat/impls/aij/seq/aij.h"
4 
5 #undef __FUNCT__
6 #define __FUNCT__ "MatToSymmetricIJ_SeqAIJ"
7 /*
8   MatToSymmetricIJ_SeqAIJ - Convert a (generally nonsymmetric) sparse AIJ matrix
9            to IJ format (ignore the "A" part) Allocates the space needed. Uses only
10            the lower triangular part of the matrix.
11 
12     Description:
13     Take the data in the row-oriented sparse storage and build the
14     IJ data for the Matrix.  Return 0 on success,row + 1 on failure
15     at that row. Produces the ij for a symmetric matrix by only using
16     the lower triangular part of the matrix.
17 
18     Input Parameters:
19 .   Matrix - matrix to convert
20 .   shiftin - the shift for the original matrix (0 or 1)
21 .   shiftout - the shift required for the ordering routine (0 or 1)
22 
23     Output Parameters:
24 .   ia     - ia part of IJ representation (row information)
25 .   ja     - ja part (column indices)
26 
27     Notes:
28     Both ia and ja may be freed with PetscFree();
29     This routine is provided for ordering routines that require a
30     symmetric structure.  It is required since those routines call
31     SparsePak routines that expect a symmetric  matrix.
32 */
33 PetscErrorCode MatToSymmetricIJ_SeqAIJ(PetscInt m,PetscInt *ai,PetscInt *aj,PetscInt shiftin,PetscInt shiftout,PetscInt **iia,PetscInt **jja)
34 {
35   PetscErrorCode ierr;
36   PetscInt       *work,*ia,*ja,*j,i,nz,row,col;
37 
38   PetscFunctionBegin;
39   /* allocate space for row pointers */
40   ierr = PetscMalloc((m+1)*sizeof(PetscInt),&ia);CHKERRQ(ierr);
41   *iia = ia;
42   ierr = PetscMemzero(ia,(m+1)*sizeof(PetscInt));CHKERRQ(ierr);
43   ierr = PetscMalloc((m+1)*sizeof(PetscInt),&work);CHKERRQ(ierr);
44 
45   /* determine the number of columns in each row */
46   ia[0] = shiftout;
47   for (row = 0; row < m; row++) {
48     nz = ai[row+1] - ai[row];
49     j  = aj + ai[row] + shiftin;
50     while (nz--) {
51        col = *j++ + shiftin;
52        if (col > row) { break;}
53        if (col != row) ia[row+1]++;
54        ia[col+1]++;
55     }
56   }
57 
58   /* shiftin ia[i] to point to next row */
59   for (i=1; i<m+1; i++) {
60     row       = ia[i-1];
61     ia[i]     += row;
62     work[i-1] = row - shiftout;
63   }
64 
65   /* allocate space for column pointers */
66   nz   = ia[m] + (!shiftin);
67   ierr = PetscMalloc(nz*sizeof(PetscInt),&ja);CHKERRQ(ierr);
68   *jja = 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