github.com/gopherd/gonum@v0.0.4/mathext/internal/amos/amos_fortran_test.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  //go:build fortran
     6  // +build fortran
     7  
     8  package amos
     9  
    10  import (
    11  	"flag"
    12  	"runtime"
    13  	"testing"
    14  
    15  	"math/rand"
    16  
    17  	"github.com/gopherd/gonum/mathext/internal/amos/amoslib"
    18  )
    19  
    20  // BUG(kortschak): Some tests here comparing the direct Go translation
    21  // of the Fortran code fail. Do not delete these tests or this file until
    22  // https://github.com/gonum/gonum/issues/1322 has been satisfactorily
    23  // resolved.
    24  var runFailing = flag.Bool("failing", false, "run known failing cases")
    25  
    26  func TestAiryFortran(t *testing.T) {
    27  	rnd := rand.New(rand.NewSource(1))
    28  	for i := 0; i < nInputs; i++ {
    29  		in := randInput(rnd)
    30  		zairytestFort(t, in.x, in.kode, in.id)
    31  	}
    32  }
    33  
    34  func TestZacaiFortran(t *testing.T) {
    35  	if !*runFailing {
    36  		t.Skip("fails")
    37  	}
    38  
    39  	switch runtime.GOARCH {
    40  	case "arm64":
    41  		t.Skipf("skipping on GOARCH=%s", runtime.GOARCH)
    42  	}
    43  	rnd := rand.New(rand.NewSource(1))
    44  	for i := 0; i < nInputs; i++ {
    45  		in := randInput(rnd)
    46  		zacaitestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi, in.kode)
    47  	}
    48  }
    49  
    50  func TestZbknuFortran(t *testing.T) {
    51  	if !*runFailing {
    52  		t.Skip("fails")
    53  	}
    54  
    55  	rnd := rand.New(rand.NewSource(1))
    56  	for i := 0; i < nInputs; i++ {
    57  		in := randInput(rnd)
    58  		zbknutestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi, in.kode)
    59  	}
    60  }
    61  
    62  func TestZasyiFortran(t *testing.T) {
    63  	if !*runFailing {
    64  		t.Skip("fails")
    65  	}
    66  
    67  	rnd := rand.New(rand.NewSource(1))
    68  	for i := 0; i < nInputs; i++ {
    69  		in := randInput(rnd)
    70  		zasyitestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi, in.kode)
    71  	}
    72  }
    73  
    74  func TestZseriFortran(t *testing.T) {
    75  	switch runtime.GOARCH {
    76  	case "arm64":
    77  		t.Skipf("skipping on GOARCH=%s", runtime.GOARCH)
    78  	}
    79  	rnd := rand.New(rand.NewSource(1))
    80  	for i := 0; i < nInputs; i++ {
    81  		in := randInput(rnd)
    82  		zseritestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi, in.kode)
    83  	}
    84  }
    85  
    86  func TestZmlriFortran(t *testing.T) {
    87  	if !*runFailing {
    88  		t.Skip("fails")
    89  	}
    90  
    91  	rnd := rand.New(rand.NewSource(1))
    92  	for i := 0; i < nInputs; i++ {
    93  		in := randInput(rnd)
    94  		zmlritestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi, in.kode)
    95  	}
    96  }
    97  
    98  func TestZksclFortran(t *testing.T) {
    99  	if !*runFailing {
   100  		t.Skip("fails")
   101  	}
   102  
   103  	rnd := rand.New(rand.NewSource(1))
   104  	for i := 0; i < nInputs; i++ {
   105  		in := randInput(rnd)
   106  		zkscltestFort(t, in.x, in.is, in.tol, in.n, in.yr, in.yi)
   107  	}
   108  }
   109  
   110  func TestZuchkFortran(t *testing.T) {
   111  	rnd := rand.New(rand.NewSource(1))
   112  	for i := 0; i < nInputs; i++ {
   113  		in := randInput(rnd)
   114  		zuchktestFort(t, in.x, in.is, in.tol)
   115  	}
   116  }
   117  
   118  func TestZs1s2Fortran(t *testing.T) {
   119  	rnd := rand.New(rand.NewSource(1))
   120  	for i := 0; i < nInputs; i++ {
   121  		in := randInput(rnd)
   122  		zs1s2testFort(t, in.x, in.is)
   123  	}
   124  }
   125  
   126  func zs1s2testFort(t *testing.T, x []float64, is []int) {
   127  	const tol = 1e-11
   128  
   129  	type data struct {
   130  		ZRR, ZRI, S1R, S1I, S2R, S2I float64
   131  		NZ                           int
   132  		ASCLE, ALIM                  float64
   133  		IUF                          int
   134  	}
   135  
   136  	input := data{
   137  		x[0], x[1], x[2], x[3], x[4], x[5],
   138  		is[0],
   139  		x[6], x[7],
   140  		is[1],
   141  	}
   142  
   143  	impl := func(input data) data {
   144  		zrr, zri, s1r, s1i, s2r, s2i, nz, ascle, alim, iuf :=
   145  			amoslib.Zs1s2Fort(input.ZRR, input.ZRI, input.S1R, input.S1I, input.S2R, input.S2I, input.NZ, input.ASCLE, input.ALIM, input.IUF)
   146  		return data{zrr, zri, s1r, s1i, s2r, s2i, nz, ascle, alim, iuf}
   147  	}
   148  
   149  	comp := func(input data) data {
   150  		zrr, zri, s1r, s1i, s2r, s2i, nz, ascle, alim, iuf :=
   151  			zs1s2Orig(input.ZRR, input.ZRI, input.S1R, input.S1I, input.S2R, input.S2I, input.NZ, input.ASCLE, input.ALIM, input.IUF)
   152  		return data{zrr, zri, s1r, s1i, s2r, s2i, nz, ascle, alim, iuf}
   153  	}
   154  
   155  	oi := impl(input)
   156  	oc := comp(input)
   157  
   158  	sameF64Approx(t, "zs1s2 zrr", oc.ZRR, oi.ZRR, tol)
   159  	sameF64Approx(t, "zs1s2 zri", oc.ZRI, oi.ZRI, tol)
   160  	sameF64Approx(t, "zs1s2 s1r", oc.S1R, oi.S1R, tol)
   161  	sameF64Approx(t, "zs1s2 s1i", oc.S1I, oi.S1I, tol)
   162  	sameF64Approx(t, "zs1s2 s2r", oc.S2R, oi.S2R, tol)
   163  	sameF64Approx(t, "zs1s2 s2i", oc.S2I, oi.S2I, tol)
   164  	sameF64Approx(t, "zs1s2 ascle", oc.ASCLE, oi.ASCLE, tol)
   165  	sameF64Approx(t, "zs1s2 alim", oc.ALIM, oi.ALIM, tol)
   166  	sameInt(t, "iuf", oc.IUF, oi.IUF)
   167  	sameInt(t, "nz", oc.NZ, oi.NZ)
   168  }
   169  
   170  func zuchktestFort(t *testing.T, x []float64, is []int, tol float64) {
   171  	t.Helper()
   172  
   173  	YR := x[0]
   174  	YI := x[1]
   175  	NZ := is[0]
   176  	ASCLE := x[2]
   177  	TOL := tol
   178  
   179  	YRfort, YIfort, NZfort, ASCLEfort, TOLfort := zuchkOrig(YR, YI, NZ, ASCLE, TOL)
   180  	YRamoslib, YIamoslib, NZamoslib, ASCLEamoslib, TOLamoslib := amoslib.ZuchkFort(YR, YI, NZ, ASCLE, TOL)
   181  
   182  	sameF64(t, "zuchk yr", YRfort, YRamoslib)
   183  	sameF64(t, "zuchk yi", YIfort, YIamoslib)
   184  	sameInt(t, "zuchk nz", NZfort, NZamoslib)
   185  	sameF64(t, "zuchk ascle", ASCLEfort, ASCLEamoslib)
   186  	sameF64(t, "zuchk tol", TOLfort, TOLamoslib)
   187  }
   188  
   189  func zkscltestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64) {
   190  	t.Helper()
   191  
   192  	ZRR := x[0]
   193  	ZRI := x[1]
   194  	FNU := x[2]
   195  	NZ := is[1]
   196  	ELIM := x[3]
   197  	ASCLE := x[4]
   198  	RZR := x[6]
   199  	RZI := x[7]
   200  
   201  	yrfort := make([]float64, len(yr))
   202  	copy(yrfort, yr)
   203  	yifort := make([]float64, len(yi))
   204  	copy(yifort, yi)
   205  	ZRRfort, ZRIfort, FNUfort, Nfort, YRfort, YIfort, NZfort, RZRfort, RZIfort, ASCLEfort, TOLfort, ELIMfort :=
   206  		zksclOrig(ZRR, ZRI, FNU, n, yrfort, yifort, NZ, RZR, RZI, ASCLE, tol, ELIM)
   207  
   208  	yramos := make([]float64, len(yr))
   209  	copy(yramos, yr)
   210  	yiamos := make([]float64, len(yi))
   211  	copy(yiamos, yi)
   212  	ZRRamoslib, ZRIamoslib, FNUamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, RZRamoslib, RZIamoslib, ASCLEamoslib, TOLamoslib, ELIMamoslib :=
   213  		amoslib.ZksclFort(ZRR, ZRI, FNU, n, yramos, yiamos, NZ, RZR, RZI, ASCLE, tol, ELIM)
   214  
   215  	sameF64(t, "zkscl zrr", ZRRfort, ZRRamoslib)
   216  	sameF64(t, "zkscl zri", ZRIfort, ZRIamoslib)
   217  	sameF64(t, "zkscl fnu", FNUfort, FNUamoslib)
   218  	sameInt(t, "zkscl n", Nfort, Namoslib)
   219  	sameInt(t, "zkscl nz", NZfort, NZamoslib)
   220  	sameF64(t, "zkscl rzr", RZRfort, RZRamoslib)
   221  	sameF64(t, "zkscl rzi", RZIfort, RZIamoslib)
   222  	sameF64(t, "zkscl ascle", ASCLEfort, ASCLEamoslib)
   223  	sameF64(t, "zkscl tol", TOLfort, TOLamoslib)
   224  	sameF64(t, "zkscl elim", ELIMfort, ELIMamoslib)
   225  
   226  	sameF64SApprox(t, "zkscl yr", YRfort, YRamoslib, 1e-14)
   227  	sameF64SApprox(t, "zkscl yi", YIfort, YIamoslib, 1e-14)
   228  }
   229  
   230  func zmlritestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64, kode int) {
   231  	t.Helper()
   232  
   233  	ZR := x[0]
   234  	ZI := x[1]
   235  	FNU := x[2]
   236  	KODE := kode
   237  	NZ := is[1]
   238  
   239  	yrfort := make([]float64, len(yr))
   240  	copy(yrfort, yr)
   241  	yifort := make([]float64, len(yi))
   242  	copy(yifort, yi)
   243  	ZRfort, ZIfort, FNUfort, KODEfort, Nfort, YRfort, YIfort, NZfort, TOLfort :=
   244  		zmlriOrig(ZR, ZI, FNU, KODE, n, yrfort, yifort, NZ, tol)
   245  
   246  	yramos := make([]float64, len(yr))
   247  	copy(yramos, yr)
   248  	yiamos := make([]float64, len(yi))
   249  	copy(yiamos, yi)
   250  	ZRamoslib, ZIamoslib, FNUamoslib, KODEamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, TOLamoslib :=
   251  		amoslib.ZmlriFort(ZR, ZI, FNU, KODE, n, yramos, yiamos, NZ, tol)
   252  
   253  	sameF64(t, "zmlri zr", ZRfort, ZRamoslib)
   254  	sameF64(t, "zmlri zi", ZIfort, ZIamoslib)
   255  	sameF64(t, "zmlri fnu", FNUfort, FNUamoslib)
   256  	sameInt(t, "zmlri kode", KODEfort, KODEamoslib)
   257  	sameInt(t, "zmlri n", Nfort, Namoslib)
   258  	sameInt(t, "zmlri nz", NZfort, NZamoslib)
   259  	sameF64(t, "zmlri tol", TOLfort, TOLamoslib)
   260  
   261  	sameF64S(t, "zmlri yr", YRfort, YRamoslib)
   262  	sameF64S(t, "zmlri yi", YIfort, YIamoslib)
   263  }
   264  
   265  func zseritestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64, kode int) {
   266  	t.Helper()
   267  
   268  	ZR := x[0]
   269  	ZI := x[1]
   270  	FNU := x[2]
   271  	KODE := kode
   272  	NZ := is[1]
   273  	ELIM := x[3]
   274  	ALIM := x[4]
   275  
   276  	yrfort := make([]float64, len(yr))
   277  	copy(yrfort, yr)
   278  	yifort := make([]float64, len(yi))
   279  	copy(yifort, yi)
   280  	ZRfort, ZIfort, FNUfort, KODEfort, Nfort, YRfort, YIfort, NZfort, TOLfort, ELIMfort, ALIMfort :=
   281  		zseriOrig(ZR, ZI, FNU, KODE, n, yrfort, yifort, NZ, tol, ELIM, ALIM)
   282  
   283  	yramos := make([]float64, len(yr))
   284  	copy(yramos, yr)
   285  	yiamos := make([]float64, len(yi))
   286  	copy(yiamos, yi)
   287  	y := make([]complex128, len(yramos))
   288  	for i, v := range yramos {
   289  		y[i] = complex(v, yiamos[i])
   290  	}
   291  
   292  	ZRamoslib, ZIamoslib, FNUamoslib, KODEamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, TOLamoslib, ELIMamoslib, ALIMamoslib :=
   293  		amoslib.ZseriFort(ZR, ZI, FNU, KODE, n, yrfort, yifort, NZ, tol, ELIM, ALIM)
   294  
   295  	sameF64(t, "zseri zr", ZRfort, ZRamoslib)
   296  	sameF64(t, "zseri zi", ZIfort, ZIamoslib)
   297  	sameF64(t, "zseri fnu", FNUfort, FNUamoslib)
   298  	sameInt(t, "zseri kode", KODEfort, KODEamoslib)
   299  	sameInt(t, "zseri n", Nfort, Namoslib)
   300  	if *runFailing {
   301  		sameInt(t, "zseri nz", NZfort, NZamoslib)
   302  	}
   303  	sameF64(t, "zseri tol", TOLfort, TOLamoslib)
   304  	sameF64(t, "zseri elim", ELIMfort, ELIMamoslib)
   305  	sameF64(t, "zseri elim", ALIMfort, ALIMamoslib)
   306  
   307  	sameF64SApprox(t, "zseri yr", YRfort, YRamoslib, 1e-9)
   308  	sameF64SApprox(t, "zseri yi", YIfort, YIamoslib, 1e-10)
   309  }
   310  
   311  func zasyitestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64, kode int) {
   312  	t.Helper()
   313  
   314  	ZR := x[0]
   315  	ZI := x[1]
   316  	FNU := x[2]
   317  	KODE := kode
   318  	NZ := is[1]
   319  	ELIM := x[3]
   320  	ALIM := x[4]
   321  	RL := x[5]
   322  
   323  	yrfort := make([]float64, len(yr))
   324  	copy(yrfort, yr)
   325  	yifort := make([]float64, len(yi))
   326  	copy(yifort, yi)
   327  	ZRfort, ZIfort, FNUfort, KODEfort, Nfort, YRfort, YIfort, NZfort, RLfort, TOLfort, ELIMfort, ALIMfort :=
   328  		zasyiOrig(ZR, ZI, FNU, KODE, n, yrfort, yifort, NZ, RL, tol, ELIM, ALIM)
   329  
   330  	yramos := make([]float64, len(yr))
   331  	copy(yramos, yr)
   332  	yiamos := make([]float64, len(yi))
   333  	copy(yiamos, yi)
   334  	ZRamoslib, ZIamoslib, FNUamoslib, KODEamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, RLamoslib, TOLamoslib, ELIMamoslib, ALIMamoslib :=
   335  		amoslib.ZasyiFort(ZR, ZI, FNU, KODE, n, yramos, yiamos, NZ, RL, tol, ELIM, ALIM)
   336  
   337  	sameF64(t, "zasyi zr", ZRfort, ZRamoslib)
   338  	sameF64(t, "zasyi zr", ZIfort, ZIamoslib)
   339  	sameF64(t, "zasyi fnu", FNUfort, FNUamoslib)
   340  	sameInt(t, "zasyi kode", KODEfort, KODEamoslib)
   341  	sameInt(t, "zasyi n", Nfort, Namoslib)
   342  	sameInt(t, "zasyi nz", NZfort, NZamoslib)
   343  	sameF64(t, "zasyi rl", RLfort, RLamoslib)
   344  	sameF64(t, "zasyi tol", TOLfort, TOLamoslib)
   345  	sameF64(t, "zasyi elim", ELIMfort, ELIMamoslib)
   346  	sameF64(t, "zasyi alim", ALIMfort, ALIMamoslib)
   347  
   348  	sameF64SApprox(t, "zasyi yr", YRfort, YRamoslib, 1e-12)
   349  	sameF64SApprox(t, "zasyi yi", YIfort, YIamoslib, 1e-12)
   350  }
   351  
   352  func zbknutestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64, kode int) {
   353  	t.Helper()
   354  
   355  	ZR := x[0]
   356  	ZI := x[1]
   357  	FNU := x[2]
   358  	KODE := kode
   359  	NZ := is[1]
   360  	ELIM := x[3]
   361  	ALIM := x[4]
   362  
   363  	yrfort := make([]float64, len(yr))
   364  	copy(yrfort, yr)
   365  	yifort := make([]float64, len(yi))
   366  	copy(yifort, yi)
   367  	ZRfort, ZIfort, FNUfort, KODEfort, Nfort, YRfort, YIfort, NZfort, TOLfort, ELIMfort, ALIMfort :=
   368  		zbknuOrig(ZR, ZI, FNU, KODE, n, yrfort, yifort, NZ, tol, ELIM, ALIM)
   369  
   370  	yramos := make([]float64, len(yr))
   371  	copy(yramos, yr)
   372  	yiamos := make([]float64, len(yi))
   373  	copy(yiamos, yi)
   374  	ZRamoslib, ZIamoslib, FNUamoslib, KODEamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, TOLamoslib, ELIMamoslib, ALIMamoslib :=
   375  		amoslib.ZbknuFort(ZR, ZI, FNU, KODE, n, yramos, yiamos, NZ, tol, ELIM, ALIM)
   376  
   377  	sameF64(t, "zbknu zr", ZRfort, ZRamoslib)
   378  	sameF64(t, "zbknu zr", ZIfort, ZIamoslib)
   379  	sameF64(t, "zbknu fnu", FNUfort, FNUamoslib)
   380  	sameInt(t, "zbknu kode", KODEfort, KODEamoslib)
   381  	sameInt(t, "zbknu n", Nfort, Namoslib)
   382  	sameInt(t, "zbknu nz", NZfort, NZamoslib)
   383  	sameF64(t, "zbknu tol", TOLfort, TOLamoslib)
   384  	sameF64(t, "zbknu elim", ELIMfort, ELIMamoslib)
   385  	sameF64(t, "zbknu alim", ALIMfort, ALIMamoslib)
   386  
   387  	sameF64SApprox(t, "zbknu yr", YRfort, YRamoslib, 1e-12)
   388  	sameF64SApprox(t, "zbknu yi", YIfort, YIamoslib, 1e-12)
   389  }
   390  
   391  func zairytestFort(t *testing.T, x []float64, kode, id int) {
   392  	const tol = 1e-8
   393  	t.Helper()
   394  
   395  	ZR := x[0]
   396  	ZI := x[1]
   397  	KODE := kode
   398  	ID := id
   399  
   400  	AIRfort, AIIfort, NZfort, IERRfort := zairyOrig(ZR, ZI, ID, KODE)
   401  	AIRamos, AIIamos, NZamos, IERRamos := amoslib.ZairyFort(ZR, ZI, ID, KODE)
   402  
   403  	sameF64Approx(t, "zairy air", AIRfort, AIRamos, tol)
   404  	sameF64Approx(t, "zairy aii", AIIfort, AIIamos, tol)
   405  	sameInt(t, "zairy nz", NZfort, NZamos)
   406  	sameInt(t, "zairy ierr", IERRfort, IERRamos)
   407  }
   408  
   409  func zacaitestFort(t *testing.T, x []float64, is []int, tol float64, n int, yr, yi []float64, kode int) {
   410  	t.Helper()
   411  
   412  	ZR := x[0]
   413  	ZI := x[1]
   414  	FNU := x[2]
   415  	KODE := kode
   416  	NZ := is[1]
   417  	MR := is[2]
   418  	ELIM := x[3]
   419  	ALIM := x[4]
   420  	RL := x[5]
   421  
   422  	yrfort := make([]float64, len(yr))
   423  	copy(yrfort, yr)
   424  	yifort := make([]float64, len(yi))
   425  	copy(yifort, yi)
   426  	ZRfort, ZIfort, FNUfort, KODEfort, MRfort, Nfort, YRfort, YIfort, NZfort, RLfort, TOLfort, ELIMfort, ALIMfort :=
   427  		zacaiOrig(ZR, ZI, FNU, KODE, MR, n, yrfort, yifort, NZ, RL, tol, ELIM, ALIM)
   428  
   429  	yramos := make([]float64, len(yr))
   430  	copy(yramos, yr)
   431  	yiamos := make([]float64, len(yi))
   432  	copy(yiamos, yi)
   433  	ZRamoslib, ZIamoslib, FNUamoslib, KODEamoslib, MRamoslib, Namoslib, YRamoslib, YIamoslib, NZamoslib, RLamoslib, TOLamoslib, ELIMamoslib, ALIMamoslib :=
   434  		amoslib.ZacaiFort(ZR, ZI, FNU, KODE, MR, n, yramos, yiamos, NZ, RL, tol, ELIM, ALIM)
   435  
   436  	sameF64(t, "zacai zr", ZRfort, ZRamoslib)
   437  	sameF64(t, "zacai zi", ZIfort, ZIamoslib)
   438  	sameF64(t, "zacai fnu", FNUfort, FNUamoslib)
   439  	sameInt(t, "zacai kode", KODEfort, KODEamoslib)
   440  	sameInt(t, "zacai mr", MRfort, MRamoslib)
   441  	sameInt(t, "zacai n", Nfort, Namoslib)
   442  	sameInt(t, "zacai nz", NZfort, NZamoslib)
   443  	sameF64(t, "zacai rl", RLfort, RLamoslib)
   444  	sameF64(t, "zacai tol", TOLfort, TOLamoslib)
   445  	sameF64(t, "zacai elim", ELIMfort, ELIMamoslib)
   446  	sameF64(t, "zacai elim", ALIMfort, ALIMamoslib)
   447  
   448  	sameF64SApprox(t, "zacai yr", YRfort, YRamoslib, 1e-12)
   449  	sameF64SApprox(t, "zacai yi", YIfort, YIamoslib, 1e-12)
   450  }