octave_packages/gsl-1.0.8/test_ellint.m

   1 ## Copyright (C) 2008   Raymond E. Rogers
   2 ##
   3 ## This program is free software; you can redistribute it and/or modify
   4 ## it under the terms of the GNU General Public License as published by
   5 ## the Free Software Foundation; either version 3 of the License, or
   6 ##  any later version.
   7 ##
   8 ## This program is distributed in the hope that it will be useful,
   9 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  11 ## GNU General Public License for more details.
  12 ##
  13 ## You should have received a copy of the GNU General Public License
  14 ## along with this program; if not, see .
  15
  16 #
  17 #
  18 #
  19 # R. Rogers v1, testing gsl_sf_ellint_Kcomp_e, gsl_sf_ellint_Ecomp_e
  20 #       implementation in octave as ellint_Kcomp, ellint_Ecomp
  21 #       against the gsl tests
  22 #       Some array verification is done; since I messed it up once.
  23
  24
  25 global GSL_DBL_EPSILON=        2.2204460492503131e-16
  26 global GSL_SQRT_DBL_EPSILON=   1.4901161193847656e-08
  27 global TEST_TOL0= (2.0*GSL_DBL_EPSILON)
  28 global TEST_TOL1=  (16.0*GSL_DBL_EPSILON)
  29 global TEST_TOL2=  (256.0*GSL_DBL_EPSILON)
  30 global TEST_TOL3=  (2048.0*GSL_DBL_EPSILON)
  31 global TEST_TOL4=  (16384.0*GSL_DBL_EPSILON)
  32 global TEST_TOL5=  (131072.0*GSL_DBL_EPSILON)
  33 global TEST_TOL6=  (1048576.0*GSL_DBL_EPSILON)
  34 global TEST_SQRT_TOL0= (2.0*GSL_SQRT_DBL_EPSILON)
  35 global TEST_SNGL=  (1.0e-06)
  36
  37 global TEST_SF_INCONS=  1
  38 global TEST_SF_ERRNEG=  2
  39 global TEST_SF_TOLBAD=  4
  40 global TEST_SF_RETBAD=  8
  41 global TEST_SF_ERRBAD=  16
  42 global TEST_SF_ERRBIG=  32
  43 global TEST_SF_EXPBAD=  64
  44 global TEST_FACTOR=     100
  45
  46
  47 function res=READ_TEST_SF_ellint(input_name)
  48 global GSL_DBL_EPSILON
  49 global GSL_SQRT_DBL_EPSILON
  50 global TEST_TOL0
  51 global TEST_TOL1
  52 global TEST_TOL2
  53 global TEST_TOL3
  54 global TEST_TOL4
  55 global TEST_TOL5
  56 global TEST_TOL6
  57 global TEST_SQRT_TOL0
  58 global TEST_SNGL
  59
  60 global TEST_SF_INCONS
  61 global TEST_SF_ERRNEG
  62 global TEST_SF_TOLBAD
  63 global TEST_SF_RETBAD
  64 global TEST_SF_ERRBAD
  65 global TEST_SF_ERRBIG
  66 global TEST_SF_EXPBAD
  67
  68 global TEST_FACTOR
  69         [source_id,source_msg]=fopen(input_name,"r","native")
  70
  71
  72                 while (! feof(source_id))
  73                         do
  74                                 input_line=fgetl(source_id);
  75                         until( index(input_line,"//") == 0);
  76
  77                         str_p=index(input_line,"gsl_sf_ellint_Kcomp_e");
  78                         if (str_p != 0)
  79                                 # Take it apart
  80                                 string_split=char(strsplit(input_line,","));
  81                                 arg1=str2double(substr(string_split(3,:),3));
  82                                 arg2=str2double(string_split(4,:));
  83                                 arg3=(string_split(5,:));
  84                                 val=str2double(string_split(6,:));
  85                                 tol=eval(string_split(7,:));
  86                                 [ret,err]=ellint_Kcomp(arg1,arg2);
  87                                 # This is to prevent extanious stops on some errors
  88                                 if ret==NaN
  89                                         ret=Inf
  90                                 endif
  91                                 if (abs((ret-val)/val)<tol*TEST_FACTOR)
  92                                         printf("\n Passed ellint_Kcomp: %e  ",arg1)
  93                                 else
  94                                         printf("\n Failed ellint_Kcomp: %s\n value=%e, computed=%e, tol=%e, returned error=%e ",input_line,val,ret,tol,err)
  95                                         printf("\n error %e", abs((ret-val)/val))
  96                                 endif
  97                         endif
  98                         str_p=index(input_line,"gsl_sf_ellint_Ecomp_e");
  99                         if (str_p != 0)
 100                                 # Take it apart
 101                                 string_split=char(strsplit(input_line,","));
 102                                 arg1=str2double(substr(string_split(3,:),3));
 103                                 arg2=str2double(string_split(4,:));
 104                                 arg3=(string_split(5,:));
 105                                 val=str2double(string_split(6,:));
 106                                 tol=eval(string_split(7,:));
 107                                 [ret,err]=ellint_Ecomp(arg1,arg2);
 108                                 # This is to prevent extanious stops on some errors
 109                                 if ret==NaN
 110                                         ret=Inf
 111                                 endif
 112                                 if (abs((ret-val)/val)<tol*TEST_FACTOR)
 113                                         printf("\n Passed ellint_Ecomp: %e  ",arg1)
 114                                 else
 115                                         printf("\n Failed ellint_Ecomp: %s\n value=%e, computed=%e, tol=%e, returned error=%e ",input_line,val,ret,tol,err)
 116                                         printf("\n error %e", abs((ret-val)/val))
 117                                 endif
 118                         endif
 119                 endwhile
 120         fclose(source_id);
 121         #lets do some array tests
 122         disp("\n array tests")
 123         disp("\n Argument array")
 124         mat1=[.1,.2; .3,.5]
 125         disp("\n ellint_Kcomp(mat1)")
 126         ellint_Kcomp(mat1,0)
 127         disp("\n ellint_Ecomp(mat1)")
 128         ellint_Ecomp(mat1,0)
 129
 130         res="";
 131 end
 132 #file_name=input("file name: ","s")
 133 file_name="test_sf.c"
 134 READ_TEST_SF_ellint(file_name);
 135
 136