Actual source code: ex1f90.F90
1: #include <petsc/finclude/petscvec.h>
2: program main
3: use petscvec
4: implicit none
6: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
7: ! Variable declarations
8: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: !
10: ! Variables:
11: ! x, y, w - vectors
12: ! z - array of vectors
13: !
14: Vec x, y, w
15: Vec, pointer :: z(:)
16: PetscInt, pointer :: ranges(:)
17: PetscReal norm, v, v1, v2
18: PetscInt n
19: PetscErrorCode ierr
20: PetscMPIInt rank
21: PetscBool flg
22: PetscScalar, parameter :: one = 1.0, two = 2.0, three = 3.0
23: PetscScalar dots(3), dot
24: PetscReal nfloat
26: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
27: ! Beginning of program
28: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
30: PetscCallA(PetscInitialize(ierr))
32: n = 20
33: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))
34: nfloat = real(n, PETSC_REAL_KIND)
35: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
37: ! Create a vector, specifying only its global dimension.
38: ! When using VecCreate(), VecSetSizes() and VecSetFromOptions(),
39: ! the vector format (currently parallel
40: ! or sequential) is determined at runtime. Also, the parallel
41: ! partitioning of the vector is determined by PETSc at runtime.
43: PetscCallA(VecCreate(PETSC_COMM_WORLD, x, ierr))
44: PetscCallA(VecSetSizes(x, PETSC_DECIDE, n, ierr))
45: PetscCallA(VecSetFromOptions(x, ierr))
47: ! Duplicate some work vectors (of the same format and
48: ! partitioning as the initial vector).
50: PetscCallA(VecDuplicate(x, y, ierr))
51: PetscCallA(VecDuplicate(x, w, ierr))
53: ! Duplicate more work vectors (of the same format and
54: ! partitioning as the initial vector). Here we duplicate
55: ! an array of vectors, which is often more convenient than
56: ! duplicating individual ones.
58: PetscCallA(VecDuplicateVecs(x, 3_PETSC_INT_KIND, z, ierr))
60: ! Set the vectors to entries to a constant value.
62: PetscCallA(VecSet(x, one, ierr))
63: PetscCallA(VecSet(y, two, ierr))
64: PetscCallA(VecSet(z(1), one, ierr))
65: PetscCallA(VecSet(z(2), two, ierr))
66: PetscCallA(VecSet(z(3), three, ierr))
68: ! Demonstrate various basic vector routines.
70: PetscCallA(VecDot(x, x, dot, ierr))
71: PetscCallA(VecMDot(x, 3_PETSC_INT_KIND, z, dots, ierr))
73: ! Note: If using a complex numbers version of PETSc, then
74: ! PETSC_USE_COMPLEX is defined in the makefiles; otherwise,
75: ! (when using real numbers) it is undefined.
77: if (rank == 0) then
78: #if defined(PETSC_USE_COMPLEX)
79: write (6, 100) int(PetscRealPart(dot))
80: write (6, 110) int(PetscRealPart(dots(1))), int(PetscRealPart(dots(2))), int(PetscRealPart(dots(3)))
81: #else
82: write (6, 100) int(dot)
83: write (6, 110) int(dots(1)), int(dots(2)), int(dots(3))
84: #endif
85: write (6, 120)
86: end if
87: 100 format('Vector length ', i6)
88: 110 format('Vector length ', 3(i6))
89: 120 format('All other values should be near zero')
91: PetscCallA(VecScale(x, two, ierr))
92: PetscCallA(VecNorm(x, NORM_2, norm, ierr))
93: v = abs(norm - 2.0*sqrt(nfloat))
94: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
95: if (rank == 0) write (6, 130) v
96: 130 format('VecScale ', 1pe9.2)
98: PetscCallA(VecCopy(x, w, ierr))
99: PetscCallA(VecNorm(w, NORM_2, norm, ierr))
100: v = abs(norm - 2.0*sqrt(nfloat))
101: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
102: if (rank == 0) write (6, 140) v
103: 140 format('VecCopy ', 1pe9.2)
105: PetscCallA(VecAXPY(y, three, x, ierr))
106: PetscCallA(VecNorm(y, NORM_2, norm, ierr))
107: v = abs(norm - 8.0*sqrt(nfloat))
108: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
109: if (rank == 0) write (6, 150) v
110: 150 format('VecAXPY ', 1pe9.2)
112: PetscCallA(VecAYPX(y, two, x, ierr))
113: PetscCallA(VecNorm(y, NORM_2, norm, ierr))
114: v = abs(norm - 18.0*sqrt(nfloat))
115: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
116: if (rank == 0) write (6, 160) v
117: 160 format('VecAYXP ', 1pe9.2)
119: PetscCallA(VecSwap(x, y, ierr))
120: PetscCallA(VecNorm(y, NORM_2, norm, ierr))
121: v = abs(norm - 2.0*sqrt(nfloat))
122: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
123: if (rank == 0) write (6, 170) v
124: 170 format('VecSwap ', 1pe9.2)
126: PetscCallA(VecNorm(x, NORM_2, norm, ierr))
127: v = abs(norm - 18.0*sqrt(nfloat))
128: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
129: if (rank == 0) write (6, 180) v
130: 180 format('VecSwap ', 1pe9.2)
132: PetscCallA(VecWAXPY(w, two, x, y, ierr))
133: PetscCallA(VecNorm(w, NORM_2, norm, ierr))
134: v = abs(norm - 38.0*sqrt(nfloat))
135: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
136: if (rank == 0) write (6, 190) v
137: 190 format('VecWAXPY ', 1pe9.2)
139: PetscCallA(VecPointwiseMult(w, y, x, ierr))
140: PetscCallA(VecNorm(w, NORM_2, norm, ierr))
141: v = abs(norm - 36.0*sqrt(nfloat))
142: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
143: if (rank == 0) write (6, 200) v
144: 200 format('VecPointwiseMult ', 1pe9.2)
146: PetscCallA(VecPointwiseDivide(w, x, y, ierr))
147: PetscCallA(VecNorm(w, NORM_2, norm, ierr))
148: v = abs(norm - 9.0*sqrt(nfloat))
149: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
150: if (rank == 0) write (6, 210) v
151: 210 format('VecPointwiseDivide ', 1pe9.2)
153: dots(1) = one
154: dots(2) = three
155: dots(3) = two
156: PetscCallA(VecSet(x, one, ierr))
157: PetscCallA(VecMAXPY(x, 3_PETSC_INT_KIND, dots, z, ierr))
158: PetscCallA(VecNorm(z(1), NORM_2, norm, ierr))
159: v = abs(norm - sqrt(nfloat))
160: if (v > -1.d-10 .and. v < 1.d-10) v = 0.0
161: PetscCallA(VecNorm(z(2), NORM_2, norm, ierr))
162: v1 = abs(norm - 2.0*sqrt(nfloat))
163: if (v1 > -1.d-10 .and. v1 < 1.d-10) v1 = 0.0
164: PetscCallA(VecNorm(z(3), NORM_2, norm, ierr))
165: v2 = abs(norm - 3.0*sqrt(nfloat))
166: if (v2 > -1.d-10 .and. v2 < 1.d-10) v2 = 0.0
167: if (rank == 0) write (6, 220) v, v1, v2
168: 220 format('VecMAXPY ', 3(1pe9.2))
170: PetscCallA(VecGetOwnershipRanges(x, ranges, ierr))
171: ! PetscCallA(VecRestoreOwnershipRanges(x, ranges, ierr))
173: ! Free work space. All PETSc objects should be destroyed when they
174: ! are no longer needed.
176: PetscCallA(VecDestroy(x, ierr))
177: PetscCallA(VecDestroy(y, ierr))
178: PetscCallA(VecDestroy(w, ierr))
179: PetscCallA(VecDestroyVecs(3_PETSC_INT_KIND, z, ierr))
180: PetscCallA(PetscFinalize(ierr))
181: end
183: !
184: !/*TEST
185: !
186: ! test:
187: ! nsize: 2
188: !
189: !TEST*/