Skip to content

Commit 35df9b8

Browse files
andreasnoacktkelman
authored andcommitted
Remove redundant uplo argument in chol family. When using Hermitian (#17909)
and Symmetric, it is redundant to also have an uplo argument and occasionally, it also gave the wrong result. This can therefore be considered a bugfix.
1 parent 8615732 commit 35df9b8

File tree

4 files changed

+66
-56
lines changed

4 files changed

+66
-56
lines changed

base/deprecated.jl

+3
Original file line numberDiff line numberDiff line change
@@ -775,6 +775,9 @@ function transpose(x)
775775
return x
776776
end
777777

778+
@deprecate cholfact!(A::Base.LinAlg.HermOrSym, uplo::Symbol, ::Type{Val{false}}) cholfact!(A, Val{false})
779+
@deprecate cholfact!(A::Base.LinAlg.HermOrSym, uplo::Symbol = :U) cholfact!(A)
780+
778781
# During the 0.5 development cycle, do not add any deprecations below this line
779782
# To be deprecated in 0.6
780783

base/linalg/cholesky.jl

+52-45
Original file line numberDiff line numberDiff line change
@@ -121,8 +121,10 @@ non_hermitian_error(f) = throw(ArgumentError("matrix is not symmetric/" *
121121

122122
# chol!. Destructive methods for computing Cholesky factor of real symmetric or Hermitian
123123
# matrix
124-
chol!(A::Hermitian) = _chol!(A.data, UpperTriangular)
125-
chol!{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}) = _chol!(A.data, UpperTriangular)
124+
chol!(A::Hermitian) =
125+
_chol!(A.uplo == 'U' ? A.data : LinAlg.copytri!(A.data, 'L', true), UpperTriangular)
126+
chol!{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}) =
127+
_chol!(A.uplo == 'U' ? A.data : LinAlg.copytri!(A.data, 'L', true), UpperTriangular)
126128
function chol!(A::StridedMatrix)
127129
ishermitian(A) || non_hermitian_error("chol!")
128130
return _chol!(A, UpperTriangular)
@@ -135,14 +137,22 @@ end
135137
function chol(A::Hermitian)
136138
T = promote_type(typeof(chol(one(eltype(A)))), Float32)
137139
AA = similar(A, T, size(A))
138-
copy!(AA, A.data)
139-
chol!(Hermitian(AA))
140+
if A.uplo == 'U'
141+
copy!(AA, A.data)
142+
else
143+
Base.ccopy!(AA, A.data)
144+
end
145+
chol!(Hermitian(AA, :U))
140146
end
141147
function chol{T<:Real,S<:AbstractMatrix}(A::Symmetric{T,S})
142148
TT = promote_type(typeof(chol(one(T))), Float32)
143149
AA = similar(A, TT, size(A))
144-
copy!(AA, A.data)
145-
chol!(Hermitian(AA))
150+
if A.uplo == 'U'
151+
copy!(AA, A.data)
152+
else
153+
Base.ccopy!(AA, A.data)
154+
end
155+
chol!(Hermitian(AA, :U))
146156
end
147157

148158
## for StridedMatrices, check that matrix is symmetric/Hermitian
@@ -170,15 +180,15 @@ chol(x::Number, args...) = _chol!(x, nothing)
170180
# cholfact!. Destructive methods for computing Cholesky factorization of real symmetric
171181
# or Hermitian matrix
172182
## No pivoting
173-
function cholfact!(A::Hermitian, uplo::Symbol, ::Type{Val{false}})
174-
if uplo == :U
183+
function cholfact!(A::Hermitian, ::Type{Val{false}})
184+
if A.uplo == :U
175185
Cholesky(_chol!(A.data, UpperTriangular).data, 'U')
176186
else
177187
Cholesky(_chol!(A.data, LowerTriangular).data, 'L')
178188
end
179189
end
180-
function cholfact!{T<:Real,S}(A::Symmetric{T,S}, uplo::Symbol, ::Type{Val{false}})
181-
if uplo == :U
190+
function cholfact!{T<:Real,S}(A::Symmetric{T,S}, ::Type{Val{false}})
191+
if A.uplo == :U
182192
Cholesky(_chol!(A.data, UpperTriangular).data, 'U')
183193
else
184194
Cholesky(_chol!(A.data, LowerTriangular).data, 'L')
@@ -187,7 +197,7 @@ end
187197

188198
### for StridedMatrices, check that matrix is symmetric/Hermitian
189199
"""
190-
cholfact!(A, uplo::Symbol, Val{false}) -> Cholesky
200+
cholfact!(A, [uplo::Symbol,] Val{false}) -> Cholesky
191201
192202
The same as `cholfact`, but saves space by overwriting the input `A`, instead
193203
of creating a copy. An `InexactError` exception is thrown if the factorisation
@@ -196,37 +206,36 @@ integer types.
196206
"""
197207
function cholfact!(A::StridedMatrix, uplo::Symbol, ::Type{Val{false}})
198208
ishermitian(A) || non_hermitian_error("cholfact!")
199-
return cholfact!(Hermitian(A), uplo, Val{false})
209+
return cholfact!(Hermitian(A, uplo), Val{false})
200210
end
201211

202212
### Default to no pivoting (and storing of upper factor) when not explicit
203-
cholfact!(A::Hermitian, uplo::Symbol = :U) = cholfact!(A, uplo, Val{false})
204-
cholfact!{T<:Real,S}(A::Symmetric{T,S}, uplo::Symbol = :U) = cholfact!(A, uplo, Val{false})
213+
cholfact!(A::Hermitian) = cholfact!(A, Val{false})
214+
cholfact!{T<:Real,S}(A::Symmetric{T,S}) = cholfact!(A, Val{false})
205215
#### for StridedMatrices, check that matrix is symmetric/Hermitian
206216
function cholfact!(A::StridedMatrix, uplo::Symbol = :U)
207217
ishermitian(A) || non_hermitian_error("cholfact!")
208-
return cholfact!(Hermitian(A), uplo)
218+
return cholfact!(Hermitian(A, uplo))
209219
end
210220

211221

212222
## With pivoting
213223
### BLAS/LAPACK element types
214224
function cholfact!{T<:BlasReal,S<:StridedMatrix}(A::RealHermSymComplexHerm{T,S},
215-
uplo::Symbol, ::Type{Val{true}}; tol = 0.0)
216-
uplochar = char_uplo(uplo)
217-
AA, piv, rank, info = LAPACK.pstrf!(uplochar, A.data, tol)
218-
return CholeskyPivoted{eltype(AA),typeof(AA)}(AA, uplochar, piv, rank, tol, info)
225+
::Type{Val{true}}; tol = 0.0)
226+
AA, piv, rank, info = LAPACK.pstrf!(A.uplo, A.data, tol)
227+
return CholeskyPivoted{eltype(AA),typeof(AA)}(AA, A.uplo, piv, rank, tol, info)
219228
end
220229

