GCC Code Coverage Report

Directory: ./
File: src/svd/tutorials/ex53.c
Date: 2026-02-22 03:58:10
Exec Total Coverage
Lines: 441 451 97.8%
Functions: 21 21 100.0%
Branches: 866 1458 59.4%
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 96 static PetscErrorCode XMatCreate(MPI_Comm comm,XMat *A,PetscInt Istart,PetscInt Iend,PetscInt M,PetscInt N)
82 {
83 96 PetscInt i;
84 96 PetscMPIInt np;
85
86
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 96 PetscFunctionBeginUser;
87
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscNew(A));
88 96 (*A)->M = M;
89 96 (*A)->N = N;
90 96 (*A)->Istart = Istart;
91 96 (*A)->Iend = Iend;
92 96 (*A)->comm = comm;
93
14/28
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 96 PetscCallMPI(MPI_Comm_size(comm,&np));
94
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscMalloc4(Iend-Istart,&(*A)->nzr,Iend-Istart,&(*A)->rcol,Iend-Istart,&(*A)->rval,np+1,&(*A)->ranges));
95
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 for (i=0; i<Iend-Istart; i++) {
96 12160 (*A)->nzr[i] = 0;
97 12160 (*A)->rcol[i] = NULL;
98 12160 (*A)->rval[i] = NULL;
99 }
100
15/30
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 6 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.
 192 PetscCallMPI(MPI_Allgather(&Istart,1,MPIU_INT,(*A)->ranges,1,MPIU_INT,comm));
101 96 (*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.
 96 PetscFunctionReturn(PETSC_SUCCESS);
103 }
104
105 96 static PetscErrorCode XMatDestroy(XMat *A)
106 {
107 96 PetscInt m,i;
108
109
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 96 PetscFunctionBeginUser;
110 96 m = (*A)->Iend-(*A)->Istart;
111
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 for (i=0; i<m; i++) {
112
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12160 if ((*A)->nzr[i]>0) {
113
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 12160 PetscCall(PetscFree2((*A)->rcol[i],(*A)->rval[i]));
114 }
115 }
116
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscFree4((*A)->nzr,(*A)->rcol,(*A)->rval,(*A)->ranges));
117
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 96 PetscCall(PetscFree(*A));
118 96 *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.
 96 PetscFunctionReturn(PETSC_SUCCESS);
120 }
121
122 10112 static PetscErrorCode XMatGetSize(XMat A,PetscInt *M,PetscInt *N)
123 {
124
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 10112 PetscFunctionBeginUser;
125 10112 *M = A->M;
126 10112 *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.
 10112 PetscFunctionReturn(PETSC_SUCCESS);
128 }
129
130 14208 static PetscErrorCode XMatGetOwnershipRange(XMat A,PetscInt *Istart,PetscInt *Iend)
131 {
132
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 14208 PetscFunctionBeginUser;
133 14208 *Istart = A->Istart;
134 14208 *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.
 14208 PetscFunctionReturn(PETSC_SUCCESS);
136 }
137
138 1280 static PetscErrorCode XMatGetOwnershipRanges(XMat A,PetscInt **ranges)
139 {
140
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 1280 PetscFunctionBeginUser;
141 1280 *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.
 1280 PetscFunctionReturn(PETSC_SUCCESS);
143 }
144
145 18944 static PetscErrorCode XMatGetRow(XMat A,PetscInt i,PetscInt *nc,PetscInt **vj,PetscScalar **va)
146 {
147 18944 PetscInt iloc=i-A->Istart;
148
149
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 18944 PetscFunctionBeginUser;
150 18944 *nc = A->nzr[iloc];
151 18944 *vj = A->rcol[iloc];
152 18944 *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.
 18944 PetscFunctionReturn(PETSC_SUCCESS);
154 }
155
156 18944 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.
 18944 PetscFunctionBeginUser;
159 18944 *nc = 0;
160 18944 *vj = NULL;
161 18944 *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.
 18944 PetscFunctionReturn(PETSC_SUCCESS);
