Actual source code: ex48.c
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
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 GSVD problem with matrices loaded from a file.\n"
12: "The command line options are:\n"
13: " -f1 <filename>, PETSc binary file containing matrix A.\n"
14: " -f2 <filename>, PETSc binary file containing matrix B (optional).\n\n";
16: #include <slepcsvd.h>
18: int main(int argc,char **argv)
19: {
20: Mat A,B; /* matrices */
21: SVD svd; /* singular value problem solver context */
22: PetscInt i,m,n,p,Istart,Iend,col[2];
23: PetscScalar vals[] = { -1, 1 };
24: char filename[PETSC_MAX_PATH_LEN];
25: PetscViewer viewer;
26: PetscBool flg,terse;
28: SlepcInitialize(&argc,&argv,(char*)0,help);
30: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
31: Load matrices that define the generalized singular value problem
32: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
34: PetscPrintf(PETSC_COMM_WORLD,"\nGeneralized singular value problem stored in file.\n\n");
35: PetscOptionsGetString(NULL,NULL,"-f1",filename,sizeof(filename),&flg);
38: #if defined(PETSC_USE_COMPLEX)
39: PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrix from a binary file...\n");
40: #else
41: PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrix from a binary file...\n");
42: #endif
43: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
44: MatCreate(PETSC_COMM_WORLD,&A);
45: MatSetFromOptions(A);
46: MatLoad(A,viewer);
47: PetscViewerDestroy(&viewer);
49: MatGetSize(A,&m,&n);
51: PetscOptionsGetString(NULL,NULL,"-f2",filename,sizeof(filename),&flg);
52: if (flg) {
53: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
54: MatCreate(PETSC_COMM_WORLD,&B);
55: MatSetFromOptions(B);
56: MatLoad(B,viewer);
57: PetscViewerDestroy(&viewer);
58: } else {
59: p = n+1;
60: PetscOptionsGetInt(NULL,NULL,"-p",&p,&flg);
61: PetscPrintf(PETSC_COMM_WORLD," Matrix B was not provided, setting B=bidiag(1,-1) with p=%" PetscInt_FMT "\n\n",p);
62: MatCreate(PETSC_COMM_WORLD,&B);
63: MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,p,n);
64: MatSetFromOptions(B);
65: MatSetUp(B);
66: MatGetOwnershipRange(B,&Istart,&Iend);
67: for (i=Istart;i<Iend;i++) {
68: col[0]=i-1; col[1]=i;
69: if (i==0) MatSetValue(B,i,col[1],vals[1],INSERT_VALUES);
70: else if (i<n) MatSetValues(B,1,&i,2,col,vals,INSERT_VALUES);
71: else if (i==n) MatSetValue(B,i,col[0],vals[0],INSERT_VALUES);
72: }
73: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
74: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
75: }
77: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
78: Create the singular value solver and set various options
79: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
81: /*
82: Create singular value solver context
83: */
84: SVDCreate(PETSC_COMM_WORLD,&svd);
86: /*
87: Set operators of GSVD problem
88: */
89: SVDSetOperators(svd,A,B);
90: SVDSetProblemType(svd,SVD_GENERALIZED);
92: /*
93: Set solver parameters at runtime
94: */
95: SVDSetFromOptions(svd);
97: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
98: Solve the problem and print solution
99: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
101: SVDSolve(svd);
103: /* show detailed info unless -terse option is given by user */
104: PetscOptionsHasName(NULL,NULL,"-terse",&terse);
105: if (terse) SVDErrorView(svd,SVD_ERROR_NORM,NULL);
106: else {
107: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);
108: SVDConvergedReasonView(svd,PETSC_VIEWER_STDOUT_WORLD);
109: SVDErrorView(svd,SVD_ERROR_NORM,PETSC_VIEWER_STDOUT_WORLD);
110: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
111: }
112: SVDDestroy(&svd);
113: MatDestroy(&A);
114: MatDestroy(&B);
115: SlepcFinalize();
116: return 0;
117: }
118: /*TEST
120: testset:
121: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
122: args: -f1 ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62a.petsc -f2 ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62b.petsc -svd_nsv 3 -terse
123: output_file: output/ex48_1.out
124: test:
125: suffix: 1
126: args: -svd_type trlanczos -svd_trlanczos_explicitmatrix {{0 1}}
127: TODO: does not work for largest singular values
128: test:
129: suffix: 1_spqr
130: args: -svd_type trlanczos -svd_trlanczos_explicitmatrix -svd_trlanczos_pc_type qr
131: requires: suitesparse
132: TODO: does not work for largest singular values
133: test:
134: suffix: 1_cross
135: args: -svd_type cross -svd_cross_explicitmatrix
136: test:
137: suffix: 1_cyclic
138: args: -svd_type cyclic -svd_cyclic_explicitmatrix
140: testset:
141: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
142: args: -f1 ${DATAFILESPATH}/matrices/complex/qc324.petsc -svd_nsv 3 -terse
143: output_file: output/ex48_2.out
144: filter: sed -e "s/30749/30748/"
145: timeoutfactor: 2
146: test:
147: suffix: 2
148: args: -svd_type trlanczos -svd_trlanczos_explicitmatrix {{0 1}} -svd_trlanczos_ksp_rtol 1e-10
149: requires: !defined(PETSCTEST_VALGRIND)
150: test:
151: suffix: 2_spqr
152: args: -svd_type trlanczos -svd_trlanczos_explicitmatrix -svd_trlanczos_pc_type qr -svd_trlanczos_ksp_rtol 1e-10
153: requires: suitesparse
154: test:
155: suffix: 2_cross
156: args: -svd_type cross -svd_cross_explicitmatrix
157: test:
158: suffix: 2_cyclic
159: args: -svd_type cyclic -svd_cyclic_explicitmatrix
161: test:
162: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES) !defined(PETSCTEST_VALGRIND)
163: args: -f1 ${DATAFILESPATH}/matrices/complex/qc324.petsc -p 320 -svd_nsv 3 -svd_type trlanczos -svd_trlanczos_ksp_rtol 1e-14 -terse
164: timeoutfactor: 2
165: suffix: 3
167: TEST*/