| Line | Branch | Exec | Source |
|---|---|---|---|
| 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[] = "Singular value decomposition of the Lauchli matrix.\n\n" | ||
| 12 | "The command line options are:\n" | ||
| 13 | " -n <n>, where <n> = matrix dimension.\n" | ||
| 14 | " -mu <mu>, where <mu> = subdiagonal value.\n\n"; | ||
| 15 | |||
| 16 | #include <slepcsvd.h> | ||
| 17 | |||
| 18 | 22 | int main(int argc,char **argv) | |
| 19 | { | ||
| 20 | 22 | Mat A; /* operator matrix */ | |
| 21 | 22 | Vec u,v; /* left and right singular vectors */ | |
| 22 | 22 | SVD svd; /* singular value problem solver context */ | |
| 23 | 22 | SVDType type; | |
| 24 | 22 | PetscReal error,tol,sigma,mu=PETSC_SQRT_MACHINE_EPSILON; | |
| 25 | 22 | PetscInt n=100,i,j,Istart,Iend,nsv,maxit,its,nconv; | |
| 26 | |||
| 27 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
22 | PetscFunctionBeginUser; |
| 28 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 29 | |||
| 30 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
| 31 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscOptionsGetReal(NULL,NULL,"-mu",&mu,NULL)); |
| 32 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nLauchli singular value decomposition, (%" PetscInt_FMT " x %" PetscInt_FMT ") mu=%g\n\n",n+1,n,(double)mu)); |
| 33 | |||
| 34 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 35 | Build the Lauchli matrix | ||
| 36 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 37 | |||
| 38 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); |
| 39 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n+1,n)); |
| 40 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatSetFromOptions(A)); |
| 41 | |||
| 42 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 43 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1840 | for (i=Istart;i<Iend;i++) { |
| 44 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1818 | if (i == 0) { |
| 45 |
7/8✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 2 times.
|
1818 | for (j=0;j<n;j++) PetscCall(MatSetValue(A,0,j,1.0,INSERT_VALUES)); |
| 46 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1818 | } else PetscCall(MatSetValue(A,i,i-1,mu,INSERT_VALUES)); |
| 47 | } | ||
| 48 | |||
| 49 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 50 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 51 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatCreateVecs(A,&v,&u)); |
| 52 | |||
| 53 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 54 | Create the singular value solver and set various options | ||
| 55 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 56 | |||
| 57 | /* | ||
| 58 | Create singular value solver context | ||
| 59 | */ | ||
| 60 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd)); |
| 61 | |||
| 62 | /* | ||
| 63 | Set operators and problem type | ||
| 64 | */ | ||
| 65 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDSetOperators(svd,A,NULL)); |
| 66 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDSetProblemType(svd,SVD_STANDARD)); |
| 67 | |||
| 68 | /* | ||
| 69 | Use thick-restart Lanczos as default solver | ||
| 70 | */ | ||
| 71 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDSetType(svd,SVDTRLANCZOS)); |
| 72 | |||
| 73 | /* | ||
| 74 | Set solver parameters at runtime | ||
| 75 | */ | ||
| 76 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDSetFromOptions(svd)); |
| 77 | |||
| 78 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 79 | Solve the singular value system | ||
| 80 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 81 | |||
| 82 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDSolve(svd)); |
| 83 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDGetIterationNumber(svd,&its)); |
| 84 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of iterations of the method: %" PetscInt_FMT "\n",its)); |
| 85 | |||
| 86 | /* | ||
| 87 | Optional: Get some information from the solver and display it | ||
| 88 | */ | ||
| 89 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDGetType(svd,&type)); |
| 90 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type)); |
| 91 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDGetDimensions(svd,&nsv,NULL,NULL)); |
| 92 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested singular values: %" PetscInt_FMT "\n",nsv)); |
| 93 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDGetTolerances(svd,&tol,&maxit)); |
| 94 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Stopping condition: tol=%.4g, maxit=%" PetscInt_FMT "\n",(double)tol,maxit)); |
| 95 | |||
| 96 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 97 | Display solution and clean up | ||
| 98 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 99 | |||
| 100 | /* | ||
| 101 | Get number of converged singular triplets | ||
| 102 | */ | ||
| 103 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDGetConverged(svd,&nconv)); |
| 104 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of converged approximate singular triplets: %" PetscInt_FMT "\n\n",nconv)); |
| 105 | |||
| 106 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
22 | if (nconv>0) { |
| 107 | /* | ||
| 108 | Display singular values and relative errors | ||
| 109 | */ | ||
| 110 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD, |
| 111 | " sigma relative error\n" | ||
| 112 | " --------------------- ------------------\n")); | ||
| 113 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1322 | for (i=0;i<nconv;i++) { |
| 114 | /* | ||
| 115 | Get converged singular triplets: i-th singular value is stored in sigma | ||
| 116 | */ | ||
| 117 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1300 | PetscCall(SVDGetSingularTriplet(svd,i,&sigma,u,v)); |
| 118 | |||
| 119 | /* | ||
| 120 | Compute the error associated to each singular triplet | ||
| 121 | */ | ||
| 122 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1300 | PetscCall(SVDComputeError(svd,i,SVD_ERROR_RELATIVE,&error)); |
| 123 | |||
| 124 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1300 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," % 6f ",(double)sigma)); |
| 125 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1300 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," % 12g\n",(double)error)); |
| 126 | } | ||
| 127 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n")); |
| 128 | } | ||
| 129 | |||
| 130 | /* | ||
| 131 | Free work space | ||
| 132 | */ | ||
| 133 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(SVDDestroy(&svd)); |
| 134 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(MatDestroy(&A)); |
| 135 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(VecDestroy(&u)); |
| 136 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
22 | PetscCall(VecDestroy(&v)); |
| 137 |
3/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
22 | PetscCall(SlepcFinalize()); |
| 138 | return 0; | ||
| 139 | } | ||
| 140 | |||
| 141 | /*TEST | ||
| 142 | |||
| 143 | testset: | ||
| 144 | filter: sed -e "s/[0-9]\.[0-9]*e[+-]\([0-9]*\)/removed/g" | ||
| 145 | requires: double | ||
| 146 | test: | ||
| 147 | suffix: 1 | ||
| 148 | test: | ||
| 149 | suffix: 1_scalapack | ||
| 150 | nsize: {{1 2}} | ||
| 151 | args: -svd_type scalapack | ||
| 152 | requires: scalapack !__float128 | ||
| 153 | |||
| 154 | TEST*/ | ||
| 155 |