Actual source code: ex110.c
petsc-3.3-p6 2013-02-11
1: static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";
3: #include <petscmat.h>
4: #include <../src/mat/impls/aij/mpi/mpiaij.h>
8: int main(int argc,char **argv) {
9: Mat A,B;
10: PetscInt i,j,column;
11: PetscInt *di,*dj,*oi,*oj;
12: PetscScalar *oa,*da,value;
13: PetscRandom rctx;
15: PetscBool equal;
16: Mat_SeqAIJ *daij,*oaij;
17: Mat_MPIAIJ *Ampiaij;
18: PetscMPIInt size,rank;
20: PetscInitialize(&argc,&argv,(char *)0,help);
21: MPI_Comm_size(PETSC_COMM_WORLD,&size);
22: if (size == 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_ARG_OUTOFRANGE,"Must run with 2 or more processes");
23: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
25: /* Create a mpiaij matrix for checking */
26: MatCreateAIJ(PETSC_COMM_WORLD,5,5,PETSC_DECIDE,PETSC_DECIDE,0,PETSC_NULL,0,PETSC_NULL,&A);
27: MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
28: MatSetUp(A);
29: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
30: PetscRandomSetFromOptions(rctx);
31:
32: for (i=5*rank; i<5*rank+5; i++) {
33: for (j=0; j<5*size; j++){
34: PetscRandomGetValue(rctx,&value);
35: column = (PetscInt) (5*size*PetscRealPart(value));
36: PetscRandomGetValue(rctx,&value);
37: MatSetValues(A,1,&i,1,&column,&value,INSERT_VALUES);
38: }
39: }
40: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
41: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
43: Ampiaij = (Mat_MPIAIJ*) A->data;
44: daij = (Mat_SeqAIJ*) Ampiaij->A->data;
45: oaij = (Mat_SeqAIJ*) Ampiaij->B->data;
46:
47: di = daij->i;
48: dj = daij->j;
49: da = daij->a;
51: oi = oaij->i;
52: oa = oaij->a;
53: PetscMalloc(oi[5]*sizeof(PetscInt),&oj);
54: PetscMemcpy(oj,oaij->j,oi[5]*sizeof(PetscInt));
55: /* modify the column entries in the non-diagonal portion back to global numbering */
56: for (i=0; i<oi[5]; i++) {
57: oj[i] = Ampiaij->garray[oj[i]];
58: }
60: MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD,5,5,PETSC_DETERMINE,PETSC_DETERMINE,di,dj,da,oi,oj,oa,&B);
61: MatSetUp(B);
62: MatEqual(A,B,&equal);
64: if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Likely a bug in MatCreateMPIAIJWithSplitArrays()");
66: /* Free spaces */
67: PetscRandomDestroy(&rctx);
68: MatDestroy(&A);
69: MatDestroy(&B);
70: PetscFree(oj);
71: PetscFinalize();
72: return(0);
73: }