github.com/jingcheng-WU/gonum@v0.9.1-0.20210323123734-f1a2a11a8f7b/lapack/testlapack/dpbtf2.go (about)

     1  // Copyright ©2017 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.
     4  
     5  package testlapack
     6  
     7  import (
     8  	"fmt"
     9  	"testing"
    10  
    11  	"golang.org/x/exp/rand"
    12  
    13  	"github.com/jingcheng-WU/gonum/blas"
    14  )
    15  
    16  type Dpbtf2er interface {
    17  	Dpbtf2(uplo blas.Uplo, n, kd int, ab []float64, ldab int) (ok bool)
    18  }
    19  
    20  // Dpbtf2Test tests Dpbtf2 on random symmetric positive definite band matrices
    21  // by checking that the Cholesky factors multiply back to the original matrix.
    22  func Dpbtf2Test(t *testing.T, impl Dpbtf2er) {
    23  	// TODO(vladimir-ch): include expected-failure test case.
    24  	rnd := rand.New(rand.NewSource(1))
    25  	for _, n := range []int{0, 1, 2, 3, 4, 5, 10, 20} {
    26  		for _, kd := range []int{0, (n + 1) / 4, (3*n - 1) / 4, (5*n + 1) / 4} {
    27  			for _, uplo := range []blas.Uplo{blas.Upper, blas.Lower} {
    28  				for _, ldab := range []int{kd + 1, kd + 1 + 7} {
    29  					dpbtf2Test(t, impl, rnd, uplo, n, kd, ldab)
    30  				}
    31  			}
    32  		}
    33  	}
    34  }
    35  
    36  func dpbtf2Test(t *testing.T, impl Dpbtf2er, rnd *rand.Rand, uplo blas.Uplo, n, kd int, ldab int) {
    37  	const tol = 1e-12
    38  
    39  	name := fmt.Sprintf("uplo=%v,n=%v,kd=%v,ldab=%v", string(uplo), n, kd, ldab)
    40  
    41  	// Generate a random symmetric positive definite band matrix.
    42  	ab := randSymBand(uplo, n, kd, ldab, rnd)
    43  
    44  	// Compute the Cholesky decomposition of A.
    45  	abFac := make([]float64, len(ab))
    46  	copy(abFac, ab)
    47  	ok := impl.Dpbtf2(uplo, n, kd, abFac, ldab)
    48  	if !ok {
    49  		t.Fatalf("%v: bad test matrix, Dpbtf2 failed", name)
    50  	}
    51  
    52  	// Reconstruct an symmetric band matrix from the Uᵀ*U or L*Lᵀ factorization, overwriting abFac.
    53  	dsbmm(uplo, n, kd, abFac, ldab)
    54  
    55  	// Compute and check the max-norm distance between the reconstructed and original matrix A.
    56  	dist := distSymBand(uplo, n, kd, abFac, ldab, ab, ldab)
    57  	if dist > tol {
    58  		t.Errorf("%v: unexpected result, diff=%v", name, dist)
    59  	}
    60  }