Actual source code: ex1f90.F90
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: DMLabel, pointer :: nolabel(:) => NULL()
16: PetscInt, target, dimension(3) :: numComp
17: PetscInt, pointer :: pNumComp(:)
18: PetscInt, target, dimension(12) :: numDof
19: PetscInt, pointer :: pNumDof(:)
20: PetscInt, target, dimension(1) :: bcField
21: PetscInt, pointer :: pBcField(:)
22: PetscInt, parameter :: zero = 0, one = 1, two = 2, eight = 8
23: PetscMPIInt :: size
24: IS, target, dimension(1) :: bcCompIS
25: IS, target, dimension(1) :: bcPointIS
26: IS, pointer :: pBcCompIS(:)
27: IS, pointer :: pBcPointIS(:)
28: PetscErrorCode :: ierr
30: call PetscInitialize(PETSC_NULL_CHARACTER, ierr)
31: if (ierr .ne. 0) then
32: print*,'Unable to initialize PETSc'
33: stop
34: endif
35: call MPI_Comm_size(PETSC_COMM_WORLD, size, ierr);CHKERRA(ierr)
36: ! Create a mesh
37: call DMCreate(PETSC_COMM_WORLD, dm, ierr);CHKERRA(ierr)
38: call DMSetType(dm, DMPLEX, ierr);CHKERRA(ierr)
39: call DMSetFromOptions(dm, ierr);CHKERRA(ierr)
40: call DMViewFromOptions(dm, PETSC_NULL_VEC, '-dm_view', ierr);CHKERRA(ierr)
41: call DMGetDimension(dm, dim, ierr);CHKERRA(ierr)
42: ! Create a scalar field u, a vector field v, and a surface vector field w
43: numFields = 3
44: numComp(1) = 1
45: numComp(2) = dim
46: numComp(3) = dim-1
47: pNumComp => numComp
48: do i = 1, numFields*(dim+1)
49: numDof(i) = 0
50: end do
51: ! Let u be defined on vertices
52: numDof(0*(dim+1)+1) = 1
53: ! Let v be defined on cells
54: numDof(1*(dim+1)+dim+1) = dim
55: ! Let v be defined on faces
56: numDof(2*(dim+1)+dim) = dim-1
57: pNumDof => numDof
58: ! Setup boundary conditions
59: numBC = 1
60: ! Test label retrieval
61: call DMGetLabel(dm, 'marker', label, ierr);CHKERRA(ierr)
62: call DMLabelGetValue(label, zero, val, ierr);CHKERRA(ierr)
63: if (size .eq. 1 .and. val .ne. -1) then
64: SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library')
65: endif
66: call DMLabelGetValue(label, eight, val, ierr);CHKERRA(ierr)
67: if (size .eq. 1 .and. val .ne. 1) then
68: SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library')
69: endif
70: ! Prescribe a Dirichlet condition on u on the boundary
71: ! Label "marker" is made by the mesh creation routine
72: bcField(1) = 0
73: pBcField => bcField
74: call ISCreateStride(PETSC_COMM_WORLD, one, zero, one, bcCompIS(1), ierr);CHKERRA(ierr)
75: pBcCompIS => bcCompIS
76: call DMGetStratumIS(dm, 'marker', one, bcPointIS(1),ierr);CHKERRA(ierr)
77: pBcPointIS => bcPointIS
78: ! Create a PetscSection with this data layout
79: call DMSetNumFields(dm, numFields,ierr);CHKERRA(ierr)
80: call DMPlexCreateSection(dm,nolabel,pNumComp,pNumDof,numBC,pBcField,pBcCompIS,pBcPointIS,PETSC_NULL_IS,section,ierr)
81: CHKERRA(ierr)
82: call ISDestroy(bcCompIS(1), ierr);CHKERRA(ierr)
83: call ISDestroy(bcPointIS(1), ierr);CHKERRA(ierr)
84: ! Name the Field variables
85: call PetscSectionSetFieldName(section, zero, 'u', ierr);CHKERRA(ierr)
86: call PetscSectionSetFieldName(section, one, 'v', ierr);CHKERRA(ierr)
87: call PetscSectionSetFieldName(section, two, 'w', ierr);CHKERRA(ierr)
88: if (size .eq. 1) then
89: call PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD, ierr);CHKERRA(ierr)
90: endif
91: ! Tell the DM to use this data layout
92: call DMSetLocalSection(dm, section, ierr);CHKERRA(ierr)
93: ! Create a Vec with this layout and view it
94: call DMGetGlobalVector(dm, u, ierr);CHKERRA(ierr)
95: call PetscViewerCreate(PETSC_COMM_WORLD, viewer, ierr);CHKERRA(ierr)
96: call PetscViewerSetType(viewer, PETSCVIEWERVTK, ierr);CHKERRA(ierr)
97: call PetscViewerFileSetName(viewer, 'sol.vtu', ierr);CHKERRA(ierr)
98: call VecView(u, viewer, ierr);CHKERRA(ierr)
99: call PetscViewerDestroy(viewer, ierr);CHKERRA(ierr)
100: call DMRestoreGlobalVector(dm, u, ierr);CHKERRA(ierr)
101: ! Cleanup
102: call PetscSectionDestroy(section, ierr);CHKERRA(ierr)
103: call DMDestroy(dm, ierr);CHKERRA(ierr)
105: call PetscFinalize(ierr)
106: end program DMPlexTestField
108: !/*TEST
109: ! build:
110: ! requires: defined(PETSC_USING_F90FREEFORM)
111: !
112: ! test:
113: ! suffix: 0
114: ! requires: triangle
115: ! args: -info :~sys,mat:
116: !
117: ! test:
118: ! suffix: 0_2
119: ! requires: triangle
120: ! nsize: 2
121: ! args: -info :~sys,mat,partitioner:
122: !
123: ! test:
124: ! suffix: 1
125: ! requires: ctetgen
126: ! args: -dm_plex_dim 3 -info :~sys,mat:
127: !
128: !TEST*/