Actual source code: ex16.c

  2: static char help[] = "Tests DMComposite routines.\n\n";

  4: #include <petscdmda.h>
  5: #include <petscdmcomposite.h>
  6: #include <petscpf.h>

 10: int main(int argc,char **argv)
 11: {
 13:   PetscInt       nredundant1 = 5,nredundant2 = 2,i;
 14:   ISLocalToGlobalMapping *ltog;
 15:   PetscMPIInt    rank,size;
 16:   PetscScalar    *redundant1,*redundant2;
 17:   DM             packer;
 18:   Vec            global,local1,local2;
 19:   PF             pf;
 20:   DM             da1,da2;
 21:   PetscViewer    sviewer;
 22:   PetscBool      gather_add = PETSC_FALSE;

 24:   PetscInitialize(&argc,&argv,(char*)0,help);
 25:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 26:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 28:   PetscOptionsGetBool(PETSC_NULL,"-gather_add",&gather_add,PETSC_NULL);

 30:   DMCompositeCreate(PETSC_COMM_WORLD,&packer);

 32:   PetscMalloc(nredundant1*sizeof(PetscScalar),&redundant1);
 33:   DMCompositeAddArray(packer,0,nredundant1);

 35:   DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,8,1,1,PETSC_NULL,&da1);
 36:   DMCreateLocalVector(da1,&local1);
 37:   DMCompositeAddDM(packer,(DM)da1);

 39:   PetscMalloc(nredundant2*sizeof(PetscScalar),&redundant2);
 40:   DMCompositeAddArray(packer,1%size,nredundant2);

 42:   DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,6,1,1,PETSC_NULL,&da2);
 43:   DMCreateLocalVector(da2,&local2);
 44:   DMCompositeAddDM(packer,(DM)da2);

 46:   DMCreateGlobalVector(packer,&global);
 47:   PFCreate(PETSC_COMM_WORLD,1,1,&pf);
 48:   PFSetType(pf,PFIDENTITY,PETSC_NULL);
 49:   PFApplyVec(pf,PETSC_NULL,global);
 50:   PFDestroy(&pf);
 51:   VecView(global,PETSC_VIEWER_STDOUT_WORLD);

 53:   DMCompositeScatter(packer,global,redundant1,local1,redundant2,local2);
 54:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
 55:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant1 array\n",rank);
 56:   PetscScalarView(nredundant1,redundant1,PETSC_VIEWER_STDOUT_WORLD);
 57:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
 58:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da1 vector\n",rank);
 59:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
 60:   VecView(local1,sviewer);
 61:   PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
 62:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
 63:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant2 array\n",rank);
 64:   PetscScalarView(nredundant2,redundant2,PETSC_VIEWER_STDOUT_WORLD);
 65:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
 66:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da2 vector\n",rank);
 67:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
 68:   VecView(local2,sviewer);
 69:   PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);

 71:   for (i=0; i<nredundant1; i++) redundant1[i] = (rank+2)*i;
 72:   for (i=0; i<nredundant2; i++) redundant2[i] = (rank+10)*i;

 74:   DMCompositeGather(packer,global,gather_add?ADD_VALUES:INSERT_VALUES,redundant1,local1,redundant2,local2);
 75:   VecView(global,PETSC_VIEWER_STDOUT_WORLD);

 77:   /* get the global numbering for each subvector/array element */
 78:   DMCompositeGetISLocalToGlobalMappings(packer,&ltog);

 80:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant1 array\n");
 81:   ISLocalToGlobalMappingView(ltog[0],PETSC_VIEWER_STDOUT_WORLD);
 82:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local1 vector\n");
 83:   ISLocalToGlobalMappingView(ltog[1],PETSC_VIEWER_STDOUT_WORLD);
 84:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant2 array\n");
 85:   ISLocalToGlobalMappingView(ltog[2],PETSC_VIEWER_STDOUT_WORLD);
 86:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local2 vector\n");
 87:   ISLocalToGlobalMappingView(ltog[3],PETSC_VIEWER_STDOUT_WORLD);

 89:   for (i=0; i<4; i++) {ISLocalToGlobalMappingDestroy(&ltog[i]);}
 90:   PetscFree(ltog);

 92:   DMDestroy(&da1);
 93:   DMDestroy(&da2);
 94:   VecDestroy(&local1);
 95:   VecDestroy(&local2);
 96:   VecDestroy(&global);
 97:   DMDestroy(&packer);
 98:   PetscFree(redundant1);
 99:   PetscFree(redundant2);
100:   PetscFinalize();
101:   return 0;
102: }
103: