1program ex47f90 2#include "petsc/finclude/petsc.h" 3#include "petsc/finclude/petscvec.h" 4 use petsc 5 use petscvec 6 implicit none 7 8 Type(tDM) :: dm 9 Type(tPetscSection) :: section 10 Character(len=PETSC_MAX_PATH_LEN) :: IOBuffer 11 PetscInt :: dof, p, pStart, pEnd, d 12 Type(tVec) :: v 13 PetscInt :: zero = 0 14 PetscInt :: one = 1 15 PetscInt :: two = 2 16 PetscScalar, Dimension(:), Pointer :: val 17 PetscScalar, pointer :: x(:) 18 PetscErrorCode :: ierr 19 20 PetscCallA(PetscInitialize(ierr)) 21 22 PetscCallA(DMCreate(PETSC_COMM_WORLD, dm, ierr)) 23 PetscCallA(DMSetType(dm, DMPLEX, ierr)) 24 PetscCallA(DMSetFromOptions(dm, ierr)) 25 PetscCallA(DMViewFromOptions(dm, PETSC_NULL_OBJECT, '-d_view', ierr)) 26 27 PetscCallA(PetscSectionCreate(PETSC_COMM_WORLD, section, ierr)) 28 PetscCallA(DMPlexGetChart(dm, pStart, pEnd, ierr)) 29 PetscCallA(PetscSectionSetChart(section, pStart, pEnd, ierr)) 30 PetscCallA(DMPlexGetHeightStratum(dm, zero, pStart, pEnd, ierr)) 31 Do p = pStart, pEnd - 1 32 PetscCallA(PetscSectionSetDof(section, p, one, ierr)) 33 End Do 34 PetscCallA(DMPlexGetDepthStratum(dm, zero, pStart, pEnd, ierr)) 35 Do p = pStart, pEnd - 1 36 PetscCallA(PetscSectionSetDof(section, p, two, ierr)) 37 End Do 38 PetscCallA(PetscSectionSetUp(section, ierr)) 39 PetscCallA(DMSetLocalSection(dm, section, ierr)) 40 PetscCallA(PetscSectionViewFromOptions(section, PETSC_NULL_OBJECT, '-s_view', ierr)) 41 42 PetscCallA(DMCreateGlobalVector(dm, v, ierr)) 43 44 PetscCallA(DMPlexGetChart(dm, pStart, pEnd, ierr)) 45 Do p = pStart, pEnd - 1 46 PetscCallA(PetscSectionGetDof(section, p, dof, ierr)) 47 Allocate (val(dof)) 48 Do d = 1, dof 49 val(d) = 100*p + d - 1 50 End Do 51 PetscCallA(VecSetValuesSection(v, section, p, val, INSERT_VALUES, ierr)) 52 DeAllocate (val) 53 End Do 54 PetscCallA(VecView(v, PETSC_VIEWER_STDOUT_WORLD, ierr)) 55 56 Do p = pStart, pEnd - 1 57 PetscCallA(PetscSectionGetDof(section, p, dof, ierr)) 58 PetscCallA(VecGetValuesSection(v, section, p, x, ierr)) 59 write (IOBuffer, *) 'Point ', p, ' dof ', dof, '\n' 60 PetscCallA(PetscPrintf(PETSC_COMM_SELF, IOBuffer, ierr)) 61 PetscCallA(VecRestoreValuesSection(v, section, p, x, ierr)) 62 End Do 63 64 PetscCallA(PetscSectionDestroy(section, ierr)) 65 PetscCallA(VecDestroy(v, ierr)) 66 PetscCallA(DMDestroy(dm, ierr)) 67 PetscCallA(PetscFinalize(ierr)) 68end program ex47f90 69 70!/*TEST 71! 72! test: 73! suffix: 0 74! args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh 75! 76!TEST*/ 77