github.com/gonum/lapack@v0.0.0-20181123203213-e4cdc5a0bff9/testlapack/dlansy.go (about)

     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.
     4  
     5  package testlapack
     6  
     7  import (
     8  	"math"
     9  	"math/rand"
    10  	"testing"
    11  
    12  	"github.com/gonum/blas"
    13  	"github.com/gonum/lapack"
    14  )
    15  
    16  type Dlansyer interface {
    17  	Dlanger
    18  	Dlansy(norm lapack.MatrixNorm, uplo blas.Uplo, n int, a []float64, lda int, work []float64) float64
    19  }
    20  
    21  func DlansyTest(t *testing.T, impl Dlansyer) {
    22  	rnd := rand.New(rand.NewSource(1))
    23  	for _, norm := range []lapack.MatrixNorm{lapack.MaxAbs, lapack.MaxColumnSum, lapack.MaxRowSum, lapack.NormFrob} {
    24  		for _, uplo := range []blas.Uplo{blas.Lower, blas.Upper} {
    25  			for _, test := range []struct {
    26  				n, lda int
    27  			}{
    28  				{1, 0},
    29  				{3, 0},
    30  
    31  				{1, 10},
    32  				{3, 10},
    33  			} {
    34  				for trial := 0; trial < 100; trial++ {
    35  					n := test.n
    36  					lda := test.lda
    37  					if lda == 0 {
    38  						lda = n
    39  					}
    40  					a := make([]float64, lda*n)
    41  					if trial == 0 {
    42  						for i := range a {
    43  							a[i] = float64(i)
    44  						}
    45  					} else {
    46  						for i := range a {
    47  							a[i] = rnd.NormFloat64()
    48  						}
    49  					}
    50  
    51  					aDense := make([]float64, n*n)
    52  					if uplo == blas.Upper {
    53  						for i := 0; i < n; i++ {
    54  							for j := i; j < n; j++ {
    55  								v := a[i*lda+j]
    56  								aDense[i*n+j] = v
    57  								aDense[j*n+i] = v
    58  							}
    59  						}
    60  					} else {
    61  						for i := 0; i < n; i++ {
    62  							for j := 0; j <= i; j++ {
    63  								v := a[i*lda+j]
    64  								aDense[i*n+j] = v
    65  								aDense[j*n+i] = v
    66  							}
    67  						}
    68  					}
    69  					work := make([]float64, n)
    70  					got := impl.Dlansy(norm, uplo, n, a, lda, work)
    71  					want := impl.Dlange(norm, n, n, aDense, n, work)
    72  					if math.Abs(want-got) > 1e-14 {
    73  						t.Errorf("Norm mismatch. norm = %c, upper = %v, n = %v, lda = %v, want %v, got %v.",
    74  							norm, uplo == blas.Upper, n, lda, got, want)
    75  					}
    76  				}
    77  			}
    78  		}
    79  	}
    80  }