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