Actual source code: ex11f90.F90
1: ! Tests DMDAGetVecGetArray()
2: #include <petsc/finclude/petscdmda.h>
3: program main
4: use petscdmda
5: implicit none
7: type(tVec) g
8: type(tDM) ada
10: PetscScalar, pointer :: x1(:), x2(:, :)
11: PetscScalar, pointer :: x3(:, :, :), x4(:, :, :, :)
12: PetscErrorCode ierr
13: PetscInt, parameter :: m = 5, n = 6, p = 4, s = 1, sw = 1
14: PetscInt i, j, k, dof
15: PetscInt xs, xl, ys, yl, zs, zl
17: PetscInt nen, nel
18: PetscInt, pointer :: elements(:)
20: PetscInt nfields
21: character(80), pointer :: namefields(:)
22: IS, pointer :: isfields(:)
23: DM, pointer :: dmfields(:)
25: PetscCallA(PetscInitialize(ierr))
27: dof = 1
28: PetscCallA(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, m, dof, sw, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
29: PetscCallA(DMSetUp(ada, ierr))
30: PetscCallA(DMGetGlobalVector(ada, g, ierr))
31: PetscCallA(DMDAGetCorners(ada, xs, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, xl, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, ierr))
32: PetscCallA(DMDAVecGetArray(ada, g, x1, ierr))
33: do i = xs, xs + xl - 1
34: x1(i) = i
35: end do
36: PetscCallA(DMDAVecRestoreArray(ada, g, x1, ierr))
37: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
38: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
39: PetscCallA(DMDestroy(ada, ierr))
41: PetscCallA(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, m, n, PETSC_DECIDE, PETSC_DECIDE, dof, s, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
42: PetscCallA(DMSetUp(ada, ierr))
43: PetscCallA(DMGetGlobalVector(ada, g, ierr))
44: PetscCallA(DMDAGetCorners(ada, xs, ys, PETSC_NULL_INTEGER, xl, yl, PETSC_NULL_INTEGER, ierr))
45: PetscCallA(DMDAVecGetArray(ada, g, x2, ierr))
46: do i = xs, xs + xl - 1
47: do j = ys, ys + yl - 1
48: x2(i, j) = i + j
49: end do
50: end do
51: PetscCallA(DMDAVecRestoreArray(ada, g, x2, ierr))
52: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
53: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
55: PetscCallA(DMDAGetElements(ada, nen, nel, elements, ierr))
56: do i = 1, nen*nel
57: PetscCheckA(elements(i) >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Error getting DMDA elements')
58: end do
59: PetscCallA(DMDARestoreElements(ada, nen, nel, elements, ierr))
60: PetscCallA(DMDestroy(ada, ierr))
62: PetscCallA(DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, m, n, p, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, s, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
63: PetscCallA(DMSetUp(ada, ierr))
64: PetscCallA(DMGetGlobalVector(ada, g, ierr))
65: PetscCallA(DMDAGetCorners(ada, xs, ys, zs, xl, yl, zl, ierr))
66: PetscCallA(DMDAVecGetArray(ada, g, x3, ierr))
67: do i = xs, xs + xl - 1
68: do j = ys, ys + yl - 1
69: do k = zs, zs + zl - 1
70: x3(i, j, k) = i + j + k
71: end do
72: end do
73: end do
74: PetscCallA(DMDAVecRestoreArray(ada, g, x3, ierr))
75: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
76: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
77: PetscCallA(DMDestroy(ada, ierr))
79: !
80: ! Same tests but now with DOF > 1, so dimensions of array are one higher
81: !
82: dof = 2
83: PetscCallA(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, m, dof, sw, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
84: PetscCallA(DMSetUp(ada, ierr))
85: PetscCallA(DMGetGlobalVector(ada, g, ierr))
86: PetscCallA(DMDAGetCorners(ada, xs, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, xl, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, ierr))
87: PetscCallA(DMDAVecGetArray(ada, g, x2, ierr))
88: do i = xs, xs + xl - 1
89: x2(0, i) = i
90: x2(1, i) = -i
91: end do
92: PetscCallA(DMDAVecRestoreArray(ada, g, x1, ierr))
93: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
94: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
96: ! some testing unrelated to the example
97: PetscCallA(DMDASetFieldName(ada, 0_PETSC_INT_KIND, 'Field 0', ierr))
98: PetscCallA(DMDASetFieldName(ada, 1_PETSC_INT_KIND, 'Field 1', ierr))
99: PetscCallA(DMCreateFieldDecomposition(ada, nfields, namefields, PETSC_NULL_IS_POINTER, PETSC_NULL_DM_POINTER, ierr))
100: ! print*,nfields,trim(namefields(1)),trim(namefields(2))
101: PetscCallA(DMDestroyFieldDecomposition(ada, nfields, namefields, PETSC_NULL_IS_POINTER, PETSC_NULL_DM_POINTER, ierr))
102: PetscCallA(DMCreateFieldDecomposition(ada, nfields, namefields, isfields, dmfields, ierr))
103: PetscCallA(DMDestroyFieldDecomposition(ada, nfields, namefields, isfields, dmfields, ierr))
105: PetscCallA(DMDestroy(ada, ierr))
107: PetscCallA(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, m, n, PETSC_DECIDE, PETSC_DECIDE, dof, s, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
108: PetscCallA(DMSetUp(ada, ierr))
109: PetscCallA(DMGetGlobalVector(ada, g, ierr))
110: PetscCallA(DMDAGetCorners(ada, xs, ys, PETSC_NULL_INTEGER, xl, yl, PETSC_NULL_INTEGER, ierr))
111: PetscCallA(DMDAVecGetArray(ada, g, x3, ierr))
112: do i = xs, xs + xl - 1
113: do j = ys, ys + yl - 1
114: x3(0, i, j) = i + j
115: x3(1, i, j) = -(i + j)
116: end do
117: end do
118: PetscCallA(DMDAVecRestoreArray(ada, g, x3, ierr))
119: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
120: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
121: PetscCallA(DMDestroy(ada, ierr))
123: dof = 3
124: PetscCallA(DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, m, n, p, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, s, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, PETSC_NULL_INTEGER_ARRAY, ada, ierr))
125: PetscCallA(DMSetUp(ada, ierr))
126: PetscCallA(DMGetGlobalVector(ada, g, ierr))
127: PetscCallA(DMDAGetCorners(ada, xs, ys, zs, xl, yl, zl, ierr))
128: PetscCallA(DMDAVecGetArray(ada, g, x4, ierr))
129: do i = xs, xs + xl - 1
130: do j = ys, ys + yl - 1
131: do k = zs, zs + zl - 1
132: x4(0, i, j, k) = i + j + k
133: x4(1, i, j, k) = -(i + j + k)
134: x4(2, i, j, k) = i + j + k
135: end do
136: end do
137: end do
138: PetscCallA(DMDAVecRestoreArray(ada, g, x4, ierr))
139: PetscCallA(VecView(g, PETSC_VIEWER_STDOUT_WORLD, ierr))
140: PetscCallA(DMRestoreGlobalVector(ada, g, ierr))
141: PetscCallA(DMDestroy(ada, ierr))
143: PetscCallA(PetscFinalize(ierr))
144: end program
146: !
147: !/*TEST
148: !
149: ! build:
150: ! requires: !complex
151: !
152: ! test:
153: ! filter: Error: grep -v "Vec Object" | grep -v "Warning: ieee_inexact is signaling"
154: !
155: !TEST*/