Actual source code: ex16.c
petsc-3.4.2 2013-07-02

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

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

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

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

     29:   PetscOptionsGetBool(NULL,"-gather_add",&gather_add,NULL);

     31:   DMCompositeCreate(PETSC_COMM_WORLD,&packer);

     33:   DMRedundantCreate(PETSC_COMM_WORLD,0,nredundant1,&dmred1);
     34:   DMCreateLocalVector(dmred1,&redundant1);
     35:   DMCompositeAddDM(packer,dmred1);

     37:   DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,8,1,1,NULL,&da1);
     38:   DMCreateLocalVector(da1,&local1);
     39:   DMCompositeAddDM(packer,da1);

     41:   DMRedundantCreate(PETSC_COMM_WORLD,1%size,nredundant2,&dmred2);
     42:   DMCreateLocalVector(dmred2,&redundant2);
     43:   DMCompositeAddDM(packer,dmred2);

     45:   DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,6,1,1,NULL,&da2);
     46:   DMCreateLocalVector(da2,&local2);
     47:   DMCompositeAddDM(packer,da2);

     49:   DMCreateGlobalVector(packer,&global);
     50:   PFCreate(PETSC_COMM_WORLD,1,1,&pf);
     51:   PFSetType(pf,PFIDENTITY,NULL);
     52:   PFApplyVec(pf,NULL,global);
     53:   PFDestroy(&pf);
     54:   VecView(global,PETSC_VIEWER_STDOUT_WORLD);

     56:   DMCompositeScatter(packer,global,redundant1,local1,redundant2,local2);
     57:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
     58:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant1 vector\n",rank);
     59:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     60:   VecView(redundant1,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 da1 vector\n",rank);
     64:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     65:   VecView(local1,sviewer);
     66:   PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     67:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
     68:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant2 vector\n",rank);
     69:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     70:   VecView(redundant2,sviewer);
     71:   PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     72:   PetscViewerASCIISynchronizedAllow(PETSC_VIEWER_STDOUT_WORLD,PETSC_TRUE);
     73:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da2 vector\n",rank);
     74:   PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
     75:   VecView(local2,sviewer);
     76:   PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);

     78:   VecGetArray(redundant1,&redundant1a);
     79:   VecGetArray(redundant2,&redundant2a);
     80:   for (i=0; i<nredundant1; i++) redundant1a[i] = (rank+2)*i;
     81:   for (i=0; i<nredundant2; i++) redundant2a[i] = (rank+10)*i;
     82:   VecRestoreArray(redundant1,&redundant1a);
     83:   VecRestoreArray(redundant2,&redundant2a);

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

     88:   /* get the global numbering for each subvector element */
     89:   DMCompositeGetISLocalToGlobalMappings(packer,&ltog);

     91:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant1 vector\n");
     92:   ISLocalToGlobalMappingView(ltog[0],PETSC_VIEWER_STDOUT_WORLD);
     93:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local1 vector\n");
     94:   ISLocalToGlobalMappingView(ltog[1],PETSC_VIEWER_STDOUT_WORLD);
     95:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant2 vector\n");
     96:   ISLocalToGlobalMappingView(ltog[2],PETSC_VIEWER_STDOUT_WORLD);
     97:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local2 vector\n");
     98:   ISLocalToGlobalMappingView(ltog[3],PETSC_VIEWER_STDOUT_WORLD);

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

    103:   DMDestroy(&da1);
    104:   DMDestroy(&dmred1);
    105:   DMDestroy(&dmred2);
    106:   DMDestroy(&da2);
    107:   VecDestroy(&redundant1);
    108:   VecDestroy(&redundant2);
    109:   VecDestroy(&local1);
    110:   VecDestroy(&local2);
    111:   VecDestroy(&global);
    112:   DMDestroy(&packer);
    113:   PetscFinalize();
    114:   return 0;
    115: }