Actual source code: ex16.c
2: static char help[] = "Tests DMComposite routines.\n\n";
4: #include petscda.h
5: #include petscpf.h
9: int main(int argc,char **argv)
10: {
12: PetscInt nredundant1 = 5,nredundant2 = 2,i,*ridx1,*ridx2,*lidx1,*lidx2,nlocal;
13: PetscMPIInt rank;
14: PetscScalar *redundant1,*redundant2;
15: DMComposite packer;
16: Vec global,local1,local2;
17: PF pf;
18: DA da1,da2;
19: PetscViewer sviewer;
21: PetscInitialize(&argc,&argv,(char*)0,help);
22: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
24: DMCompositeCreate(PETSC_COMM_WORLD,&packer);
26: PetscMalloc(nredundant1*sizeof(PetscScalar),&redundant1);
27: DMCompositeAddArray(packer,0,nredundant1);
29: DACreate1d(PETSC_COMM_WORLD,DA_NONPERIODIC,8,1,1,PETSC_NULL,&da1);
30: DACreateLocalVector(da1,&local1);
31: DMCompositeAddDM(packer,(DM)da1);
33: PetscMalloc(nredundant2*sizeof(PetscScalar),&redundant2);
34: DMCompositeAddArray(packer,0,nredundant2);
36: DACreate1d(PETSC_COMM_WORLD,DA_NONPERIODIC,6,1,1,PETSC_NULL,&da2);
37: DACreateLocalVector(da2,&local2);
38: DMCompositeAddDM(packer,(DM)da2);
40: DMCompositeCreateGlobalVector(packer,&global);
41: PFCreate(PETSC_COMM_WORLD,1,1,&pf);
42: PFSetType(pf,PFIDENTITY,PETSC_NULL);
43: PFApplyVec(pf,PETSC_NULL,global);
44: PFDestroy(pf);
45: VecView(global,PETSC_VIEWER_STDOUT_WORLD);
47: DMCompositeScatter(packer,global,redundant1,local1,redundant2,local2);
48: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant1 array\n",rank);
49: PetscScalarView(nredundant1,redundant1,PETSC_VIEWER_STDOUT_WORLD);
50: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da1 vector\n",rank);
51: PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
52: VecView(local1,sviewer);
53: PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
54: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant2 array\n",rank);
55: PetscScalarView(nredundant2,redundant2,PETSC_VIEWER_STDOUT_WORLD);
56: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da2 vector\n",rank);
57: PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
58: VecView(local2,sviewer);
59: PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
61: for (i=0; i<nredundant1; i++) redundant1[i] = (rank+2)*i;
62: for (i=0; i<nredundant2; i++) redundant2[i] = (rank+10)*i;
64: DMCompositeGather(packer,global,redundant1,local1,redundant2,local2);
65: VecView(global,PETSC_VIEWER_STDOUT_WORLD);
67: /* get the global numbering for each subvector/array element */
68: DMCompositeGetGlobalIndices(packer,&ridx1,&lidx1,&ridx2,&lidx2);
69:
70: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] Global numbering of redundant1 array\n",rank);
71: PetscIntView(nredundant1,ridx1,PETSC_VIEWER_STDOUT_WORLD);
72: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] Global numbering of local1 vector\n",rank);
73: VecGetSize(local1,&nlocal);
74: PetscIntView(nlocal,lidx1,PETSC_VIEWER_STDOUT_WORLD);
75: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] Global numbering of redundant2 array\n",rank);
76: PetscIntView(nredundant2,ridx2,PETSC_VIEWER_STDOUT_WORLD);
77: PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] Global numbering of local2 vector\n",rank);
78: VecGetSize(local2,&nlocal);
79: PetscIntView(nlocal,lidx2,PETSC_VIEWER_STDOUT_WORLD);
80:
81: PetscFree(ridx1);
82: PetscFree(lidx1);
83: PetscFree(ridx2);
84: PetscFree(lidx2);
86: DADestroy(da1);
87: DADestroy(da2);
88: VecDestroy(local1);
89: VecDestroy(local2);
90: VecDestroy(global);
91: DMCompositeDestroy(packer);
92: PetscFree(redundant1);
93: PetscFree(redundant2);
94: PetscFinalize();
95: return 0;
96: }
97: