Actual source code: ex11f90.F90
petsc-3.14.6 2021-03-30
1: program main
2: !-----------------------------------------------------------------------
3: !
4: ! Tests DMDAGetVecGetArray()
5: !-----------------------------------------------------------------------
6: !
8: #include <petsc/finclude/petscdm.h>
9: use petsc
10: implicit none
12: Type(tVec) g
13: Type(tDM) ada
15: PetscScalar,pointer :: x1(:),x2(:,:)
16: PetscScalar,pointer :: x3(:,:,:),x4(:,:,:,:)
17: PetscErrorCode ierr
18: PetscInt m,n,p,dof,s,i,j,k,xs,xl
19: PetscInt ys,yl
20: PetscInt zs,zl,sw
22: m = 5
23: n = 6
24: p = 4;
25: s = 1
26: dof = 1
27: sw = 1
28: CALL PetscInitialize(PETSC_NULL_CHARACTER,ierr)
29: if (ierr .ne. 0) then
30: print*,'Unable to initialize PETSc'
31: stop
32: endif
33: call DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,sw,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
34: call DMSetUp(ada,ierr);CHKERRA(ierr)
35: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
36: call DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
37: call DMDAVecGetArrayF90(ada,g,x1,ierr);CHKERRA(ierr)
38: do i=xs,xs+xl-1
39: ! CHKMEMQ
40: x1(i) = i
41: ! CHKMEMQ
42: enddo
43: call DMDAVecRestoreArrayF90(ada,g,x1,ierr);CHKERRA(ierr)
44: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
45: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
46: call DMDestroy(ada,ierr);CHKERRA(ierr)
48: call DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s, &
49: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
50: call DMSetUp(ada,ierr);CHKERRA(ierr)
51: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
52: call DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER,xl,yl,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
53: call DMDAVecGetArrayF90(ada,g,x2,ierr);CHKERRA(ierr)
54: do i=xs,xs+xl-1
55: do j=ys,ys+yl-1
56: ! CHKMEMQ
57: x2(i,j) = i + j
58: ! CHKMEMQ
59: enddo
60: enddo
61: call DMDAVecRestoreArrayF90(ada,g,x2,ierr);CHKERRA(ierr)
62: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
63: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
64: call DMDestroy(ada,ierr);CHKERRA(ierr)
66: call DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX, m,n,p,PETSC_DECIDE,PETSC_DECIDE, &
67: & PETSC_DECIDE,dof,s,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
68: call DMSetUp(ada,ierr);CHKERRA(ierr)
69: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
70: call DMDAGetCorners(ada,xs,ys,zs,xl,yl,zl,ierr);CHKERRA(ierr)
71: call DMDAVecGetArrayF90(ada,g,x3,ierr);CHKERRA(ierr)
72: do i=xs,xs+xl-1
73: do j=ys,ys+yl-1
74: do k=zs,zs+zl-1
75: ! CHKMEMQ
76: x3(i,j,k) = i + j + k
77: ! CHKMEMQ
78: enddo
79: enddo
80: enddo
81: call DMDAVecRestoreArrayF90(ada,g,x3,ierr);CHKERRA(ierr)
82: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
83: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
84: call DMDestroy(ada,ierr);CHKERRA(ierr)
86: !
87: ! Same tests but now with DOF > 1, so dimensions of array are one higher
88: !
89: dof = 2
90: call DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,sw,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
91: call DMSetUp(ada,ierr);CHKERRA(ierr)
92: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
93: call DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
94: call DMDAVecGetArrayF90(ada,g,x2,ierr);CHKERRA(ierr)
95: do i=xs,xs+xl-1
96: ! CHKMEMQ
97: x2(0,i) = i
98: x2(1,i) = -i
99: ! CHKMEMQ
100: enddo
101: call DMDAVecRestoreArrayF90(ada,g,x1,ierr);CHKERRA(ierr)
102: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
103: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
104: call DMDestroy(ada,ierr);CHKERRA(ierr)
106: dof = 2
107: call DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s, &
108: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
109: call DMSetUp(ada,ierr);CHKERRA(ierr)
110: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
111: call DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER,xl,yl,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
112: call DMDAVecGetArrayF90(ada,g,x3,ierr);CHKERRA(ierr)
113: do i=xs,xs+xl-1
114: do j=ys,ys+yl-1
115: ! CHKMEMQ
116: x3(0,i,j) = i + j
117: x3(1,i,j) = -(i + j)
118: ! CHKMEMQ
119: enddo
120: enddo
121: call DMDAVecRestoreArrayF90(ada,g,x3,ierr);CHKERRA(ierr)
122: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
123: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
124: call DMDestroy(ada,ierr);CHKERRA(ierr)
126: dof = 3
127: call 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, &
128: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr);CHKERRA(ierr)
129: call DMSetUp(ada,ierr);CHKERRA(ierr)
130: call DMGetGlobalVector(ada,g,ierr);CHKERRA(ierr)
131: call DMDAGetCorners(ada,xs,ys,zs,xl,yl,zl,ierr);CHKERRA(ierr)
132: call DMDAVecGetArrayF90(ada,g,x4,ierr);CHKERRA(ierr)
133: do i=xs,xs+xl-1
134: do j=ys,ys+yl-1
135: do k=zs,zs+zl-1
136: ! CHKMEMQ
137: x4(0,i,j,k) = i + j + k
138: x4(1,i,j,k) = -(i + j + k)
139: x4(2,i,j,k) = i + j + k
140: ! CHKMEMQ
141: enddo
142: enddo
143: enddo
144: call DMDAVecRestoreArrayF90(ada,g,x4,ierr);CHKERRA(ierr)
145: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr)
146: call DMRestoreGlobalVector(ada,g,ierr);CHKERRA(ierr)
147: call DMDestroy(ada,ierr);CHKERRA(ierr)
149: CALL PetscFinalize(ierr)
150: END PROGRAM
152: !
153: !/*TEST
154: !
155: ! build:
156: ! requires: !complex
157: !
158: ! test:
159: ! filter: Error: grep -v "Vec Object" | grep -v "Warning: ieee_inexact is signaling"
160: !
161: !TEST*/