GCC Code Coverage Report

Directory: ./
File: src/sys/classes/fn/impls/invsqrt/fninvsqrt.c
Date: 2026-10-22 02:08:06
Exec Total Coverage
Lines: 134 182 73.6%
Functions: 10 12 83.3%
Branches: 254 746 34.0%
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 Inverse square root function x^(-1/2)
12 */
13
14 #include <slepc/private/fnimpl.h> /*I "slepcfn.h" I*/
15 #include <slepcblaslapack.h>
16
17 72 static PetscErrorCode FNEvaluateFunction_Invsqrt(FN fn,PetscScalar x,PetscScalar *y)
18 {
19
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 72 PetscFunctionBegin;
20
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 72 PetscCheck(x!=0.0,PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Function not defined in the requested value");
21 #if !PetscDefined(USE_COMPLEX)
22
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 36 PetscCheck(x>0.0,PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Function not defined in the requested value");
23 #endif
24 72 *y = 1.0/PetscSqrtScalar(x);
25
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 72 PetscFunctionReturn(PETSC_SUCCESS);
26 }
27
28 72 static PetscErrorCode FNEvaluateDerivative_Invsqrt(FN fn,PetscScalar x,PetscScalar *y)
29 {
30
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 72 PetscFunctionBegin;
31
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 72 PetscCheck(x!=0.0,PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Derivative not defined in the requested value");
32 #if !PetscDefined(USE_COMPLEX)
33
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 36 PetscCheck(x>0.0,PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Derivative not defined in the requested value");
34 #endif
35
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 72 *y = -1.0/(2.0*PetscPowScalarReal(x,1.5));
36
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 72 PetscFunctionReturn(PETSC_SUCCESS);
37 }
38
39 ✗ static PetscErrorCode FNEvaluateFunctionMat_Invsqrt_Schur(FN fn,Mat A,Mat B)
40 {
41 ✗ PetscBLASInt n=0,ld,*ipiv;
42 ✗ PetscScalar *Ba,*Wa;
43 ✗ PetscInt m;
44 ✗ Mat W;
45
46 ✗ PetscFunctionBegin;
47 ✗ PetscCall(FN_AllocateWorkMat(fn,A,&W));
48 ✗ if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
49 ✗ PetscCall(MatDenseGetArray(B,&Ba));
50 ✗ PetscCall(MatDenseGetArray(W,&Wa));
51 /* compute B = sqrtm(A) */
52 ✗ PetscCall(MatGetSize(A,&m,NULL));
53 ✗ PetscCall(PetscBLASIntCast(m,&n));
54 ✗ ld = n;
55 ✗ PetscCall(FNSqrtmSchur(fn,n,Ba,n,PETSC_FALSE));
56 /* compute B = A\B */
57 ✗ PetscCall(PetscMalloc1(ld,&ipiv));
58 ✗ PetscCallLAPACKInfo("LAPACKgesv",LAPACKgesv_(&n,&n,Wa,&ld,ipiv,Ba,&ld,&info));
59 ✗ PetscCall(PetscLogFlops(2.0*n*n*n/3.0+2.0*n*n*n));
60 ✗ PetscCall(PetscFree(ipiv));
61 ✗ PetscCall(MatDenseRestoreArray(W,&Wa));
62 ✗ PetscCall(MatDenseRestoreArray(B,&Ba));
63 ✗ PetscCall(FN_FreeWorkMat(fn,&W));
64 ✗ PetscFunctionReturn(PETSC_SUCCESS);
65 }
66
67 ✗ static PetscErrorCode FNEvaluateFunctionMatVec_Invsqrt_Schur(FN fn,Mat A,Vec v)
68 {
69 ✗ PetscBLASInt n=0,ld,*ipiv,one=1;
70 ✗ PetscScalar *Ba,*Wa;
71 ✗ PetscInt m;
72 ✗ Mat B,W;
73
74 ✗ PetscFunctionBegin;
75 ✗ PetscCall(FN_AllocateWorkMat(fn,A,&B));
76 ✗ PetscCall(FN_AllocateWorkMat(fn,A,&W));
77 ✗ PetscCall(MatDenseGetArray(B,&Ba));
78 ✗ PetscCall(MatDenseGetArray(W,&Wa));
79 /* compute B_1 = sqrtm(A)*e_1 */
80 ✗ PetscCall(MatGetSize(A,&m,NULL));
81 ✗ PetscCall(PetscBLASIntCast(m,&n));
82 ✗ ld = n;
83 ✗ PetscCall(FNSqrtmSchur(fn,n,Ba,n,PETSC_TRUE));
84 /* compute B_1 = A\B_1 */
85 ✗ PetscCall(PetscMalloc1(ld,&ipiv));
86 ✗ PetscCallLAPACKInfo("LAPACKgesv",LAPACKgesv_(&n,&one,Wa,&ld,ipiv,Ba,&ld,&info));
87 ✗ PetscCall(PetscFree(ipiv));
88 ✗ PetscCall(MatDenseRestoreArray(W,&Wa));
89 ✗ PetscCall(MatDenseRestoreArray(B,&Ba));
90 ✗ PetscCall(MatGetColumnVector(B,v,0));
91 ✗ PetscCall(FN_FreeWorkMat(fn,&W));
92 ✗ PetscCall(FN_FreeWorkMat(fn,&B));
93 ✗ PetscFunctionReturn(PETSC_SUCCESS);
94 }
95
96 120 static PetscErrorCode FNEvaluateFunctionMat_Invsqrt_DBP(FN fn,Mat A,Mat B)
97 {
98 120 PetscBLASInt n=0;
99 120 PetscScalar *T;
100 120 PetscInt m;
101
102
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 120 PetscFunctionBegin;
103
5/8
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 120 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
104
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseGetArray(B,&T));
105
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatGetSize(A,&m,NULL));
106
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(PetscBLASIntCast(m,&n));
107
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(FNSqrtmDenmanBeavers(fn,n,T,n,PETSC_TRUE));
108
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseRestoreArray(B,&T));
109
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 24 PetscFunctionReturn(PETSC_SUCCESS);
110 }
111
112 120 static PetscErrorCode FNEvaluateFunctionMat_Invsqrt_NS(FN fn,Mat A,Mat B)
113 {
114 120 PetscBLASInt n=0;
115 120 PetscScalar *T;
116 120 PetscInt m;
117
118
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 120 PetscFunctionBegin;
119
5/8
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 120 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
120
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseGetArray(B,&T));
121
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatGetSize(A,&m,NULL));
122
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(PetscBLASIntCast(m,&n));
123
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(FNSqrtmNewtonSchulz(fn,n,T,n,PETSC_TRUE));
124
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseRestoreArray(B,&T));
125
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 24 PetscFunctionReturn(PETSC_SUCCESS);
126 }
127
128 120 static PetscErrorCode FNEvaluateFunctionMat_Invsqrt_Sadeghi(FN fn,Mat A,Mat B)
129 {
130 120 PetscBLASInt n=0,ld,*ipiv;
131 120 PetscScalar *Ba,*Wa;
132 120 PetscInt m;
133 120 Mat W;
134
135
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 120 PetscFunctionBegin;
136
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(FN_AllocateWorkMat(fn,A,&W));
137
5/8
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 120 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
138
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseGetArray(B,&Ba));
139
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseGetArray(W,&Wa));
140 /* compute B = sqrtm(A) */
141
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatGetSize(A,&m,NULL));
142
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(PetscBLASIntCast(m,&n));
143 120 ld = n;
144
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(FNSqrtmSadeghi(fn,n,Ba,n));
145 /* compute B = A\B */
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.
 120 PetscCall(PetscMalloc1(ld,&ipiv));
