Float16, Float32, Float64, BigFloat

Author: jbrea

src/owent.jl

@@ -99,25 +99,7 @@ function owent(h::T, a::T) where {T <: Real}
     # calculate Owen's T
 
     if a ≤ T(0.999999)
-
-
-        x = (-T(0.9987710072524261), -T(0.9935301722663508), -T(0.9841245837228269), -T(0.9705915925462473), -T(0.9529877031604309), -T(0.9313866907065543), -T(0.9058791367155696)
-        , -T(0.8765720202742479), -T(0.8435882616243935), -T(0.8070662040294426), -T(0.7671590325157404), -T(0.7240341309238146), -T(0.6778723796326639), -T(0.6288673967765136)
-        , -T(0.5772247260839727), -T(0.5231609747222331), -T(0.4669029047509584), -T(0.4086864819907167), -T(0.34875588629216075), -T(0.28736248735545555), -T(0.22476379039468905)
-        , -T(0.1612223560688917), -T(0.0970046992094627), -T(0.03238017096286937), T(0.03238017096286937), T(0.0970046992094627), T(0.1612223560688917), T(0.22476379039468905)
-        , T(0.28736248735545555), T(0.34875588629216075), T(0.4086864819907167), T(0.4669029047509584), T(0.5231609747222331), T(0.5772247260839727), T(0.6288673967765136)
-        , T(0.6778723796326639), T(0.7240341309238146), T(0.7671590325157404), T(0.8070662040294426), T(0.8435882616243935), T(0.8765720202742479), T(0.9058791367155696)
-        , T(0.9313866907065543), T(0.9529877031604309), T(0.9705915925462473), T(0.9841245837228269), T(0.9935301722663508), T(0.9987710072524261))
-
-        w = (T(0.0031533460523059122), T(0.0073275539012762885), T(0.011477234579234613), T(0.015579315722943824), T(0.01961616045735561), T(0.023570760839324363)
-        , T(0.027426509708356944), T(0.031167227832798003), T(0.03477722256477052), T(0.038241351065830737), T(0.04154508294346467), T(0.0446745608566943), T(0.04761665849249045)
-        , T(0.05035903555385445), T(0.05289018948519363), T(0.055199503699984116), T(0.05727729210040322), T(0.05911483969839564), T(0.06070443916589387), T(0.06203942315989268)
-        , T(0.06311419228625402), T(0.06392423858464813), T(0.06446616443594998), T(0.06473769681268386), T(0.06473769681268386), T(0.06446616443594998), T(0.06392423858464813)
-        , T(0.06311419228625402), T(0.06203942315989268), T(0.06070443916589387), T(0.05911483969839564), T(0.05727729210040322), T(0.055199503699984116)
-        , T(0.05289018948519363), T(0.05035903555385445), T(0.04761665849249045), T(0.0446745608566943), T(0.04154508294346467), T(0.038241351065830737), T(0.03477722256477052)
-        , T(0.031167227832798003), T(0.027426509708356944), T(0.023570760839324363), T(0.01961616045735561), T(0.015579315722943824), T(0.011477234579234613), T(0.0073275539012762885)
-        , T(0.0031533460523059122))
-
+        x, w = gauss_legendre(T)
         return sum(w .* t2.(h, a, x))
     else
         # a > 0.999999, T6 from paper (quadrature using QuadGK would also work, but be slower)
@@ -133,3 +115,19 @@ end
 t2(h::T, a, x) where T = T(inv4π)*a*exp(-T(0.5)*(h^2)*(one(h)+(a*x)^2))/(one(h)+(a*x)^2)
 
 owent(h::Real, a::Real) = owent(promote(h,a)...)
