vendor/gonum.org/v1/gonum/blas/testblas/common.go

   1 // Copyright ©2014 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 testblas
   6
   7 import (
   8         "math"
   9         "math/cmplx"
  10         "testing"
  11
  12         "gonum.org/v1/gonum/blas"
  13 )
  14
  15 // throwPanic will throw unexpected panics if true, or will just report them as errors if false
  16 const throwPanic = true
  17
  18 var znan = cmplx.NaN()
  19
  20 func dTolEqual(a, b float64) bool {
  21         if math.IsNaN(a) && math.IsNaN(b) {
  22                 return true
  23         }
  24         if a == b {
  25                 return true
  26         }
  27         m := math.Max(math.Abs(a), math.Abs(b))
  28         if m > 1 {
  29                 a /= m
  30                 b /= m
  31         }
  32         if math.Abs(a-b) < 1e-14 {
  33                 return true
  34         }
  35         return false
  36 }
  37
  38 func dSliceTolEqual(a, b []float64) bool {
  39         if len(a) != len(b) {
  40                 return false
  41         }
  42         for i := range a {
  43                 if !dTolEqual(a[i], b[i]) {
  44                         return false
  45                 }
  46         }
  47         return true
  48 }
  49
  50 func dStridedSliceTolEqual(n int, a []float64, inca int, b []float64, incb int) bool {
  51         ia := 0
  52         ib := 0
  53         if inca <= 0 {
  54                 ia = -(n - 1) * inca
  55         }
  56         if incb <= 0 {
  57                 ib = -(n - 1) * incb
  58         }
  59         for i := 0; i < n; i++ {
  60                 if !dTolEqual(a[ia], b[ib]) {
  61                         return false
  62                 }
  63                 ia += inca
  64                 ib += incb
  65         }
  66         return true
  67 }
  68
  69 func dSliceEqual(a, b []float64) bool {
  70         if len(a) != len(b) {
  71                 return false
  72         }
  73         for i := range a {
  74                 if !dTolEqual(a[i], b[i]) {
  75                         return false
  76                 }
  77         }
  78         return true
  79 }
  80
  81 func dCopyTwoTmp(x, xTmp, y, yTmp []float64) {
  82         if len(x) != len(xTmp) {
  83                 panic("x size mismatch")
  84         }
  85         if len(y) != len(yTmp) {
  86                 panic("y size mismatch")
  87         }
  88         copy(xTmp, x)
  89         copy(yTmp, y)
  90 }
  91
  92 // returns true if the function panics
  93 func panics(f func()) (b bool) {
  94         defer func() {
  95                 err := recover()
  96                 if err != nil {
  97                         b = true
  98                 }
  99         }()
 100         f()
 101         return
 102 }
 103
 104 func testpanics(f func(), name string, t *testing.T) {
 105         b := panics(f)
 106         if !b {
 107                 t.Errorf("%v should panic and does not", name)
 108         }
 109 }
 110
 111 func sliceOfSliceCopy(a [][]float64) [][]float64 {
 112         n := make([][]float64, len(a))
 113         for i := range a {
 114                 n[i] = make([]float64, len(a[i]))
 115                 copy(n[i], a[i])
 116         }
 117         return n
 118 }
 119
 120 func sliceCopy(a []float64) []float64 {
 121         n := make([]float64, len(a))
 122         copy(n, a)
 123         return n
 124 }
 125
 126 func flatten(a [][]float64) []float64 {
 127         if len(a) == 0 {
 128                 return nil
 129         }
 130         m := len(a)
 131         n := len(a[0])
 132         s := make([]float64, m*n)
 133         for i := 0; i < m; i++ {
 134                 for j := 0; j < n; j++ {
 135                         s[i*n+j] = a[i][j]
 136                 }
 137         }
 138         return s
 139 }
 140
 141 func unflatten(a []float64, m, n int) [][]float64 {
 142         s := make([][]float64, m)
 143         for i := 0; i < m; i++ {
 144                 s[i] = make([]float64, n)
 145                 for j := 0; j < n; j++ {
 146                         s[i][j] = a[i*n+j]
 147                 }
 148         }
 149         return s
 150 }
 151
 152 // flattenTriangular turns the upper or lower triangle of a dense slice of slice
 153 // into a single slice with packed storage. a must be a square matrix.
 154 func flattenTriangular(a [][]float64, ul blas.Uplo) []float64 {
 155         m := len(a)
 156         aFlat := make([]float64, m*(m+1)/2)
 157         var k int
 158         if ul == blas.Upper {
 159                 for i := 0; i < m; i++ {
 160                         k += copy(aFlat[k:], a[i][i:])
 161                 }
 162                 return aFlat
 163         }
 164         for i := 0; i < m; i++ {
 165                 k += copy(aFlat[k:], a[i][:i+1])
 166         }
 167         return aFlat
 168 }
 169
 170 // flattenBanded turns a dense banded slice of slice into the compact banded matrix format
 171 func flattenBanded(a [][]float64, ku, kl int) []float64 {
 172         m := len(a)
 173         n := len(a[0])
 174         if ku < 0 || kl < 0 {
 175                 panic("testblas: negative band length")
 176         }
 177         nRows := m
 178         nCols := (ku + kl + 1)
 179         aflat := make([]float64, nRows*nCols)
 180         for i := range aflat {
 181                 aflat[i] = math.NaN()
 182         }
 183         // loop over the rows, and then the bands
 184         // elements in the ith row stay in the ith row
 185         // order in bands is kept
 186         for i := 0; i < nRows; i++ {
 187                 min := -kl
 188                 if i-kl < 0 {
 189                         min = -i
 190                 }
 191                 max := ku
 192                 if i+ku >= n {
 193                         max = n - i - 1
 194                 }
 195                 for j := min; j <= max; j++ {
 196                         col := kl + j
 197                         aflat[i*nCols+col] = a[i][i+j]
 198                 }
 199         }
 200         return aflat
 201 }
 202
 203 // makeIncremented takes a slice with inc == 1 and makes an incremented version
 204 // and adds extra values on the end
 205 func makeIncremented(x []float64, inc int, extra int) []float64 {
 206         if inc == 0 {
 207                 panic("zero inc")
 208         }
 209         absinc := inc
 210         if absinc < 0 {
 211                 absinc = -inc
 212         }
 213         xcopy := make([]float64, len(x))
 214         if inc > 0 {
 215                 copy(xcopy, x)
 216         } else {
 217                 for i := 0; i < len(x); i++ {
 218                         xcopy[i] = x[len(x)-i-1]
 219                 }
 220         }
 221
 222         // don't use NaN because it makes comparison hard
 223         // Do use a weird unique value for easier debugging
 224         counter := 100.0
 225         var xnew []float64
 226         for i, v := range xcopy {
 227                 xnew = append(xnew, v)
 228                 if i != len(x)-1 {
 229                         for j := 0; j < absinc-1; j++ {
 230                                 xnew = append(xnew, counter)
 231                                 counter++
 232                         }
 233                 }
 234         }
 235         for i := 0; i < extra; i++ {
 236                 xnew = append(xnew, counter)
 237                 counter++
 238         }
 239         return xnew
 240 }
 241
 242 func abs(x int) int {
 243         if x < 0 {
 244                 return -x
 245         }
 246         return x
 247 }
 248
 249 func allPairs(x, y []int) [][2]int {
 250         var p [][2]int
 251         for _, v0 := range x {
 252                 for _, v1 := range y {
 253                         p = append(p, [2]int{v0, v1})
 254                 }
 255         }
 256         return p
 257 }
 258
 259 func sameFloat64(a, b float64) bool {
 260         return a == b || math.IsNaN(a) && math.IsNaN(b)
 261 }
 262
 263 func sameComplex128(x, y complex128) bool {
 264         return sameFloat64(real(x), real(y)) && sameFloat64(imag(x), imag(y))
 265 }
 266
 267 func zsame(x, y []complex128) bool {
 268         if len(x) != len(y) {
 269                 return false
 270         }
 271         for i, v := range x {
 272                 w := y[i]
 273                 if !sameComplex128(v, w) {
 274                         return false
 275                 }
 276         }
 277         return true
 278 }
 279
 280 // zSameAtNonstrided returns whether elements at non-stride positions of vectors
 281 // x and y are same.
 282 func zSameAtNonstrided(x, y []complex128, inc int) bool {
 283         if len(x) != len(y) {
 284                 return false
 285         }
 286         if inc < 0 {
 287                 inc = -inc
 288         }
 289         for i, v := range x {
 290                 if i%inc == 0 {
 291                         continue
 292                 }
 293                 w := y[i]
 294                 if !sameComplex128(v, w) {
 295                         return false
 296                 }
 297         }
 298         return true
 299 }
 300
 301 // zEqualApproxAtStrided returns whether elements at stride positions of vectors
 302 // x and y are approximately equal within tol.
 303 func zEqualApproxAtStrided(x, y []complex128, inc int, tol float64) bool {
 304         if len(x) != len(y) {
 305                 return false
 306         }
 307         if inc < 0 {
 308                 inc = -inc
 309         }
 310         for i := 0; i < len(x); i += inc {
 311                 v := x[i]
 312                 w := y[i]
 313                 if !(cmplx.Abs(v-w) <= tol) {
 314                         return false
 315                 }
 316         }
 317         return true
 318 }
 319
 320 func makeZVector(data []complex128, inc int) []complex128 {
 321         if inc == 0 {
 322                 panic("bad test")
 323         }
 324         if len(data) == 0 {
 325                 return nil
 326         }
 327         inc = abs(inc)
 328         x := make([]complex128, (len(data)-1)*inc+1)
 329         for i := range x {
 330                 x[i] = znan
 331         }
 332         for i, v := range data {
 333                 x[i*inc] = v
 334         }
 335         return x
 336 }
 337
 338 func makeZGeneral(data []complex128, m, n int, ld int) []complex128 {
 339         if m < 0 || n < 0 {
 340                 panic("bad test")
 341         }
 342         if data != nil && len(data) != m*n {
 343                 panic("bad test")
 344         }
 345         if ld < max(1, n) {
 346                 panic("bad test")
 347         }
 348         if m == 0 || n == 0 {
 349                 return nil
 350         }
 351         a := make([]complex128, (m-1)*ld+n)
 352         for i := range a {
 353                 a[i] = znan
 354         }
 355         if data != nil {
 356                 for i := 0; i < m; i++ {
 357                         copy(a[i*ld:i*ld+n], data[i*n:i*n+n])
 358                 }
 359         }
 360         return a
 361 }
 362
 363 func max(a, b int) int {
 364         if a < b {
 365                 return b
 366         }
 367         return a
 368 }
 369
 370 // zPack returns the uplo triangle of an n×n matrix A in packed format.
 371 func zPack(uplo blas.Uplo, n int, a []complex128, lda int) []complex128 {
 372         if n == 0 {
 373                 return nil
 374         }
 375         ap := make([]complex128, n*(n+1)/2)
 376         var ii int
 377         if uplo == blas.Upper {
 378                 for i := 0; i < n; i++ {
 379                         for j := i; j < n; j++ {
 380                                 ap[ii] = a[i*lda+j]
 381                                 ii++
 382                         }
 383                 }
 384         } else {
 385                 for i := 0; i < n; i++ {
 386                         for j := 0; j <= i; j++ {
 387                                 ap[ii] = a[i*lda+j]
 388                                 ii++
 389                         }
 390                 }
 391         }
 392         return ap
 393 }
 394
 395 // zUnpackAsHermitian returns an n×n general Hermitian matrix (with stride n)
 396 // whose packed uplo triangle is stored on entry in ap.
 397 func zUnpackAsHermitian(uplo blas.Uplo, n int, ap []complex128) []complex128 {
 398         if n == 0 {
 399                 return nil
 400         }
 401         a := make([]complex128, n*n)
 402         lda := n
 403         var ii int
 404         if uplo == blas.Upper {
 405                 for i := 0; i < n; i++ {
 406                         for j := i; j < n; j++ {
 407                                 a[i*lda+j] = ap[ii]
 408                                 if i != j {
 409                                         a[j*lda+i] = cmplx.Conj(ap[ii])
 410                                 }
 411                                 ii++
 412                         }
 413                 }
 414         } else {
 415                 for i := 0; i < n; i++ {
 416                         for j := 0; j <= i; j++ {
 417                                 a[i*lda+j] = ap[ii]
 418                                 if i != j {
 419                                         a[j*lda+i] = cmplx.Conj(ap[ii])
 420                                 }
 421                                 ii++
 422                         }
 423                 }
 424         }
 425         return a
 426 }