147
16/32
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 8 times.
✗ 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 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.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
 120 PetscCallLAPACKInfo("LAPACKgesv",LAPACKgesv_(&n,&n,Wa,&ld,ipiv,Ba,&ld,&info));
148
3/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 taken 2 times.
 120 PetscCall(PetscLogFlops(2.0*n*n*n/3.0+2.0*n*n*n));
149
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 2 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.
 120 PetscCall(PetscFree(ipiv));
150
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseRestoreArray(W,&Wa));
151
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 120 PetscCall(MatDenseRestoreArray(B,&Ba));
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.
 120 PetscCall(FN_FreeWorkMat(fn,&W));
153
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 24 PetscFunctionReturn(PETSC_SUCCESS);
154 }
155
156 #if PetscDefined(HAVE_CUDA)
157 24 PetscErrorCode FNEvaluateFunctionMat_Invsqrt_NS_CUDA(FN fn,Mat A,Mat B)
158 {
159 24 PetscBLASInt n=0;
160 24 PetscScalar *Ba;
161 24 PetscInt m;
162
163 24 PetscFunctionBegin;
164
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
 24 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
165
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDAGetArray(B,&Ba));
166
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatGetSize(A,&m,NULL));
167
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(PetscBLASIntCast(m,&n));
168
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(FNSqrtmNewtonSchulz_CUDA(fn,n,Ba,n,PETSC_TRUE));
169
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDARestoreArray(B,&Ba));
170 PetscFunctionReturn(PETSC_SUCCESS);
171 }
172
173 #if PetscDefined(HAVE_MAGMA)
174 #include <slepcmagma.h>
175
176 24 PetscErrorCode FNEvaluateFunctionMat_Invsqrt_DBP_CUDAm(FN fn,Mat A,Mat B)
177 {
178 24 PetscBLASInt n=0;
179 24 PetscScalar *T;
180 24 PetscInt m;
181
182 24 PetscFunctionBegin;
183
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
 24 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
184
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDAGetArray(B,&T));
185
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatGetSize(A,&m,NULL));
186
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(PetscBLASIntCast(m,&n));
187
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(FNSqrtmDenmanBeavers_CUDAm(fn,n,T,n,PETSC_TRUE));
188
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDARestoreArray(B,&T));
189 PetscFunctionReturn(PETSC_SUCCESS);
190 }
191
192 24 PetscErrorCode FNEvaluateFunctionMat_Invsqrt_Sadeghi_CUDAm(FN fn,Mat A,Mat B)
193 {
194 24 PetscBLASInt n=0,ld,*ipiv;
195 24 PetscScalar *Ba,*Wa;
196 24 PetscInt m;
197 24 Mat W;
198
199 24 PetscFunctionBegin;
200
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(FN_AllocateWorkMat(fn,A,&W));
201
2/4
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
 24 if (A!=B) PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN));
