Actual source code: ex14.c
slepc-main 2024-12-17
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: static char help[] = "Solves a singular value problem with the matrix loaded from a file.\n"
12: "This example works for both real and complex numbers.\n\n"
13: "The command line options are:\n"
14: " -file <filename>, where <filename> = matrix file in PETSc binary form.\n\n";
16: #include <slepcsvd.h>
18: int main(int argc,char **argv)
19: {
20: Mat A; /* operator matrix */
21: SVD svd; /* singular value problem solver context */
22: SVDType type;
23: SVDStop stop;
24: PetscReal tol,thres;
25: PetscInt nsv,maxit,its;
26: char filename[PETSC_MAX_PATH_LEN];
27: PetscViewer viewer;
28: PetscBool flg,terse;
30: PetscFunctionBeginUser;
31: PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
33: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
34: Load the operator matrix that defines the singular value problem
35: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
37: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nSingular value problem stored in file.\n\n"));
38: PetscCall(PetscOptionsGetString(NULL,NULL,"-file",filename,sizeof(filename),&flg));
39: PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must indicate a file name with the -file option");
41: #if defined(PETSC_USE_COMPLEX)
42: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrix from a binary file...\n"));
43: #else
44: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrix from a binary file...\n"));
45: #endif
46: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer));
47: PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
48: PetscCall(MatSetFromOptions(A));
49: PetscCall(MatLoad(A,viewer));
50: PetscCall(PetscViewerDestroy(&viewer));
52: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
53: Create the singular value solver and set various options
54: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
56: /*
57: Create singular value solver context
58: */
59: PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd));
61: /*
62: Set operator
63: */
64: PetscCall(SVDSetOperators(svd,A,NULL));
66: /*
67: Set solver parameters at runtime
68: */
69: PetscCall(SVDSetFromOptions(svd));
71: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: Solve the singular value system
73: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
75: PetscCall(SVDSolve(svd));
76: PetscCall(SVDGetIterationNumber(svd,&its));
77: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of iterations of the method: %" PetscInt_FMT "\n",its));
79: /*
80: Optional: Get some information from the solver and display it
81: */
82: PetscCall(SVDGetType(svd,&type));
83: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type));
84: PetscCall(SVDGetStoppingTest(svd,&stop));
85: if (stop!=SVD_STOP_THRESHOLD) {
86: PetscCall(SVDGetDimensions(svd,&nsv,NULL,NULL));
87: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested singular values: %" PetscInt_FMT "\n",nsv));
88: } else {
89: PetscCall(SVDGetThreshold(svd,&thres,NULL));
90: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Using threshold: %.4g\n",(double)thres));
91: }
92: PetscCall(SVDGetTolerances(svd,&tol,&maxit));
93: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Stopping condition: tol=%.4g, maxit=%" PetscInt_FMT "\n",(double)tol,maxit));
95: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
96: Display solution and clean up
97: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
99: /* show detailed info unless -terse option is given by user */
100: PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
101: if (terse) PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,NULL));
102: else {
103: PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
104: PetscCall(SVDConvergedReasonView(svd,PETSC_VIEWER_STDOUT_WORLD));
105: PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
106: PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
107: }
108: PetscCall(SVDDestroy(&svd));
109: PetscCall(MatDestroy(&A));
110: PetscCall(SlepcFinalize());
111: return 0;
112: }
113: /*TEST
115: testset:
116: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
117: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/rdb200.petsc -terse
118: test:
119: suffix: 1
120: args: -svd_nsv 4 -svd_standard -svd_ncv 12 -svd_type {{trlanczos lanczos randomized cross}}
121: filter: grep -v method
122: test:
123: suffix: 1_scalapack
124: nsize: {{1 2 3}}
125: args: -svd_nsv 4 -svd_type scalapack
126: requires: scalapack
127: test:
128: suffix: 1_elemental
129: nsize: {{1 2 3}}
130: args: -svd_nsv 4 -svd_type elemental
131: output_file: output/ex14_1_scalapack.out
132: filter: sed -e "s/elemental/scalapack/"
133: requires: elemental
134: test:
135: suffix: 1_ksvd
136: nsize: 6
137: args: -svd_nsv 4 -svd_type ksvd -svd_ksvd_eigen_method {{mrrr elpa}} -svd_ksvd_polar_method {{qdwh zolopd}}
138: output_file: output/ex14_1_scalapack.out
139: filter: sed -e "s/ksvd/scalapack/"
140: requires: elpa polar ksvd
141: test:
142: suffix: 2
143: args: -svd_nsv 2 -svd_type cyclic -svd_cyclic_explicitmatrix -svd_cyclic_st_type sinvert -svd_cyclic_eps_target 12.5 -svd_cyclic_st_ksp_type preonly -svd_cyclic_st_pc_type lu -svd_view
144: filter: grep -v tolerance
145: test:
146: suffix: 2_cross
147: args: -svd_nsv 2 -svd_type cross -svd_cross_explicitmatrix -svd_cross_st_type sinvert -svd_cross_eps_target 100.0
148: filter: grep -v tolerance
150: testset:
151: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
152: args: -file ${DATAFILESPATH}/matrices/complex/qc324.petsc -terse
153: test:
154: suffix: 1_complex
155: args: -svd_nsv 4
156: test:
157: suffix: 1_complex_scalapack
158: nsize: {{1 2 3}}
159: args: -svd_nsv 4 -svd_type scalapack
160: requires: scalapack
161: test:
162: suffix: 1_complex_elemental
163: nsize: {{1 2 3}}
164: args: -svd_nsv 4 -svd_type elemental
165: requires: elemental
166: test:
167: suffix: 2_complex
168: args: -svd_nsv 2 -svd_type cyclic -svd_cyclic_explicitmatrix -svd_cyclic_st_type sinvert -svd_cyclic_eps_target 12.0 -svd_cyclic_st_ksp_type preonly -svd_cyclic_st_pc_type lu -svd_view
169: filter: grep -v tolerance
171: testset:
172: args: -svd_nsv 5 -svd_type randomized -svd_max_it 1 -svd_conv_maxit
173: test:
174: suffix: 3
175: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/rdb200.petsc
176: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
177: test:
178: suffix: 3_complex
179: args: -file ${DATAFILESPATH}/matrices/complex/qc324.petsc
180: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
181: filter: sed -e 's/[0-9][0-9]$//'
183: testset:
184: args: -svd_type {{trlanczos lanczos cross}} -terse
185: filter: grep -v method
186: test:
187: suffix: 4
188: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/rdb200.petsc -svd_threshold_relative 0.9 -svd_ncv {{26 12}}
189: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
190: test:
191: suffix: 4_complex
192: args: -file ${DATAFILESPATH}/matrices/complex/qc324.petsc -svd_threshold_relative 0.6 -svd_ncv {{18 10}}
193: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
195: TEST*/