Actual source code: ex11f90.F
petsc-3.6.1 2015-08-06
1: !-----------------------------------------------------------------------
2: !
3: ! Tests DMDAGetVecGetArray()
4: !-----------------------------------------------------------------------
5: !
7: !#define PETSC_USE_FORTRAN_MODULES 1
8: #include <petsc/finclude/petscsysdef.h>
9: #include <petsc/finclude/petscvecdef.h>
10: #include <petsc/finclude/petscdmdef.h>
11: #if defined(PETSC_USE_FORTRAN_MODULES) || defined(PETSC_USE_FORTRAN_DATATYPES)
12: use petsc
13: #endif
14: implicit none
15: #if !defined(PETSC_USE_FORTRAN_MODULES) && !defined(PETSC_USE_FORTRAN_DATATYPES)
16: #include <petsc/finclude/petscsys.h>
17: #include <petsc/finclude/petscvec.h>
18: #include <petsc/finclude/petscdm.h>
19: #include <petsc/finclude/petscdmda.h>
20: #include <petsc/finclude/petscvec.h90>
21: #include <petsc/finclude/petscdmda.h90>
22: #include <petsc/finclude/petscviewer.h>
23: #endif
25: #if defined(PETSC_USE_FORTRAN_DATATYPES)
26: Type(Vec) g
27: Type(DM) ada
28: #else
29: Vec g
30: DM ada
31: #endif
32: PetscScalar,pointer :: x1(:),x2(:,:)
33: PetscScalar,pointer :: x3(:,:,:),x4(:,:,:,:)
34: PetscErrorCode ierr
35: PetscInt m,n,p,dof,s,i,j,k,xs,xl
36: PetscInt ys,yl
37: PetscInt zs,zl
39: m = 5
40: n = 6
41: p = 4;
42: s = 1
43: dof = 1
44: CALL PetscInitialize(PETSC_NULL_CHARACTER,ierr)
45: call DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,1, &
46: & PETSC_NULL_INTEGER,ada,ierr)
47: call DMGetGlobalVector(ada,g,ierr)
48: call DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER, &
49: & xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr)
50: call DMDAVecGetArrayF90(ada,g,x1,ierr)
51: do i=xs,xs+xl-1
52: ! CHKMEMQ
53: x1(i) = i
54: ! CHKMEMQ
55: enddo
56: call DMDAVecRestoreArrayF90(ada,g,x1,ierr)
57: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
58: call DMRestoreGlobalVector(ada,g,ierr)
59: call DMDestroy(ada,ierr)
61: call DMDACreate2d(PETSC_COMM_WORLD, &
62: & DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, &
63: & DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s, &
64: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)
65: call DMGetGlobalVector(ada,g,ierr)
66: call DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER, &
67: & xl,yl,PETSC_NULL_INTEGER,ierr)
68: call DMDAVecGetArrayF90(ada,g,x2,ierr)
69: do i=xs,xs+xl-1
70: do j=ys,ys+yl-1
71: ! CHKMEMQ
72: x2(i,j) = i + j
73: ! CHKMEMQ
74: enddo
75: enddo
76: call DMDAVecRestoreArrayF90(ada,g,x2,ierr)
77: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
78: call DMRestoreGlobalVector(ada,g,ierr)
79: call DMDestroy(ada,ierr)
81: call DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE, &
82: & DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, &
83: & DMDA_STENCIL_BOX, m,n,p,PETSC_DECIDE,PETSC_DECIDE, &
84: & PETSC_DECIDE,dof,s, &
85: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER, &
86: & PETSC_NULL_INTEGER,ada,ierr)
87: call DMGetGlobalVector(ada,g,ierr)
88: call DMDAGetCorners(ada,xs,ys,zs, &
89: & xl,yl,zl,ierr)
90: call DMDAVecGetArrayF90(ada,g,x3,ierr)
91: do i=xs,xs+xl-1
92: do j=ys,ys+yl-1
93: do k=zs,zs+zl-1
94: ! CHKMEMQ
95: x3(i,j,k) = i + j + k
96: ! CHKMEMQ
97: enddo
98: enddo
99: enddo
100: call DMDAVecRestoreArrayF90(ada,g,x3,ierr)
101: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
102: call DMRestoreGlobalVector(ada,g,ierr)
103: call DMDestroy(ada,ierr)
105: !
106: ! Same tests but now with DOF > 1, so dimensions of array are one higher
107: !
108: dof = 2
109: CALL PetscInitialize(PETSC_NULL_CHARACTER,ierr)
110: call DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,1, &
111: & PETSC_NULL_INTEGER,ada,ierr)
112: call DMGetGlobalVector(ada,g,ierr)
113: call DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER, &
114: & xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr)
115: call DMDAVecGetArrayF90(ada,g,x2,ierr)
116: do i=xs,xs+xl-1
117: ! CHKMEMQ
118: x2(0,i) = i
119: x2(1,i) = -i
120: ! CHKMEMQ
121: enddo
122: call DMDAVecRestoreArrayF90(ada,g,x1,ierr)
123: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
124: call DMRestoreGlobalVector(ada,g,ierr)
125: call DMDestroy(ada,ierr)
127: dof = 2
128: call DMDACreate2d(PETSC_COMM_WORLD, &
129: & DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, &
130: & DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s, &
131: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)
132: call DMGetGlobalVector(ada,g,ierr)
133: call DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER, &
134: & xl,yl,PETSC_NULL_INTEGER,ierr)
135: call DMDAVecGetArrayF90(ada,g,x3,ierr)
136: do i=xs,xs+xl-1
137: do j=ys,ys+yl-1
138: ! CHKMEMQ
139: x3(0,i,j) = i + j
140: x3(1,i,j) = -(i + j)
141: ! CHKMEMQ
142: enddo
143: enddo
144: call DMDAVecRestoreArrayF90(ada,g,x3,ierr)
145: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
146: call DMRestoreGlobalVector(ada,g,ierr)
147: call DMDestroy(ada,ierr)
149: dof = 3
150: call DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE, &
151: & DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, &
152: & DMDA_STENCIL_BOX,m,n,p,PETSC_DECIDE,PETSC_DECIDE, &
153: & PETSC_DECIDE,dof,s, &
154: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER, &
155: & PETSC_NULL_INTEGER,ada,ierr)
156: call DMGetGlobalVector(ada,g,ierr)
157: call DMDAGetCorners(ada,xs,ys,zs, &
158: & xl,yl,zl,ierr)
159: call DMDAVecGetArrayF90(ada,g,x4,ierr)
160: do i=xs,xs+xl-1
161: do j=ys,ys+yl-1
162: do k=zs,zs+zl-1
163: ! CHKMEMQ
164: x4(0,i,j,k) = i + j + k
165: x4(1,i,j,k) = -(i + j + k)
166: x4(2,i,j,k) = i + j + k
167: ! CHKMEMQ
168: enddo
169: enddo
170: enddo
171: call DMDAVecRestoreArrayF90(ada,g,x4,ierr)
172: call VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)
173: call DMRestoreGlobalVector(ada,g,ierr)
174: call DMDestroy(ada,ierr)
176: CALL PetscFinalize(ierr)
177: stop
178: END PROGRAM