Actual source code: ex1f90.F90
petsc-3.10.5 2019-03-28
1: program DMPlexTestField
2: #include "petsc/finclude/petscdmplex.h"
3: #include "petsc/finclude/petscdmlabel.h"
4: use petscdmplex
5: use petscsys
6: implicit none
8: DM :: dm
9: DMLabel :: label
10: Vec :: u
11: PetscViewer :: viewer
12: PetscSection :: section
13: PetscInt :: dim,numFields,numBC
14: PetscInt :: i,val
15: PetscInt, target, dimension(3) :: numComp
16: PetscInt, pointer :: pNumComp(:)
17: PetscInt, target, dimension(12) :: numDof
18: PetscInt, pointer :: pNumDof(:)
19: PetscInt, target, dimension(1) :: bcField
20: PetscInt, pointer :: pBcField(:)
21: PetscInt, parameter :: zero = 0, one = 1, two = 2, eight = 8
22: IS, target, dimension(1) :: bcCompIS
23: IS, target, dimension(1) :: bcPointIS
24: IS, pointer :: pBcCompIS(:)
25: IS, pointer :: pBcPointIS(:)
26: PetscBool :: interpolate,flg
27: PetscErrorCode :: ierr
29: call PetscInitialize(PETSC_NULL_CHARACTER, ierr)
30: if (ierr .ne. 0) then
31: print*,'Unable to initialize PETSc'
32: stop
33: endif
34: dim = 2
35: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-dim', dim,flg, ierr);CHKERRA(ierr)
36: interpolate = PETSC_TRUE
37: ! Create a mesh
38: call DMPlexCreateBoxMesh(PETSC_COMM_WORLD, dim, PETSC_TRUE, PETSC_NULL_INTEGER, PETSC_NULL_REAL, PETSC_NULL_REAL, PETSC_NULL_INTEGER, interpolate, dm, ierr);CHKERRA(ierr)
39: ! Create a scalar field u, a vector field v, and a surface vector field w
40: numFields = 3
41: numComp(1) = 1
42: numComp(2) = dim
43: numComp(3) = dim-1
44: pNumComp => numComp
45: do i = 1, numFields*(dim+1)
46: numDof(i) = 0
47: end do
48: ! Let u be defined on vertices
49: numDof(0*(dim+1)+1) = 1
50: ! Let v be defined on cells
51: numDof(1*(dim+1)+dim+1) = dim
52: ! Let v be defined on faces
53: numDof(2*(dim+1)+dim) = dim-1
54: pNumDof => numDof
55: ! Setup boundary conditions
56: numBC = 1
57: ! Test label retrieval
58: call DMGetLabel(dm, 'marker', label, ierr);CHKERRA(ierr)
59: call DMLabelGetValue(label, zero, val, ierr);CHKERRA(ierr)
60: if (val .ne. -1) then
61: CHKERRA(1)
62: endif
63: call DMLabelGetValue(label, eight, val, ierr);CHKERRA(ierr)
64: if (val .ne. 1) then
65: CHKERRA(1)
66: endif
67: ! Prescribe a Dirichlet condition on u on the boundary
68: ! Label "marker" is made by the mesh creation routine
69: bcField(1) = 0
70: pBcField => bcField
71: call ISCreateStride(PETSC_COMM_WORLD, one, zero, one, bcCompIS(1), ierr);CHKERRA(ierr)
72: pBcCompIS => bcCompIS
73: call DMGetStratumIS(dm, 'marker', one, bcPointIS(1),ierr);CHKERRA(ierr)
74: pBcPointIS => bcPointIS
75: ! Create a PetscSection with this data layout
76: call DMPlexCreateSection(dm,dim,numFields,pNumComp,pNumDof,numBC,pBcField,pBcCompIS,pBcPointIS,PETSC_NULL_IS,section,ierr)
77: CHKERRA(ierr)
78: call ISDestroy(bcCompIS(1), ierr);CHKERRA(ierr)
79: call ISDestroy(bcPointIS(1), ierr);CHKERRA(ierr)
80: ! Name the Field variables
81: call PetscSectionSetFieldName(section, zero, 'u', ierr);CHKERRA(ierr)
82: call PetscSectionSetFieldName(section, one, 'v', ierr);CHKERRA(ierr)
83: call PetscSectionSetFieldName(section, two, 'w', ierr);CHKERRA(ierr)
84: call PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD, ierr);CHKERRA(ierr)
85: ! Tell the DM to use this data layout
86: call DMSetSection(dm, section, ierr);CHKERRA(ierr)
87: ! Create a Vec with this layout and view it
88: call DMGetGlobalVector(dm, u, ierr);CHKERRA(ierr)
89: call PetscViewerCreate(PETSC_COMM_WORLD, viewer, ierr);CHKERRA(ierr)
90: call PetscViewerSetType(viewer, PETSCVIEWERVTK, ierr);CHKERRA(ierr)
91: call PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_VTK, ierr);CHKERRA(ierr)
92: call PetscViewerFileSetName(viewer, 'sol.vtk', ierr);CHKERRA(ierr)
93: call VecView(u, viewer, ierr);CHKERRA(ierr)
94: call PetscViewerDestroy(viewer, ierr);CHKERRA(ierr)
95: call DMRestoreGlobalVector(dm, u, ierr);CHKERRA(ierr)
96: ! Cleanup
97: call PetscSectionDestroy(section, ierr);CHKERRA(ierr)
98: call DMDestroy(dm, ierr);CHKERRA(ierr)
100: call PetscFinalize(ierr)
101: end program DMPlexTestField
103: !/*TEST
104: ! build:
105: ! requires: define(PETSC_USING_F90FREEFORM)
106: !
107: ! test:
108: ! suffix: 0
109: ! requires: triangle
110: !
111: ! test:
112: ! suffix: 1
113: ! requires: ctetgen
114: ! args: -dim 3
115: !
116: !TEST*/