vendor/gonum.org/v1/gonum/mat/lq_test.go

   1 // Copyright ©2013 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 mat
   6
   7 import (
   8         "testing"
   9
  10         "golang.org/x/exp/rand"
  11 )
  12
  13 func TestLQ(t *testing.T) {
  14         for _, test := range []struct {
  15                 m, n int
  16         }{
  17                 {5, 5},
  18                 {5, 10},
  19         } {
  20                 m := test.m
  21                 n := test.n
  22                 a := NewDense(m, n, nil)
  23                 for i := 0; i < m; i++ {
  24                         for j := 0; j < n; j++ {
  25                                 a.Set(i, j, rand.NormFloat64())
  26                         }
  27                 }
  28                 var want Dense
  29                 want.Clone(a)
  30
  31                 var lq LQ
  32                 lq.Factorize(a)
  33                 q := lq.QTo(nil)
  34
  35                 if !isOrthonormal(q, 1e-10) {
  36                         t.Errorf("Q is not orthonormal: m = %v, n = %v", m, n)
  37                 }
  38
  39                 l := lq.LTo(nil)
  40
  41                 var got Dense
  42                 got.Mul(l, q)
  43                 if !EqualApprox(&got, &want, 1e-12) {
  44                         t.Errorf("LQ does not equal original matrix. \nWant: %v\nGot: %v", want, got)
  45                 }
  46         }
  47 }
  48
  49 func TestSolveLQ(t *testing.T) {
  50         for _, trans := range []bool{false, true} {
  51                 for _, test := range []struct {
  52                         m, n, bc int
  53                 }{
  54                         {5, 5, 1},
  55                         {5, 10, 1},
  56                         {5, 5, 3},
  57                         {5, 10, 3},
  58                 } {
  59                         m := test.m
  60                         n := test.n
  61                         bc := test.bc
  62                         a := NewDense(m, n, nil)
  63                         for i := 0; i < m; i++ {
  64                                 for j := 0; j < n; j++ {
  65                                         a.Set(i, j, rand.Float64())
  66                                 }
  67                         }
  68                         br := m
  69                         if trans {
  70                                 br = n
  71                         }
  72                         b := NewDense(br, bc, nil)
  73                         for i := 0; i < br; i++ {
  74                                 for j := 0; j < bc; j++ {
  75                                         b.Set(i, j, rand.Float64())
  76                                 }
  77                         }
  78                         var x Dense
  79                         lq := &LQ{}
  80                         lq.Factorize(a)
  81                         lq.Solve(&x, trans, b)
  82
  83                         // Test that the normal equations hold.
  84                         // A^T * A * x = A^T * b if !trans
  85                         // A * A^T * x = A * b if trans
  86                         var lhs Dense
  87                         var rhs Dense
  88                         if trans {
  89                                 var tmp Dense
  90                                 tmp.Mul(a, a.T())
  91                                 lhs.Mul(&tmp, &x)
  92                                 rhs.Mul(a, b)
  93                         } else {
  94                                 var tmp Dense
  95                                 tmp.Mul(a.T(), a)
  96                                 lhs.Mul(&tmp, &x)
  97                                 rhs.Mul(a.T(), b)
  98                         }
  99                         if !EqualApprox(&lhs, &rhs, 1e-10) {
 100                                 t.Errorf("Normal equations do not hold.\nLHS: %v\n, RHS: %v\n", lhs, rhs)
 101                         }
 102                 }
 103         }
 104         // TODO(btracey): Add in testOneInput when it exists.
 105 }
 106
 107 func TestSolveLQVec(t *testing.T) {
 108         for _, trans := range []bool{false, true} {
 109                 for _, test := range []struct {
 110                         m, n int
 111                 }{
 112                         {5, 5},
 113                         {5, 10},
 114                 } {
 115                         m := test.m
 116                         n := test.n
 117                         a := NewDense(m, n, nil)
 118                         for i := 0; i < m; i++ {
 119                                 for j := 0; j < n; j++ {
 120                                         a.Set(i, j, rand.Float64())
 121                                 }
 122                         }
 123                         br := m
 124                         if trans {
 125                                 br = n
 126                         }
 127                         b := NewVecDense(br, nil)
 128                         for i := 0; i < br; i++ {
 129                                 b.SetVec(i, rand.Float64())
 130                         }
 131                         var x VecDense
 132                         lq := &LQ{}
 133                         lq.Factorize(a)
 134                         lq.SolveVec(&x, trans, b)
 135
 136                         // Test that the normal equations hold.
 137                         // A^T * A * x = A^T * b if !trans
 138                         // A * A^T * x = A * b if trans
 139                         var lhs Dense
 140                         var rhs Dense
 141                         if trans {
 142                                 var tmp Dense
 143                                 tmp.Mul(a, a.T())
 144                                 lhs.Mul(&tmp, &x)
 145                                 rhs.Mul(a, b)
 146                         } else {
 147                                 var tmp Dense
 148                                 tmp.Mul(a.T(), a)
 149                                 lhs.Mul(&tmp, &x)
 150                                 rhs.Mul(a.T(), b)
 151                         }
 152                         if !EqualApprox(&lhs, &rhs, 1e-10) {
 153                                 t.Errorf("Normal equations do not hold.\nLHS: %v\n, RHS: %v\n", lhs, rhs)
 154                         }
 155                 }
 156         }
 157         // TODO(btracey): Add in testOneInput when it exists.
 158 }
 159
 160 func TestSolveLQCond(t *testing.T) {
 161         for _, test := range []*Dense{
 162                 NewDense(2, 2, []float64{1, 0, 0, 1e-20}),
 163                 NewDense(2, 3, []float64{1, 0, 0, 0, 1e-20, 0}),
 164         } {
 165                 m, _ := test.Dims()
 166                 var lq LQ
 167                 lq.Factorize(test)
 168                 b := NewDense(m, 2, nil)
 169                 var x Dense
 170                 if err := lq.Solve(&x, false, b); err == nil {
 171                         t.Error("No error for near-singular matrix in matrix solve.")
 172                 }
 173
 174                 bvec := NewVecDense(m, nil)
 175                 var xvec VecDense
 176                 if err := lq.SolveVec(&xvec, false, bvec); err == nil {
 177                         t.Error("No error for near-singular matrix in matrix solve.")
 178                 }
 179         }
 180 }