163 }
164
165 97280 static PetscErrorCode XMatChangeRow(XMat A,PetscInt i,PetscInt nc,PetscInt vj[],PetscScalar va[])
166 {
167 97280 PetscInt iloc = i-A->Istart;
168
169
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 97280 PetscFunctionBeginUser;
170
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 97280 if (A->nzr[iloc]>0) {
171
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96704 PetscCall(PetscFree2(A->rcol[iloc],A->rval[iloc]));
172 }
173 97280 A->nzr[iloc] = nc;
174 97280 A->rcol[iloc] = vj;
175 97280 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.
 97280 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 1811200 static PetscBool getcolpos(const PetscInt rcol[],PetscInt nz,PetscInt j,PetscInt *k)
183 {
184 1811200 PetscInt ka,kb;
185
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 1811200 if (nz==0) {
186 191616 *k = 0;
187 191616 return PETSC_FALSE;
188 }
189
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 1619584 if (j<=rcol[0]) {
190 540032 *k = 0;
191 540032 return (rcol[0]==j)? PETSC_TRUE: PETSC_FALSE;
192 }
193
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 1079552 if (rcol[nz-1]<=j) {
194 365440 *k = nz-(rcol[nz-1]==j);
195 365440 return (rcol[nz-1]==j)? PETSC_TRUE: PETSC_FALSE;
196 }
197 /* Bisection: rcol[ka] < j < rcol[kb] */
198 714112 ka = 0; kb = nz-1;
199
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 2216384 while (ka+1<kb) {
200 1560384 *k = (ka+kb)/2;
201
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 1560384 if (rcol[*k]<=j) {
202
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 784256 if (rcol[*k]==j) return PETSC_TRUE;
203 else ka = *k;
204 } else kb = *k;
205 }
206 656000 *k = kb;
207 656000 return PETSC_FALSE;
208 }
209
210 /* Only local elements can be set */
211 11520 static PetscErrorCode XMatSetValue(XMat A,PetscInt i,PetscInt j,PetscScalar x,InsertMode dummy)
212 {
213 11520 PetscInt nz,iloc,k,kx,*col1,*col2;
214 11520 PetscScalar *val1,*val2;
215
216
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 11520 PetscFunctionBeginUser;
217
4/10
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 8 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 8 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 11520 PetscCheck(i>=A->Istart && i<A->Iend && j>=0 && j<A->N,A->comm,PETSC_ERR_ARG_OUTOFRANGE,"Incorrect index of row/column");
218 11520 iloc = i-A->Istart;
219 11520 nz = A->nzr[iloc];
220 11520 col1 = A->rcol[iloc];
221 11520 val1 = A->rval[iloc];
222
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 11520 if (!getcolpos(col1,nz,j,&k)) {
223 11520 A->nzr[iloc]++;
224 /* Replace row */
225
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 11520 PetscCall(PetscMalloc2(nz+1,&col2,nz+1,&val2));
226 11520 A->rcol[iloc] = col2;
227 11520 A->rval[iloc] = val2;
228
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 8 times.
 11520 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 8 times.
 11520 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 8 times.
 11520 if (nz>0) {
237 ✗ PetscCall(PetscFree2(col1,val1));
238 }
239 11520 col1 = col2;
240 11520 val1 = val2;
241 }
242 11520 col1[k] = j;
243 11520 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.
 11520 PetscFunctionReturn(PETSC_SUCCESS);
245 }
246
247 /* Convert to PETSc Mat */
248 96 static PetscErrorCode XMatConvert(XMat XA,Mat A)
249 {
250 96 PetscInt i,iloc;
251
252
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 96 PetscFunctionBeginUser;
253
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 for (i=XA->Istart; i<XA->Iend; i++) {
254 12160 iloc = i-XA->Istart;
255
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 12160 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
258
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 96 static PetscErrorCode XMatGetNumberNonzeros(XMat A,PetscInt *nz)
264 {
265 96 PetscInt i;
266
267
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 96 PetscFunctionBeginUser;
268 96 *nz = 0;
269
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 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.
 96 PetscFunctionReturn(PETSC_SUCCESS);
271 }
272
273 7104 static PetscErrorCode XMatCreateColsNzIter(XMat A,PetscInt j1,PetscInt j2,ColsNzIter *iter)
274 {
275
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 7104 PetscFunctionBeginUser;
276
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(PetscNew(iter));
277 7104 (*iter)->A = A;
278 7104 (*iter)->j1 = j1;
279 7104 (*iter)->j2 = j2;
280 7104 (*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.
 7104 PetscFunctionReturn(PETSC_SUCCESS);
282 }
283
284 7104 static PetscErrorCode ColitDestroy(ColsNzIter *iter)
285 {
286
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 7104 PetscFunctionBeginUser;
287
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 7104 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 92480 static PetscErrorCode ColitNextRow(ColsNzIter iter,PetscInt *i,PetscScalar **paij1,PetscScalar **paij2)
292 {
293 92480 PetscInt iloc = iter->iloc,*rcol,nz,k,k1,k2,kx,jx,*rcol2;
294 92480 PetscScalar *rval,*rval2;
295 92480 PetscBool found1, found2;
296 92480 XMat A=iter->A;
297
298
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 92480 PetscFunctionBeginUser;
299 92480 *i = -1;
300 92480 *paij1 = *paij2 = NULL;
301
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 92480 if (iloc>=0) {
302 91904 rcol = A->rcol[iloc];
303 91904 nz = A->nzr[iloc];
304 1707776 while (PETSC_TRUE) {
305 899840 found1 = getcolpos(rcol,nz,iter->j1,&k1);
306
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 899840 found2 = getcolpos(PetscSafePointerPlusOffset(rcol,k1+found1),nz-k1-found1,iter->j2,&k2);
307
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 899840 if (found1 || found2) break;
308 814464 iloc++;
309
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 814464 if (iloc>=A->Iend-A->Istart) {
310 iloc = -1;
311 break;
312 }
313 807936 rcol = A->rcol[iloc];
314 807936 nz = A->nzr[iloc];
315 }
316
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 91904 if (iloc>=0) {
317 85376 k2 += k1 + 1;
318 85376 rval = A->rval[iloc];
319 85376 *i = A->Istart+iloc;
320
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 85376 if (!found1 || !found2) {
321 /* Copy row to a new one */
322 78464 rcol2 = rcol;
323 78464 rval2 = rval;
324
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 78464 PetscCall(PetscMalloc2(nz+1,&rcol,nz+1,&rval));
325
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 78464 if (found1) {
326 43776 kx = k2;
327 43776 jx = iter->j2;
328 } else {
329 34688 kx = k1;
330 34688 jx = iter->j1;
331 }
332
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 511808 for (k=0; k<kx; k++) {
333 433344 rcol[k] = rcol2[k];
334 433344 rval[k] = rval2[k];
335 }
336 78464 rcol[kx] = jx;
337 78464 rval[kx] = 0.0;
338
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 582272 for (k=kx+1; k<nz+1; k++) {
339 503808 rcol[k] = rcol2[k-1];
340 503808 rval[k] = rval2[k-1];
341 }
342
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 78464 PetscCall(XMatChangeRow(A,*i,nz+1,rcol,rval));
343 }
344 85376 *paij1 = &rval[k1];
345 85376 *paij2 = &rval[k2];
346 85376 iloc++;
347
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 85376 if (iloc>=A->Iend-A->Istart) iloc = -1;
348 }
349 91904 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 1280 static PetscErrorCode SendrecvRow(XMat A,PetscInt nc1,const PetscInt vj1[],const PetscScalar va1[],PetscInt i2,PetscInt *nc2,PetscInt vj2[],PetscScalar va2[])
355 {
356 1280 PetscInt *ranges,N=A->N;
357 1280 PetscMPIInt rank,naux,len1,len2;
358 1280 MPI_Status st;
359
360
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 1280 PetscFunctionBeginUser;
361
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 1280 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 8 times.
✓ Branch 1 taken 8 times.
 1920 while (ranges[rank+1]<=i2) rank++;
366 /* Send row i1, receive row i2 */
367
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 1280 PetscCall(PetscMPIIntCast(nc1,&len1));
368
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 1280 PetscCall(PetscMPIIntCast(N,&len2));
369
14/28
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 1280 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 6 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 1280 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 6 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 1280 PetscCallMPI(MPI_Get_count(&st,MPIU_INT,&naux));
372 1280 *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.
 1280 PetscFunctionReturn(PETSC_SUCCESS);
374 }
375
376 10112 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 10112 PetscInt k1,k2,k,n1=nc1+nc2,*vjjaux=iwork;
379 10112 PetscScalar *vaa1aux=swork,*vaa2aux=swork+n1;
380
381
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 10112 PetscFunctionBeginUser;
382 10112 *nc=0;
383
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 84480 for (k1=k2=k=0; k1<nc1 && k2<nc2; ) {
384
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 74368 if (vj1[k1]<vj2[k2]) {
385 /* Take next element from first row */
386 34592 vaa1aux[k] = va1[k1];
387 34592 vaa2aux[k] = 0.0;
388 34592 vjjaux[k++] = vj1[k1++];
389
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 39776 } else if (vj1[k1]>vj2[k2]) {
390 /* Take next element from second row */
391 34496 vaa1aux[k] = 0.0;
392 34496 vaa2aux[k] = va2[k2];
393 34496 vjjaux[k++] = vj2[k2++];
394 } else {
395 /* Take next element from both rows */
396 5280 vaa1aux[k] = va1[k1];
397 5280 vaa2aux[k] = va2[k2++];
398 5280 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 8 times.
✓ Branch 1 taken 8 times.
 20448 while (k1<nc1) {
404 /* Take next element from first row */
405 10336 vaa1aux[k] = va1[k1];
406 10336 vaa2aux[k] = 0.0;
407 10336 vjjaux[k++] = vj1[k1++];
408 }
409
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 20640 while (k2<nc2) {
410 /* Take next element from second row */
411 10528 vaa1aux[k] = 0.0;
412 10528 vaa2aux[k] = va2[k2];
413 10528 vjjaux[k++] = vj2[k2++];
414 }
415 10112 *nc=k;
416
417 /* Copy rows (each row must be allocated separately)*/
418
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PetscMalloc2(k,vjj1,k,vaa1));
419
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PetscMalloc2(k,vjj2,k,vaa2));
420
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PetscArraycpy(*vjj1,vjjaux,k));
421
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PetscArraycpy(*vjj2,vjjaux,k));
422
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PetscArraycpy(*vaa1,vaa1aux,k));
423
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 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 1280 static inline void MyAxpby(PetscInt n,PetscScalar a,PetscScalar x[],PetscScalar b,const PetscScalar y[])
428 {
429 1280 PetscInt i;
430
431
4/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 8 times.
 10240 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 14208 static PetscErrorCode RotateRows(XMat A,PetscInt i1,PetscInt i2,PetscReal c,PetscReal s,PetscInt *nzloc,PetscInt *iwork,PetscScalar *swork)
436 {
437 14208 PetscInt M,N,Istart,Iend,nc,nc1,nc2,*vj1=NULL,*vj2=NULL,*vjj1=NULL,*vjj2=NULL,*iworkx=iwork;
438 14208 PetscBLASInt nc_,one=1;
439 14208 PetscScalar *va1=NULL,*va2=NULL,*vaa1=NULL,*vaa2=NULL,*sworkx=swork;
440 14208 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.
 14208 PetscFunctionBeginUser;
446
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 14208 PetscCall(XMatGetOwnershipRange(A,&Istart,&Iend));
447
4/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 8 times.
 14208 i1mine = (Istart<=i1 && i1<Iend)? PETSC_TRUE: PETSC_FALSE;
448
4/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 8 times.
 14208 i2mine = (Istart<=i2 && i2<Iend)? PETSC_TRUE: PETSC_FALSE;
449
450
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 14208 if (i1mine || i2mine) {
451
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(XMatGetSize(A,&M,&N));
452
453 /* Get local row(s) (at least 1, possibly 2) */
454 10112 time_now(t0);
455
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 10112 if (i1mine) PetscCall(XMatGetRow(A,i1,&nc1,&vj1,&va1));
456
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 10112 if (i2mine) PetscCall(XMatGetRow(A,i2,&nc2,&vj2,&va2));
457 10112 time_diff(tt_getr,t0,t,t0);
458 /* Get remote row (at most 1, possibly none) */
459
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
 10112 if (!i2mine) {
460 640 vj2 = iworkx;
461 640 iworkx += N;
462 640 va2 = sworkx;
463 640 sworkx += N;
464
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 640 PetscCall(SendrecvRow(A,nc1,vj1,va1,i2,&nc2,vj2,va2));
465
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 9472 } else if (!i1mine) {
466 640 vj1 = iworkx;
467 640 iworkx += N;
468 640 va1 = sworkx;
469 640 sworkx += N;
470
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 640 PetscCall(SendrecvRow(A,nc2,vj2,va2,i1,&nc1,vj1,va1));
471 }
472
473
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10112 PetscCall(PadRows(nc1,vj1,va1,nc2,vj2,va2,&nc,&vjj1,&vaa1,&vjj2,&vaa2,iworkx,sworkx));
474
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 10112 if (i1mine) {
475 9472 *nzloc += nc-nc1;
476
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 9472 PetscCall(XMatRestoreRow(A,i1,&nc1,&vj1,&va1));
477 }
478
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 10112 if (i2mine) {
479 9472 *nzloc += nc-nc2;
480
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 9472 PetscCall(XMatRestoreRow(A,i2,&nc2,&vj2,&va2));
481 }
482
483 /* Compute rotation and update matrix */
484
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 10112 if (nc) {
485
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 10048 PetscCall(PetscBLASIntCast(nc,&nc_));
486
12/22
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 10048 if (i1mine && i2mine) PetscCallBLAS("BLASrot",BLASMIXEDrot_(&nc_,vaa1,&one,vaa2,&one,&c,&s));
487
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 1280 else if (i1mine) {
488 640 MyAxpby(nc,c,vaa1,s,vaa2);
489
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 640 PetscCall(PetscFree2(vjj2,vaa2));
490 } else {
491 640 MyAxpby(nc,c,vaa2,-s,vaa1);
492
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 640 PetscCall(PetscFree2(vjj1,vaa1));
493 }
494 10048 time_now(t0);
495
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 10048 if (i1mine) PetscCall(XMatChangeRow(A,i1,nc,vjj1,vaa1));
496
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 10048 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 7104 static PetscErrorCode RotateCols(XMat A,PetscInt j1,PetscInt j2,PetscReal c,PetscReal s,PetscInt *nzloc)
505 {
506 7104 PetscInt i;
507 7104 PetscScalar aux1,aux2,*paij1,*paij2;
508 7104 ColsNzIter iter;
509
510
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 7104 PetscFunctionBeginUser;
511
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(XMatCreateColsNzIter(A,j1,j2,&iter));
512 177856 while (PETSC_TRUE) {
513
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 92480 PetscCall(ColitNextRow(iter,&i,&paij1,&paij2));
514
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 92480 if (i<0) break;
515
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 85376 if (*paij1**paij2==0.0) (*nzloc)++;
516 85376 aux1 = *paij1*c+*paij2*s;
517 85376 aux2 = -*paij1*s+*paij2*c;
518 85376 *paij1 = aux1;
519 85376 *paij2 = aux2;
520 }
521
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 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 21312 static PetscErrorCode GetIndexPair(PetscRandom rand,PetscInt n,PetscInt *i1,PetscInt *i2)
527 {
528 21312 PetscInt aux;
529 21312 PetscReal x;
530
531
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 21312 PetscFunctionBeginUser;
532
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 21312 PetscCall(PetscRandomGetValueReal(rand,&x));
533 21312 *i1 = (PetscInt)(x*n);
534
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 21312 PetscCall(PetscRandomGetValueReal(rand,&x));
535 21312 *i2 = (PetscInt)(x*(n-1));
536
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 21312 if (*i2>=*i1) (*i2)++;
537
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 21312 if (*i1>*i2) {
538 10656 aux = *i1;
539 10656 *i1 = *i2;
540 10656 *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 96 int main(int argc,char **argv)
546 {
547 96 Mat A,Omega; /* operator matrix, signature matrix */
548 96 XMat XA;
549 96 SVD svd; /* singular value problem solver context */
550 96 Vec u,v,vomega,*U,*V;
551 96 MatType Atype;
552 96 PetscReal tol,lev1=0.0,lev2=0.0,k=1e3,dens=0.1,log0,loginc,sq,x,angle,c,s;
553 96 PetscScalar *swork;
554 96 PetscInt P=110,Q=80,N=100,i,j,Istart,Iend,nconv,nsv,i1,i2,nroth,nrotv,*iwork,nnzwanted,nz,nzloc;
555 96 PetscBool flg,flgp,flgq,flgn,terse,skiporth=PETSC_FALSE,progress=PETSC_FALSE;
556 96 PetscRandom rand;
557 96 MatInfo Ainfo;
558 #ifdef TIMING
559 PetscReal t,t0;
560 #endif
561
562
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 96 PetscFunctionBeginUser;
563
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
564
565
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetInt(NULL,NULL,"-p",&P,&flgp));
566
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetInt(NULL,NULL,"-q",&Q,&flgq));
567
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&N,&flgn));
568
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetReal(NULL,NULL,"-k",&k,&flg));
569
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetReal(NULL,NULL,"-d",&dens,&flg));
570
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetBool(NULL,NULL,"-prog",&progress,&flg));
571
2/4
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 8 times.
 96 if (!flgp && flgn) P=N;
572
2/4
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 8 times.
 96 if (!flgq && flgn) Q=N;
573
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 96 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 8 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 96 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 8 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nHyperbolic singular value problem.\n\n"));
582
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
586
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,P+Q,N));
587
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatSetFromOptions(A));
588
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
592
593
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(XMatCreate(MPI_COMM_WORLD,&XA,Istart,Iend,P+Q,N));
594 96 log0 = PetscLogReal(1/k);
595 96 loginc = -log0/(N-1);
596 96 sq = PetscSqrtReal(5.0/4.0);
597
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 for (i=Istart;i<Iend;i++) {
598
4/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 8 times.
 12160 if (i<P && i<N) {
599
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 6400 if (i<Q) x = sq*PetscExpReal(log0+i*loginc);
600 1280 else x = PetscExpReal(log0+i*loginc);
601
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 6400 PetscCall(XMatSetValue(XA,i,i,x,INSERT_VALUES));
602
3/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 5760 } else if (i>=P && i-P<N) {
603 5120 j = i-P;
604
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 12160 PetscCall(XMatSetValue(XA,i,j,PetscExpReal(log0+j*loginc),INSERT_VALUES));
605 }
606 }
607
608 /* Apply plane rotations */
609 96 nnzwanted = dens*(P+Q)*N;
610
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(XMatGetNumberNonzeros(XA,&nzloc));
611
15/30
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 6 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.
 192 PetscCallMPI(MPI_Allreduce(&nzloc,&nz,1,MPIU_INT,MPI_SUM,MPI_COMM_WORLD));
612
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscRandomCreate(PETSC_COMM_SELF,&rand));
613
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscRandomSetFromOptions(rand));
614 96 nroth = nrotv = 0;
615
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscMalloc2(2*N,&iwork,3*N,&swork));
616 time_now(t0);
617 7200 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(PetscRandomGetValueReal(rand,&angle));
620 7104 c=PetscCosReal(angle);
621 7104 s=PetscSinReal(angle);
622
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(GetIndexPair(rand,P,&i1,&i2));
623
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(RotateRows(XA,i1,i2,c,s,&nzloc,iwork,swork));
624 7104 time_diff(tt_rotr,t0,t,t0);
625 7104 nroth++;
626 /* Plane rotation on 2 of the bottom Q rows*/
627
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 7104 if (Q>1) {
628
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(PetscRandomGetValueReal(rand,&angle));
629 7104 c=PetscCosReal(angle);
630 7104 s=PetscSinReal(angle);
631
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(GetIndexPair(rand,Q,&i1,&i2));
632
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(RotateRows(XA,P+i1,P+i2,c,s,&nzloc,iwork,swork));
633 7104 time_diff(tt_rotr,t0,t,t0);
634 7104 nroth++;
635 }
636 /* Plane rotation on 2 columns */
637
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(PetscRandomGetValueReal(rand,&angle));
638 7104 c=PetscCosReal(angle);
639 7104 s=PetscSinReal(angle);
640
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(GetIndexPair(rand,N,&i1,&i2));
641
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 7104 PetscCall(RotateCols(XA,i1,i2,c,s,&nzloc));
642 7104 time_diff(tt_rotc,t0,t,t0);
643 7104 nrotv++;
644 /* Update total number of non-zeros */
645
15/30
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 6 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.
 14208 PetscCallMPI(MPI_Allreduce(&nzloc,&nz,1,MPIU_INT,MPI_SUM,MPI_COMM_WORLD));
