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 | #include <slepc/private/dsimpl.h> | ||
12 | #include <slepcblaslapack.h> | ||
13 | |||
14 | struct HRtr | ||
15 | { | ||
16 | PetscScalar *data; | ||
17 | PetscInt m; | ||
18 | PetscInt idx[2]; | ||
19 | PetscInt n[2]; | ||
20 | PetscScalar tau[2]; | ||
21 | PetscReal alpha; | ||
22 | PetscReal cs; | ||
23 | PetscReal sn; | ||
24 | PetscInt type; | ||
25 | }; | ||
26 | |||
27 | /* | ||
28 | Generates a hyperbolic rotation | ||
29 | if x1*x1 - x2*x2 != 0 | ||
30 | r = sqrt(|x1*x1 - x2*x2|) | ||
31 | c = x1/r s = x2/r | ||
32 | |||
33 | | c -s||x1| |d*r| | ||
34 | |-s c||x2| = | 0 | | ||
35 | where d = 1 for type==1 and -1 for type==2 | ||
36 | Returns the condition number of the reduction | ||
37 | */ | ||
38 | 219908 | static PetscErrorCode HRGen(PetscReal x1,PetscReal x2,PetscInt *type,PetscReal *c,PetscReal *s,PetscReal *r,PetscReal *cond) | |
39 | { | ||
40 | 219908 | PetscReal t,n2,xa,xb; | |
41 | 219908 | PetscInt type_; | |
42 | |||
43 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
219908 | PetscFunctionBegin; |
44 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
219908 | if (x2==0.0) { |
45 | ✗ | *r = PetscAbsReal(x1); *c = (x1>=0.0)?1.0:-1.0; *s = 0.0; | |
46 | ✗ | if (type) *type = 1; | |
47 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
48 | } | ||
49 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
219908 | if (PetscAbsReal(x1) == PetscAbsReal(x2)) { |
50 | /* hyperbolic rotation doesn't exist */ | ||
51 | ✗ | *c = *s = *r = 0.0; | |
52 | ✗ | if (type) *type = 0; | |
53 | ✗ | *cond = PETSC_MAX_REAL; | |
54 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
55 | } | ||
56 | |||
57 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
219908 | if (PetscAbsReal(x1)>PetscAbsReal(x2)) { |
58 | xa = x1; xb = x2; type_ = 1; | ||
59 | } else { | ||
60 | 39317 | xa = x2; xb = x1; type_ = 2; | |
61 | } | ||
62 | 219908 | t = xb/xa; | |
63 | 219908 | n2 = PetscAbsReal(1 - t*t); | |
64 | 219908 | *r = PetscSqrtReal(n2)*PetscAbsReal(xa); | |
65 | 219908 | *c = x1/(*r); | |
66 | 219908 | *s = x2/(*r); | |
67 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
219908 | if (type_ == 2) *r *= -1; |
68 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
219908 | if (type) *type = type_; |
69 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
219908 | if (cond) *cond = (PetscAbsReal(*c) + PetscAbsReal(*s))/PetscAbsReal(PetscAbsReal(*c) - PetscAbsReal(*s)); |
70 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
45433 | PetscFunctionReturn(PETSC_SUCCESS); |
71 | } | ||
72 | |||
73 | /* | ||
74 | |c s| | ||
75 | Applies an hyperbolic rotator |s c| | ||
76 | |c s| | ||
77 | [x1 x2]|s c| | ||
78 | */ | ||
79 | 532453 | static PetscErrorCode HRApply(PetscInt n,PetscScalar *x1,PetscInt inc1,PetscScalar *x2,PetscInt inc2,PetscReal c,PetscReal s) | |
80 | { | ||
81 | 532453 | PetscInt i; | |
82 | 532453 | PetscReal t; | |
83 | 532453 | PetscScalar tmp; | |
84 | |||
85 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
532453 | PetscFunctionBegin; |
86 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
532453 | if (PetscAbsReal(c)>PetscAbsReal(s)) { /* Type I */ |
87 | 455680 | t = s/c; | |
88 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
14047996 | for (i=0;i<n;i++) { |
89 | 13592316 | x1[i*inc1] = c*x1[i*inc1] + s*x2[i*inc2]; | |
90 | 13592316 | x2[i*inc2] = t*x1[i*inc1] + x2[i*inc2]/c; | |
91 | } | ||
92 | } else { /* Type II */ | ||
93 | 76773 | t = c/s; | |
94 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1695381 | for (i=0;i<n;i++) { |
95 | 1618608 | tmp = x1[i*inc1]; | |
96 | 1618608 | x1[i*inc1] = c*x1[i*inc1] + s*x2[i*inc2]; | |
97 | 1618608 | x2[i*inc2] = t*x1[i*inc1] + tmp/s; | |
98 | } | ||
99 | } | ||
100 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
532453 | PetscFunctionReturn(PETSC_SUCCESS); |
101 | } | ||
102 | |||
103 | /* | ||
104 | Reduction to tridiagonal-diagonal form (see F. Tisseur, SIMAX 26(1), 2004). | ||
105 | |||
106 | Input: | ||
107 | A symmetric (only lower triangular part is referred) | ||
108 | s vector +1 and -1 (signature matrix) | ||
109 | Output: | ||
110 | d,e | ||
111 | s | ||
112 | Q s-orthogonal matrix with Q^T*A*Q = T (symmetric tridiagonal matrix) | ||
113 | */ | ||
114 | 10001 | static PetscErrorCode TridiagDiag_HHR(PetscInt n,PetscScalar *A,PetscInt lda,PetscReal *s,PetscScalar* Q,PetscInt ldq,PetscBool flip,PetscReal *d,PetscReal *e,PetscInt *perm_,PetscScalar *work,PetscReal *rwork,PetscBLASInt *iwork) | |
115 | { | ||
116 | 10001 | PetscInt i,j,k,*ii,*jj,i0=0,ik=0,type; | |
117 | 10001 | PetscInt nwu=0; | |
118 | 10001 | PetscReal *ss,cond=1.0,cs,sn,r; | |
119 | 10001 | PetscScalar tau,t,*AA; | |
120 | 10001 | PetscBLASInt n0,n1,ni,inc=1,m,n_,lda_,ldq_,*perm,tmp; | |
121 | 10001 | PetscBool breakdown = PETSC_TRUE; | |
122 | |||
123 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10001 | PetscFunctionBegin; |
124 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
10001 | if (n<3) { |
125 | ✗ | if (n==1) Q[0]=1; | |
126 | ✗ | if (n==2) { | |
127 | ✗ | Q[0] = Q[1+ldq] = 1; | |
128 | ✗ | Q[1] = Q[ldq] = 0; | |
129 | } | ||
130 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
131 | } | ||
132 |
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.
|
10001 | PetscCall(PetscBLASIntCast(lda,&lda_)); |
133 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10001 | PetscCall(PetscBLASIntCast(n,&n_)); |
134 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10001 | PetscCall(PetscBLASIntCast(ldq,&ldq_)); |
135 | 10001 | ss = rwork; | |
136 | 10001 | perm = iwork; | |
137 | 10001 | AA = work; | |
138 |
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.
|
229968 | for (i=0;i<n;i++) PetscCall(PetscArraycpy(AA+i*n,A+i*lda,n)); |
139 | 10001 | nwu += n*n; | |
140 | 10001 | k=0; | |
141 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
20002 | while (breakdown && k<n) { |
142 | 10001 | breakdown = PETSC_FALSE; | |
143 | /* Classify (and flip) A and s according to sign */ | ||
144 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10001 | if (flip) { |
145 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
229528 | for (i=0;i<n;i++) { |
146 |
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.
|
219567 | PetscCall(PetscBLASIntCast(n-1-perm_[i],&perm[i])); |
147 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
219567 | if (perm[i]==0) i0 = i; |
148 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
219567 | if (perm[i]==k) ik = i; |
149 | } | ||
150 | } else { | ||
151 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
440 | for (i=0;i<n;i++) { |
152 |
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.
|
400 | PetscCall(PetscBLASIntCast(perm_[i],&perm[i])); |
153 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
400 | if (perm[i]==0) i0 = i; |
154 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
400 | if (perm[i]==k) ik = i; |
155 | } | ||
156 | } | ||
157 | 10001 | perm[ik] = 0; | |
158 |
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.
|
10001 | PetscCall(PetscBLASIntCast(k,&perm[i0])); |
159 | 10001 | i=1; | |
160 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
198755 | while (i<n-1 && s[perm[i-1]]==s[perm[0]]) { |
161 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
188754 | if (s[perm[i]]!=s[perm[0]]) { |
162 | 162033 | j=i+1; | |
163 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
177855 | while (j<n-1 && s[perm[j]]!=s[perm[0]])j++; |
164 | 162033 | tmp = perm[i]; perm[i] = perm[j]; perm[j] = tmp; | |
165 | } | ||
166 | 188754 | i++; | |
167 | } | ||
168 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
229968 | for (i=0;i<n;i++) { |
169 | 219967 | ss[i] = s[perm[i]]; | |
170 | } | ||
171 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10001 | if (flip) { |
172 | ii = &j; | ||
173 | jj = &i; | ||
174 | } else { | ||
175 | 40 | ii = &i; | |
176 | 40 | jj = &j; | |
177 | } | ||
178 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
229968 | for (i=0;i<n;i++) |
179 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
5840550 | for (j=0;j<n;j++) |
180 | 5620583 | A[i+j*lda] = AA[perm[*ii]+perm[*jj]*n]; | |
181 | /* Initialize Q */ | ||
182 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
229968 | for (i=0;i<n;i++) { |
183 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
219967 | PetscCall(PetscArrayzero(Q+i*ldq,n)); |
184 | 219967 | Q[perm[i]+i*ldq] = 1.0; | |
185 | } | ||
186 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
199038 | for (ni=1;ni<n && ss[ni]==ss[0]; ni++); |
187 | 10001 | n0 = ni-1; | |
188 | 10001 | n1 = n_-ni; | |
189 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
209966 | for (j=0;j<n-2;j++) { |
190 |
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.
|
199965 | PetscCall(PetscBLASIntCast(n-j-1,&m)); |
191 | /* Forming and applying reflectors */ | ||
192 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
199965 | if (n0 > 1) { |
193 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
187336 | PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n0,A+ni-n0+j*lda,A+ni-n0+j*lda+1,&inc,&tau)); |
194 | /* Apply reflector */ | ||
195 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
187336 | if (PetscAbsScalar(tau) != 0.0) { |
196 | 187336 | t=*(A+ni-n0+j*lda); *(A+ni-n0+j*lda)=1.0; | |
197 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
187336 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&m,&n0,A+ni-n0+j*lda,&inc,&tau,A+j+1+(j+1)*lda,&lda_,work+nwu)); |
198 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
187336 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n0,&m,A+ni-n0+j*lda,&inc,&tau,A+j+1+(j+1)*lda,&lda_,work+nwu)); |
199 | /* Update Q */ | ||
200 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
187336 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n0,A+ni-n0+j*lda,&inc,&tau,Q+(j+1)*ldq,&ldq_,work+nwu)); |
201 | 187336 | *(A+ni-n0+j*lda) = t; | |
202 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2490743 | for (i=1;i<n0;i++) { |
203 | 2303407 | *(A+ni-n0+j*lda+i) = 0.0; *(A+j+(ni-n0+i)*lda) = 0.0; | |
204 | } | ||
205 | 187336 | *(A+j+(ni-n0)*lda) = *(A+ni-n0+j*lda); | |
206 | } | ||
207 | } | ||
208 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
199965 | if (n1 > 1) { |
209 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
14014 | PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n1,A+n-n1+j*lda,A+n-n1+j*lda+1,&inc,&tau)); |
210 | /* Apply reflector */ | ||
211 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
14014 | if (PetscAbsScalar(tau) != 0.0) { |
212 | 14014 | t=*(A+n-n1+j*lda); *(A+n-n1+j*lda)=1.0; | |
213 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
14014 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&m,&n1,A+n-n1+j*lda,&inc,&tau,A+j+1+(n-n1)*lda,&lda_,work+nwu)); |
214 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
14014 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n1,&m,A+n-n1+j*lda,&inc,&tau,A+n-n1+(j+1)*lda,&lda_,work+nwu)); |
215 | /* Update Q */ | ||
216 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
14014 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n1,A+n-n1+j*lda,&inc,&tau,Q+(n-n1)*ldq,&ldq_,work+nwu)); |
217 | 14014 | *(A+n-n1+j*lda) = t; | |
218 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
96072 | for (i=1;i<n1;i++) { |
219 | 82058 | *(A+n-n1+i+j*lda) = 0.0; *(A+j+(n-n1+i)*lda) = 0.0; | |
220 | } | ||
221 | 14014 | *(A+j+(n-n1)*lda) = *(A+n-n1+j*lda); | |
222 | } | ||
223 | } | ||
224 | /* Hyperbolic rotation */ | ||
225 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
199965 | if (n0 > 0 && n1 > 0) { |
226 |
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.
|
104877 | PetscCall(HRGen(PetscRealPart(A[ni-n0+j*lda]),PetscRealPart(A[n-n1+j*lda]),&type,&cs,&sn,&r,&cond)); |
227 | /* Check condition number */ | ||
228 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
104877 | if (cond > 1.0/(10*PETSC_SQRT_MACHINE_EPSILON)) { |
229 | ✗ | breakdown = PETSC_TRUE; | |
230 | ✗ | k++; | |
231 | ✗ | PetscCheck(k<n && !flip,PETSC_COMM_SELF,PETSC_ERR_SUP,"Breakdown in construction of hyperbolic transformation"); | |
232 | break; | ||
233 | } | ||
234 | 104877 | A[ni-n0+j*lda] = r; A[n-n1+j*lda] = 0.0; | |
235 | 104877 | A[j+(ni-n0)*lda] = r; A[j+(n-n1)*lda] = 0.0; | |
236 | /* Apply to A */ | ||
237 |
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.
|
104877 | PetscCall(HRApply(m,A+j+1+(ni-n0)*lda,1,A+j+1+(n-n1)*lda,1,cs,-sn)); |
238 |
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.
|
104877 | PetscCall(HRApply(m,A+ni-n0+(j+1)*lda,lda,A+n-n1+(j+1)*lda,lda,cs,-sn)); |
239 | |||
240 | /* Update Q */ | ||
241 |
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.
|
104877 | PetscCall(HRApply(n,Q+(ni-n0)*ldq,1,Q+(n-n1)*ldq,1,cs,-sn)); |
242 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
104877 | if (type==2) { |
243 | 8472 | ss[ni-n0] = -ss[ni-n0]; ss[n-n1] = -ss[n-n1]; | |
244 | 8472 | n0++;ni++;n1--; | |
245 | } | ||
246 | } | ||
247 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
199965 | if (n0>0) n0--; |
248 | 11164 | else n1--; | |
249 | } | ||
250 | } | ||
251 | |||
252 | /* flip matrices */ | ||
253 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10001 | if (flip) { |
254 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
219567 | for (i=0;i<n-1;i++) { |
255 | 209606 | d[i] = PetscRealPart(A[n-i-1+(n-i-1)*lda]); | |
256 | 209606 | e[i] = PetscRealPart(A[n-i-1+(n-i-2)*lda]); | |
257 | 209606 | s[i] = ss[n-i-1]; | |
258 | } | ||
259 | 9961 | s[n-1] = ss[0]; | |
260 | 9961 | d[n-1] = PetscRealPart(A[0]); | |
261 |
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.
|
229528 | for (i=0;i<n;i++) PetscCall(PetscArraycpy(work+i*n,Q+i*ldq,n)); |
262 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
229528 | for (i=0;i<n;i++) |
263 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
5836150 | for (j=0;j<n;j++) |
264 | 5616583 | Q[i+j*ldq] = work[i+(n-j-1)*n]; | |
265 | } else { | ||
266 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
400 | for (i=0;i<n-1;i++) { |
267 | 360 | d[i] = PetscRealPart(A[i+i*lda]); | |
268 | 360 | e[i] = PetscRealPart(A[i+1+i*lda]); | |
269 | 360 | s[i] = ss[i]; | |
270 | } | ||
271 | 40 | s[n-1] = ss[n-1]; | |
272 | 40 | d[n-1] = PetscRealPart(A[n-1 + (n-1)*lda]); | |
273 | } | ||
274 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
1983 | PetscFunctionReturn(PETSC_SUCCESS); |
275 | } | ||
276 | |||
277 | 335846 | static PetscErrorCode MadeHRtr(PetscInt sz,PetscInt n,PetscInt idx0,PetscInt n0,PetscInt idx1,PetscInt n1,struct HRtr *tr1,struct HRtr *tr2,PetscReal *ncond,PetscScalar *work) | |
278 | { | ||
279 | 335846 | PetscScalar *x,*y; | |
280 | 335846 | PetscReal ncond2=1.0; | |
281 | 335846 | PetscBLASInt n0_,n1_,inc=1; | |
282 | |||
283 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
335846 | PetscFunctionBegin; |
284 | /* Hyperbolic transformation to make zeros in x */ | ||
285 | 335846 | x = tr1->data; | |
286 | 335846 | tr1->n[0] = n0; | |
287 | 335846 | tr1->n[1] = n1; | |
288 | 335846 | tr1->idx[0] = idx0; | |
289 | 335846 | tr1->idx[1] = idx1; | |
290 |
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.
|
335846 | PetscCall(PetscBLASIntCast(tr1->n[0],&n0_)); |
291 |
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.
|
335846 | PetscCall(PetscBLASIntCast(tr1->n[1],&n1_)); |
292 |
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.
|
335846 | if (tr1->n[0] > 1) PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n0_,x+tr1->idx[0],x+tr1->idx[0]+1,&inc,tr1->tau)); |
293 |
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.
|
335846 | if (tr1->n[1]> 1) PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n1_,x+tr1->idx[1],x+tr1->idx[1]+1,&inc,tr1->tau+1)); |
294 |
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.
|
335846 | if (tr1->idx[0]<tr1->idx[1]) PetscCall(HRGen(PetscRealPart(x[tr1->idx[0]]),PetscRealPart(x[tr1->idx[1]]),&tr1->type,&tr1->cs,&tr1->sn,&tr1->alpha,ncond)); |
295 | else { | ||
296 | 235037 | tr1->alpha = PetscRealPart(x[tr1->idx[0]]); | |
297 | 235037 | *ncond = 1.0; | |
298 | } | ||
299 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
335846 | if (sz==2) { |
300 | 21826 | y = tr2->data; | |
301 | /* Apply first transformation to second column */ | ||
302 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
21826 | if (tr1->n[0] > 1 && PetscAbsScalar(tr1->tau[0])!=0.0) { |
303 | 20728 | x[tr1->idx[0]] = 1.0; | |
304 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
20728 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n0_,&inc,x+tr1->idx[0],&inc,tr1->tau,y+tr1->idx[0],&n0_,work)); |
305 | } | ||
306 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
21826 | if (tr1->n[1] > 1 && PetscAbsScalar(tr1->tau[1])!=0.0) { |
307 | 7716 | x[tr1->idx[1]] = 1.0; | |
308 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
7716 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n1_,&inc,x+tr1->idx[1],&inc,tr1->tau+1,y+tr1->idx[1],&n1_,work)); |
309 | } | ||
310 |
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.
|
21826 | if (tr1->idx[0]<tr1->idx[1]) PetscCall(HRApply(1,y+tr1->idx[0],1,y+tr1->idx[1],1,tr1->cs,-tr1->sn)); |
311 | 21826 | tr2->n[0] = tr1->n[0]; | |
312 | 21826 | tr2->n[1] = tr1->n[1]; | |
313 | 21826 | tr2->idx[0] = tr1->idx[0]; | |
314 | 21826 | tr2->idx[1] = tr1->idx[1]; | |
315 |
3/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
21826 | if (tr1->idx[0]<tr1->idx[1] && tr1->type==2) { |
316 | 10913 | tr2->idx[1]++; tr2->n[1]--; tr2->n[0]++; | |
317 | } | ||
318 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
21826 | if (tr2->n[0]>0) { |
319 | 21826 | tr2->n[0]--; tr2->idx[0]++; | |
320 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
21826 | if (tr2->n[1]==0) tr2->idx[1] = tr2->idx[0]; |
321 | } else { | ||
322 | ✗ | tr2->n[1]--; tr2->idx[1]++; tr2->idx[0] = tr2->idx[1]; | |
323 | } | ||
324 | /* Hyperbolic transformation to make zeros in y */ | ||
325 |
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.
|
21826 | PetscCall(PetscBLASIntCast(tr2->n[0],&n0_)); |
326 |
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.
|
21826 | PetscCall(PetscBLASIntCast(tr2->n[1],&n1_)); |
327 |
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.
|
21826 | if (tr2->n[0] > 1) PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n0_,y+tr2->idx[0],y+tr2->idx[0]+1,&inc,tr2->tau)); |
328 |
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.
|
21826 | if (tr2->n[1]> 1) PetscCallBLAS("LAPACKlarfg",LAPACKlarfg_(&n1_,y+tr2->idx[1],y+tr2->idx[1]+1,&inc,tr2->tau+1)); |
329 |
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.
|
21826 | if (tr2->idx[0]<tr2->idx[1]) PetscCall(HRGen(PetscRealPart(y[tr2->idx[0]]),PetscRealPart(y[tr2->idx[1]]),&tr2->type,&tr2->cs,&tr2->sn,&tr2->alpha,&ncond2)); |
330 | else { | ||
331 | 7604 | tr2->alpha = PetscRealPart(y[tr2->idx[0]]); | |
332 | 7604 | ncond2 = 1.0; | |
333 | } | ||
334 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
21826 | if (ncond2>*ncond) *ncond = ncond2; |
335 | } | ||
336 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
67713 | PetscFunctionReturn(PETSC_SUCCESS); |
337 | } | ||
338 | |||
339 | /* | ||
340 | Auxiliary function to try perform one iteration of hr routine, | ||
341 | checking condition number. If it is < tolD, apply the | ||
342 | transformation to H and R, if not, ok=false and it do nothing | ||
343 | tolE, tolerance to exchange complex pairs to improve conditioning | ||
344 | */ | ||
345 | 324933 | static PetscErrorCode TryHRIt(PetscInt n,PetscInt j,PetscInt sz,PetscScalar *H,PetscInt ldh,PetscScalar *R,PetscInt ldr,PetscReal *s,PetscBool *exg,PetscBool *ok,PetscInt *n0,PetscInt *n1,PetscInt *idx0,PetscInt *idx1,PetscReal *cond,PetscScalar *work) | |
346 | { | ||
347 | 324933 | struct HRtr *tr1,*tr2,tr1_t,tr2_t,tr1_te,tr2_te; | |
348 | 324933 | PetscScalar *x,*y; | |
349 | 324933 | PetscReal ncond=0,ncond_e; | |
350 | 324933 | PetscInt nwu=0,i,d=1; | |
351 | 324933 | PetscBLASInt n0_,n1_,inc=1,mh,mr,n_,ldr_,ldh_; | |
352 | 324933 | PetscReal tolD = 1e+5; | |
353 | |||
354 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
324933 | PetscFunctionBegin; |
355 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
324933 | if (cond) *cond = 1.0; |
356 |
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.
|
324933 | PetscCall(PetscBLASIntCast(n,&n_)); |
357 |
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.
|
324933 | PetscCall(PetscBLASIntCast(ldr,&ldr_)); |
358 |
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.
|
324933 | PetscCall(PetscBLASIntCast(ldh,&ldh_)); |
359 | 324933 | x = work+nwu; | |
360 | 324933 | nwu += n; | |
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.
|
324933 | PetscCall(PetscArraycpy(x,R+j*ldr,n)); |
362 | 324933 | *exg = PETSC_FALSE; | |
363 | 324933 | *ok = PETSC_TRUE; | |
364 | 324933 | tr1_t.data = x; | |
365 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (sz==1) { |
366 | /* Hyperbolic transformation to make zeros in x */ | ||
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.
|
314020 | PetscCall(MadeHRtr(sz,n,*idx0,*n0,*idx1,*n1,&tr1_t,NULL,&ncond,work+nwu)); |
368 | /* Check condition number to single column*/ | ||
369 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
314020 | if (ncond>tolD) *ok = PETSC_FALSE; |
370 | tr1 = &tr1_t; | ||
371 | tr2 = &tr2_t; | ||
372 | } else { | ||
373 | 10913 | y = work+nwu; | |
374 | 10913 | nwu += n; | |
375 |
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.
|
10913 | PetscCall(PetscArraycpy(y,R+(j+1)*ldr,n)); |
376 | 10913 | tr2_t.data = y; | |
377 |
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.
|
10913 | PetscCall(MadeHRtr(sz,n,*idx0,*n0,*idx1,*n1,&tr1_t,&tr2_t,&ncond,work+nwu)); |
378 | /* Computing hyperbolic transformations also for exchanged vectors */ | ||
379 | 10913 | tr1_te.data = work+nwu; | |
380 | 10913 | nwu += n; | |
381 |
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.
|
10913 | PetscCall(PetscArraycpy(tr1_te.data,R+(j+1)*ldr,n)); |
382 | 10913 | tr2_te.data = work+nwu; | |
383 | 10913 | nwu += n; | |
384 |
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.
|
10913 | PetscCall(PetscArraycpy(tr2_te.data,R+j*ldr,n)); |
385 |
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.
|
10913 | PetscCall(MadeHRtr(sz,n,*idx0,*n0,*idx1,*n1,&tr1_te,&tr2_te,&ncond_e,work+nwu)); |
386 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10913 | if (ncond > d*ncond_e) { |
387 | 6204 | *exg = PETSC_TRUE; | |
388 | 6204 | tr1 = &tr1_te; | |
389 | 6204 | tr2 = &tr2_te; | |
390 | 6204 | ncond = ncond_e; | |
391 | } else { | ||
392 | tr1 = &tr1_t; | ||
393 | tr2 = &tr2_t; | ||
394 | } | ||
395 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
10913 | if (ncond>tolD) *ok = PETSC_FALSE; |
396 | } | ||
397 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
324933 | if (*ok) { |
398 | /* Everything is OK, apply transformations to R and H */ | ||
399 | /* First column */ | ||
400 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
324933 | if (cond && *cond<ncond) *cond = ncond; |
401 | 324933 | x = tr1->data; | |
402 |
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.
|
324933 | PetscCall(PetscBLASIntCast(tr1->n[0],&n0_)); |
403 |
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.
|
324933 | PetscCall(PetscBLASIntCast(tr1->n[1],&n1_)); |
404 |
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.
|
324933 | PetscCall(PetscBLASIntCast(n-j-sz,&mr)); |
405 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
324933 | if (tr1->n[0] > 1 && PetscAbsScalar(tr1->tau[0])!=0.0) { |
406 | 302394 | x[tr1->idx[0]] = 1.0; | |
407 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
302394 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n0_,&mr,x+tr1->idx[0],&inc,tr1->tau,R+(j+sz)*ldr+tr1->idx[0],&ldr_,work+nwu)); |
408 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
302394 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n0_,x+tr1->idx[0],&inc,tr1->tau,H+(tr1->idx[0])*ldh,&ldh_,work+nwu)); |
409 | } | ||
410 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
324933 | if (tr1->n[1] > 1 && PetscAbsScalar(tr1->tau[1])!=0.0) { |
411 | 50965 | x[tr1->idx[1]] = 1.0; | |
412 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
50965 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n1_,&mr,x+tr1->idx[1],&inc,tr1->tau+1,R+(j+sz)*ldr+tr1->idx[1],&ldr_,work+nwu)); |
413 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
50965 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n1_,x+tr1->idx[1],&inc,tr1->tau+1,H+(tr1->idx[1])*ldh,&ldh_,work+nwu)); |
414 | } | ||
415 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (tr1->idx[0]<tr1->idx[1]) { |
416 |
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.
|
89896 | PetscCall(HRApply(mr,R+(j+sz)*ldr+tr1->idx[0],ldr,R+(j+sz)*ldr+tr1->idx[1],ldr,tr1->cs,-tr1->sn)); |
417 |
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.
|
89896 | if (tr1->type==1) PetscCall(HRApply(n,H+(tr1->idx[0])*ldh,1,H+(tr1->idx[1])*ldh,1,tr1->cs,tr1->sn)); |
418 | else { | ||
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.
|
13879 | PetscCall(HRApply(n,H+(tr1->idx[0])*ldh,1,H+(tr1->idx[1])*ldh,1,-tr1->cs,-tr1->sn)); |
420 | 13879 | s[tr1->idx[0]] = -s[tr1->idx[0]]; | |
421 | 13879 | s[tr1->idx[1]] = -s[tr1->idx[1]]; | |
422 | } | ||
423 | } | ||
424 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6703128 | for (i=0;i<tr1->idx[0];i++) *(R+j*ldr+i) = x[i]; |
425 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6567374 | for (i=tr1->idx[0]+1;i<n;i++) *(R+j*ldr+i) = 0.0; |
426 | 324933 | *(R+j*ldr+tr1->idx[0]) = tr1->alpha; | |
427 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (sz==2) { |
428 | 10913 | y = tr2->data; | |
429 | /* Second column */ | ||
430 |
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.
|
10913 | PetscCall(PetscBLASIntCast(tr2->n[0],&n0_)); |
431 |
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.
|
10913 | PetscCall(PetscBLASIntCast(tr2->n[1],&n1_)); |
432 |
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.
|
10913 | PetscCall(PetscBLASIntCast(n-j-sz,&mr)); |
433 |
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.
|
10913 | PetscCall(PetscBLASIntCast(n-tr2->idx[0],&mh)); |
434 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
10913 | if (tr2->n[0] > 1 && PetscAbsScalar(tr2->tau[0])!=0.0) { |
435 | 10101 | y[tr2->idx[0]] = 1.0; | |
436 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
10101 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n0_,&mr,y+tr2->idx[0],&inc,tr2->tau,R+(j+2)*ldr+tr2->idx[0],&ldr_,work+nwu)); |
437 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
10101 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n0_,y+tr2->idx[0],&inc,tr2->tau,H+(tr2->idx[0])*ldh,&ldh_,work+nwu)); |
438 | } | ||
439 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
10913 | if (tr2->n[1] > 1 && PetscAbsScalar(tr2->tau[1])!=0.0) { |
440 | 3508 | y[tr2->idx[1]] = 1.0; | |
441 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
3508 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("L",&n1_,&mr,y+tr2->idx[1],&inc,tr2->tau+1,R+(j+2)*ldr+tr2->idx[1],&ldr_,work+nwu)); |
442 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
3508 | PetscCallBLAS("LAPACKlarf",LAPACKlarf_("R",&n_,&n1_,y+tr2->idx[1],&inc,tr2->tau+1,H+(tr2->idx[1])*ldh,&ldh_,work+nwu)); |
443 | } | ||
444 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10913 | if (tr2->idx[0]<tr2->idx[1]) { |
445 |
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.
|
8102 | PetscCall(HRApply(mr,R+(j+2)*ldr+tr2->idx[0],ldr,R+(j+2)*ldr+tr2->idx[1],ldr,tr2->cs,-tr2->sn)); |
446 |
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.
|
8102 | if (tr2->type==1) PetscCall(HRApply(n,H+(tr2->idx[0])*ldh,1,H+(tr2->idx[1])*ldh,1,tr2->cs,tr2->sn)); |
447 | else { | ||
448 |
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.
|
6343 | PetscCall(HRApply(n,H+(tr2->idx[0])*ldh,1,H+(tr2->idx[1])*ldh,1,-tr2->cs,-tr2->sn)); |
449 | 6343 | s[tr2->idx[0]] = -s[tr2->idx[0]]; | |
450 | 6343 | s[tr2->idx[1]] = -s[tr2->idx[1]]; | |
451 | } | ||
452 | } | ||
453 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
122246 | for (i=0;i<tr2->idx[0]-1;i++) *(R+(j+1)*ldr+i) = y[i]; |
454 | 10913 | *(R+(j+1)*ldr+tr2->idx[0]-1) = y[tr2->idx[0]-1]; | |
455 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
268913 | for (i=tr2->idx[0]+1;i<n;i++) *(R+(j+1)*ldr+i) = 0.0; |
456 | 10913 | *(R+(j+1)*ldr+tr2->idx[0]) = tr2->alpha; | |
457 | 10913 | *n0 = tr2->n[0]; | |
458 | 10913 | *n1 = tr2->n[1]; | |
459 | 10913 | *idx0 = tr2->idx[0]; | |
460 | 10913 | *idx1 = tr2->idx[1]; | |
461 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
10913 | if (tr2->idx[0]<tr2->idx[1] && tr2->type==2) { |
462 | 6343 | (*idx1)++; (*n1)--; (*n0)++; | |
463 | } | ||
464 | } else { | ||
465 | 314020 | *n0 = tr1->n[0]; | |
466 | 314020 | *n1 = tr1->n[1]; | |
467 | 314020 | *idx0 = tr1->idx[0]; | |
468 | 314020 | *idx1 = tr1->idx[1]; | |
469 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
314020 | if (tr1->idx[0]<tr1->idx[1] && tr1->type==2) { |
470 | 9568 | (*idx1)++; (*n1)--; (*n0)++; | |
471 | } | ||
472 | } | ||
473 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (*n0>0) { |
474 | 313560 | (*n0)--; (*idx0)++; | |
475 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
313560 | if (*n1==0) *idx1 = *idx0; |
476 | } else { | ||
477 | 11373 | (*n1)--; (*idx1)++; *idx0 = *idx1; | |
478 | } | ||
479 | } | ||
480 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
65465 | PetscFunctionReturn(PETSC_SUCCESS); |
481 | } | ||
482 | |||
483 | /* | ||
484 | compute V = HR whit H s-orthogonal and R upper triangular | ||
485 | (work: work space of minimum size 6*nv) | ||
486 | */ | ||
487 | 10520 | PetscErrorCode DSPseudoOrthog_HR(PetscInt *nv,PetscScalar *V,PetscInt ldv,PetscReal *s,PetscScalar *R,PetscInt ldr,PetscBLASInt *perm,PetscBLASInt *cmplxEig,PetscBool *breakdown,PetscScalar *work) | |
488 | { | ||
489 | 10520 | PetscInt i,j,n,n0,n1,np,idx0,idx1,sz=1,k=0,nwu=0; | |
490 | 10520 | PetscBLASInt t1,t2; | |
491 | 10520 | PetscScalar *col1,*col2; | |
492 | 10520 | PetscBool exg=PETSC_FALSE,ok=PETSC_FALSE; | |
493 | |||
494 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10520 | PetscFunctionBegin; |
495 | 10520 | n = *nv; | |
496 | 10520 | col1 = work+nwu; | |
497 | 10520 | nwu += n; | |
498 | 10520 | col2 = work+nwu; | |
499 | 10520 | nwu += n; | |
500 | /* Sort R and s according to sing(s) */ | ||
501 | 10520 | np = 0; | |
502 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
346366 | for (i=0;i<n;i++) if (s[i]>0) np++; |
503 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10520 | if (s[0]>0) n1 = np; |
504 | 994 | else n1 = n-np; | |
505 | 10520 | n0 = 0; | |
506 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
346366 | for (i=0;i<n;i++) { |
507 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
335846 | if (s[i]==s[0]) { |
508 | 304251 | s[n0] = s[0]; | |
509 |
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.
|
304251 | PetscCall(PetscBLASIntCast(i,&perm[n0++])); |
510 |
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.
|
335846 | } else PetscCall(PetscBLASIntCast(i,&perm[n1++])); |
511 | } | ||
512 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
42115 | for (i=n0;i<n;i++) s[i] = -s[0]; |
513 | 10520 | n1 -= n0; | |
514 | 10520 | idx0 = 0; | |
515 | 10520 | idx1 = n0; | |
516 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10520 | if (idx1==n) idx1=idx0; |
517 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
346366 | for (i=0;i<n;i++) { |
518 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
13672574 | for (j=0;j<n;j++) R[j*ldr+i] = V[j*ldv+perm[i]]; |
519 | } | ||
520 | /* Initialize H */ | ||
521 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
346366 | for (i=0;i<n;i++) { |
522 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
335846 | PetscCall(PetscArrayzero(V+i*ldv,n)); |
523 | 335846 | V[perm[i]+i*ldv] = 1.0; | |
524 | } | ||
525 |
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.
|
346366 | for (i=0;i<n;i++) PetscCall(PetscBLASIntCast(i,&perm[i])); |
526 | j = 0; | ||
527 | 335453 | while (j<n-k) { | |
528 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (cmplxEig[j]==0) sz=1; |
529 | 10913 | else sz=2; | |
530 |
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.
|
324933 | PetscCall(TryHRIt(n,j,sz,V,ldv,R,ldr,s,&exg,&ok,&n0,&n1,&idx0,&idx1,NULL,work+nwu)); |
531 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
324933 | if (ok) { |
532 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324933 | if (exg) cmplxEig[j] = -cmplxEig[j]; |
533 | 324933 | j = j+sz; | |
534 | } else { /* to be discarded */ | ||
535 | ✗ | k = k+1; | |
536 | ✗ | if (cmplxEig[j]==0) { | |
537 | ✗ | if (j<n) { | |
538 | ✗ | t1 = perm[j]; | |
539 | ✗ | for (i=j;i<n-1;i++) perm[i] = perm[i+1]; | |
540 | ✗ | perm[n-1] = t1; | |
541 | ✗ | t1 = cmplxEig[j]; | |
542 | ✗ | for (i=j;i<n-1;i++) cmplxEig[i] = cmplxEig[i+1]; | |
543 | ✗ | cmplxEig[n-1] = t1; | |
544 | ✗ | PetscCall(PetscArraycpy(col1,R+j*ldr,n)); | |
545 | ✗ | for (i=j;i<n-1;i++) PetscCall(PetscArraycpy(R+i*ldr,R+(i+1)*ldr,n)); | |
546 | ✗ | PetscCall(PetscArraycpy(R+(n-1)*ldr,col1,n)); | |
547 | } | ||
548 | } else { | ||
549 | ✗ | k = k+1; | |
550 | ✗ | if (j<n-1) { | |
551 | ✗ | t1 = perm[j]; t2 = perm[j+1]; | |
552 | ✗ | for (i=j;i<n-2;i++) perm[i] = perm[i+2]; | |
553 | ✗ | perm[n-2] = t1; perm[n-1] = t2; | |
554 | ✗ | t1 = cmplxEig[j]; t2 = cmplxEig[j+1]; | |
555 | ✗ | for (i=j;i<n-2;i++) cmplxEig[i] = cmplxEig[i+2]; | |
556 | ✗ | cmplxEig[n-2] = t1; cmplxEig[n-1] = t2; | |
557 | ✗ | PetscCall(PetscArraycpy(col1,R+j*ldr,n)); | |
558 | ✗ | PetscCall(PetscArraycpy(col2,R+(j+1)*ldr,n)); | |
559 | ✗ | for (i=j;i<n-2;i++) PetscCall(PetscArraycpy(R+i*ldr,R+(i+2)*ldr,n)); | |
560 | ✗ | PetscCall(PetscArraycpy(R+(n-2)*ldr,col1,n)); | |
561 |
4/8✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 8 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 2 times.
|
335453 | PetscCall(PetscArraycpy(R+(n-1)*ldr,col2,n)); |
562 | } | ||
563 | } | ||
564 | } | ||
565 | } | ||
566 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
10520 | if (k!=0) { |
567 | ✗ | if (breakdown) *breakdown = PETSC_TRUE; | |
568 | ✗ | *nv = n-k; | |
569 | } | ||
570 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
2092 | PetscFunctionReturn(PETSC_SUCCESS); |
571 | } | ||
572 | |||
573 | 10500 | PetscErrorCode DSGHIEPOrthogEigenv(DS ds,DSMatType mat,PetscScalar *wr,PetscScalar *wi,PetscBool accum) | |
574 | { | ||
575 | 10500 | PetscInt lws,nwus=0,nwui=0,lwi,off,n,nv,ld,i,ldr,l; | |
576 | 10500 | const PetscScalar *B,*W; | |
577 | 10500 | PetscScalar *Q,*X,*R,*ts,szero=0.0,sone=1.0; | |
578 | 10500 | PetscReal *s,vi,vr,tr,*d,*e; | |
579 | 10500 | PetscBLASInt ld_,n_,nv_,*perm,*cmplxEig; | |
580 | |||
581 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10500 | PetscFunctionBegin; |
582 | 10500 | l = ds->l; | |
583 | 10500 | n = ds->n-l; | |
584 |
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.
|
10500 | PetscCall(PetscBLASIntCast(n,&n_)); |
585 | 10500 | ld = ds->ld; | |
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.
|
10500 | PetscCall(PetscBLASIntCast(ld,&ld_)); |
587 | 10500 | off = l*ld+l; | |
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.
|
10500 | PetscCall(DSGetArrayReal(ds,DS_MAT_D,&s)); |
589 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10500 | if (!ds->compact) { |
590 |
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.
|
40 | PetscCall(MatDenseGetArrayRead(ds->omat[DS_MAT_B],&B)); |
591 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
440 | for (i=l;i<ds->n;i++) s[i] = PetscRealPart(B[i*ld+i]); |
592 |
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.
|
40 | PetscCall(MatDenseRestoreArrayRead(ds->omat[DS_MAT_B],&B)); |
593 | } | ||
594 | 10500 | lws = n*n+7*n; | |
595 | 10500 | lwi = 2*n; | |
596 |
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.
|
10500 | PetscCall(DSAllocateWork_Private(ds,lws,0,lwi)); |
597 | 10500 | R = ds->work+nwus; | |
598 | 10500 | nwus += n*n; | |
599 | 10500 | ldr = n; | |
600 | 10500 | perm = ds->iwork + nwui; | |
601 | 10500 | nwui += n; | |
602 | 10500 | cmplxEig = ds->iwork+nwui; | |
603 |
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.
|
10500 | PetscCall(MatDenseGetArray(ds->omat[mat],&X)); |
604 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
335273 | for (i=0;i<n;i++) { |
605 | #if defined(PETSC_USE_COMPLEX) | ||
606 | 47006 | vi = PetscImaginaryPart(wr[l+i]); | |
607 | #else | ||
608 |
1/2✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
|
277767 | vi = wi?PetscRealPart(wi[l+i]):0.0; |
609 | #endif | ||
610 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
324773 | if (vi!=0) { |
611 | 10913 | cmplxEig[i] = 1; | |
612 | 10913 | cmplxEig[i+1] = 2; | |
613 | 10913 | i++; | |
614 | 313860 | } else cmplxEig[i] = 0; | |
615 | } | ||
616 | 10500 | nv = n; | |
617 | |||
618 | /* Perform HR decomposition */ | ||
619 | /* Hyperbolic rotators */ | ||
620 |
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.
|
10500 | PetscCall(DSPseudoOrthog_HR(&nv,X+off,ld,s+l,R,ldr,perm,cmplxEig,NULL,ds->work+nwus)); |
621 | /* Sort wr,wi perm */ | ||
622 | 10500 | ts = ds->work+nwus; | |
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.
|
10500 | PetscCall(PetscArraycpy(ts,wr+l,n)); |
624 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
346186 | for (i=0;i<n;i++) wr[i+l] = ts[perm[i]]; |
625 | #if !defined(PETSC_USE_COMPLEX) | ||
626 |
1/2✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
|
8714 | if (wi) { |
627 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
8714 | PetscCall(PetscArraycpy(ts,wi+l,n)); |
628 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
296220 | for (i=0;i<n;i++) wi[i+l] = ts[perm[i]]; |
629 | } | ||
630 | #endif | ||
631 | /* Projected Matrix */ | ||
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.
|
10500 | PetscCall(DSGetArrayReal(ds,DS_MAT_T,&d)); |
633 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10500 | PetscCall(PetscArrayzero(d+2*ld,ld)); |
634 | 10500 | e = d+ld; | |
635 | 10500 | d[l+nv-1] = PetscRealPart(wr[l+nv-1]*s[l+nv-1]); | |
636 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
325110 | for (i=0;i<nv-1;i++) { |
637 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
314610 | if (cmplxEig[i]==0) { /* Real */ |
638 | 303697 | d[l+i] = PetscRealPart(wr[l+i]*s[l+i]); | |
639 | 303697 | e[l+i] = 0.0; | |
640 | } else { | ||
641 | 10913 | vr = PetscRealPart(wr[l+i]); | |
642 | #if defined(PETSC_USE_COMPLEX) | ||
643 | 1174 | vi = PetscImaginaryPart(wr[l+i]); | |
644 | #else | ||
645 |
1/2✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
|
9739 | vi = wi?PetscRealPart(wi[l+i]):0.0; |
646 | #endif | ||
647 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10913 | if (cmplxEig[i]==-1) vi = -vi; |
648 | 10913 | tr = PetscRealPart((R[i+(i+1)*ldr]/R[i+i*ldr]))*vi; | |
649 | 10913 | d[l+i] = (vr-tr)*s[l+i]; | |
650 | 10913 | d[l+i+1] = (vr+tr)*s[l+i+1]; | |
651 | 10913 | e[l+i] = PetscRealPart(s[l+i]*(R[(i+1)+(i+1)*ldr]/R[i+i*ldr])*vi); | |
652 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10913 | if (i<nv-2) e[l+i+1] = 0.0; |
653 | i++; | ||
654 | } | ||
655 | } | ||
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.
|
10500 | PetscCall(DSRestoreArrayReal(ds,DS_MAT_T,&d)); |
657 |
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.
|
10500 | PetscCall(DSRestoreArrayReal(ds,DS_MAT_D,&s)); |
658 | /* accumulate previous Q */ | ||
659 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10500 | if (accum) { |
660 |
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.
|
10470 | PetscCall(MatDenseGetArray(ds->omat[DS_MAT_Q],&Q)); |
661 |
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.
|
10470 | PetscCall(PetscBLASIntCast(nv,&nv_)); |
662 |
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.
|
10470 | PetscCall(DSAllocateMat_Private(ds,DS_MAT_W)); |
663 |
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.
|
10470 | PetscCall(MatCopy(ds->omat[DS_MAT_Q],ds->omat[DS_MAT_W],SAME_NONZERO_PATTERN)); |
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.
|
10470 | PetscCall(MatDenseGetArrayRead(ds->omat[DS_MAT_W],&W)); |
665 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
10470 | PetscCallBLAS("BLASgemm",BLASgemm_("N","N",&n_,&nv_,&n_,&sone,W+off,&ld_,X+off,&ld_,&szero,Q+off,&ld_)); |
666 |
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.
|
10470 | PetscCall(MatDenseRestoreArrayRead(ds->omat[DS_MAT_W],&W)); |
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.
|
10470 | PetscCall(MatDenseRestoreArray(ds->omat[DS_MAT_Q],&Q)); |
668 | } | ||
669 | 10500 | ds->t = nv+l; | |
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.
|
10500 | PetscCall(MatDenseRestoreArray(ds->omat[mat],&X)); |
671 |
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.
|
10500 | if (!ds->compact) PetscCall(DSSwitchFormat_GHIEP(ds,PETSC_FALSE)); |
672 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
2088 | PetscFunctionReturn(PETSC_SUCCESS); |
673 | } | ||
674 | |||
675 | /* | ||
676 | Reduce to tridiagonal-diagonal pair by means of TridiagDiag_HHR. | ||
677 | */ | ||
678 | 10500 | PetscErrorCode DSIntermediate_GHIEP(DS ds) | |
679 | { | ||
680 | 10500 | PetscInt i,ld,off; | |
681 | 10500 | PetscInt nwall,nwallr,nwalli; | |
682 | 10500 | PetscScalar *A,*B,*Q; | |
683 | 10500 | PetscReal *d,*e,*s; | |
684 | |||
685 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10500 | PetscFunctionBegin; |
686 | 10500 | ld = ds->ld; | |
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.
|
10500 | PetscCall(MatDenseGetArray(ds->omat[DS_MAT_A],&A)); |
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.
|
10500 | PetscCall(MatDenseGetArray(ds->omat[DS_MAT_B],&B)); |
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.
|
10500 | PetscCall(MatDenseGetArray(ds->omat[DS_MAT_Q],&Q)); |
690 |
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.
|
10500 | PetscCall(DSGetArrayReal(ds,DS_MAT_T,&d)); |
691 |
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.
|
10500 | PetscCall(DSGetArrayReal(ds,DS_MAT_D,&s)); |
692 | 10500 | e = d+ld; | |
693 | 10500 | off = ds->l+ds->l*ld; | |
694 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10500 | PetscCall(PetscArrayzero(Q,ld*ld)); |
695 | 10500 | nwall = ld*ld+ld; | |
696 | 10500 | nwallr = ld; | |
697 | 10500 | nwalli = ld; | |
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.
|
10500 | PetscCall(DSAllocateWork_Private(ds,nwall,nwallr,nwalli)); |
699 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
357842 | for (i=0;i<ds->n;i++) Q[i+i*ld]=1.0; |
700 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
346186 | for (i=0;i<ds->n-ds->l;i++) *(ds->perm+i)=i; |
701 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10500 | if (ds->compact) { |
702 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10460 | if (ds->state < DS_STATE_INTERMEDIATE) { |
703 |
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.
|
9961 | PetscCall(DSSwitchFormat_GHIEP(ds,PETSC_FALSE)); |
704 |
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.
|
9961 | PetscCall(TridiagDiag_HHR(ds->k-ds->l+1,A+off,ld,s+ds->l,Q+off,ld,PETSC_TRUE,d+ds->l,e+ds->l,ds->perm,ds->work,ds->rwork,ds->iwork)); |
705 | 9961 | ds->k = ds->l; | |
706 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9961 | PetscCall(PetscArrayzero(d+2*ld+ds->l,ds->n-ds->l)); |
707 | } | ||
708 | } else { | ||
709 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
40 | if (ds->state < DS_STATE_INTERMEDIATE) { |
710 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
440 | for (i=0;i<ds->n;i++) s[i] = PetscRealPart(B[i+i*ld]); |
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.
|
40 | PetscCall(TridiagDiag_HHR(ds->n-ds->l,A+off,ld,s+ds->l,Q+off,ld,PETSC_FALSE,d+ds->l,e+ds->l,ds->perm,ds->work,ds->rwork,ds->iwork)); |
712 |
4/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(PetscArrayzero(d+2*ld,ds->n)); |
713 | 40 | ds->k = ds->l; | |
714 |
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.
|
40 | PetscCall(DSSwitchFormat_GHIEP(ds,PETSC_FALSE)); |
715 | ✗ | } else PetscCall(DSSwitchFormat_GHIEP(ds,PETSC_TRUE)); | |
716 | } | ||
717 |
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.
|
10500 | PetscCall(MatDenseRestoreArray(ds->omat[DS_MAT_A],&A)); |
718 |
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.
|
10500 | PetscCall(MatDenseRestoreArray(ds->omat[DS_MAT_B],&B)); |
719 |
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.
|
10500 | PetscCall(MatDenseRestoreArray(ds->omat[DS_MAT_Q],&Q)); |
720 |
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.
|
10500 | PetscCall(DSRestoreArrayReal(ds,DS_MAT_T,&d)); |
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.
|
10500 | PetscCall(DSRestoreArrayReal(ds,DS_MAT_D,&s)); |
722 |
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
2088 | PetscFunctionReturn(PETSC_SUCCESS); |
723 | } | ||
724 |