gonum.org/v1/gonum@v0.14.0/lapack/testlapack/dlascl.go (about)

     1  // Copyright ©2016 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  	"math"
    10  	"testing"
    11  
    12  	"golang.org/x/exp/rand"
    13  
    14  	"gonum.org/v1/gonum/lapack"
    15  )
    16  
    17  type Dlascler interface {
    18  	Dlascl(kind lapack.MatrixType, kl, ku int, cfrom, cto float64, m, n int, a []float64, lda int)
    19  }
    20  
    21  func DlasclTest(t *testing.T, impl Dlascler) {
    22  	const tol = 1e-15
    23  
    24  	rnd := rand.New(rand.NewSource(1))
    25  	for ti, test := range []struct {
    26  		m, n int
    27  	}{
    28  		{0, 0},
    29  		{1, 1},
    30  		{1, 10},
    31  		{10, 1},
    32  		{2, 2},
    33  		{2, 11},
    34  		{11, 2},
    35  		{3, 3},
    36  		{3, 11},
    37  		{11, 3},
    38  		{11, 11},
    39  		{11, 100},
    40  		{100, 11},
    41  	} {
    42  		m := test.m
    43  		n := test.n
    44  		for _, extra := range []int{0, 11} {
    45  			for _, kind := range []lapack.MatrixType{lapack.General, lapack.UpperTri, lapack.LowerTri} {
    46  				a := randomGeneral(m, n, n+extra, rnd)
    47  				aCopy := cloneGeneral(a)
    48  				cfrom := rnd.NormFloat64()
    49  				cto := rnd.NormFloat64()
    50  				scale := cto / cfrom
    51  
    52  				impl.Dlascl(kind, -1, -1, cfrom, cto, m, n, a.Data, a.Stride)
    53  
    54  				prefix := fmt.Sprintf("Case #%v: kind=%v,m=%v,n=%v,extra=%v", ti, kind, m, n, extra)
    55  				if !generalOutsideAllNaN(a) {
    56  					t.Errorf("%v: out-of-range write to A", prefix)
    57  				}
    58  				switch kind {
    59  				case lapack.UpperTri:
    60  					var mod bool
    61  				loopLower:
    62  					for i := 0; i < m; i++ {
    63  						for j := 0; j < min(i, n); j++ {
    64  							if a.Data[i*a.Stride+j] != aCopy.Data[i*aCopy.Stride+j] {
    65  								mod = true
    66  								break loopLower
    67  							}
    68  						}
    69  					}
    70  					if mod {
    71  						t.Errorf("%v: unexpected modification in lower triangle of A", prefix)
    72  					}
    73  				case lapack.LowerTri:
    74  					var mod bool
    75  				loopUpper:
    76  					for i := 0; i < m; i++ {
    77  						for j := i + 1; j < n; j++ {
    78  							if a.Data[i*a.Stride+j] != aCopy.Data[i*aCopy.Stride+j] {
    79  								mod = true
    80  								break loopUpper
    81  							}
    82  						}
    83  					}
    84  					if mod {
    85  						t.Errorf("%v: unexpected modification in upper triangle of A", prefix)
    86  					}
    87  				}
    88  
    89  				var resid float64
    90  				switch kind {
    91  				case lapack.General:
    92  					for i := 0; i < m; i++ {
    93  						for j := 0; j < n; j++ {
    94  							want := scale * aCopy.Data[i*aCopy.Stride+j]
    95  							got := a.Data[i*a.Stride+j]
    96  							resid = math.Max(resid, math.Abs(want-got))
    97  						}
    98  					}
    99  				case lapack.UpperTri:
   100  					for i := 0; i < m; i++ {
   101  						for j := i; j < n; j++ {
   102  							want := scale * aCopy.Data[i*aCopy.Stride+j]
   103  							got := a.Data[i*a.Stride+j]
   104  							resid = math.Max(resid, math.Abs(want-got))
   105  						}
   106  					}
   107  				case lapack.LowerTri:
   108  					for i := 0; i < m; i++ {
   109  						for j := 0; j <= min(i, n-1); j++ {
   110  							want := scale * aCopy.Data[i*aCopy.Stride+j]
   111  							got := a.Data[i*a.Stride+j]
   112  							resid = math.Max(resid, math.Abs(want-got))
   113  						}
   114  					}
   115  				}
   116  				if resid > tol*float64(max(m, n)) {
   117  					t.Errorf("%v: unexpected result; residual=%v, want<=%v", prefix, resid, tol*float64(max(m, n)))
   118  				}
   119  			}
   120  		}
   121  	}
   122  }