Line data Source code
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 :
11 : static char help[] = "SVD problem with user-defined stopping test.\n\n"
12 : "The command line options are:\n"
13 : " -m <m>, where <m> = matrix rows.\n"
14 : " -n <n>, where <n> = matrix columns (defaults to m+2).\n\n";
15 :
16 : #include <slepcsvd.h>
17 : #include <petsctime.h>
18 :
19 : /*
20 : This example computes the singular values of a rectangular bidiagonal matrix
21 :
22 : | 1 2 |
23 : | 1 2 |
24 : | 1 2 |
25 : A = | . . |
26 : | . . |
27 : | 1 2 |
28 : | 1 2 |
29 : */
30 :
31 : /*
32 : Function for user-defined stopping test.
33 :
34 : Checks that the computing time has not exceeded the deadline.
35 : */
36 282 : PetscErrorCode MyStoppingTest(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,SVDConvergedReason *reason,void *ctx)
37 : {
38 282 : PetscLogDouble now,deadline = *(PetscLogDouble*)ctx;
39 :
40 282 : PetscFunctionBeginUser;
41 : /* check if usual termination conditions are met */
42 282 : PetscCall(SVDStoppingBasic(svd,its,max_it,nconv,nev,reason,NULL));
43 282 : if (*reason==SVD_CONVERGED_ITERATING) {
44 : /* check if deadline has expired */
45 282 : PetscCall(PetscTime(&now));
46 282 : if (now>deadline) *reason = SVD_CONVERGED_USER;
47 : }
48 282 : PetscFunctionReturn(PETSC_SUCCESS);
49 : }
50 :
51 1 : int main(int argc,char **argv)
52 : {
53 1 : Mat A;
54 1 : SVD svd;
55 1 : SVDConvergedReason reason;
56 1 : PetscInt m=20,n,Istart,Iend,i,col[2],nconv;
57 1 : PetscReal seconds=2.5; /* maximum time allowed for computation */
58 1 : PetscLogDouble deadline; /* time to abort computation */
59 1 : PetscScalar value[] = { 1, 2 };
60 1 : PetscBool terse,flg;
61 1 : PetscViewer viewer;
62 :
63 1 : PetscFunctionBeginUser;
64 1 : PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
65 :
66 1 : PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL));
67 1 : PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,&flg));
68 1 : if (!flg) n=m+2;
69 1 : PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nRectangular bidiagonal matrix, m=%" PetscInt_FMT " n=%" PetscInt_FMT "\n\n",m,n));
70 1 : PetscCall(PetscOptionsGetReal(NULL,NULL,"-seconds",&seconds,NULL));
71 1 : deadline = seconds;
72 1 : PetscCall(PetscTimeAdd(&deadline));
73 :
74 : /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
75 : Generate the matrix
76 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
77 :
78 1 : PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
79 1 : PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n));
80 1 : PetscCall(MatSetFromOptions(A));
81 1 : PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
82 751 : for (i=Istart;i<Iend;i++) {
83 750 : col[0]=i; col[1]=i+1;
84 750 : if (i<n-1) PetscCall(MatSetValues(A,1,&i,2,col,value,INSERT_VALUES));
85 750 : else if (i==n-1) PetscCall(MatSetValue(A,i,col[0],value[0],INSERT_VALUES));
86 : }
87 1 : PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
88 1 : PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
89 :
90 : /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
91 : Compute singular values
92 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
93 :
94 1 : PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd));
95 1 : PetscCall(SVDSetOperators(svd,A,NULL));
96 1 : PetscCall(SVDSetWhichSingularTriplets(svd,SVD_SMALLEST));
97 1 : PetscCall(SVDSetTolerances(svd,PETSC_CURRENT,1000));
98 1 : PetscCall(SVDSetType(svd,SVDTRLANCZOS));
99 1 : PetscCall(SVDSetStoppingTestFunction(svd,MyStoppingTest,&deadline,NULL));
100 1 : PetscCall(SVDSetFromOptions(svd));
101 :
102 : /* call the solver */
103 1 : PetscCall(SVDSolve(svd));
104 :
105 : /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106 : Display solution and clean up
107 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
108 :
109 : /* show detailed info unless -terse option is given by user */
110 1 : PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
111 1 : if (terse) PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
112 : else {
113 1 : PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer));
114 1 : PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO_DETAIL));
115 1 : PetscCall(SVDGetConvergedReason(svd,&reason));
116 1 : if (reason!=SVD_CONVERGED_USER) {
117 0 : PetscCall(SVDConvergedReasonView(svd,viewer));
118 0 : PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer));
119 : } else {
120 1 : PetscCall(SVDGetConverged(svd,&nconv));
121 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"SVD solve finished with %" PetscInt_FMT " converged eigenpairs; reason=%s\n",nconv,SVDConvergedReasons[reason]));
122 : }
123 1 : PetscCall(PetscViewerPopFormat(viewer));
124 : }
125 1 : PetscCall(SVDDestroy(&svd));
126 1 : PetscCall(MatDestroy(&A));
127 1 : PetscCall(SlepcFinalize());
128 : return 0;
129 : }
130 :
131 : /*TEST
132 :
133 : test:
134 : suffix: 1
135 : args: -m 750 -seconds 0.1 -svd_max_it 10000
136 : requires: !single
137 :
138 : TEST*/
|