221230
### Non BLAS/LAPACK element types (generic). Since generic fallback for pivoted Cholesky
222231
### is not implemented yet we throw an error
223-
cholfact!{T<:Real,S}(A::RealHermSymComplexHerm{T,S}, uplo::Symbol, ::Type{Val{true}};
232+
cholfact!{T<:Real,S}(A::RealHermSymComplexHerm{T,S}, ::Type{Val{true}};
224233
tol = 0.0) =
225234
throw(ArgumentError("generic pivoted Cholesky factorization is not implemented yet"))
226235

227236
### for StridedMatrices, check that matrix is symmetric/Hermitian
228237
"""
229-
cholfact!(A, uplo::Symbol, Val{true}; tol = 0.0) -> CholeskyPivoted
238+
cholfact!(A, [uplo::Symbol,] Val{true}; tol = 0.0) -> CholeskyPivoted
230239
231240
The same as `cholfact`, but saves space by overwriting the input `A`, instead
232241
of creating a copy. An `InexactError` exception is thrown if the
@@ -235,63 +244,61 @@ e.g. for integer types.
235244
"""
236245
function cholfact!(A::StridedMatrix, uplo::Symbol, ::Type{Val{true}}; tol = 0.0)
237246
ishermitian(A) || non_hermitian_error("cholfact!")
238-
return cholfact!(Hermitian(A), uplo, Val{true}; tol = tol)
247+
return cholfact!(Hermitian(A, uplo), Val{true}; tol = tol)
239248
end
240249

241-
242-
243250
# cholfact. Non-destructive methods for computing Cholesky factorization of real symmetric
244251
# or Hermitian matrix
245252
## No pivoting
246-
cholfact(A::Hermitian, uplo::Symbol, ::Type{Val{false}}) =
247-
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)), uplo, Val{false})
248-
cholfact{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}, uplo::Symbol, ::Type{Val{false}}) =
249-
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)), uplo, Val{false})
253+
cholfact(A::Hermitian, ::Type{Val{false}}) =
254+
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)), Val{false})
255+
cholfact{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}, ::Type{Val{false}}) =
256+
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)), Val{false})
250257

