Actual source code: test4f.F
slepc-3.17.1 2022-04-11
1: !
2: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
4: ! Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5: !
6: ! This file is part of SLEPc.
7: ! SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: !
10: ! Program usage: mpiexec -n <np> ./test4f [-help] [-n <n>] [-m <m>] [all SLEPc options]
11: !
12: ! Description: Singular value decomposition of a bidiagonal matrix.
13: !
14: ! | 1 2 |
15: ! | 1 2 |
16: ! | 1 2 |
17: ! A = | . . |
18: ! | . . |
19: ! | 1 2 |
20: ! | 1 2 |
21: !
22: ! The command line options are:
23: ! -m <m>, where <m> = matrix rows.
24: ! -n <n>, where <n> = matrix columns (defaults to m+2).
25: !
26: ! ----------------------------------------------------------------------
27: !
28: program main
29: #include <slepc/finclude/slepcsvd.h>
30: use slepcsvd
31: implicit none
33: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
34: ! Declarations
35: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
36: !
37: Mat A, B
38: SVD svd
39: SVDConv conv;
40: SVDStop stp;
41: SVDWhich which;
42: SVDConvergedReason reason;
43: PetscInt m, n, i, Istart
44: PetscInt col(2), its, Iend
45: PetscScalar val(2)
46: SVDProblemType ptype
47: PetscMPIInt rank
48: PetscErrorCode ierr
49: PetscBool flg, tmode
50: PetscViewerAndFormat vf
52: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
53: ! Beginning of program
54: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
56: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
57: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
58: m = 20
59: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER, &
60: & '-m',m,flg,ierr)
61: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER, &
62: & '-n',n,flg,ierr)
63: if (.not. flg) n = m+2
65: if (rank .eq. 0) then
66: write(*,100) m, n
67: endif
68: 100 format (/'Bidiagonal matrix, m =',I3,', n=',I3,' (Fortran)')
70: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71: ! Build the Lauchli matrix
72: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: call MatCreate(PETSC_COMM_WORLD,A,ierr)
75: call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n,ierr)
76: call MatSetFromOptions(A,ierr)
77: call MatSetUp(A,ierr)
79: call MatGetOwnershipRange(A,Istart,Iend,ierr)
80: val(1) = 1.0
81: val(2) = 2.0
82: do i=Istart,Iend-1
83: col(1) = i
84: col(2) = i+1
85: if (i .le. n-1) then
86: call MatSetValue(A,i,col(1),val(1),INSERT_VALUES,ierr)
87: end if
88: if (i .lt. n-1) then
89: call MatSetValue(A,i,col(2),val(2),INSERT_VALUES,ierr)
90: end if
91: enddo
93: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
94: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
96: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
97: ! Compute singular values
98: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
100: call SVDCreate(PETSC_COMM_WORLD,svd,ierr)
101: call SVDSetOperators(svd,A,PETSC_NULL_MAT,ierr)
103: ! ** test some interface functions
104: call SVDGetOperators(svd,B,PETSC_NULL_MAT,ierr)
105: call MatView(B,PETSC_VIEWER_STDOUT_WORLD,ierr)
106: call SVDSetConvergenceTest(svd,SVD_CONV_ABS,ierr)
107: call SVDSetStoppingTest(svd,SVD_STOP_BASIC,ierr)
109: ! ** query properties and print them
110: call SVDGetProblemType(svd,ptype,ierr)
111: if (rank .eq. 0) then
112: write(*,105) ptype
113: endif
114: 105 format (/' Problem type = ',I2)
115: call SVDIsGeneralized(svd,flg,ierr)
116: if (flg .and. rank .eq. 0) then
117: write(*,*) 'generalized'
118: endif
119: call SVDGetImplicitTranspose(svd,tmode,ierr)
120: if (rank .eq. 0) then
121: if (tmode) then
122: write(*,110) 'implicit'
123: else
124: write(*,110) 'explicit'
125: endif
126: endif
127: 110 format (' Transpose mode is',A9)
128: call SVDGetConvergenceTest(svd,conv,ierr)
129: if (rank .eq. 0) then
130: write(*,120) conv
131: endif
132: 120 format (' Convergence test is',I2)
133: call SVDGetStoppingTest(svd,stp,ierr)
134: if (rank .eq. 0) then
135: write(*,130) stp
136: endif
137: 130 format (' Stopping test is',I2)
138: call SVDGetWhichSingularTriplets(svd,which,ierr)
139: if (rank .eq. 0) then
140: if (which .eq. SVD_LARGEST) then
141: write(*,140) 'largest'
142: else
143: write(*,140) 'smallest'
144: endif
145: endif
146: 140 format (' Which =',A9)
148: call PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, &
149: & PETSC_VIEWER_DEFAULT,vf,ierr)
150: call SVDMonitorSet(svd,SVDMONITORFIRST,vf, &
151: & PetscViewerAndFormatDestroy,ierr)
152: call SVDMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD, &
153: & PETSC_VIEWER_DEFAULT,PETSC_NULL_VEC,vf,ierr)
154: call SVDMonitorSet(svd,SVDMONITORCONVERGED,vf, &
155: & SVDMonitorConvergedDestroy,ierr)
156: call SVDMonitorCancel(svd,ierr)
158: ! ** call the solver
159: call SVDSetFromOptions(svd,ierr)
160: call SVDSolve(svd,ierr)
161: call SVDGetConvergedReason(svd,reason,ierr)
162: if (rank .eq. 0) then
163: write(*,150) reason
164: endif
165: 150 format (' Converged reason:',I2)
166: call SVDGetIterationNumber(svd,its,ierr)
167: ! if (rank .eq. 0) then
168: ! write(*,160) its
169: ! endif
170: !160 format (' Number of iterations of the method:',I4)
172: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
173: ! Display solution and clean up
174: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
176: call SVDErrorView(svd,SVD_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr)
177: call SVDDestroy(svd,ierr)
178: call MatDestroy(A,ierr)
180: call SlepcFinalize(ierr)
181: end
183: !/*TEST
184: !
185: ! test:
186: ! suffix: 1
187: ! args: -svd_type {{lanczos trlanczos cross cyclic randomized}}
188: ! filter: sed -e 's/2.99255/2.99254/'
189: !
190: !TEST*/