gonum.org/v1/gonum@v0.14.0/lapack/testlapack/matgen_test.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  	"testing"
     9  
    10  	"golang.org/x/exp/rand"
    11  
    12  	"gonum.org/v1/gonum/blas/blas64"
    13  )
    14  
    15  func TestDlagsy(t *testing.T) {
    16  	const tol = 1e-14
    17  	rnd := rand.New(rand.NewSource(1))
    18  	for _, n := range []int{0, 1, 2, 3, 4, 5, 10, 50} {
    19  		for _, lda := range []int{0, 2*n + 1} {
    20  			if lda == 0 {
    21  				lda = max(1, n)
    22  			}
    23  			// D is the identity matrix I.
    24  			d := make([]float64, n)
    25  			for i := range d {
    26  				d[i] = 1
    27  			}
    28  			// Allocate an n×n symmetric matrix A and fill it with NaNs.
    29  			a := nanSlice(n * lda)
    30  			work := make([]float64, 2*n)
    31  			// Compute A = U * D * Uᵀ where U is a random orthogonal matrix.
    32  			Dlagsy(n, 0, d, a, lda, rnd, work)
    33  			// A should be the identity matrix because
    34  			//  A = U * D * Uᵀ = U * I * Uᵀ = U * Uᵀ = I.
    35  			dist := distFromIdentity(n, a, lda)
    36  			if dist > tol {
    37  				t.Errorf("Case n=%v,lda=%v: |A-I|=%v is too large", n, lda, dist)
    38  			}
    39  		}
    40  	}
    41  }
    42  
    43  func TestDlagge(t *testing.T) {
    44  	const tol = 1e-14
    45  	rnd := rand.New(rand.NewSource(1))
    46  	for _, n := range []int{0, 1, 2, 3, 4, 5, 10, 50} {
    47  		for _, lda := range []int{0, 2*n + 1} {
    48  			if lda == 0 {
    49  				lda = max(1, n)
    50  			}
    51  			d := make([]float64, n)
    52  			for i := range d {
    53  				d[i] = 1
    54  			}
    55  			a := blas64.General{
    56  				Rows:   n,
    57  				Cols:   n,
    58  				Stride: lda,
    59  				Data:   nanSlice(n * lda),
    60  			}
    61  			work := make([]float64, a.Rows+a.Cols)
    62  
    63  			Dlagge(a.Rows, a.Cols, 0, 0, d, a.Data, a.Stride, rnd, work)
    64  
    65  			if resid := residualOrthogonal(a, false); resid > tol {
    66  				t.Errorf("Case n=%v,lda=%v: unexpected result", n, lda)
    67  			}
    68  		}
    69  	}
    70  }
    71  
    72  func TestRandomOrthogonal(t *testing.T) {
    73  	const tol = 1e-14
    74  	rnd := rand.New(rand.NewSource(1))
    75  	for n := 1; n <= 20; n++ {
    76  		q := randomOrthogonal(n, rnd)
    77  		if resid := residualOrthogonal(q, false); resid > tol {
    78  			t.Errorf("Case n=%v: Q not orthogonal; resid=%v, want<=%v", n, resid, tol)
    79  		}
    80  	}
    81  }