1 static const char help[] = "Test freeing of MPI types in PetscSF\n\n"; 2 3 #include <petscvec.h> 4 #include <petscsf.h> 5 #include <petscviewer.h> 6 7 int main(int argc, char **argv) 8 { 9 PetscInt ierr; 10 PetscSF sf; 11 Vec A,Aout; 12 PetscScalar *bufA; 13 PetscScalar *bufAout; 14 PetscMPIInt rank, size; 15 PetscInt nroots, nleaves; 16 PetscInt i; 17 PetscInt *ilocal; 18 PetscSFNode *iremote; 19 PetscBool test_dupped_type; 20 MPI_Datatype contig; 21 22 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 23 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr); 24 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRMPI(ierr); 25 26 if (size != 1) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for one MPI process"); 27 28 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,"","PetscSF type freeing options","none");CHKERRQ(ierr); 29 test_dupped_type = PETSC_FALSE; 30 ierr = PetscOptionsBool("-test_dupped_type", "Test dupped input type","",test_dupped_type,&test_dupped_type,NULL);CHKERRQ(ierr); 31 ierr = PetscOptionsEnd();CHKERRQ(ierr); 32 33 ierr = PetscSFCreate(PETSC_COMM_WORLD,&sf);CHKERRQ(ierr); 34 ierr = PetscSFSetFromOptions(sf);CHKERRQ(ierr); 35 36 nleaves = 1; 37 nroots = 1; 38 ierr = PetscMalloc1(nleaves,&ilocal);CHKERRQ(ierr); 39 40 for (i = 0; i<nleaves; i++) { 41 ilocal[i] = i; 42 } 43 44 ierr = PetscMalloc1(nleaves,&iremote);CHKERRQ(ierr); 45 iremote[0].rank = 0; 46 iremote[0].index = 0; 47 ierr = PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_OWN_POINTER,iremote,PETSC_OWN_POINTER);CHKERRQ(ierr); 48 ierr = PetscSFSetUp(sf);CHKERRQ(ierr); 49 ierr = PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 50 ierr = VecCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 51 ierr = VecSetSizes(A,4,PETSC_DETERMINE);CHKERRQ(ierr); 52 ierr = VecSetFromOptions(A);CHKERRQ(ierr); 53 ierr = VecSetUp(A);CHKERRQ(ierr); 54 55 ierr = VecDuplicate(A,&Aout);CHKERRQ(ierr); 56 ierr = VecGetArray(A,&bufA);CHKERRQ(ierr); 57 for (i=0; i<4; i++) { 58 bufA[i] = (PetscScalar)i; 59 } 60 ierr = VecRestoreArray(A,&bufA);CHKERRQ(ierr); 61 62 ierr = VecGetArrayRead(A,(const PetscScalar**)&bufA);CHKERRQ(ierr); 63 ierr = VecGetArray(Aout,&bufAout);CHKERRQ(ierr); 64 65 ierr = MPI_Type_contiguous(4, MPIU_SCALAR, &contig);CHKERRMPI(ierr); 66 ierr = MPI_Type_commit(&contig);CHKERRMPI(ierr); 67 68 if (test_dupped_type) { 69 MPI_Datatype tmp; 70 ierr = MPI_Type_dup(contig, &tmp);CHKERRMPI(ierr); 71 ierr = MPI_Type_free(&contig);CHKERRMPI(ierr); 72 contig = tmp; 73 } 74 for (i=0;i<10000;i++) { 75 ierr = PetscSFBcastBegin(sf,contig,bufA,bufAout,MPI_REPLACE);CHKERRQ(ierr); 76 ierr = PetscSFBcastEnd(sf,contig,bufA,bufAout,MPI_REPLACE);CHKERRQ(ierr); 77 } 78 ierr = VecRestoreArrayRead(A,(const PetscScalar**)&bufA);CHKERRQ(ierr); 79 ierr = VecRestoreArray(Aout,&bufAout);CHKERRQ(ierr); 80 81 ierr = VecView(Aout,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 82 ierr = VecDestroy(&A);CHKERRQ(ierr); 83 ierr = VecDestroy(&Aout);CHKERRQ(ierr); 84 ierr = PetscSFDestroy(&sf);CHKERRQ(ierr); 85 ierr = MPI_Type_free(&contig);CHKERRMPI(ierr); 86 ierr = PetscFinalize(); 87 return ierr; 88 } 89 90 /*TEST 91 92 test: 93 suffix: basic 94 args: -sf_type basic 95 96 test: 97 suffix: basic_dupped 98 args: -test_dupped_type -sf_type basic 99 100 test: 101 suffix: window 102 filter: grep -v "type" | grep -v "sort" 103 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create allocate dynamic}} 104 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW) 105 106 test: 107 suffix: window_dupped 108 filter: grep -v "type" | grep -v "sort" 109 args: -test_dupped_type -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create allocate dynamic}} 110 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW) 111 112 test: 113 suffix: window_shared 114 output_file: output/ex3_window.out 115 filter: grep -v "type" | grep -v "sort" 116 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared 117 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) defined(PETSC_HAVE_MPI_ONE_SIDED) !defined(PETSC_HAVE_I_MPI_NUMVERSION) 118 119 test: 120 suffix: window_dupped_shared 121 output_file: output/ex3_window_dupped.out 122 filter: grep -v "type" | grep -v "sort" 123 args: -test_dupped_type -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared 124 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) defined(PETSC_HAVE_MPI_ONE_SIDED) !defined(PETSC_HAVE_I_MPI_NUMVERSION) 125 126 TEST*/ 127