+
+# 48-point Gauss-Legendre quadrature (Arblib.hypgeom_legendre_p_ui_root!)
+gauss_legendre(::Type{Float64}) =
+    (0.9987710072524261, 0.9935301722663508, 0.9841245837228269, 0.9705915925462473, 0.9529877031604309, 0.9313866907065543, 0.9058791367155696, 0.8765720202742479, 0.8435882616243935, 0.8070662040294426, 0.7671590325157404, 0.7240341309238146, 0.6778723796326639, 0.6288673967765136, 0.5772247260839727, 0.523160974722233, 0.4669029047509584, 0.4086864819907167, 0.34875588629216075, 0.28736248735545555, 0.22476379039468905, 0.1612223560688917, 0.0970046992094627, 0.03238017096286936, -0.03238017096286936, -0.0970046992094627, -0.1612223560688917, -0.22476379039468905, -0.28736248735545555, -0.34875588629216075, -0.4086864819907167, -0.4669029047509584, -0.523160974722233, -0.5772247260839727, -0.6288673967765136, -0.6778723796326639, -0.7240341309238146, -0.7671590325157404, -0.8070662040294426, -0.8435882616243935, -0.8765720202742479, -0.9058791367155696, -0.9313866907065543, -0.9529877031604309, -0.9705915925462473, -0.9841245837228269, -0.9935301722663508, -0.9987710072524261),
+    (0.0031533460523058385, 0.0073275539012762625, 0.01147723457923454, 0.015579315722943849, 0.01961616045735553, 0.02357076083932438, 0.027426509708356948, 0.03116722783279809, 0.03477722256477044, 0.03824135106583071, 0.04154508294346475, 0.04467456085669428, 0.04761665849249048, 0.05035903555385447, 0.05289018948519367, 0.055199503699984165, 0.057277292100403214, 0.059114839698395635, 0.06070443916589388, 0.062039423159892665, 0.06311419228625402, 0.06392423858464819, 0.06446616443595009, 0.06473769681268392, 0.06473769681268392, 0.06446616443595009, 0.06392423858464819, 0.06311419228625402, 0.062039423159892665, 0.06070443916589388, 0.059114839698395635, 0.057277292100403214, 0.055199503699984165, 0.05289018948519367, 0.05035903555385447, 0.04761665849249048, 0.04467456085669428, 0.04154508294346475, 0.03824135106583071, 0.03477722256477044, 0.03116722783279809, 0.027426509708356948, 0.02357076083932438, 0.01961616045735553, 0.015579315722943849, 0.01147723457923454, 0.0073275539012762625, 0.0031533460523058385)
+# 24-point Gauss-Legendre quadrature (Arblib.hypgeom_legendre_p_ui_root!)
+gauss_legendre(::Type{Float32}) =
+    (0.9951872f0, 0.9747286f0, 0.93827456f0, 0.88641554f0, 0.82000196f0, 0.74012417f0, 0.64809364f0, 0.5454215f0, 0.43379351f0, 0.31504267f0, 0.19111887f0, 0.064056896f0, -0.064056896f0, -0.19111887f0, -0.31504267f0, -0.43379351f0, -0.5454215f0, -0.64809364f0, -0.74012417f0, -0.82000196f0, -0.88641554f0, -0.93827456f0, -0.9747286f0, -0.9951872f0),
+    (0.01234123f0, 0.02853139f0, 0.044277437f0, 0.059298586f0, 0.07334648f0, 0.086190164f0, 0.097618654f0, 0.10744427f0, 0.115505666f0, 0.12167047f0, 0.12583746f0, 0.1279382f0, 0.1279382f0, 0.12583746f0, 0.12167047f0, 0.115505666f0, 0.10744427f0, 0.097618654f0, 0.086190164f0, 0.07334648f0, 0.059298586f0, 0.044277437f0, 0.02853139f0, 0.01234123f0)
+gauss_legendre(::Type{Float16}) =
+    (Float16(0.9814), Float16(0.9043), Float16(0.77), Float16(0.5874), Float16(0.368), Float16(0.1252), Float16(-0.1252), Float16(-0.368), Float16(-0.5874), Float16(-0.77), Float16(-0.9043), Float16(-0.9814)),
+    (Float16(0.04718), Float16(0.10693), Float16(0.16), Float16(0.2031), Float16(0.2335), Float16(0.2491), Float16(0.2491), Float16(0.2335), Float16(0.2031), Float16(0.16), Float16(0.10693), Float16(0.04718))
+# 48-point Gauss-Legendre quadrature (Arblib.hypgeom_legendre_p_ui_root!)
+gauss_legendre(::Type{BigFloat}) =
+    (big"0.9987710072524261186005414915631136400889376502767210386129404813754588436074878", big"0.9935301722663507575479287508490741183566147495946719296171518380987546182067713", big"0.9841245837228268577445836000265988305892392234173847299576501679855297780009794", big"0.9705915925462472504614119838006600573024339116308837060283723521653233091284874", big"0.9529877031604308607229606660257183432085413318239187368639476034939458705853333", big"0.9313866907065543331141743801016012677199970856189504298706048642530730422171056", big"0.9058791367155696728220748356710117883122621998274108453524854254710168231209838", big"0.8765720202742478859056935548050967545616485337299619927478757518746727101403824", big"0.8435882616243935307110898445196560498708870117375524015149131998988410546898503", big"0.8070662040294426270825530430245384459730130294604153865758629418121821540044232", big"0.7671590325157403392538554375229690536226423308482073722351285886640508368078524", big"0.7240341309238146546744822334936652465850928122807223627293663025733514606200864", big"0.6778723796326639052118512806759090588499546790260486130710406429754946468798164", big"0.6288673967765136239951649330699946520249089997901617709817329945195319139770715", big"0.5772247260839727038178092385404787728539972861401955280523973994277369963343583", big"0.5231609747222330336782258691375085262891876218118841075802295472194144547473473", big"0.4669029047509584045449288616507985092368121042585169441818691951347943934426029", big"0.4086864819907167299162254958146332864599228429948880647711509833256205384841253", big"0.3487558862921607381598179372704079161343096499683925760321229677812815940686729", big"0.2873624873554555767358864613167976878515583058010397789085000321689998442687597", big"0.224763790394689061224865440174692277438561804041654806164742641045181941897513", big"0.161222356068891718056437390783497694774374379741895117703242637556516342099581", big"0.09700469920946269893005395585362452015273622930093698643058076594480403626262214", big"0.03238017096286936203332224315213444204596280236151809242500322001737781920338223", big"-0.03238017096286936203332224315213444204596280236151809242500322001737781920338223", big"-0.09700469920946269893005395585362452015273622930093698643058076594480403626262214", big"-0.161222356068891718056437390783497694774374379741895117703242637556516342099581", big"-0.224763790394689061224865440174692277438561804041654806164742641045181941897513", big"-0.2873624873554555767358864613167976878515583058010397789085000321689998442687597", big"-0.3487558862921607381598179372704079161343096499683925760321229677812815940686729", big"-0.4086864819907167299162254958146332864599228429948880647711509833256205384841253", big"-0.4669029047509584045449288616507985092368121042585169441818691951347943934426029", big"-0.5231609747222330336782258691375085262891876218118841075802295472194144547473473", big"-0.5772247260839727038178092385404787728539972861401955280523973994277369963343583", big"-0.6288673967765136239951649330699946520249089997901617709817329945195319139770715", big"-0.6778723796326639052118512806759090588499546790260486130710406429754946468798164", big"-0.7240341309238146546744822334936652465850928122807223627293663025733514606200864", big"-0.7671590325157403392538554375229690536226423308482073722351285886640508368078524", big"-0.8070662040294426270825530430245384459730130294604153865758629418121821540044232", big"-0.8435882616243935307110898445196560498708870117375524015149131998988410546898503", big"-0.8765720202742478859056935548050967545616485337299619927478757518746727101403824", big"-0.9058791367155696728220748356710117883122621998274108453524854254710168231209838", big"-0.9313866907065543331141743801016012677199970856189504298706048642530730422171056", big"-0.9529877031604308607229606660257183432085413318239187368639476034939458705853333", big"-0.9705915925462472504614119838006600573024339116308837060283723521653233091284874", big"-0.9841245837228268577445836000265988305892392234173847299576501679855297780009794", big"-0.9935301722663507575479287508490741183566147495946719296171518380987546182067713", big"-0.9987710072524261186005414915631136400889376502767210386129404813754588436074878"),
+    (big"0.003153346052305838632677311543891487578283938831693622295209493250319586438316842", big"0.007327553901276262102383979621786550058707902559201353274881829548806980072502799", big"0.01147723457923453948959266760909162808642050630874764065376681674103503658508731", big"0.01557931572294384872817695583446031397637626899155246951309343105269243335619984", big"0.01961616045735552781446071965221270969581303773413223918112083050740924629812146", big"0.02357076083932437914051930137844923022172973852218859873423906486456506379639118", big"0.0274265097083569482000738362625058204511841551616509759972809374993765019410236", big"0.03116722783279808890206575684635441945428534148356953550954371886143141262424302", big"0.03477722256477043889254858596380241059728139690706809871800663617967672335903626", big"0.03824135106583070631721725652371561786382396835498228892925819103405053922410909", big"0.04154508294346474921405882236106479775347282603403806308273482122272582562965843", big"0.04467456085669428041944858712585039498846278686250200843292144633919149051230188", big"0.04761665849249047482590662347892983015799806674344968539676989627880988507905503", big"0.05035903555385447495780761908786560603299409302590633069379205724693441466024811", big"0.05289018948519366709550505626469891466172648563310918638649123384829276249063296", big"0.05519950369998416286820349519163543900445092560756100054805625793058523675145725", big"0.05727729210040321570515023468470057624152712300411207753884993747681745421856388", big"0.05911483969839563574647481743351991065965560255705499855629113348583514270048057", big"0.06070443916589388005296923202782047788526086425647775511151144466063789427975123", big"0.06203942315989266390419778413759851830638339966509146156903781450273903590161649", big"0.06311419228625402565712602275023331812741364337110079121114724790803811921086588", big"0.06392423858464818662390620182551540891897408498264299989087420749955378258611148", big"0.06446616443595008220650419365770506572569192445553030876055845653739235337295456", big"0.06473769681268392250302493873659155355208191894663651001456309552308307891126462", big"0.06473769681268392250302493873659155355208191894663651001456309552308307891126462", big"0.06446616443595008220650419365770506572569192445553030876055845653739235337295456", big"0.06392423858464818662390620182551540891897408498264299989087420749955378258611148", big"0.06311419228625402565712602275023331812741364337110079121114724790803811921086588", big"0.06203942315989266390419778413759851830638339966509146156903781450273903590161649", big"0.06070443916589388005296923202782047788526086425647775511151144466063789427975123", big"0.05911483969839563574647481743351991065965560255705499855629113348583514270048057", big"0.05727729210040321570515023468470057624152712300411207753884993747681745421856388", big"0.05519950369998416286820349519163543900445092560756100054805625793058523675145725", big"0.05289018948519366709550505626469891466172648563310918638649123384829276249063296", big"0.05035903555385447495780761908786560603299409302590633069379205724693441466024811", big"0.04761665849249047482590662347892983015799806674344968539676989627880988507905503", big"0.04467456085669428041944858712585039498846278686250200843292144633919149051230188", big"0.04154508294346474921405882236106479775347282603403806308273482122272582562965843", big"0.03824135106583070631721725652371561786382396835498228892925819103405053922410909", big"0.03477722256477043889254858596380241059728139690706809871800663617967672335903626", big"0.03116722783279808890206575684635441945428534148356953550954371886143141262424302", big"0.0274265097083569482000738362625058204511841551616509759972809374993765019410236", big"0.02357076083932437914051930137844923022172973852218859873423906486456506379639118", big"0.01961616045735552781446071965221270969581303773413223918112083050740924629812146", big"0.01557931572294384872817695583446031397637626899155246951309343105269243335619984", big"0.01147723457923453948959266760909162808642050630874764065376681674103503658508731", big"0.007327553901276262102383979621786550058707902559201353274881829548806980072502799", big"0.003153346052305838632677311543891487578283938831693622295209493250319586438316842")