vendor/gonum.org/v1/gonum/lapack/testlapack/dlartg.go

   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         "testing"
  10
  11         "golang.org/x/exp/rand"
  12
  13         "gonum.org/v1/gonum/floats"
  14 )
  15
  16 type Dlartger interface {
  17         Dlartg(f, g float64) (cs, sn, r float64)
  18 }
  19
  20 func DlartgTest(t *testing.T, impl Dlartger) {
  21         const tol = 1e-14
  22         // safmn2 and safmx2 are copied from native.Dlartg.
  23         safmn2 := math.Pow(dlamchB, math.Trunc(math.Log(dlamchS/dlamchE)/math.Log(dlamchB)/2))
  24         safmx2 := 1 / safmn2
  25         rnd := rand.New(rand.NewSource(1))
  26         for i := 0; i < 1000; i++ {
  27                 var f float64
  28                 var fHuge bool
  29                 switch rnd.Intn(3) {
  30                 case 0:
  31                         // Huge f.
  32                         fHuge = true
  33                         f = math.Pow(10, 10-20*rnd.Float64()) * safmx2
  34                 case 1:
  35                         // Tiny f.
  36                         f = math.Pow(10, 10-20*rnd.Float64()) * safmn2
  37                 default:
  38                         f = rnd.NormFloat64()
  39                 }
  40                 if rnd.Intn(2) == 0 {
  41                         f *= -1
  42                 }
  43
  44                 var g float64
  45                 var gHuge bool
  46                 switch rnd.Intn(3) {
  47                 case 0:
  48                         // Huge g.
  49                         gHuge = true
  50                         g = math.Pow(10, 10-20*rnd.Float64()) * safmx2
  51                 case 1:
  52                         // Tiny g.
  53                         g = math.Pow(10, 10-20*rnd.Float64()) * safmn2
  54                 default:
  55                         g = rnd.NormFloat64()
  56                 }
  57                 if rnd.Intn(2) == 0 {
  58                         g *= -1
  59                 }
  60
  61                 cs, sn, r := impl.Dlartg(f, g)
  62
  63                 rWant := cs*f + sn*g
  64                 if !floats.EqualWithinAbsOrRel(math.Abs(rWant), math.Abs(r), tol, tol) {
  65                         t.Errorf("Case f=%v,g=%v: unexpected r. Want %v, got %v", f, g, rWant, r)
  66                 }
  67                 oneTest := cs*cs + sn*sn
  68                 if math.Abs(oneTest-1) > tol {
  69                         t.Errorf("Case f=%v,g=%v: expected cs^2+sn^2==1, got %v", f, g, oneTest)
  70                 }
  71                 if !fHuge && !gHuge {
  72                         zeroTest := -sn*f + cs*g
  73                         if math.Abs(zeroTest) > tol {
  74                                 t.Errorf("Case f=%v,g=%v: expected zero, got %v", f, g, zeroTest)
  75                         }
  76                 }
  77                 if math.Abs(f) > math.Abs(g) && cs < 0 {
  78                         t.Errorf("Case f=%v,g=%v: unexpected negative cs %v", f, g, cs)
  79                 }
  80         }
  81         for i := 0; i < 100; i++ {
  82                 cs, sn, _ := impl.Dlartg(rnd.NormFloat64(), 0)
  83                 if cs != 1 {
  84                         t.Errorf("Unexpected cs for g=0. Want 1, got %v", cs)
  85                 }
  86                 if sn != 0 {
  87                         t.Errorf("Unexpected sn for g=0. Want 0, got %v", sn)
  88                 }
  89         }
  90         for i := 0; i < 100; i++ {
  91                 cs, sn, _ := impl.Dlartg(0, rnd.NormFloat64())
  92                 if cs != 0 {
  93                         t.Errorf("Unexpected cs for f=0. Want 0, got %v", cs)
  94                 }
  95                 if sn != 1 {
  96                         t.Errorf("Unexpected sn for f=0. Want 1, got %v", sn)
  97                 }
  98         }
  99 }