251258
### for StridedMatrices, check that matrix is symmetric/Hermitian
252259
"""
253-
cholfact(A, uplo::Symbol, Val{false}) -> Cholesky
260+
cholfact(A, [uplo::Symbol,] Val{false}) -> Cholesky
254261
255262
Compute the Cholesky factorization of a dense symmetric positive definite matrix `A`
256-
and return a `Cholesky` factorization.
257-
The `uplo` argument may be `:L` for using the lower part or `:U` for the upper part of `A`.
263+
and return a `Cholesky` factorization. The matrix `A` can either be a `Symmetric` or `Hermitian`
264+
`StridedMatrix` or a *perfectly* symmetric or Hermitian `StridedMatrix`. In the latter case,
265+
the optional argument `uplo` may be `:L` for using the lower part or `:U` for the upper part of `A`.
258266
The default is to use `:U`.
259267
The triangular Cholesky factor can be obtained from the factorization `F` with: `F[:L]` and `F[:U]`.
260268
The following functions are available for `Cholesky` objects: `size`, `\\`, `inv`, `det`.
261269
A `PosDefException` exception is thrown in case the matrix is not positive definite.
262270
"""
263271
function cholfact(A::StridedMatrix, uplo::Symbol, ::Type{Val{false}})
264272
ishermitian(A) || non_hermitian_error("cholfact")
265-
return cholfact(Hermitian(A), uplo, Val{false})
273+
return cholfact(Hermitian(A, uplo), Val{false})
266274
end
267275

268276
### Default to no pivoting (and storing of upper factor) when not explicit
269-
cholfact(A::Hermitian, uplo::Symbol = :U) = cholfact(A, uplo, Val{false})
270-
cholfact{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}, uplo::Symbol = :U) =
271-
cholfact(A, uplo, Val{false})
277+
cholfact(A::Hermitian) = cholfact(A, Val{false})
278+
cholfact{T<:Real,S<:StridedMatrix}(A::Symmetric{T,S}) = cholfact(A, Val{false})
272279
#### for StridedMatrices, check that matrix is symmetric/Hermitian
273280
function cholfact(A::StridedMatrix, uplo::Symbol = :U)
274281
ishermitian(A) || non_hermitian_error("cholfact")
275-
return cholfact(Hermitian(A), uplo)
282+
return cholfact(Hermitian(A, uplo))
276283
end
277284

278285

279286
## With pivoting
280-
cholfact(A::Hermitian, uplo::Symbol, ::Type{Val{true}}; tol = 0.0) =
287+
cholfact(A::Hermitian, ::Type{Val{true}}; tol = 0.0) =
281288
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)),
282-
uplo, Val{true}; tol = tol)
283-
cholfact{T<:Real,S<:StridedMatrix}(A::RealHermSymComplexHerm{T,S}, uplo::Symbol,
284-
::Type{Val{true}}; tol = 0.0) =
289+
Val{true}; tol = tol)
290+
cholfact{T<:Real,S<:StridedMatrix}(A::RealHermSymComplexHerm{T,S}, ::Type{Val{true}}; tol = 0.0) =
285291
cholfact!(copy_oftype(A, promote_type(typeof(chol(one(eltype(A)))),Float32)),
286-
uplo, Val{true}; tol = tol)
292+
Val{true}; tol = tol)
287293