646
5/10
✗ Branch 0 not taken.
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 6 times.
✓ Branch 5 taken 6 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✓ Branch 9 taken 2 times.
 14304 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 8 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.
 96 if (progress) PetscCall(PetscFPrintf(MPI_COMM_WORLD,stderr,"\r"));
649
650
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscFree2(iwork,swork));
651
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscRandomDestroy(&rand));
652
653
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(XMatConvert(XA,A));
654
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(XMatDestroy(&XA));
655
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatGetInfo(A,MAT_GLOBAL_SUM,&Ainfo));
656
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatCreateVecs(A,NULL,&vomega));
665
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecSet(vomega,1.0));
666
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 5216 for (i=PetscMax(P,Istart);i<Iend;i++) {
667
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 5120 PetscCall(VecSetValue(vomega,i,0.5,INSERT_VALUES));
668 }
669
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecAssemblyBegin(vomega));
670
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecAssemblyEnd(vomega));
671
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatDiagonalScale(A,vomega,NULL));
672
673 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
674 Create the signature
675 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
676
677
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecSet(vomega,1.0));
678
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecGetOwnershipRange(vomega,&Istart,&Iend));
679
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 12256 for (i=Istart;i<Iend;i++) {
680
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 12160 if (i>=P) PetscCall(VecSetValue(vomega,i,-1.0,INSERT_VALUES));
681 }
682
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecAssemblyBegin(vomega));
683
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecAssemblyEnd(vomega));
684
685
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatGetType(A,&Atype));
686
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatCreate(PETSC_COMM_WORLD,&Omega));
687
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatSetSizes(Omega,PETSC_DECIDE,PETSC_DECIDE,P+Q,P+Q));
688
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatSetType(Omega,Atype));
689
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd));
696
697
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDSetOperators(svd,A,NULL));
698
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDSetSignature(svd,vomega));
699
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDSetProblemType(svd,SVD_HYPERBOLIC));
700
701
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDSetFromOptions(svd));
702
703 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
704 Solve the problem, display solution
705 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
706
707
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatCreateVecs(A,&v,&u));
708
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
712
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(PetscOptionsGetBool(NULL,NULL,"-skiporth",&skiporth,NULL));
722
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDGetConverged(svd,&nconv));
723
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDGetDimensions(svd,&nsv,NULL,NULL));
724 96 nconv = PetscMin(nconv,nsv);
725
2/4
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 96 if (nconv>0 && !skiporth) {
726
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDGetTolerances(svd,&tol,NULL));
727
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDuplicateVecs(u,nconv,&U));
728
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDuplicateVecs(v,nconv,&V));
729
7/8
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 6 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 2 times.
 576 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecCheckOrthonormality(U,nconv,NULL,nconv,Omega,NULL,&lev1));
731
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecCheckOrthonormality(V,nconv,NULL,nconv,NULL,NULL,&lev2));
732
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 96 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDestroyVecs(nconv,&U));
735
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDestroyVecs(nconv,&V));
736 }
737
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDestroy(&u));
738
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDestroy(&v));
739
740 /* free work space */
741
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(SVDDestroy(&svd));
742
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatDestroy(&Omega));
743
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(MatDestroy(&A));
744
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 96 PetscCall(VecDestroy(&vomega));
745
2/6
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 96 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