Actual source code: ex12.c

petsc-3.5.4 2015-05-23
Report Typos and Errors
```  2: static char help[] = "Reads a PETSc matrix and vector from a file appends the vector the matrix\n\n";

4: /*T
5:    Concepts: Mat^ordering a matrix - loading a binary matrix and vector;
9:    Processors: 1
10: T*/

12: /*
13:   Include "petscmat.h" so that we can use matrices.
14:   automatically includes:
15:      petscsys.h       - base PETSc routines   petscvec.h    - vectors
16:      petscmat.h    - matrices
17:      petscis.h     - index sets            petscviewer.h - viewers
18: */
19: #include <petscmat.h>
20: #include <../src/mat/impls/aij/seq/aij.h>

24: PetscErrorCode PadMatrix(Mat A,Vec v,PetscScalar c,Mat *B)
25: {
27:   PetscInt       n    = A->rmap->n,i,*cnt,*indices;
28:   Mat_SeqAIJ     *aij = (Mat_SeqAIJ*)A->data;
29:   PetscScalar    *vv;

32:   VecGetArray(v,&vv);
33:   PetscMalloc1(n,&indices);
34:   for (i=0; i<n; i++) indices[i] = i;

36:   /* determine number of nonzeros per row in the new matrix */
37:   PetscMalloc1((n+1),&cnt);
38:   for (i=0; i<n; i++) {
39:     cnt[i] = aij->i[i+1] - aij->i[i] + (vv[i] != 0.0);
40:   }
41:   cnt[n] = 1;
42:   for (i=0; i<n; i++) {
43:     cnt[n] += (vv[i] != 0.0);
44:   }
45:   MatCreateSeqAIJ(PETSC_COMM_SELF,n+1,n+1,0,cnt,B);
46:   MatSetOption(*B,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE);

48:   /* copy over the matrix entries from the matrix and then the vector */
49:   for (i=0; i<n; i++) {
50:     MatSetValues(*B,1,&i,aij->i[i+1] - aij->i[i],aij->j + aij->i[i],aij->a + aij->i[i],INSERT_VALUES);
51:   }
52:   MatSetValues(*B,1,&n,n,indices,vv,INSERT_VALUES);
53:   MatSetValues(*B,n,indices,1,&n,vv,INSERT_VALUES);
54:   MatSetValues(*B,1,&n,1,&n,&c,INSERT_VALUES);

56:   MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
57:   MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
58:   VecRestoreArray(v,&vv);
59:   PetscFree(cnt);
60:   PetscFree(indices);
61:   return(0);
62: }

67: int main(int argc,char **args)
68: {
69:   Mat            A,B;
70:   PetscViewer    fd;                        /* viewer */
71:   char           file[PETSC_MAX_PATH_LEN];  /* input file name */
73:   PetscBool      flg;
74:   Vec            v;

76:   PetscInitialize(&argc,&args,(char*)0,help);

79:   /*
80:      Determine files from which we read the two linear systems
81:      (matrix and right-hand-side vector).
82:   */
83:   PetscOptionsGetString(NULL,"-f0",file,PETSC_MAX_PATH_LEN,&flg);
84:   if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate binary file with the -f0 option");

88:   MatCreate(PETSC_COMM_WORLD,&A);
89:   MatSetType(A,MATSEQAIJ);