288294
### for StridedMatrices, check that matrix is symmetric/Hermitian
289295
"""
290-
cholfact(A, uplo::Symbol, Val{true}; tol = 0.0) -> CholeskyPivoted
296+
cholfact(A, [uplo::Symbol,] Val{true}; tol = 0.0) -> CholeskyPivoted
291297
292298
Compute the pivoted Cholesky factorization of a dense symmetric positive semi-definite matrix `A`
293-
and return a `CholeskyPivoted` factorization.
294-
The `uplo` argument may be `:L` for using the lower part or `:U` for the upper part of `A`.
299+
and return a `CholeskyPivoted` factorization. The matrix `A` can either be a `Symmetric` or `Hermitian`
300+
`StridedMatrix` or a *perfectly* symmetric or Hermitian `StridedMatrix`. In the latter case,
301+
the optional argument `uplo` may be `:L` for using the lower part or `:U` for the upper part of `A`.
295302
The default is to use `:U`.
296303
The triangular Cholesky factor can be obtained from the factorization `F` with: `F[:L]` and `F[:U]`.
297304
The following functions are available for `PivotedCholesky` objects: `size`, `\\`, `inv`, `det`, and `rank`.
@@ -300,7 +307,7 @@ For negative values, the tolerance is the machine precision.
300307
"""
301308
function cholfact(A::StridedMatrix, uplo::Symbol, ::Type{Val{true}}; tol = 0.0)
302309
ishermitian(A) || non_hermitian_error("cholfact")
303-
return cholfact(Hermitian(A), uplo, Val{true}; tol = tol)
310+
return cholfact(Hermitian(A, uplo), Val{true}; tol = tol)
304311
end
305312

306313
## Number

doc/stdlib/linalg.rst

+6-6
Original file line numberDiff line numberDiff line change
@@ -308,17 +308,17 @@ Linear algebra functions in Julia are largely implemented by calling functions f
308308
309309
Compute the square root of a non-negative number ``x``\ .
310310

311-
.. function:: cholfact(A, uplo::Symbol, Val{false}) -> Cholesky
311+
.. function:: cholfact(A, [uplo::Symbol,] Val{false}) -> Cholesky
312312

313313
.. Docstring generated from Julia source
314314
315-
Compute the Cholesky factorization of a dense symmetric positive definite matrix ``A`` and return a ``Cholesky`` factorization. The ``uplo`` argument may be ``:L`` for using the lower part or ``:U`` for the upper part of ``A``\ . The default is to use ``:U``\ . The triangular Cholesky factor can be obtained from the factorization ``F`` with: ``F[:L]`` and ``F[:U]``\ . The following functions are available for ``Cholesky`` objects: ``size``\ , ``\``\ , ``inv``\ , ``det``\ . A ``PosDefException`` exception is thrown in case the matrix is not positive definite.
315+
Compute the Cholesky factorization of a dense symmetric positive definite matrix ``A`` and return a ``Cholesky`` factorization. The matrix ``A`` can either be a ``Symmetric`` or ``Hermitian`` ``StridedMatrix`` or a *perfectly* symmetric or Hermitian ``StridedMatrix``\ . In the latter case, the optional argument ``uplo`` may be ``:L`` for using the lower part or ``:U`` for the upper part of ``A``\ . The default is to use ``:U``\ . The triangular Cholesky factor can be obtained from the factorization ``F`` with: ``F[:L]`` and ``F[:U]``\ . The following functions are available for ``Cholesky`` objects: ``size``\ , ``\``\ , ``inv``\ , ``det``\ . A ``PosDefException`` exception is thrown in case the matrix is not positive definite.
316316

317-
.. function:: cholfact(A, uplo::Symbol, Val{true}; tol = 0.0) -> CholeskyPivoted
317+
.. function:: cholfact(A, [uplo::Symbol,] Val{true}; tol = 0.0) -> CholeskyPivoted
318318

