Actual source code: ex1f90.F90
1: program main
2: #include <petsc/finclude/petscvec.h>
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: PetscReal norm,v,v1,v2
17: PetscInt n,ithree
18: PetscErrorCode ierr
19: PetscMPIInt rank
20: PetscBool flg
21: PetscScalar one,two,three
22: PetscScalar dots(3),dot
23: PetscReal nfloat
25: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
26: ! Beginning of program
27: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
29: PetscCallA(PetscInitialize(ierr))
30: one = 1.0
31: two = 2.0
32: three = 3.0
33: n = 20
34: ithree = 3
36: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr))
37: nfloat = n
38: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))
40: ! Create a vector, specifying only its global dimension.
41: ! When using VecCreate(), VecSetSizes() and VecSetFromOptions(),
42: ! the vector format (currently parallel
43: ! or sequential) is determined at runtime. Also, the parallel
44: ! partitioning of the vector is determined by PETSc at runtime.
46: PetscCallA(VecCreate(PETSC_COMM_WORLD,x,ierr))
47: PetscCallA(VecSetSizes(x,PETSC_DECIDE,n,ierr))
48: PetscCallA(VecSetFromOptions(x,ierr))
50: ! Duplicate some work vectors (of the same format and
51: ! partitioning as the initial vector).
53: PetscCallA(VecDuplicate(x,y,ierr))
54: PetscCallA(VecDuplicate(x,w,ierr))
56: ! Duplicate more work vectors (of the same format and
57: ! partitioning as the initial vector). Here we duplicate
58: ! an array of vectors, which is often more convenient than
59: ! duplicating individual ones.
61: PetscCallA(VecDuplicateVecsF90(x,ithree,z,ierr))
63: ! Set the vectors to entries to a constant value.
65: PetscCallA(VecSet(x,one,ierr))
66: PetscCallA(VecSet(y,two,ierr))
67: PetscCallA(VecSet(z(1),one,ierr))
68: PetscCallA(VecSet(z(2),two,ierr))
69: PetscCallA(VecSet(z(3),three,ierr))
71: ! Demonstrate various basic vector routines.
73: PetscCallA(VecDot(x,x,dot,ierr))
74: PetscCallA(VecMDot(x,ithree,z,dots,ierr))
76: ! Note: If using a complex numbers version of PETSc, then
77: ! PETSC_USE_COMPLEX is defined in the makefiles; otherwise,
78: ! (when using real numbers) it is undefined.
80: if (rank .eq. 0) then
81: #if defined(PETSC_USE_COMPLEX)
82: write(6,100) int(PetscRealPart(dot))
83: write(6,110) int(PetscRealPart(dots(1))),int(PetscRealPart(dots(2))),int(PetscRealPart(dots(3)))
84: #else
85: write(6,100) int(dot)
86: write(6,110) int(dots(1)),int(dots(2)),int(dots(3))
87: #endif
88: write(6,120)
89: endif
90: 100 format ('Vector length ',i6)
91: 110 format ('Vector length ',3(i6))
92: 120 format ('All other values should be near zero')
94: PetscCallA(VecScale(x,two,ierr))
95: PetscCallA(VecNorm(x,NORM_2,norm,ierr))
96: v = abs(norm-2.0*sqrt(nfloat))
97: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
98: if (rank .eq. 0) write(6,130) v
99: 130 format ('VecScale ',1pe9.2)
101: PetscCallA(VecCopy(x,w,ierr))
102: PetscCallA(VecNorm(w,NORM_2,norm,ierr))
103: v = abs(norm-2.0*sqrt(nfloat))
104: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
105: if (rank .eq. 0) write(6,140) v
106: 140 format ('VecCopy ',1pe9.2)
108: PetscCallA(VecAXPY(y,three,x,ierr))
109: PetscCallA(VecNorm(y,NORM_2,norm,ierr))
110: v = abs(norm-8.0*sqrt(nfloat))
111: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
112: if (rank .eq. 0) write(6,150) v
113: 150 format ('VecAXPY ',1pe9.2)
115: PetscCallA(VecAYPX(y,two,x,ierr))
116: PetscCallA(VecNorm(y,NORM_2,norm,ierr))
117: v = abs(norm-18.0*sqrt(nfloat))
118: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
119: if (rank .eq. 0) write(6,160) v
120: 160 format ('VecAYXP ',1pe9.2)
122: PetscCallA(VecSwap(x,y,ierr))
123: PetscCallA(VecNorm(y,NORM_2,norm,ierr))
124: v = abs(norm-2.0*sqrt(nfloat))
125: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
126: if (rank .eq. 0) write(6,170) v
127: 170 format ('VecSwap ',1pe9.2)
129: PetscCallA(VecNorm(x,NORM_2,norm,ierr))
130: v = abs(norm-18.0*sqrt(nfloat))
131: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
132: if (rank .eq. 0) write(6,180) v
133: 180 format ('VecSwap ',1pe9.2)
135: PetscCallA(VecWAXPY(w,two,x,y,ierr))
136: PetscCallA(VecNorm(w,NORM_2,norm,ierr))
137: v = abs(norm-38.0*sqrt(nfloat))
138: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
139: if (rank .eq. 0) write(6,190) v
140: 190 format ('VecWAXPY ',1pe9.2)
142: PetscCallA(VecPointwiseMult(w,y,x,ierr))
143: PetscCallA(VecNorm(w,NORM_2,norm,ierr))
144: v = abs(norm-36.0*sqrt(nfloat))
145: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
146: if (rank .eq. 0) write(6,200) v
147: 200 format ('VecPointwiseMult ',1pe9.2)
149: PetscCallA(VecPointwiseDivide(w,x,y,ierr))
150: PetscCallA(VecNorm(w,NORM_2,norm,ierr))
151: v = abs(norm-9.0*sqrt(nfloat))
152: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
153: if (rank .eq. 0) write(6,210) v
154: 210 format ('VecPointwiseDivide ',1pe9.2)
156: dots(1) = one
157: dots(2) = three
158: dots(3) = two
159: PetscCallA(VecSet(x,one,ierr))
160: PetscCallA(VecMAXPY(x,ithree,dots,z,ierr))
161: PetscCallA(VecNorm(z(1),NORM_2,norm,ierr))
162: v = abs(norm-sqrt(nfloat))
163: if (v .gt. -1.d-10 .and. v .lt. 1.d-10) v = 0.0
164: PetscCallA(VecNorm(z(2),NORM_2,norm,ierr))
165: v1 = abs(norm-2.0*sqrt(nfloat))
166: if (v1 .gt. -1.d-10 .and. v1 .lt. 1.d-10) v1 = 0.0
167: PetscCallA(VecNorm(z(3),NORM_2,norm,ierr))
168: v2 = abs(norm-3.0*sqrt(nfloat))
169: if (v2 .gt. -1.d-10 .and. v2 .lt. 1.d-10) v2 = 0.0
170: if (rank .eq. 0) write(6,220) v,v1,v2
171: 220 format ('VecMAXPY ',3(1pe9.2))
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(VecDestroyVecsF90(ithree,z,ierr))
180: PetscCallA(PetscFinalize(ierr))
182: end
184: !
185: !/*TEST
186: !
187: ! test:
188: ! nsize: 2
189: !
190: !TEST*/