mpmath_special_functions_test_generator.py 3.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108
  1. #!/usr/bin/env python3
  2. import mpmath as mp
  3. import mpmath_riccati_bessel as mrb
  4. import mpmath_input_arguments as mia
  5. class test_data:
  6. def __init__(self, list_to_parse):
  7. self.evaluated_data = []
  8. self.comment = list_to_parse[0]
  9. if self.comment[:2] != '//': raise ValueError('Not a comment')
  10. self.ending = list_to_parse[-1]
  11. if self.ending != '};': raise ValueError('For C++ we expect closing };')
  12. self.list_to_parse = list_to_parse
  13. class update_special_functions_evaluations:
  14. def __init__(self, filename='default_out.h', complex_arguments = []):
  15. self.evaluated_data = []
  16. self.test_setup = []
  17. self.filename = filename
  18. self.complex_arguments = complex_arguments
  19. def read_evaluated_data(self):
  20. with open(self.filename, 'r') as in_file:
  21. content = in_file.readlines()
  22. content = [x.strip() for x in content]
  23. while '' in content:
  24. record_end_index = content.index('')
  25. new_record = content[:record_end_index]
  26. self.evaluated_data.append(test_data(new_record))
  27. content = content[record_end_index+1:]
  28. self.evaluated_data.append(test_data(new_record))
  29. def get_n_list(z, max_number_of_elements = 10):
  30. nmax = mrb.LeRu_cutoff(z)
  31. factor = nmax**(1/(max_number_of_elements-2))
  32. n_list = [int(factor**i) for i in range(max_number_of_elements-1) ]
  33. n_list.append(0)
  34. n_set = set(n_list)
  35. return sorted(n_set)
  36. def get_test_data_nlist(z_record, output_dps, n):
  37. x = str(z_record[0])
  38. mr = str(z_record[1][0])
  39. mi = str(z_record[1][1])
  40. z_str = ''
  41. try:
  42. z = mp.mpf(x)*mp.mpc(mr, mi)
  43. D1nz = mrb.D1(n,z)
  44. z_str=('{{'+
  45. mp.nstr(z.real, output_dps*2)+','+
  46. mp.nstr(z.imag, output_dps*2)+'},'+
  47. str(n)+',{'+
  48. mp.nstr(D1nz.real, output_dps)+','+
  49. mp.nstr(D1nz.imag, output_dps)+'},'+
  50. mp.nstr(mp.fabs(D1nz.real* 10**-output_dps),2)+',' +
  51. mp.nstr(mp.fabs(D1nz.imag* 10**-output_dps),2)+
  52. '},\n')
  53. except:
  54. pass
  55. return z_str
  56. def get_test_data(Du_test, output_dps, max_num_elements_of_n_list):
  57. output_str = ('// complex(z), n, complex(D1(n,z)), abs_err_real, abs_err_imag\n'+
  58. 'std::vector< std::tuple< std::complex<double>, int, std::complex<double>, double, double > >'+
  59. 'D1_test_'+str(output_dps)+'digits = {\n')
  60. for z_record in Du_test:
  61. x = str(z_record[0])
  62. mr = str(z_record[1][0])
  63. mi = str(z_record[1][1])
  64. mp.mp.dps = 20
  65. z = mp.mpf(x)*mp.mpc(mr, mi)
  66. n_list = get_n_list(z, max_num_elements_of_n_list)
  67. print(z, n_list)
  68. for n in n_list:
  69. mp.mp.dps = 20
  70. old_z_string = get_test_data_nlist(z_record, output_dps, n)
  71. mp.mp.dps = 37
  72. new_z_string = get_test_data_nlist(z_record, output_dps, n)
  73. while old_z_string != new_z_string:
  74. new_dps = int(mp.mp.dps * 1.41)
  75. if new_dps > 300: break
  76. mp.mp.dps = new_dps
  77. print("New dps = ", mp.mp.dps)
  78. old_z_string = new_z_string
  79. new_z_string = get_test_data_nlist(z_record, output_dps, n)
  80. output_str += new_z_string
  81. output_str += '};\n'
  82. return output_str
  83. def main ():
  84. sf_evals = update_special_functions_evaluations(filename='test_spec_functions_data.hpp',
  85. complex_arguments = mia.complex_arguments)
  86. sf_evals.read_evaluated_data()
  87. print(sf_evals.evaluated_data[1].list_to_parse)
  88. output_dps = 16
  89. max_num_elements_of_nlist = 51
  90. # out_filename = 'test_spec_functions_data.hpp'
  91. # output_str = get_test_data(mia.complex_arguments, output_dps, max_num_elements_of_nlist)
  92. # with open(out_filename, 'w') as out_file:
  93. # out_file.write(output_str)
  94. main()