Actual source code: test4f.F90
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: #include <slepc/finclude/slepcsvd.h>
29: program test4f
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: PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, ierr))
57: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
58: m = 20
59: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-m', m, flg, ierr))
60: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))
61: if (.not. flg) n = m + 2
63: if (rank == 0) then
64: write (*, '(/a,i3,a,i3,a)') 'Bidiagonal matrix, m =', m, ', n=', n, ' (Fortran)'
65: end if
67: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
68: ! Build the Lauchli matrix
69: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71: PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
72: PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m, n, ierr))
73: PetscCallA(MatSetFromOptions(A, ierr))
75: PetscCallA(MatGetOwnershipRange(A, Istart, Iend, ierr))
76: val(1) = 1.0
77: val(2) = 2.0
78: do i = Istart, Iend - 1
79: col(1) = i
80: col(2) = i + 1
81: if (i < n) then
82: PetscCallA(MatSetValue(A, i, col(1), val(1), INSERT_VALUES, ierr))
83: end if
84: if (i < n - 1) then
85: PetscCallA(MatSetValue(A, i, col(2), val(2), INSERT_VALUES, ierr))
86: end if
87: end do
89: PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr))
90: PetscCallA(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY, ierr))
92: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
93: ! Compute singular values
94: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
96: PetscCallA(SVDCreate(PETSC_COMM_WORLD, svd, ierr))
97: PetscCallA(SVDSetOperators(svd, A, PETSC_NULL_MAT, ierr))
99: ! ** test some interface functions
100: PetscCallA(SVDGetOperators(svd, B, PETSC_NULL_MAT, ierr))
101: PetscCallA(MatView(B, PETSC_VIEWER_STDOUT_WORLD, ierr))
102: PetscCallA(SVDSetConvergenceTest(svd, SVD_CONV_ABS, ierr))
103: PetscCallA(SVDSetStoppingTest(svd, SVD_STOP_BASIC, ierr))
105: ! ** query properties and print them
106: PetscCallA(SVDGetProblemType(svd, ptype, ierr))
107: if (rank == 0) then
108: write (*, '(/a,i2)') ' Problem type = ', ptype
109: end if
110: PetscCallA(SVDIsGeneralized(svd, flg, ierr))
111: if (flg .and. rank == 0) then
112: write (*, *) 'generalized'
113: end if
114: PetscCallA(SVDGetImplicitTranspose(svd, tmode, ierr))
115: if (rank == 0) then
116: if (tmode) then
117: write (*, *) ' Transpose mode is implicit'
118: else
119: write (*, *) ' Transpose mode is explicit'
120: end if
121: end if
122: PetscCallA(SVDGetConvergenceTest(svd, conv, ierr))
123: if (rank == 0) then
124: write (*, '(a,i2)') ' Convergence test is', conv
125: end if
126: PetscCallA(SVDGetStoppingTest(svd, stp, ierr))
127: if (rank == 0) then
128: write (*, '(a,i2)') ' Stopping test is', stp
129: end if
130: PetscCallA(SVDGetWhichSingularTriplets(svd, which, ierr))
131: if (rank == 0) then
132: if (which == SVD_LARGEST) then
133: write (*, *) ' Which = largest'
134: else
135: write (*, *) ' Which = smallest'
136: end if
137: end if
139: PetscCallA(PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_DEFAULT, vf, ierr))
140: PetscCallA(SVDMonitorSet(svd, SVDMONITORFIRST, vf, PetscViewerAndFormatDestroy, ierr))
141: PetscCallA(SVDMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_DEFAULT, PETSC_NULL, vf, ierr))
142: PetscCallA(SVDMonitorSet(svd, SVDMONITORCONVERGED, vf, PetscViewerAndFormatDestroy, ierr))
143: PetscCallA(SVDMonitorCancel(svd, ierr))
145: ! ** call the solver
146: PetscCallA(SVDSetFromOptions(svd, ierr))
147: PetscCallA(SVDSolve(svd, ierr))
148: PetscCallA(SVDGetConvergedReason(svd, reason, ierr))
149: if (rank == 0) then
150: write (*, '(a,i2)') ' Converged reason:', reason
151: end if
152: PetscCallA(SVDGetIterationNumber(svd, its, ierr))
153: ! if (rank==0) then
154: ! write(*,'(a,i4)') ' Number of iterations of the method:', its
155: ! end if
157: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
158: ! Display solution and clean up
159: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
161: PetscCallA(SVDErrorView(svd, SVD_ERROR_RELATIVE, PETSC_NULL_VIEWER, ierr))
162: PetscCallA(SVDDestroy(svd, ierr))
163: PetscCallA(MatDestroy(A, ierr))
165: PetscCallA(SlepcFinalize(ierr))
166: end program test4f
168: !/*TEST
169: !
170: ! test:
171: ! suffix: 1
172: ! args: -svd_type {{lanczos trlanczos cross cyclic randomized}}
173: ! filter: sed -e 's/2.99255/2.99254/'
174: !
175: !TEST*/