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