202
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDAGetArray(B,&Ba));
203
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDAGetArray(W,&Wa));
204 /* compute B = sqrtm(A) */
205
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatGetSize(A,&m,NULL));
206
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(PetscBLASIntCast(m,&n));
207 24 ld = n;
208
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(FNSqrtmSadeghi_CUDAm(fn,n,Ba,n));
209 /* compute B = A\B */
210
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(SlepcMagmaInit());
211
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(PetscMalloc1(ld,&ipiv));
212
2/6
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 24 PetscCallMAGMA(magma_xgesv_gpu,n,n,Wa,ld,ipiv,Ba,ld);
213 24 PetscCall(PetscLogFlops(2.0*n*n*n/3.0+2.0*n*n*n));
214
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(PetscFree(ipiv));
215
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDARestoreArray(W,&Wa));
216
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(MatDenseCUDARestoreArray(B,&Ba));
217
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
 24 PetscCall(FN_FreeWorkMat(fn,&W));
218 PetscFunctionReturn(PETSC_SUCCESS);
219 }
220 #endif /* PETSC_HAVE_MAGMA */
221 #endif /* PETSC_HAVE_CUDA */
222
223 72 static PetscErrorCode FNView_Invsqrt(FN fn,PetscViewer viewer)
224 {
225 72 PetscBool isascii;
226 72 char str[50];
227 72 const char *methodname[] = {
228 "Schur method for inv(A)*sqrtm(A)",
229 "Denman-Beavers (product form)",
230 "Newton-Schulz iteration",
231 "Sadeghi iteration"
232 };
233 72 const int nmeth=PETSC_STATIC_ARRAY_LENGTH(methodname);
234
235
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 72 PetscFunctionBegin;
236
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.
 72 PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
237
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
 72 if (isascii) {
238
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
 72 if (fn->beta==(PetscScalar)1.0) {
239 ✗ if (fn->alpha==(PetscScalar)1.0) PetscCall(PetscViewerASCIIPrintf(viewer," inverse square root: x^(-1/2)\n"));
240 else {
241 ✗ PetscCall(SlepcSNPrintfScalar(str,sizeof(str),fn->alpha,PETSC_TRUE));
242 ✗ PetscCall(PetscViewerASCIIPrintf(viewer," inverse square root: (%s*x)^(-1/2)\n",str));
243 }
244 } else {
245
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.
 72 PetscCall(SlepcSNPrintfScalar(str,sizeof(str),fn->beta,PETSC_TRUE));
246
1/8
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 72 if (fn->alpha==(PetscScalar)1.0) PetscCall(PetscViewerASCIIPrintf(viewer," inverse square root: %s*x^(-1/2)\n",str));
247 else {
248
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.
 72 PetscCall(PetscViewerASCIIPrintf(viewer," inverse square root: %s",str));
249
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.
 72 PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
250
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.
 72 PetscCall(SlepcSNPrintfScalar(str,sizeof(str),fn->alpha,PETSC_TRUE));
251
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.
 72 PetscCall(PetscViewerASCIIPrintf(viewer,"*(%s*x)^(-1/2)\n",str));
252
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.
 72 PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
253 }
254 }
255
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.
 72 if (fn->method<nmeth) PetscCall(PetscViewerASCIIPrintf(viewer," computing matrix functions with: %s\n",methodname[fn->method]));
256 }
257
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 12 PetscFunctionReturn(PETSC_SUCCESS);
258 }
259
260 /*MC
261 FNINVSQRT - FNINVSQRT = "invsqrt" - The inverse square root function $f(x)=x^{-\frac{1}{2}}$.
262
263 Level: beginner
264
265 .seealso: [](sec:fn), `FN`, `FNType`, `FNSetType()`
266 M*/
267
268 72 SLEPC_EXTERN PetscErrorCode FNCreate_Invsqrt(FN fn)
269 {
270
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 72 PetscFunctionBegin;
271 72 fn->ops->evaluatefunction = FNEvaluateFunction_Invsqrt;
272 72 fn->ops->evaluatederivative = FNEvaluateDerivative_Invsqrt;
273 72 fn->ops->evaluatefunctionmat[0] = FNEvaluateFunctionMat_Invsqrt_Schur;
274 72 fn->ops->evaluatefunctionmat[1] = FNEvaluateFunctionMat_Invsqrt_DBP;
275 72 fn->ops->evaluatefunctionmat[2] = FNEvaluateFunctionMat_Invsqrt_NS;
276 72 fn->ops->evaluatefunctionmat[3] = FNEvaluateFunctionMat_Invsqrt_Sadeghi;
277 #if PetscDefined(HAVE_CUDA)
278 24 fn->ops->evaluatefunctionmatcuda[2] = FNEvaluateFunctionMat_Invsqrt_NS_CUDA;
279 #if PetscDefined(HAVE_MAGMA)
280 24 fn->ops->evaluatefunctionmatcuda[1] = FNEvaluateFunctionMat_Invsqrt_DBP_CUDAm;
281 24 fn->ops->evaluatefunctionmatcuda[3] = FNEvaluateFunctionMat_Invsqrt_Sadeghi_CUDAm;
282 #endif /* PETSC_HAVE_MAGMA */
283 #endif /* PETSC_HAVE_CUDA */
284 72 fn->ops->evaluatefunctionmatvec[0] = FNEvaluateFunctionMatVec_Invsqrt_Schur;
285 72 fn->ops->view = FNView_Invsqrt;
286
7/14
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
 72 PetscFunctionReturn(PETSC_SUCCESS);
287 }
288