319319
.. Docstring generated from Julia source
320320
321-
Compute the pivoted Cholesky factorization of a dense symmetric positive semi-definite matrix ``A`` and return a ``CholeskyPivoted`` factorization. The ``uplo`` argument may be ``:L`` for using the lower part or ``:U`` for the upper part of ``A``\ . The default is to use ``:U``\ . The triangular Cholesky factor can be obtained from the factorization ``F`` with: ``F[:L]`` and ``F[:U]``\ . The following functions are available for ``PivotedCholesky`` objects: ``size``\ , ``\``\ , ``inv``\ , ``det``\ , and ``rank``\ . The argument ``tol`` determines the tolerance for determining the rank. For negative values, the tolerance is the machine precision.
321+
Compute the pivoted Cholesky factorization of a dense symmetric positive semi-definite matrix ``A`` and return a ``CholeskyPivoted`` factorization. The matrix ``A`` can either be a ``Symmetric`` or ``Hermitian`` ``StridedMatrix`` or a *perfectly* symmetric or Hermitian ``StridedMatrix``\ . In the latter case, the optional argument ``uplo`` may be ``:L`` for using the lower part or ``:U`` for the upper part of ``A``\ . The default is to use ``:U``\ . The triangular Cholesky factor can be obtained from the factorization ``F`` with: ``F[:L]`` and ``F[:U]``\ . The following functions are available for ``PivotedCholesky`` objects: ``size``\ , ``\``\ , ``inv``\ , ``det``\ , and ``rank``\ . The argument ``tol`` determines the tolerance for determining the rank. For negative values, the tolerance is the machine precision.
322322

323323
.. function:: cholfact(A; shift = 0.0, perm = Int[]) -> CHOLMOD.Factor
324324

@@ -344,13 +344,13 @@ Linear algebra functions in Julia are largely implemented by calling functions f
344344
This method uses the CHOLMOD library from SuiteSparse, which only supports doubles or complex doubles. Input matrices not of those element types will be converted to ``SparseMatrixCSC{Float64}`` or ``SparseMatrixCSC{Complex128}`` as appropriate.
345345

346346

347-
.. function:: cholfact!(A, uplo::Symbol, Val{false}) -> Cholesky
347+
.. function:: cholfact!(A, [uplo::Symbol,] Val{false}) -> Cholesky
348348

349349
.. Docstring generated from Julia source
350350
351351
The same as ``cholfact``\ , but saves space by overwriting the input ``A``\ , instead of creating a copy. An ``InexactError`` exception is thrown if the factorisation produces a number not representable by the element type of ``A``\ , e.g. for integer types.
352352

353-
.. function:: cholfact!(A, uplo::Symbol, Val{true}; tol = 0.0) -> CholeskyPivoted
353+
.. function:: cholfact!(A, [uplo::Symbol,] Val{true}; tol = 0.0) -> CholeskyPivoted
354354

355355
.. Docstring generated from Julia source
356356

test/linalg/cholesky.jl

+5-5
Original file line numberDiff line numberDiff line change
@@ -74,15 +74,15 @@ for eltya in (Float32, Float64, Complex64, Complex128, BigFloat, Int)
7474
@test full(lapd) apd
7575
l = lapd[:L]
7676
@test l*l' apd
77-
@test triu(capd.factors) lapd[:U]
78-
@test tril(lapd.factors) capd[:L]
77+
@test capd[:U] lapd[:U]
78+
@test lapd[:L] capd[:L]
7979
if eltya <: Real
8080
capds = cholfact(apds)
81-
lapds = cholfact(apds, :L)
81+
lapds = cholfact(full(apds), :L)
8282
ls = lapds[:L]
8383
@test ls*ls' apd
84-
@test triu(capds.factors) lapds[:U]
85-
@test tril(lapds.factors) capds[:L]
84+
@test capds[:U] lapds[:U]
85+
@test lapds[:L] capds[:L]
8686
end
8787

8888
#pivoted upper Cholesky

0 commit comments

Comments
 (0)