1 static char help[] = "Tests for VecGetValuesSection / VecSetValuesSection \n\n"; 2 3 #include <petscdmplex.h> 4 5 int main(int argc, char **argv) { 6 DM dm; 7 Vec v; 8 PetscSection section; 9 PetscScalar val[2]; 10 PetscInt pStart, pEnd, p; 11 12 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 13 14 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm)); 15 PetscCall(DMSetType(dm, DMPLEX)); 16 PetscCall(DMSetFromOptions(dm)); 17 PetscCall(DMViewFromOptions(dm, NULL, "-d_view")); 18 19 PetscCall(PetscSectionCreate(PetscObjectComm((PetscObject)dm), §ion)); 20 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd)); 21 PetscCall(PetscSectionSetChart(section, pStart, pEnd)); 22 PetscCall(DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd)); 23 for (p = pStart; p < pEnd; p++) PetscCall(PetscSectionSetDof(section, p, 1)); 24 PetscCall(DMPlexGetDepthStratum(dm, 0, &pStart, &pEnd)); 25 for (p = pStart; p < pEnd; p++) PetscCall(PetscSectionSetDof(section, p, 2)); 26 PetscCall(PetscSectionSetUp(section)); 27 PetscCall(DMSetLocalSection(dm, section)); 28 PetscCall(PetscSectionViewFromOptions(section, NULL, "-s_view")); 29 30 PetscCall(DMCreateGlobalVector(dm, &v)); 31 PetscCall(VecViewFromOptions(v, NULL, "-v_view")); 32 33 /* look through all cells and change "cell values" */ 34 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd)); 35 for (p = pStart; p < pEnd; ++p) { 36 PetscInt dof; 37 38 PetscCall(PetscSectionGetDof(section, p, &dof)); 39 for (PetscInt d = 0; d < dof; ++d) val[d] = 100 * p + d; 40 PetscCall(VecSetValuesSection(v, section, p, val, INSERT_VALUES)); 41 } 42 PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD)); 43 44 for (p = pStart; p < pEnd; ++p) { 45 PetscScalar *x; 46 PetscInt dof; 47 48 PetscCall(PetscSectionGetDof(section, p, &dof)); 49 PetscCall(VecGetValuesSection(v, section, p, &x)); 50 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Point #%" PetscInt_FMT " %" PetscInt_FMT " dof\n", p, dof)); 51 } 52 53 PetscCall(VecDestroy(&v)); 54 PetscCall(PetscSectionDestroy(§ion)); 55 PetscCall(DMDestroy(&dm)); 56 PetscCall(PetscFinalize()); 57 return 0; 58 } 59 60 /*TEST 61 62 test: 63 suffix: 0 64 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh 65 66 TEST*/ 67