Actual source code: ex14.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7:
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: static char help[] = "Solves a singular value problem with the matrix loaded from a file.\n"
23: "This example works for both real and complex numbers.\n\n"
24: "The command line options are:\n"
25: " -file <filename>, where <filename> = matrix file in PETSc binary form.\n\n";
27: #include slepcsvd.h
31: int main( int argc, char **argv )
32: {
33: Mat A; /* operator matrix */
34: SVD svd; /* singular value problem solver context */
35: const SVDType type;
36: PetscReal error, tol, sigma;
38: PetscInt nsv, maxit, i, its, nconv;
39: char filename[256];
40: PetscViewer viewer;
41: PetscTruth flg;
44: SlepcInitialize(&argc,&argv,(char*)0,help);
46: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: Load the operator matrix that defines the singular value problem
48: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
50: PetscPrintf(PETSC_COMM_WORLD,"\nSingular value problem stored in file.\n\n");
51: PetscOptionsGetString(PETSC_NULL,"-file",filename,256,&flg);
52: if (!flg) {
53: SETERRQ(1,"Must indicate a file name with the -file option.");
54: }
56: #if defined(PETSC_USE_COMPLEX)
57: PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrix from a binary file...\n");
58: #else
59: PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrix from a binary file...\n");
60: #endif
61: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
62: MatLoad(viewer,MATAIJ,&A);
63: PetscViewerDestroy(viewer);
65: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
66: Create the singular value solver and set various options
67: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
69: /*
70: Create singular value solver context
71: */
72: SVDCreate(PETSC_COMM_WORLD,&svd);
74: /*
75: Set operator
76: */
77: SVDSetOperator(svd,A);
79: /*
80: Set solver parameters at runtime
81: */
82: SVDSetFromOptions(svd);
84: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
85: Solve the singular value system
86: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
88: SVDSolve(svd);
89: SVDGetIterationNumber(svd, &its);
90: PetscPrintf(PETSC_COMM_WORLD," Number of iterations of the method: %d\n",its);
92: /*
93: Optional: Get some information from the solver and display it
94: */
95: SVDGetType(svd,&type);
96: PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type);
97: SVDGetDimensions(svd,&nsv,PETSC_NULL,PETSC_NULL);
98: PetscPrintf(PETSC_COMM_WORLD," Number of requested singular values: %d\n",nsv);
99: SVDGetTolerances(svd,&tol,&maxit);
100: PetscPrintf(PETSC_COMM_WORLD," Stopping condition: tol=%.4g, maxit=%d\n",tol,maxit);
102: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103: Display solution and clean up
104: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
106: /*
107: Get number of converged singular triplets
108: */
109: SVDGetConverged(svd,&nconv);
110: PetscPrintf(PETSC_COMM_WORLD," Number of converged approximate singular triplets: %d\n\n",nconv);
112: if (nconv>0) {
113: /*
114: Display singular values and relative errors
115: */
116: PetscPrintf(PETSC_COMM_WORLD,
117: " sigma residual norm\n"
118: " --------------------- ------------------\n" );
119: for( i=0; i<nconv; i++ ) {
120: /*
121: Get converged singular triplets: i-th singular value is stored in sigma
122: */
123: SVDGetSingularTriplet(svd,i,&sigma,PETSC_NULL,PETSC_NULL);
125: /*
126: Compute the error associated to each singular triplet
127: */
128: SVDComputeRelativeError(svd,i,&error);
129: PetscPrintf(PETSC_COMM_WORLD," % 6f ",sigma);
130: PetscPrintf(PETSC_COMM_WORLD," % 12g\n",error);
131: }
132: PetscPrintf(PETSC_COMM_WORLD,"\n" );
133: }
134:
135: /*
136: Free work space
137: */
138: SVDDestroy(svd);
139: MatDestroy(A);
140: SlepcFinalize();
141: return 0;
142: }