1 // Copyright ©2015 The Gonum Authors. All rights reserved.
2 // Use of this source code is governed by a BSD-style
3 // license that can be found in the LICENSE file.
10 "golang.org/x/exp/rand"
12 "gonum.org/v1/gonum/blas"
13 "gonum.org/v1/gonum/blas/blas64"
14 "gonum.org/v1/gonum/floats"
17 type Dpotrfer interface {
18 Dpotrf(ul blas.Uplo, n int, a []float64, lda int) (ok bool)
21 func DpotrfTest(t *testing.T, impl Dpotrfer) {
23 rnd := rand.New(rand.NewSource(1))
24 bi := blas64.Implementation()
25 for _, uplo := range []blas.Uplo{blas.Upper, blas.Lower} {
26 for tc, test := range []struct {
53 // Random diagonal matrix D with positive entries.
54 d := make([]float64, n)
55 Dlatm1(d, 4, 10000, false, 1, rnd)
57 // Construct a positive definite matrix A as
59 // where U is a random orthogonal matrix.
64 a := make([]float64, n*lda)
65 Dlagsy(n, 0, d, a, lda, rnd, make([]float64, 2*n))
67 aCopy := make([]float64, len(a))
70 ok := impl.Dpotrf(uplo, n, a, lda)
72 t.Errorf("Case %v: unexpected failure for positive definite matrix", tc)
78 for i := 0; i < n; i++ {
79 for j := 0; j < i; j++ {
84 for i := 0; i < n; i++ {
85 for j := i + 1; j < n; j++ {
93 ans := make([]float64, len(a))
97 bi.Dsyrk(uplo, blas.Trans, n, n, 1, a, lda, 0, ans, lda)
100 bi.Dsyrk(uplo, blas.NoTrans, n, n, 1, a, lda, 0, ans, lda)
106 for i := 0; i < n; i++ {
107 for j := i; j < n; j++ {
108 if !floats.EqualWithinAbsOrRel(ans[i*lda+j], aCopy[i*lda+j], tol, tol) {
114 for i := 0; i < n; i++ {
115 for j := 0; j <= i; j++ {
116 if !floats.EqualWithinAbsOrRel(ans[i*lda+j], aCopy[i*lda+j], tol, tol) {
123 t.Errorf("Case %v (uplo=%v,n=%v,lda=%v): unexpected result", tc, uplo, n, lda)
126 // Make one element of D negative so that A is not
127 // positive definite, and check that Dpotrf fails.
129 Dlagsy(n, 0, d, a, lda, rnd, make([]float64, 2*n))
130 ok = impl.Dpotrf(uplo, n, a, lda)
132 t.Errorf("Case %v: unexpected success for not positive definite matrix", tc)