Actual source code: ex213.c
2: static char help[] = "Tests MatMPIBAIJSetPreallocationCSR()\n\n";
4: /*T
5: Concepts: partitioning
6: Processors: 4
7: T*/
9: /*
10: Include "petscmat.h" so that we can use matrices. Note that this file
11: automatically includes:
12: petscsys.h - base PETSc routines petscvec.h - vectors
13: petscmat.h - matrices
14: petscis.h - index sets
15: petscviewer.h - viewers
16: */
17: #include <petscmat.h>
19: int main(int argc,char **args)
20: {
21: Mat A;
22: PetscInt *ia,*ja, bs = 2;
23: PetscInt N = 9, n;
24: PetscInt rstart, rend, row, col;
25: PetscInt i;
26: PetscMPIInt rank,size;
27: Vec v;
29: PetscInitialize(&argc,&args,(char*)0,help);
30: MPI_Comm_size(PETSC_COMM_WORLD,&size);
32: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
34: /* Get a partition range based on the vector size */
35: VecCreateMPI(PETSC_COMM_WORLD, PETSC_DECIDE, N, &v);
36: VecGetLocalSize(v, &n);
37: VecGetOwnershipRange(v, &rstart, &rend);
38: VecDestroy(&v);
40: PetscMalloc1(n+1,&ia);
41: PetscMalloc1(3*n,&ja);
43: /* Construct a tri-diagonal CSR indexing */
44: i = 1;
45: ia[0] = 0;
46: for (row = rstart; row < rend; row++)
47: {
48: ia[i] = ia[i-1];
50: /* diagonal */
51: col = row;
52: {
53: ja[ia[i]] = col;
54: ia[i]++;
55: }
57: /* lower diagonal */
58: col = row-1;
59: if (col >= 0)
60: {
61: ja[ia[i]] = col;
62: ia[i]++;
63: }
65: /* upper diagonal */
66: col = row+1;
67: if (col < N)
68: {
69: ja[ia[i]] = col;
70: ia[i]++;
71: }
72: i++;
73: }
75: MatCreate(PETSC_COMM_WORLD, &A);
76: MatSetSizes(A, n, n, PETSC_DETERMINE, PETSC_DETERMINE);
77: MatSetType(A,MATMPIAIJ);
78: MatMPIAIJSetPreallocationCSR(A, ia, ja, NULL);
79: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
80: MatDestroy(&A);
82: MatCreate(PETSC_COMM_WORLD, &A);
83: MatSetSizes(A, bs*n, bs*n, PETSC_DETERMINE, PETSC_DETERMINE);
84: MatSetType(A,MATMPIBAIJ);
85: MatMPIBAIJSetPreallocationCSR(A, bs, ia, ja, NULL);
86: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
87: MatDestroy(&A);
89: PetscFree(ia);
90: PetscFree(ja);
91: PetscFinalize();
92: return 0;
93: }
95: /*TEST
97: test:
98: nsize: 4
100: TEST*/