Actual source code: test12.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[] = "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: PetscFunctionReturn(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;
63: SlepcInitialize(&argc,&argv,(char*)0,help);
65: PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
66: PetscOptionsGetInt(NULL,NULL,"-n",&n,&flg);
67: if (!flg) n=m+2;
68: PetscPrintf(PETSC_COMM_WORLD,"\nRectangular bidiagonal matrix, m=%" PetscInt_FMT " n=%" PetscInt_FMT "\n\n",m,n);
69: PetscOptionsGetReal(NULL,NULL,"-seconds",&seconds,NULL);
70: deadline = seconds;
71: PetscTimeAdd(&deadline);
73: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: Generate the matrix
75: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
77: MatCreate(PETSC_COMM_WORLD,&A);
78: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);
79: MatSetFromOptions(A);
80: MatSetUp(A);
81: MatGetOwnershipRange(A,&Istart,&Iend);
82: for (i=Istart;i<Iend;i++) {
83: col[0]=i; col[1]=i+1;
84: if (i<n-1) MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
85: else if (i==n-1) MatSetValue(A,i,col[0],value[0],INSERT_VALUES);
86: }
87: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
88: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
90: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
91: Compute singular values
92: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
94: SVDCreate(PETSC_COMM_WORLD,&svd);
95: SVDSetOperators(svd,A,NULL);
96: SVDSetWhichSingularTriplets(svd,SVD_SMALLEST);
97: SVDSetTolerances(svd,PETSC_DEFAULT,1000);
98: SVDSetType(svd,SVDTRLANCZOS);
99: SVDSetStoppingTestFunction(svd,MyStoppingTest,&deadline,NULL);
100: SVDSetFromOptions(svd);
102: /* call the solver */
103: SVDSolve(svd);
105: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106: Display solution and clean up
107: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
109: /* show detailed info unless -terse option is given by user */
110: PetscOptionsHasName(NULL,NULL,"-terse",&terse);
111: if (terse) SVDErrorView(svd,SVD_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD);
112: else {
113: PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer);
114: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO_DETAIL);
115: SVDGetConvergedReason(svd,&reason);
116: if (reason!=SVD_CONVERGED_USER) {
117: SVDConvergedReasonView(svd,viewer);
118: SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer);
119: } else {
120: SVDGetConverged(svd,&nconv);
121: PetscViewerASCIIPrintf(viewer,"SVD solve finished with %" PetscInt_FMT " converged eigenpairs; reason=%s\n",nconv,SVDConvergedReasons[reason]);
122: }
123: PetscViewerPopFormat(viewer);
124: }
125: SVDDestroy(&svd);
126: MatDestroy(&A);
127: SlepcFinalize();
128: return 0;
129: }
131: /*TEST
133: test:
134: suffix: 1
135: args: -m 750 -seconds 0.1 -svd_max_it 10000
137: TEST*/