| 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[] = "Partial hyperbolic singular value decomposition (HSVD) of matrix generated by plane rotations.\n\n" | ||
| 12 | "The command line options are:\n" | ||
| 13 | " -p <p>, number of rows in upper part of matrix A.\n" | ||
| 14 | " -q <q>, number of rows in lower part of matrix A.\n" | ||
| 15 | " -n <n>, number of columns of matrix A.\n" | ||
| 16 | " -k <k>, condition number.\n" | ||
| 17 | " -d <d>, density.\n" | ||
| 18 | " -prog, show progress of matrix generation\n"; | ||
| 19 | |||
| 20 | #include <slepcsvd.h> | ||
| 21 | #include <slepcblaslapack.h> | ||
| 22 | |||
| 23 | /* The problem matrix is the (p+q)-by-n matrix | ||
| 24 | |||
| 25 | A = | Q1 * D1 * U | | ||
| 26 | | 1/2*Q2 * D2 * U | | ||
| 27 | |||
| 28 | with p>=n, where Q1 (p-by-p), Q2 (q-by-q) and U (n-by-n) are orthonormal | ||
| 29 | matrices, D1 is a diagonal p-by-n matrix with elements | ||
| 30 | D1(i,i) = sqrt(5/4)*d(i), if i<=q | ||
| 31 | = d(i), if i>q | ||
| 32 | D2 is a diagonal q-by-n matrix with elements D2(i,i)= d(i), | ||
| 33 | and d is an n-element vector with elements distributed exponentially | ||
| 34 | from 1/k to 1. | ||
| 35 | |||
| 36 | The signature matrix is Omega = diag(I_p,-I_q) | ||
| 37 | |||
| 38 | The hyperbolic singular values of A are equal to the elements of vector d. | ||
| 39 | |||
| 40 | A is generated by random plane rotations applied to the diagonal matrices | ||
| 41 | D1 and D2 | ||
| 42 | |||
| 43 | This test case is based on the paper: | ||
| 44 | |||
| 45 | [1] A. Bojanczyk, N.J. Higham, and H. Patel, "Solving the Indefinite Least | ||
| 46 | Squares Problem by Hyperbolic QR Factorization", SIAM J. Matrix Anal. | ||
| 47 | Appl. 24(4):914-931 2003 | ||
| 48 | */ | ||
| 49 | |||
| 50 | /* Timing */ | ||
| 51 | #ifdef TIMING | ||
| 52 | static PetscReal tt_rotr,tt_rotc,tt_assmr=0.0,tt_assmc=0.0,tt_getr=0.0,tt_getc=0.0,tt_setvalr=0.0,tt_setvalc=0.0; | ||
| 53 | #define time_now(t) t=MPI_Wtime(); | ||
| 54 | #define time_diff(tacum,t0,t,t1) {t=MPI_Wtime(); tacum += t-t0; t1=t;} | ||
| 55 | #else | ||
| 56 | #define time_diff(tacum,t0,t,t1) | ||
| 57 | #define time_now(t) | ||
| 58 | #endif | ||
| 59 | |||
| 60 | struct _p_XMat { | ||
| 61 | MPI_Comm comm; | ||
| 62 | PetscInt M,N; /* global number of rows and columns */ | ||
| 63 | PetscInt Istart,Iend; /* Index of 1st row and last local rows */ | ||
| 64 | PetscInt *nzr; /* number of non-zeros in each row */ | ||
| 65 | PetscInt **rcol; /* row columns */ | ||
| 66 | PetscInt *ranges; /* proc row ranges */ | ||
| 67 | PetscScalar **rval; /* row values */ | ||
| 68 | }; | ||
| 69 | |||
| 70 | typedef struct _p_XMat *XMat; | ||
| 71 | |||
| 72 | /* Iterator over non-zero rows of 2 columns */ | ||
| 73 | struct _p_ColsNzIter { | ||
| 74 | XMat A; | ||
| 75 | PetscInt j1,j2; /* Columns considered */ | ||
| 76 | PetscInt iloc; /* current local row */ | ||
| 77 | }; | ||
| 78 | |||
| 79 | typedef struct _p_ColsNzIter *ColsNzIter; | ||
| 80 | |||
| 81 | 120 | static PetscErrorCode XMatCreate(MPI_Comm comm,XMat *A,PetscInt Istart,PetscInt Iend,PetscInt M,PetscInt N) | |
| 82 | { | ||
| 83 | 120 | PetscInt i; | |
| 84 | 120 | PetscMPIInt np; | |
| 85 | |||
| 86 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
120 | PetscFunctionBeginUser; |
| 87 |
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.
|
120 | PetscCall(PetscNew(A)); |
| 88 | 120 | (*A)->M = M; | |
| 89 | 120 | (*A)->N = N; | |
| 90 | 120 | (*A)->Istart = Istart; | |
| 91 | 120 | (*A)->Iend = Iend; | |
| 92 | 120 | (*A)->comm = comm; | |
| 93 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
120 | PetscCallMPI(MPI_Comm_size(comm,&np)); |
| 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.
|
120 | PetscCall(PetscMalloc4(Iend-Istart,&(*A)->nzr,Iend-Istart,&(*A)->rcol,Iend-Istart,&(*A)->rval,np+1,&(*A)->ranges)); |
| 95 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=0; i<Iend-Istart; i++) { |
| 96 | 15200 | (*A)->nzr[i] = 0; | |
| 97 | 15200 | (*A)->rcol[i] = NULL; | |
| 98 | 15200 | (*A)->rval[i] = NULL; | |
| 99 | } | ||
| 100 |
15/30✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 2 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✓ Branch 27 taken 2 times.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
|
240 | PetscCallMPI(MPI_Allgather(&Istart,1,MPIU_INT,(*A)->ranges,1,MPIU_INT,comm)); |
| 101 | 120 | (*A)->ranges[np]=M; | |
| 102 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
120 | PetscFunctionReturn(PETSC_SUCCESS); |
| 103 | } | ||
| 104 | |||
| 105 | 120 | static PetscErrorCode XMatDestroy(XMat *A) | |
| 106 | { | ||
| 107 | 120 | PetscInt m,i; | |
| 108 | |||
| 109 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
120 | PetscFunctionBeginUser; |
| 110 | 120 | m = (*A)->Iend-(*A)->Istart; | |
| 111 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=0; i<m; i++) { |
| 112 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15200 | if ((*A)->nzr[i]>0) { |
| 113 |
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.
|
15200 | PetscCall(PetscFree2((*A)->rcol[i],(*A)->rval[i])); |
| 114 | } | ||
| 115 | } | ||
| 116 |
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.
|
120 | PetscCall(PetscFree4((*A)->nzr,(*A)->rcol,(*A)->rval,(*A)->ranges)); |
| 117 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
120 | PetscCall(PetscFree(*A)); |
| 118 | 120 | *A = NULL; | |
| 119 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
120 | PetscFunctionReturn(PETSC_SUCCESS); |
| 120 | } | ||
| 121 | |||
| 122 | 12640 | static PetscErrorCode XMatGetSize(XMat A,PetscInt *M,PetscInt *N) | |
| 123 | { | ||
| 124 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
12640 | PetscFunctionBeginUser; |
| 125 | 12640 | *M = A->M; | |
| 126 | 12640 | *N = A->N; | |
| 127 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
12640 | PetscFunctionReturn(PETSC_SUCCESS); |
| 128 | } | ||
| 129 | |||
| 130 | 17760 | static PetscErrorCode XMatGetOwnershipRange(XMat A,PetscInt *Istart,PetscInt *Iend) | |
| 131 | { | ||
| 132 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
17760 | PetscFunctionBeginUser; |
| 133 | 17760 | *Istart = A->Istart; | |
| 134 | 17760 | *Iend = A->Iend; | |
| 135 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
17760 | PetscFunctionReturn(PETSC_SUCCESS); |
| 136 | } | ||
| 137 | |||
| 138 | 1600 | static PetscErrorCode XMatGetOwnershipRanges(XMat A,PetscInt **ranges) | |
| 139 | { | ||
| 140 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1600 | PetscFunctionBeginUser; |
| 141 | 1600 | *ranges = A->ranges; | |
| 142 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1600 | PetscFunctionReturn(PETSC_SUCCESS); |
| 143 | } | ||
| 144 | |||
| 145 | 23680 | static PetscErrorCode XMatGetRow(XMat A,PetscInt i,PetscInt *nc,PetscInt **vj,PetscScalar **va) | |
| 146 | { | ||
| 147 | 23680 | PetscInt iloc=i-A->Istart; | |
| 148 | |||
| 149 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
23680 | PetscFunctionBeginUser; |
| 150 | 23680 | *nc = A->nzr[iloc]; | |
| 151 | 23680 | *vj = A->rcol[iloc]; | |
| 152 | 23680 | *va = A->rval[iloc]; | |
| 153 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
23680 | PetscFunctionReturn(PETSC_SUCCESS); |
| 154 | } | ||
| 155 | |||
| 156 | 23680 | static PetscErrorCode XMatRestoreRow(XMat A,PetscInt i,PetscInt *nc,PetscInt **vj,PetscScalar **va) | |
| 157 | { | ||
| 158 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
23680 | PetscFunctionBeginUser; |
| 159 | 23680 | *nc = 0; | |
| 160 | 23680 | *vj = NULL; | |
| 161 | 23680 | *va = NULL; | |
| 162 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
23680 | PetscFunctionReturn(PETSC_SUCCESS); |
| 163 | } | ||
| 164 | |||
| 165 | 121600 | static PetscErrorCode XMatChangeRow(XMat A,PetscInt i,PetscInt nc,PetscInt vj[],PetscScalar va[]) | |
| 166 | { | ||
| 167 | 121600 | PetscInt iloc = i-A->Istart; | |
| 168 | |||
| 169 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
121600 | PetscFunctionBeginUser; |
| 170 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
121600 | if (A->nzr[iloc]>0) { |
| 171 |
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.
|
120880 | PetscCall(PetscFree2(A->rcol[iloc],A->rval[iloc])); |
| 172 | } | ||
| 173 | 121600 | A->nzr[iloc] = nc; | |
| 174 | 121600 | A->rcol[iloc] = vj; | |
| 175 | 121600 | A->rval[iloc] = va; | |
| 176 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
121600 | PetscFunctionReturn(PETSC_SUCCESS); |
| 177 | } | ||
| 178 | |||
| 179 | /* Given a row with non-zero elements in columns indicated by rcol, get the | ||
| 180 | position where the non-zero element of column j is, or where it should be inserted. | ||
| 181 | Returns true if there is a non-zero element for the column, false otherwise */ | ||
| 182 | 2264000 | static PetscBool getcolpos(const PetscInt rcol[],PetscInt nz,PetscInt j,PetscInt *k) | |
| 183 | { | ||
| 184 | 2264000 | PetscInt ka,kb; | |
| 185 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2264000 | if (nz==0) { |
| 186 | 239520 | *k = 0; | |
| 187 | 239520 | return PETSC_FALSE; | |
| 188 | } | ||
| 189 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2024480 | if (j<=rcol[0]) { |
| 190 | 675040 | *k = 0; | |
| 191 | 675040 | return (rcol[0]==j)? PETSC_TRUE: PETSC_FALSE; | |
| 192 | } | ||
| 193 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1349440 | if (rcol[nz-1]<=j) { |
| 194 | 456800 | *k = nz-(rcol[nz-1]==j); | |
| 195 | 456800 | return (rcol[nz-1]==j)? PETSC_TRUE: PETSC_FALSE; | |
| 196 | } | ||
| 197 | /* Bisection: rcol[ka] < j < rcol[kb] */ | ||
| 198 | 892640 | ka = 0; kb = nz-1; | |
| 199 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2770480 | while (ka+1<kb) { |
| 200 | 1950480 | *k = (ka+kb)/2; | |
| 201 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1950480 | if (rcol[*k]<=j) { |
| 202 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
980320 | if (rcol[*k]==j) return PETSC_TRUE; |
| 203 | else ka = *k; | ||
| 204 | } else kb = *k; | ||
| 205 | } | ||
| 206 | 820000 | *k = kb; | |
| 207 | 820000 | return PETSC_FALSE; | |
| 208 | } | ||
| 209 | |||
| 210 | /* Only local elements can be set */ | ||
| 211 | 14400 | static PetscErrorCode XMatSetValue(XMat A,PetscInt i,PetscInt j,PetscScalar x,InsertMode dummy) | |
| 212 | { | ||
| 213 | 14400 | PetscInt nz,iloc,k,kx,*col1,*col2; | |
| 214 | 14400 | PetscScalar *val1,*val2; | |
| 215 | |||
| 216 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
14400 | PetscFunctionBeginUser; |
| 217 |
4/10✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 10 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 10 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
|
14400 | PetscCheck(i>=A->Istart && i<A->Iend && j>=0 && j<A->N,A->comm,PETSC_ERR_ARG_OUTOFRANGE,"Incorrect index of row/column"); |
| 218 | 14400 | iloc = i-A->Istart; | |
| 219 | 14400 | nz = A->nzr[iloc]; | |
| 220 | 14400 | col1 = A->rcol[iloc]; | |
| 221 | 14400 | val1 = A->rval[iloc]; | |
| 222 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
14400 | if (!getcolpos(col1,nz,j,&k)) { |
| 223 | 14400 | A->nzr[iloc]++; | |
| 224 | /* Replace row */ | ||
| 225 |
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.
|
14400 | PetscCall(PetscMalloc2(nz+1,&col2,nz+1,&val2)); |
| 226 | 14400 | A->rcol[iloc] = col2; | |
| 227 | 14400 | A->rval[iloc] = val2; | |
| 228 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
14400 | for (kx=0; kx<k; kx++) { |
| 229 | ✗ | col2[kx] = col1[kx]; | |
| 230 | ✗ | val2[kx] = val1[kx]; | |
| 231 | } | ||
| 232 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
14400 | for (; kx<nz; kx++) { |
| 233 | ✗ | col2[kx+1] = col1[kx]; | |
| 234 | ✗ | val2[kx+1] = val1[kx]; | |
| 235 | } | ||
| 236 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
14400 | if (nz>0) { |
| 237 | ✗ | PetscCall(PetscFree2(col1,val1)); | |
| 238 | } | ||
| 239 | 14400 | col1 = col2; | |
| 240 | 14400 | val1 = val2; | |
| 241 | } | ||
| 242 | 14400 | col1[k] = j; | |
| 243 | 14400 | val1[k] = x; | |
| 244 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
14400 | PetscFunctionReturn(PETSC_SUCCESS); |
| 245 | } | ||
| 246 | |||
| 247 | /* Convert to PETSc Mat */ | ||
| 248 | 120 | static PetscErrorCode XMatConvert(XMat XA,Mat A) | |
| 249 | { | ||
| 250 | 120 | PetscInt i,iloc; | |
| 251 | |||
| 252 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
120 | PetscFunctionBeginUser; |
| 253 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=XA->Istart; i<XA->Iend; i++) { |
| 254 | 15200 | iloc = i-XA->Istart; | |
| 255 |
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.
|
15200 | PetscCall(MatSetValues(A,1,&i,XA->nzr[iloc],XA->rcol[iloc],XA->rval[iloc],INSERT_VALUES)); |
| 256 | } | ||
| 257 |
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.
|
120 | PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 258 |
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.
|
120 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 259 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
24 | PetscFunctionReturn(PETSC_SUCCESS); |
| 260 | } | ||
| 261 | |||
| 262 | /* gets the number of nonzeros in the matrix on this MPI rank */ | ||
| 263 | 120 | static PetscErrorCode XMatGetNumberNonzeros(XMat A,PetscInt *nz) | |
| 264 | { | ||
| 265 | 120 | PetscInt i; | |
| 266 | |||
| 267 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
120 | PetscFunctionBeginUser; |
| 268 | 120 | *nz = 0; | |
| 269 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=0; i<A->Iend-A->Istart; i++) *nz += A->nzr[i]; |
| 270 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
120 | PetscFunctionReturn(PETSC_SUCCESS); |
| 271 | } | ||
| 272 | |||
| 273 | 8880 | static PetscErrorCode XMatCreateColsNzIter(XMat A,PetscInt j1,PetscInt j2,ColsNzIter *iter) | |
| 274 | { | ||
| 275 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
8880 | PetscFunctionBeginUser; |
| 276 |
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.
|
8880 | PetscCall(PetscNew(iter)); |
| 277 | 8880 | (*iter)->A = A; | |
| 278 | 8880 | (*iter)->j1 = j1; | |
| 279 | 8880 | (*iter)->j2 = j2; | |
| 280 | 8880 | (*iter)->iloc = 0; | |
| 281 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
8880 | PetscFunctionReturn(PETSC_SUCCESS); |
| 282 | } | ||
| 283 | |||
| 284 | 8880 | static PetscErrorCode ColitDestroy(ColsNzIter *iter) | |
| 285 | { | ||
| 286 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
8880 | PetscFunctionBeginUser; |
| 287 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
8880 | PetscCall(PetscFree(*iter)); |
| 288 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1776 | PetscFunctionReturn(PETSC_SUCCESS); |
| 289 | } | ||
| 290 | |||
| 291 | 115600 | static PetscErrorCode ColitNextRow(ColsNzIter iter,PetscInt *i,PetscScalar **paij1,PetscScalar **paij2) | |
| 292 | { | ||
| 293 | 115600 | PetscInt iloc = iter->iloc,*rcol,nz,k,k1,k2,kx,jx,*rcol2; | |
| 294 | 115600 | PetscScalar *rval,*rval2; | |
| 295 | 115600 | PetscBool found1, found2; | |
| 296 | 115600 | XMat A=iter->A; | |
| 297 | |||
| 298 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
115600 | PetscFunctionBeginUser; |
| 299 | 115600 | *i = -1; | |
| 300 | 115600 | *paij1 = *paij2 = NULL; | |
| 301 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
115600 | if (iloc>=0) { |
| 302 | 114880 | rcol = A->rcol[iloc]; | |
| 303 | 114880 | nz = A->nzr[iloc]; | |
| 304 | 2134720 | while (PETSC_TRUE) { | |
| 305 | 1124800 | found1 = getcolpos(rcol,nz,iter->j1,&k1); | |
| 306 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1124800 | found2 = getcolpos(PetscSafePointerPlusOffset(rcol,k1+found1),nz-k1-found1,iter->j2,&k2); |
| 307 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1124800 | if (found1 || found2) break; |
| 308 | 1018080 | iloc++; | |
| 309 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1018080 | if (iloc>=A->Iend-A->Istart) { |
| 310 | iloc = -1; | ||
| 311 | break; | ||
| 312 | } | ||
| 313 | 1009920 | rcol = A->rcol[iloc]; | |
| 314 | 1009920 | nz = A->nzr[iloc]; | |
| 315 | } | ||
| 316 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
114880 | if (iloc>=0) { |
| 317 | 106720 | k2 += k1 + 1; | |
| 318 | 106720 | rval = A->rval[iloc]; | |
| 319 | 106720 | *i = A->Istart+iloc; | |
| 320 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
106720 | if (!found1 || !found2) { |
| 321 | /* Copy row to a new one */ | ||
| 322 | 98080 | rcol2 = rcol; | |
| 323 | 98080 | rval2 = rval; | |
| 324 |
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.
|
98080 | PetscCall(PetscMalloc2(nz+1,&rcol,nz+1,&rval)); |
| 325 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
98080 | if (found1) { |
| 326 | 54720 | kx = k2; | |
| 327 | 54720 | jx = iter->j2; | |
| 328 | } else { | ||
| 329 | 43360 | kx = k1; | |
| 330 | 43360 | jx = iter->j1; | |
| 331 | } | ||
| 332 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
639760 | for (k=0; k<kx; k++) { |
| 333 | 541680 | rcol[k] = rcol2[k]; | |
| 334 | 541680 | rval[k] = rval2[k]; | |
| 335 | } | ||
| 336 | 98080 | rcol[kx] = jx; | |
| 337 | 98080 | rval[kx] = 0.0; | |
| 338 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
727840 | for (k=kx+1; k<nz+1; k++) { |
| 339 | 629760 | rcol[k] = rcol2[k-1]; | |
| 340 | 629760 | rval[k] = rval2[k-1]; | |
| 341 | } | ||
| 342 |
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.
|
98080 | PetscCall(XMatChangeRow(A,*i,nz+1,rcol,rval)); |
| 343 | } | ||
| 344 | 106720 | *paij1 = &rval[k1]; | |
| 345 | 106720 | *paij2 = &rval[k2]; | |
| 346 | 106720 | iloc++; | |
| 347 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
106720 | if (iloc>=A->Iend-A->Istart) iloc = -1; |
| 348 | } | ||
| 349 | 114880 | iter->iloc = iloc; | |
| 350 | } | ||
| 351 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
23120 | PetscFunctionReturn(PETSC_SUCCESS); |
| 352 | } | ||
| 353 | |||
| 354 | 1600 | static PetscErrorCode SendrecvRow(XMat A,PetscInt nc1,const PetscInt vj1[],const PetscScalar va1[],PetscInt i2,PetscInt *nc2,PetscInt vj2[],PetscScalar va2[]) | |
| 355 | { | ||
| 356 | 1600 | PetscInt *ranges,N=A->N; | |
| 357 | 1600 | PetscMPIInt rank,naux,len1,len2; | |
| 358 | 1600 | MPI_Status st; | |
| 359 | |||
| 360 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1600 | PetscFunctionBeginUser; |
| 361 |
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.
|
1600 | PetscCall(XMatGetOwnershipRanges(A,&ranges)); |
| 362 |
2/4✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
|
320 | PetscAssertAbort(i2>=0 && i2<A->M,A->comm,PETSC_ERR_ARG_OUTOFRANGE,"Incorrect index of row"); |
| 363 | /* Find owner of row i2 */ | ||
| 364 | 320 | rank=0; | |
| 365 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2400 | while (ranges[rank+1]<=i2) rank++; |
| 366 | /* Send row i1, receive row i2 */ | ||
| 367 |
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.
|
1600 | PetscCall(PetscMPIIntCast(nc1,&len1)); |
| 368 |
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.
|
1600 | PetscCall(PetscMPIIntCast(N,&len2)); |
| 369 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
1600 | PetscCallMPI(MPI_Sendrecv(vj1,len1,MPIU_INT,rank,0,vj2,len2,MPIU_INT,rank,0,A->comm,&st)); |
| 370 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
1600 | PetscCallMPI(MPI_Sendrecv(va1,len1,MPIU_SCALAR,rank,0,va2,len2,MPIU_SCALAR,rank,0,A->comm,MPI_STATUS_IGNORE)); |
| 371 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
1600 | PetscCallMPI(MPI_Get_count(&st,MPIU_INT,&naux)); |
| 372 | 1600 | *nc2 = (PetscInt)naux; | |
| 373 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1600 | PetscFunctionReturn(PETSC_SUCCESS); |
| 374 | } | ||
| 375 | |||
| 376 | 12640 | static PetscErrorCode PadRows(PetscInt nc1,const PetscInt vj1[],const PetscScalar va1[],PetscInt nc2,const PetscInt vj2[],const PetscScalar va2[],PetscInt *nc,PetscInt **vjj1,PetscScalar **vaa1,PetscInt **vjj2,PetscScalar **vaa2,PetscInt *iwork,PetscScalar *swork) | |
| 377 | { | ||
| 378 | 12640 | PetscInt k1,k2,k,n1=nc1+nc2,*vjjaux=iwork; | |
| 379 | 12640 | PetscScalar *vaa1aux=swork,*vaa2aux=swork+n1; | |
| 380 | |||
| 381 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
12640 | PetscFunctionBeginUser; |
| 382 | 12640 | *nc=0; | |
| 383 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
105600 | for (k1=k2=k=0; k1<nc1 && k2<nc2; ) { |
| 384 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
92960 | if (vj1[k1]<vj2[k2]) { |
| 385 | /* Take next element from first row */ | ||
| 386 | 43240 | vaa1aux[k] = va1[k1]; | |
| 387 | 43240 | vaa2aux[k] = 0.0; | |
| 388 | 43240 | vjjaux[k++] = vj1[k1++]; | |
| 389 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
49720 | } else if (vj1[k1]>vj2[k2]) { |
| 390 | /* Take next element from second row */ | ||
| 391 | 43120 | vaa1aux[k] = 0.0; | |
| 392 | 43120 | vaa2aux[k] = va2[k2]; | |
| 393 | 43120 | vjjaux[k++] = vj2[k2++]; | |
| 394 | } else { | ||
| 395 | /* Take next element from both rows */ | ||
| 396 | 6600 | vaa1aux[k] = va1[k1]; | |
| 397 | 6600 | vaa2aux[k] = va2[k2++]; | |
| 398 | 6600 | vjjaux[k++] = vj1[k1++]; | |
| 399 | } | ||
| 400 | } | ||
| 401 | |||
| 402 | /* We reached the end of one of the rows. Continue with the other one */ | ||
| 403 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
25560 | while (k1<nc1) { |
| 404 | /* Take next element from first row */ | ||
| 405 | 12920 | vaa1aux[k] = va1[k1]; | |
| 406 | 12920 | vaa2aux[k] = 0.0; | |
| 407 | 12920 | vjjaux[k++] = vj1[k1++]; | |
| 408 | } | ||
| 409 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
25800 | while (k2<nc2) { |
| 410 | /* Take next element from second row */ | ||
| 411 | 13160 | vaa1aux[k] = 0.0; | |
| 412 | 13160 | vaa2aux[k] = va2[k2]; | |
| 413 | 13160 | vjjaux[k++] = vj2[k2++]; | |
| 414 | } | ||
| 415 | 12640 | *nc=k; | |
| 416 | |||
| 417 | /* Copy rows (each row must be allocated separately)*/ | ||
| 418 |
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.
|
12640 | PetscCall(PetscMalloc2(k,vjj1,k,vaa1)); |
| 419 |
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.
|
12640 | PetscCall(PetscMalloc2(k,vjj2,k,vaa2)); |
| 420 |
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.
|
12640 | PetscCall(PetscArraycpy(*vjj1,vjjaux,k)); |
| 421 |
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.
|
12640 | PetscCall(PetscArraycpy(*vjj2,vjjaux,k)); |
| 422 |
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.
|
12640 | PetscCall(PetscArraycpy(*vaa1,vaa1aux,k)); |
| 423 |
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.
|
12640 | PetscCall(PetscArraycpy(*vaa2,vaa2aux,k)); |
| 424 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
2528 | PetscFunctionReturn(PETSC_SUCCESS); |
| 425 | } | ||
| 426 | |||
| 427 | 1600 | static inline void MyAxpby(PetscInt n,PetscScalar a,PetscScalar x[],PetscScalar b,const PetscScalar y[]) | |
| 428 | { | ||
| 429 | 1600 | PetscInt i; | |
| 430 | |||
| 431 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
12800 | for (i=0; i<n; i++) x[i] = a*x[i]+b*y[i]; |
| 432 | } | ||
| 433 | |||
| 434 | /* Apply plane rotation on rows i1, i2 of A */ | ||
| 435 | 17760 | static PetscErrorCode RotateRows(XMat A,PetscInt i1,PetscInt i2,PetscReal c,PetscReal s,PetscInt *nzloc,PetscInt *iwork,PetscScalar *swork) | |
| 436 | { | ||
| 437 | 17760 | PetscInt M,N,Istart,Iend,nc,nc1,nc2,*vj1=NULL,*vj2=NULL,*vjj1=NULL,*vjj2=NULL,*iworkx=iwork; | |
| 438 | 17760 | PetscBLASInt nc_,one=1; | |
| 439 | 17760 | PetscScalar *va1=NULL,*va2=NULL,*vaa1=NULL,*vaa2=NULL,*sworkx=swork; | |
| 440 | 17760 | PetscBool i1mine, i2mine; | |
| 441 | #ifdef TIMING | ||
| 442 | PetscReal t,t0; | ||
| 443 | #endif | ||
| 444 | |||
| 445 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
17760 | PetscFunctionBeginUser; |
| 446 |
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.
|
17760 | PetscCall(XMatGetOwnershipRange(A,&Istart,&Iend)); |
| 447 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
17760 | i1mine = (Istart<=i1 && i1<Iend)? PETSC_TRUE: PETSC_FALSE; |
| 448 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
17760 | i2mine = (Istart<=i2 && i2<Iend)? PETSC_TRUE: PETSC_FALSE; |
| 449 | |||
| 450 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
17760 | if (i1mine || i2mine) { |
| 451 |
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.
|
12640 | PetscCall(XMatGetSize(A,&M,&N)); |
| 452 | |||
| 453 | /* Get local row(s) (at least 1, possibly 2) */ | ||
| 454 | 12640 | time_now(t0); | |
| 455 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
12640 | if (i1mine) PetscCall(XMatGetRow(A,i1,&nc1,&vj1,&va1)); |
| 456 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
12640 | if (i2mine) PetscCall(XMatGetRow(A,i2,&nc2,&vj2,&va2)); |
| 457 | 12640 | time_diff(tt_getr,t0,t,t0); | |
| 458 | /* Get remote row (at most 1, possibly none) */ | ||
| 459 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
|
12640 | if (!i2mine) { |
| 460 | 800 | vj2 = iworkx; | |
| 461 | 800 | iworkx += N; | |
| 462 | 800 | va2 = sworkx; | |
| 463 | 800 | sworkx += N; | |
| 464 |
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.
|
800 | PetscCall(SendrecvRow(A,nc1,vj1,va1,i2,&nc2,vj2,va2)); |
| 465 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
11840 | } else if (!i1mine) { |
| 466 | 800 | vj1 = iworkx; | |
| 467 | 800 | iworkx += N; | |
| 468 | 800 | va1 = sworkx; | |
| 469 | 800 | sworkx += N; | |
| 470 |
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.
|
800 | PetscCall(SendrecvRow(A,nc2,vj2,va2,i1,&nc1,vj1,va1)); |
| 471 | } | ||
| 472 | |||
| 473 |
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.
|
12640 | PetscCall(PadRows(nc1,vj1,va1,nc2,vj2,va2,&nc,&vjj1,&vaa1,&vjj2,&vaa2,iworkx,sworkx)); |
| 474 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
12640 | if (i1mine) { |
| 475 | 11840 | *nzloc += nc-nc1; | |
| 476 |
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.
|
11840 | PetscCall(XMatRestoreRow(A,i1,&nc1,&vj1,&va1)); |
| 477 | } | ||
| 478 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
12640 | if (i2mine) { |
| 479 | 11840 | *nzloc += nc-nc2; | |
| 480 |
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.
|
11840 | PetscCall(XMatRestoreRow(A,i2,&nc2,&vj2,&va2)); |
| 481 | } | ||
| 482 | |||
| 483 | /* Compute rotation and update matrix */ | ||
| 484 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
12640 | if (nc) { |
| 485 |
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.
|
12560 | PetscCall(PetscBLASIntCast(nc,&nc_)); |
| 486 |
12/22✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✓ Branch 19 taken 2 times.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
|
12560 | if (i1mine && i2mine) PetscCallBLAS("BLASrot",BLASMIXEDrot_(&nc_,vaa1,&one,vaa2,&one,&c,&s)); |
| 487 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1600 | else if (i1mine) { |
| 488 | 800 | MyAxpby(nc,c,vaa1,s,vaa2); | |
| 489 |
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.
|
800 | PetscCall(PetscFree2(vjj2,vaa2)); |
| 490 | } else { | ||
| 491 | 800 | MyAxpby(nc,c,vaa2,-s,vaa1); | |
| 492 |
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.
|
800 | PetscCall(PetscFree2(vjj1,vaa1)); |
| 493 | } | ||
| 494 | 12560 | time_now(t0); | |
| 495 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
12560 | if (i1mine) PetscCall(XMatChangeRow(A,i1,nc,vjj1,vaa1)); |
| 496 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
12560 | if (i2mine) PetscCall(XMatChangeRow(A,i2,nc,vjj2,vaa2)); |
| 497 | 3552 | time_diff(tt_setvalr,t0,t,t0); | |
| 498 | } | ||
| 499 | } | ||
| 500 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
3552 | PetscFunctionReturn(PETSC_SUCCESS); |
| 501 | } | ||
| 502 | |||
| 503 | /* Apply plane rotation on columns j1, j2 of A */ | ||
| 504 | 8880 | static PetscErrorCode RotateCols(XMat A,PetscInt j1,PetscInt j2,PetscReal c,PetscReal s,PetscInt *nzloc) | |
| 505 | { | ||
| 506 | 8880 | PetscInt i; | |
| 507 | 8880 | PetscScalar aux1,aux2,*paij1,*paij2; | |
| 508 | 8880 | ColsNzIter iter; | |
| 509 | |||
| 510 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
8880 | PetscFunctionBeginUser; |
| 511 |
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.
|
8880 | PetscCall(XMatCreateColsNzIter(A,j1,j2,&iter)); |
| 512 | 222320 | while (PETSC_TRUE) { | |
| 513 |
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.
|
115600 | PetscCall(ColitNextRow(iter,&i,&paij1,&paij2)); |
| 514 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
115600 | if (i<0) break; |
| 515 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
106720 | if (*paij1**paij2==0.0) (*nzloc)++; |
| 516 | 106720 | aux1 = *paij1*c+*paij2*s; | |
| 517 | 106720 | aux2 = -*paij1*s+*paij2*c; | |
| 518 | 106720 | *paij1 = aux1; | |
| 519 | 106720 | *paij2 = aux2; | |
| 520 | } | ||
| 521 |
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.
|
8880 | PetscCall(ColitDestroy(&iter)); |
| 522 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1776 | PetscFunctionReturn(PETSC_SUCCESS); |
| 523 | } | ||
| 524 | |||
| 525 | /* Generate a pair of indices i1, i2, where 0 <= i1 < i2 < n */ | ||
| 526 | 26640 | static PetscErrorCode GetIndexPair(PetscRandom rand,PetscInt n,PetscInt *i1,PetscInt *i2) | |
| 527 | { | ||
| 528 | 26640 | PetscInt aux; | |
| 529 | 26640 | PetscReal x; | |
| 530 | |||
| 531 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
26640 | PetscFunctionBeginUser; |
| 532 |
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.
|
26640 | PetscCall(PetscRandomGetValueReal(rand,&x)); |
| 533 | 26640 | *i1 = (PetscInt)(x*n); | |
| 534 |
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.
|
26640 | PetscCall(PetscRandomGetValueReal(rand,&x)); |
| 535 | 26640 | *i2 = (PetscInt)(x*(n-1)); | |
| 536 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
26640 | if (*i2>=*i1) (*i2)++; |
| 537 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
26640 | if (*i1>*i2) { |
| 538 | 13320 | aux = *i1; | |
| 539 | 13320 | *i1 = *i2; | |
| 540 | 13320 | *i2 = aux; | |
| 541 | } | ||
| 542 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
5328 | PetscFunctionReturn(PETSC_SUCCESS); |
| 543 | } | ||
| 544 | |||
| 545 | 120 | int main(int argc,char **argv) | |
| 546 | { | ||
| 547 | 120 | Mat A,Omega; /* operator matrix, signature matrix */ | |
| 548 | 120 | XMat XA; | |
| 549 | 120 | SVD svd; /* singular value problem solver context */ | |
| 550 | 120 | Vec u,v,vomega,*U,*V; | |
| 551 | 120 | MatType Atype; | |
| 552 | 120 | PetscReal tol,lev1=0.0,lev2=0.0,k=1e3,dens=0.1,log0,loginc,sq,x,angle,c,s; | |
| 553 | 120 | PetscScalar *swork; | |
| 554 | 120 | PetscInt P=110,Q=80,N=100,i,j,Istart,Iend,nconv,nsv,i1,i2,nroth,nrotv,*iwork,nnzwanted,nz,nzloc; | |
| 555 | 120 | PetscBool flg,flgp,flgq,flgn,terse,skiporth=PETSC_FALSE,progress=PETSC_FALSE; | |
| 556 | 120 | PetscRandom rand; | |
| 557 | 120 | MatInfo Ainfo; | |
| 558 | #ifdef TIMING | ||
| 559 | PetscReal t,t0; | ||
| 560 | #endif | ||
| 561 | |||
| 562 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
120 | PetscFunctionBeginUser; |
| 563 |
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.
|
120 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 564 | |||
| 565 |
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.
|
120 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-p",&P,&flgp)); |
| 566 |
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.
|
120 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-q",&Q,&flgq)); |
| 567 |
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.
|
120 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&N,&flgn)); |
| 568 |
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.
|
120 | PetscCall(PetscOptionsGetReal(NULL,NULL,"-k",&k,&flg)); |
| 569 |
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.
|
120 | PetscCall(PetscOptionsGetReal(NULL,NULL,"-d",&dens,&flg)); |
| 570 |
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.
|
120 | PetscCall(PetscOptionsGetBool(NULL,NULL,"-prog",&progress,&flg)); |
| 571 |
2/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
120 | if (!flgp && flgn) P=N; |
| 572 |
2/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
120 | if (!flgq && flgn) Q=N; |
| 573 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
120 | PetscCheck(P>=N,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Parameter p must be >= n"); |
| 574 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
120 | PetscCheck(k>=1.0,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Parameter k must be >= 1.0"); |
| 575 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
120 | PetscCheck(dens<=1.0,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Parameter d must be <= 1.0"); |
| 576 | |||
| 577 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 578 | Build matrix that defines the hyperbolic singular value problem | ||
| 579 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 580 | |||
| 581 |
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.
|
120 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nHyperbolic singular value problem.\n\n")); |
| 582 |
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.
|
120 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n Matrix dimensions: (%" PetscInt_FMT "+%" PetscInt_FMT ")x%" PetscInt_FMT,P,Q,N)); |
| 583 |
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.
|
120 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,", using signature Omega=blkdiag(I_%" PetscInt_FMT ",-I_%" PetscInt_FMT ")\n\n",P,Q)); |
| 584 | |||
| 585 |
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.
|
120 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); |
| 586 |
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.
|
120 | PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,P+Q,N)); |
| 587 |
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.
|
120 | PetscCall(MatSetFromOptions(A)); |
| 588 |
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.
|
120 | PetscCall(MatSetOption(A,MAT_IGNORE_OFF_PROC_ENTRIES,PETSC_TRUE)); |
| 589 | |||
| 590 | /* Set diagonals */ | ||
| 591 |
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.
|
120 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 592 | |||
| 593 |
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.
|
120 | PetscCall(XMatCreate(MPI_COMM_WORLD,&XA,Istart,Iend,P+Q,N)); |
| 594 | 120 | log0 = PetscLogReal(1/k); | |
| 595 | 120 | loginc = -log0/(N-1); | |
| 596 | 120 | sq = PetscSqrtReal(5.0/4.0); | |
| 597 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=Istart;i<Iend;i++) { |
| 598 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
15200 | if (i<P && i<N) { |
| 599 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
8000 | if (i<Q) x = sq*PetscExpReal(log0+i*loginc); |
| 600 | 1600 | else x = PetscExpReal(log0+i*loginc); | |
| 601 |
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.
|
8000 | PetscCall(XMatSetValue(XA,i,i,x,INSERT_VALUES)); |
| 602 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
7200 | } else if (i>=P && i-P<N) { |
| 603 | 6400 | j = i-P; | |
| 604 |
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.
|
15200 | PetscCall(XMatSetValue(XA,i,j,PetscExpReal(log0+j*loginc),INSERT_VALUES)); |
| 605 | } | ||
| 606 | } | ||
| 607 | |||
| 608 | /* Apply plane rotations */ | ||
| 609 | 120 | nnzwanted = dens*(P+Q)*N; | |
| 610 |
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.
|
120 | PetscCall(XMatGetNumberNonzeros(XA,&nzloc)); |
| 611 |
15/30✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 2 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✓ Branch 27 taken 2 times.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
|
240 | PetscCallMPI(MPI_Allreduce(&nzloc,&nz,1,MPIU_INT,MPI_SUM,MPI_COMM_WORLD)); |
| 612 |
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.
|
120 | PetscCall(PetscRandomCreate(PETSC_COMM_SELF,&rand)); |
| 613 |
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.
|
120 | PetscCall(PetscRandomSetFromOptions(rand)); |
| 614 | 120 | nroth = nrotv = 0; | |
| 615 |
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.
|
120 | PetscCall(PetscMalloc2(2*N,&iwork,3*N,&swork)); |
| 616 | time_now(t0); | ||
| 617 | 9000 | while (nz<0.95*nnzwanted) { | |
| 618 | /* Plane rotation on 2 of the top P rows*/ | ||
| 619 |
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.
|
8880 | PetscCall(PetscRandomGetValueReal(rand,&angle)); |
| 620 | 8880 | c=PetscCosReal(angle); | |
| 621 | 8880 | s=PetscSinReal(angle); | |
| 622 |
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.
|
8880 | PetscCall(GetIndexPair(rand,P,&i1,&i2)); |
| 623 |
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.
|
8880 | PetscCall(RotateRows(XA,i1,i2,c,s,&nzloc,iwork,swork)); |
| 624 | 8880 | time_diff(tt_rotr,t0,t,t0); | |
| 625 | 8880 | nroth++; | |
| 626 | /* Plane rotation on 2 of the bottom Q rows*/ | ||
| 627 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
8880 | if (Q>1) { |
| 628 |
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.
|
8880 | PetscCall(PetscRandomGetValueReal(rand,&angle)); |
| 629 | 8880 | c=PetscCosReal(angle); | |
| 630 | 8880 | s=PetscSinReal(angle); | |
| 631 |
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.
|
8880 | PetscCall(GetIndexPair(rand,Q,&i1,&i2)); |
| 632 |
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.
|
8880 | PetscCall(RotateRows(XA,P+i1,P+i2,c,s,&nzloc,iwork,swork)); |
| 633 | 8880 | time_diff(tt_rotr,t0,t,t0); | |
| 634 | 8880 | nroth++; | |
| 635 | } | ||
| 636 | /* Plane rotation on 2 columns */ | ||
| 637 |
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.
|
8880 | PetscCall(PetscRandomGetValueReal(rand,&angle)); |
| 638 | 8880 | c=PetscCosReal(angle); | |
| 639 | 8880 | s=PetscSinReal(angle); | |
| 640 |
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.
|
8880 | PetscCall(GetIndexPair(rand,N,&i1,&i2)); |
| 641 |
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.
|
8880 | PetscCall(RotateCols(XA,i1,i2,c,s,&nzloc)); |
| 642 | 8880 | time_diff(tt_rotc,t0,t,t0); | |
| 643 | 8880 | nrotv++; | |
| 644 | /* Update total number of non-zeros */ | ||
| 645 |
15/30✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 2 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✓ Branch 27 taken 2 times.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
|
17760 | PetscCallMPI(MPI_Allreduce(&nzloc,&nz,1,MPIU_INT,MPI_SUM,MPI_COMM_WORLD)); |
| 646 |
5/10✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 8 times.
✓ Branch 5 taken 8 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✓ Branch 9 taken 2 times.
|
17880 | if (progress) PetscCall(PetscFPrintf(MPI_COMM_WORLD,stderr,"\rProgress: %.2f%%",(double)nz/nnzwanted*100)); |
| 647 | } | ||
| 648 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
120 | if (progress) PetscCall(PetscFPrintf(MPI_COMM_WORLD,stderr,"\r")); |
| 649 | |||
| 650 |
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.
|
120 | PetscCall(PetscFree2(iwork,swork)); |
| 651 |
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.
|
120 | PetscCall(PetscRandomDestroy(&rand)); |
| 652 | |||
| 653 |
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.
|
120 | PetscCall(XMatConvert(XA,A)); |
| 654 |
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.
|
120 | PetscCall(XMatDestroy(&XA)); |
| 655 |
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.
|
120 | PetscCall(MatGetInfo(A,MAT_GLOBAL_SUM,&Ainfo)); |
| 656 |
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.
|
120 | PetscCall(PetscPrintf(MPI_COMM_WORLD," Matrix created. %" PetscInt_FMT " rotations applied. nnz=%.0f. Density: %.4g\n\n",nroth+nrotv,Ainfo.nz_used,(double)Ainfo.nz_used/(1.0*(P+Q)*N))); |
| 657 | |||
| 658 | #ifdef TIMING | ||
| 659 | PetscCall(PetscPrintf(MPI_COMM_WORLD,"#T rot-rows: %6.3f get-rows: %6.3f setval-rows: %6.3f assm-rows: %.3f\n",tt_rotr,tt_getr,tt_setvalr,tt_assmr)); | ||
| 660 | PetscCall(PetscPrintf(MPI_COMM_WORLD,"#T rot cols: %6.3f get-cols: %6.3f setval-cols: %6.3f assm-cols: %.3f\n",tt_rotc,tt_getc,tt_setvalc,tt_assmc)); | ||
| 661 | #endif | ||
| 662 | |||
| 663 | /* scale by 0.5 the lower Q rows of A */ | ||
| 664 |
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.
|
120 | PetscCall(MatCreateVecs(A,NULL,&vomega)); |
| 665 |
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.
|
120 | PetscCall(VecSet(vomega,1.0)); |
| 666 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6520 | for (i=PetscMax(P,Istart);i<Iend;i++) { |
| 667 |
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.
|
6400 | PetscCall(VecSetValue(vomega,i,0.5,INSERT_VALUES)); |
| 668 | } | ||
| 669 |
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.
|
120 | PetscCall(VecAssemblyBegin(vomega)); |
| 670 |
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.
|
120 | PetscCall(VecAssemblyEnd(vomega)); |
| 671 |
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.
|
120 | PetscCall(MatDiagonalScale(A,vomega,NULL)); |
| 672 | |||
| 673 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 674 | Create the signature | ||
| 675 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 676 | |||
| 677 |
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.
|
120 | PetscCall(VecSet(vomega,1.0)); |
| 678 |
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.
|
120 | PetscCall(VecGetOwnershipRange(vomega,&Istart,&Iend)); |
| 679 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
15320 | for (i=Istart;i<Iend;i++) { |
| 680 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
15200 | if (i>=P) PetscCall(VecSetValue(vomega,i,-1.0,INSERT_VALUES)); |
| 681 | } | ||
| 682 |
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.
|
120 | PetscCall(VecAssemblyBegin(vomega)); |
| 683 |
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.
|
120 | PetscCall(VecAssemblyEnd(vomega)); |
| 684 | |||
| 685 |
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.
|
120 | PetscCall(MatGetType(A,&Atype)); |
| 686 |
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.
|
120 | PetscCall(MatCreate(PETSC_COMM_WORLD,&Omega)); |
| 687 |
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.
|
120 | PetscCall(MatSetSizes(Omega,PETSC_DECIDE,PETSC_DECIDE,P+Q,P+Q)); |
| 688 |
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.
|
120 | PetscCall(MatSetType(Omega,Atype)); |
| 689 |
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.
|
120 | PetscCall(MatDiagonalSet(Omega,vomega,INSERT_VALUES)); |
| 690 | |||
| 691 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 692 | Create the singular value solver and set various options | ||
| 693 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 694 | |||
| 695 |
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.
|
120 | PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd)); |
| 696 | |||
| 697 |
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.
|
120 | PetscCall(SVDSetOperators(svd,A,NULL)); |
| 698 |
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.
|
120 | PetscCall(SVDSetSignature(svd,vomega)); |
| 699 |
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.
|
120 | PetscCall(SVDSetProblemType(svd,SVD_HYPERBOLIC)); |
| 700 | |||
| 701 |
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.
|
120 | PetscCall(SVDSetFromOptions(svd)); |
| 702 | |||
| 703 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 704 | Solve the problem, display solution | ||
| 705 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 706 | |||
| 707 |
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.
|
120 | PetscCall(MatCreateVecs(A,&v,&u)); |
| 708 |
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.
|
120 | PetscCall(SVDSolve(svd)); |
| 709 | |||
| 710 | /* show detailed info unless -terse option is given by user */ | ||
| 711 |
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.
|
120 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 712 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
120 | if (terse) PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,NULL)); |
| 713 | else { | ||
| 714 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
| 715 | ✗ | PetscCall(SVDConvergedReasonView(svd,PETSC_VIEWER_STDOUT_WORLD)); | |
| 716 | ✗ | PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,PETSC_VIEWER_STDOUT_WORLD)); | |
| 717 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
| 718 | } | ||
| 719 | |||
| 720 | /* check orthogonality */ | ||
| 721 |
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.
|
120 | PetscCall(PetscOptionsGetBool(NULL,NULL,"-skiporth",&skiporth,NULL)); |
| 722 |
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.
|
120 | PetscCall(SVDGetConverged(svd,&nconv)); |
| 723 |
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.
|
120 | PetscCall(SVDGetDimensions(svd,&nsv,NULL,NULL)); |
| 724 | 120 | nconv = PetscMin(nconv,nsv); | |
| 725 |
2/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
120 | if (nconv>0 && !skiporth) { |
| 726 |
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.
|
120 | PetscCall(SVDGetTolerances(svd,&tol,NULL)); |
| 727 |
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.
|
120 | PetscCall(VecDuplicateVecs(u,nconv,&U)); |
| 728 |
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.
|
120 | PetscCall(VecDuplicateVecs(v,nconv,&V)); |
| 729 |
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.
|
720 | for (i=0;i<nconv;i++) PetscCall(SVDGetSingularTriplet(svd,i,NULL,U[i],V[i])); |
| 730 |
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.
|
120 | PetscCall(VecCheckOrthonormality(U,nconv,NULL,nconv,Omega,NULL,&lev1)); |
| 731 |
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.
|
120 | PetscCall(VecCheckOrthonormality(V,nconv,NULL,nconv,NULL,NULL,&lev2)); |
| 732 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
120 | if (lev1+lev2<20*tol) PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Level of orthogonality below the tolerance\n")); |
| 733 | ✗ | else PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Level of orthogonality: %g (U) %g (V)\n",(double)lev1,(double)lev2)); | |
| 734 |
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.
|
120 | PetscCall(VecDestroyVecs(nconv,&U)); |
| 735 |
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.
|
120 | PetscCall(VecDestroyVecs(nconv,&V)); |
| 736 | } | ||
| 737 |
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.
|
120 | PetscCall(VecDestroy(&u)); |
| 738 |
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.
|
120 | PetscCall(VecDestroy(&v)); |
| 739 | |||
| 740 | /* free work space */ | ||
| 741 |
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.
|
120 | PetscCall(SVDDestroy(&svd)); |
| 742 |
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.
|
120 | PetscCall(MatDestroy(&Omega)); |
| 743 |
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.
|
120 | PetscCall(MatDestroy(&A)); |
| 744 |
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.
|
120 | PetscCall(VecDestroy(&vomega)); |
| 745 |
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.
|
120 | PetscCall(SlepcFinalize()); |
| 746 | return 0; | ||
| 747 | } | ||
| 748 | |||
| 749 | /*TEST | ||
| 750 | |||
| 751 | testset: | ||
| 752 | args: -svd_nsv 5 -d .15 -terse | ||
| 753 | output_file: output/ex53_1.out | ||
| 754 | nsize: {{1 2}} | ||
| 755 | test: | ||
| 756 | args: -svd_type trlanczos -ds_parallel {{redundant synchronized}} | ||
| 757 | suffix: 1_trlanczos | ||
| 758 | test: | ||
| 759 | args: -svd_type cross | ||
| 760 | suffix: 1_cross | ||
| 761 | test: | ||
| 762 | args: -svd_type cyclic | ||
| 763 | requires: !single | ||
| 764 | suffix: 1_cyclic | ||
| 765 | |||
| 766 | TEST*/ | ||
| 767 |