Actual source code: ex110.c
1: static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";
3: #include <petscmat.h>
5: int main(int argc,char **argv)
6: {
7: Mat A,B;
8: PetscInt i,j,column,*ooj;
9: PetscInt *di,*dj,*oi,*oj,nd;
10: const PetscInt *garray;
11: PetscScalar *oa,*da;
12: PetscScalar value;
13: PetscRandom rctx;
14: PetscBool equal,done;
15: Mat AA,AB;
16: PetscMPIInt size,rank;
18: PetscInitialize(&argc,&argv,(char*)0,help);
19: MPI_Comm_size(PETSC_COMM_WORLD,&size);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
23: /* Create a mpiaij matrix for checking */
24: MatCreateAIJ(PETSC_COMM_WORLD,5,5,PETSC_DECIDE,PETSC_DECIDE,0,NULL,0,NULL,&A);
25: MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
26: MatSetUp(A);
27: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
28: PetscRandomSetFromOptions(rctx);
30: for (i=5*rank; i<5*rank+5; i++) {
31: for (j=0; j<5*size; j++) {
32: PetscRandomGetValue(rctx,&value);
33: column = (PetscInt) (5*size*PetscRealPart(value));
34: PetscRandomGetValue(rctx,&value);
35: MatSetValues(A,1,&i,1,&column,&value,INSERT_VALUES);
36: }
37: }
38: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
39: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
41: MatMPIAIJGetSeqAIJ(A,&AA,&AB,&garray);
42: MatGetRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);
43: MatSeqAIJGetArray(AA,&da);
44: MatGetRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);
45: MatSeqAIJGetArray(AB,&oa);
47: PetscMalloc1(oi[5],&ooj);
48: PetscArraycpy(ooj,oj,oi[5]);
49: /* modify the column entries in the non-diagonal portion back to global numbering */
50: for (i=0; i<oi[5]; i++) {
51: ooj[i] = garray[ooj[i]];
52: }
54: MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD,5,5,PETSC_DETERMINE,PETSC_DETERMINE,di,dj,da,oi,ooj,oa,&B);
55: MatSetUp(B);
56: MatEqual(A,B,&equal);
58: MatRestoreRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);
59: MatSeqAIJRestoreArray(AA,&da);
60: MatRestoreRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);
61: MatSeqAIJRestoreArray(AB,&oa);
65: /* Free spaces */
66: PetscFree(ooj);
67: PetscRandomDestroy(&rctx);
68: MatDestroy(&A);
69: MatDestroy(&B);
70: PetscFree(oj);
71: PetscFinalize();
72: return 0;
73: }
75: /*TEST
77: test:
78: nsize: 3
80: TEST*/