github.com/gonum/lapack@v0.0.0-20181123203213-e4cdc5a0bff9/testlapack/dlanv2.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  	"math/rand"
    11  	"testing"
    12  )
    13  
    14  type Dlanv2er interface {
    15  	Dlanv2(a, b, c, d float64) (aa, bb, cc, dd float64, rt1r, rt1i, rt2r, rt2i float64, cs, sn float64)
    16  }
    17  
    18  func Dlanv2Test(t *testing.T, impl Dlanv2er) {
    19  	rnd := rand.New(rand.NewSource(1))
    20  	for i := 0; i < 1000; i++ {
    21  		a := rnd.NormFloat64()
    22  		b := rnd.NormFloat64()
    23  		c := rnd.NormFloat64()
    24  		d := rnd.NormFloat64()
    25  		aa, bb, cc, dd, rt1r, rt1i, rt2r, rt2i, cs, sn := impl.Dlanv2(a, b, c, d)
    26  
    27  		mat := fmt.Sprintf("[%v %v; %v %v]", a, b, c, d)
    28  		if cc == 0 {
    29  			if rt1i != 0 || rt2i != 0 {
    30  				t.Errorf("Unexpected complex eigenvalues for %v", mat)
    31  			}
    32  		} else {
    33  			if aa != dd {
    34  				t.Errorf("Diagonal elements not equal for %v: got [%v %v]", mat, aa, dd)
    35  			}
    36  			if bb*cc >= 0 {
    37  				t.Errorf("Non-diagonal elements have the same sign for %v: got [%v %v]", mat, bb, cc)
    38  			} else {
    39  				im := math.Sqrt(-bb * cc)
    40  				if math.Abs(rt1i-im) > 1e-14 && math.Abs(rt1i+im) > 1e-14 {
    41  					t.Errorf("Unexpected imaginary part of eigenvalue for %v: got %v, want %v or %v", mat, rt1i, im, -im)
    42  				}
    43  				if math.Abs(rt2i-im) > 1e-14 && math.Abs(rt2i+im) > 1e-14 {
    44  					t.Errorf("Unexpected imaginary part of eigenvalue for %v: got %v, want %v or %v", mat, rt2i, im, -im)
    45  				}
    46  			}
    47  		}
    48  		if rt1r != aa && rt1r != dd {
    49  			t.Errorf("Unexpected real part of eigenvalue for %v: got %v, want %v or %v", mat, rt1r, aa, dd)
    50  		}
    51  		if rt2r != aa && rt2r != dd {
    52  			t.Errorf("Unexpected real part of eigenvalue for %v: got %v, want %v or %v", mat, rt2r, aa, dd)
    53  		}
    54  		if math.Abs(math.Hypot(cs, sn)-1) > 1e-14 {
    55  			t.Errorf("Unexpected unitary matrix for %v: got cs %v, sn %v", mat, cs, sn)
    56  		}
    57  
    58  		gota := cs*(aa*cs-bb*sn) - sn*(cc*cs-dd*sn)
    59  		gotb := cs*(aa*sn+bb*cs) - sn*(cc*sn+dd*cs)
    60  		gotc := sn*(aa*cs-bb*sn) + cs*(cc*cs-dd*sn)
    61  		gotd := sn*(aa*sn+bb*cs) + cs*(cc*sn+dd*cs)
    62  		if math.Abs(gota-a) > 1e-14 ||
    63  			math.Abs(gotb-b) > 1e-14 ||
    64  			math.Abs(gotc-c) > 1e-14 ||
    65  			math.Abs(gotd-d) > 1e-14 {
    66  			t.Errorf("Unexpected factorization: got [%v %v; %v %v], want [%v %v; %v %v]", gota, gotb, gotc, gotd, a, b, c, d)
    67  		}
    68  	}
    69  }