tsai
/
mathnet-numerics
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Merge pull request #369 from Jdsleppy/erf-performance 

Make the coefficient arrays for ErfImp and ErfInvImp static fields.
pull/374/head
Christoph Ruegg 11 years ago
parent
aad5c5303e c3b9f3a3f8
commit
31147883ad
1 changed files with 242 additions and 65 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 307
      src/Numerics/SpecialFunctions/Erf.cs

307
src/Numerics/SpecialFunctions/Erf.cs
View File

@ -48,6 +48,226 @@ namespace MathNet.Numerics
/// </summary>
public partial class SpecialFunctions
{
/// <summary>
/// **************************************
/// COEFFICIENTS FOR METHOD ErfImp *
/// **************************************
/// </summary>
/// <summary> Polynomial coefficients for a numerator of ErfImp
/// calculation for Erf(x) in the interval [1e-10, 0.5].
/// </summary>
private static readonly double[] ErfImpAn = { 0.00337916709551257388990745, -0.00073695653048167948530905, -0.374732337392919607868241, 0.0817442448733587196071743, -0.0421089319936548595203468, 0.0070165709512095756344528, -0.00495091255982435110337458, 0.000871646599037922480317225 };
/// <summary> Polynomial coefficients for adenominator of ErfImp
/// calculation for Erf(x) in the interval [1e-10, 0.5].
/// </summary>
private static readonly double[] ErfImpAd = { 1, -0.218088218087924645390535, 0.412542972725442099083918, -0.0841891147873106755410271, 0.0655338856400241519690695, -0.0120019604454941768171266, 0.00408165558926174048329689, -0.000615900721557769691924509 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [0.5, 0.75].
/// </summary>
private static readonly double[] ErfImpBn = { -0.0361790390718262471360258, 0.292251883444882683221149, 0.281447041797604512774415, 0.125610208862766947294894, 0.0274135028268930549240776, 0.00250839672168065762786937 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [0.5, 0.75].
/// </summary>
private static readonly double[] ErfImpBd = { 1, 1.8545005897903486499845, 1.43575803037831418074962, 0.582827658753036572454135, 0.124810476932949746447682, 0.0113724176546353285778481 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [0.75, 1.25].
/// </summary>
private static readonly double[] ErfImpCn = { -0.0397876892611136856954425, 0.153165212467878293257683, 0.191260295600936245503129, 0.10276327061989304213645, 0.029637090615738836726027, 0.0046093486780275489468812, 0.000307607820348680180548455 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [0.75, 1.25].
/// </summary>
private static readonly double[] ErfImpCd = { 1, 1.95520072987627704987886, 1.64762317199384860109595, 0.768238607022126250082483, 0.209793185936509782784315, 0.0319569316899913392596356, 0.00213363160895785378615014 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [1.25, 2.25].
/// </summary>
private static readonly double[] ErfImpDn = { -0.0300838560557949717328341, 0.0538578829844454508530552, 0.0726211541651914182692959, 0.0367628469888049348429018, 0.00964629015572527529605267, 0.00133453480075291076745275, 0.778087599782504251917881e-4 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [1.25, 2.25].
/// </summary>
private static readonly double[] ErfImpDd = { 1, 1.75967098147167528287343, 1.32883571437961120556307, 0.552528596508757581287907, 0.133793056941332861912279, 0.0179509645176280768640766, 0.00104712440019937356634038, -0.106640381820357337177643e-7 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [2.25, 3.5].
/// </summary>
private static readonly double[] ErfImpEn = { -0.0117907570137227847827732, 0.014262132090538809896674, 0.0202234435902960820020765, 0.00930668299990432009042239, 0.00213357802422065994322516, 0.00025022987386460102395382, 0.120534912219588189822126e-4 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [2.25, 3.5].
/// </summary>
private static readonly double[] ErfImpEd = { 1, 1.50376225203620482047419, 0.965397786204462896346934, 0.339265230476796681555511, 0.0689740649541569716897427, 0.00771060262491768307365526, 0.000371421101531069302990367 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [3.5, 5.25].
/// </summary>
private static readonly double[] ErfImpFn = { -0.00546954795538729307482955, 0.00404190278731707110245394, 0.0054963369553161170521356, 0.00212616472603945399437862, 0.000394984014495083900689956, 0.365565477064442377259271e-4, 0.135485897109932323253786e-5 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [3.5, 5.25].
/// </summary>
private static readonly double[] ErfImpFd = { 1, 1.21019697773630784832251, 0.620914668221143886601045, 0.173038430661142762569515, 0.0276550813773432047594539, 0.00240625974424309709745382, 0.891811817251336577241006e-4, -0.465528836283382684461025e-11 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [5.25, 8].
/// </summary>
private static readonly double[] ErfImpGn = { -0.00270722535905778347999196, 0.0013187563425029400461378, 0.00119925933261002333923989, 0.00027849619811344664248235, 0.267822988218331849989363e-4, 0.923043672315028197865066e-6 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [5.25, 8].
/// </summary>
private static readonly double[] ErfImpGd = { 1, 0.814632808543141591118279, 0.268901665856299542168425, 0.0449877216103041118694989, 0.00381759663320248459168994, 0.000131571897888596914350697, 0.404815359675764138445257e-11 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [8, 11.5].
/// </summary>
private static readonly double[] ErfImpHn = { -0.00109946720691742196814323, 0.000406425442750422675169153, 0.000274499489416900707787024, 0.465293770646659383436343e-4, 0.320955425395767463401993e-5, 0.778286018145020892261936e-7 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [8, 11.5].
/// </summary>
private static readonly double[] ErfImpHd = { 1, 0.588173710611846046373373, 0.139363331289409746077541, 0.0166329340417083678763028, 0.00100023921310234908642639, 0.24254837521587225125068e-4 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [11.5, 17].
/// </summary>
private static readonly double[] ErfImpIn = { -0.00056907993601094962855594, 0.000169498540373762264416984, 0.518472354581100890120501e-4, 0.382819312231928859704678e-5, 0.824989931281894431781794e-7 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [11.5, 17].
/// </summary>
private static readonly double[] ErfImpId = { 1, 0.339637250051139347430323, 0.043472647870310663055044, 0.00248549335224637114641629, 0.535633305337152900549536e-4, -0.117490944405459578783846e-12 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [17, 24].
/// </summary>
private static readonly double[] ErfImpJn = { -0.000241313599483991337479091, 0.574224975202501512365975e-4, 0.115998962927383778460557e-4, 0.581762134402593739370875e-6, 0.853971555085673614607418e-8 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [17, 24].
/// </summary>
private static readonly double[] ErfImpJd = { 1, 0.233044138299687841018015, 0.0204186940546440312625597, 0.000797185647564398289151125, 0.117019281670172327758019e-4 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [24, 38].
/// </summary>
private static readonly double[] ErfImpKn = { -0.000146674699277760365803642, 0.162666552112280519955647e-4, 0.269116248509165239294897e-5, 0.979584479468091935086972e-7, 0.101994647625723465722285e-8 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [24, 38].
/// </summary>
private static readonly double[] ErfImpKd = { 1, 0.165907812944847226546036, 0.0103361716191505884359634, 0.000286593026373868366935721, 0.298401570840900340874568e-5 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [38, 60].
/// </summary>
private static readonly double[] ErfImpLn = { -0.583905797629771786720406e-4, 0.412510325105496173512992e-5, 0.431790922420250949096906e-6, 0.993365155590013193345569e-8, 0.653480510020104699270084e-10 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [38, 60].
/// </summary>
private static readonly double[] ErfImpLd = { 1, 0.105077086072039915406159, 0.00414278428675475620830226, 0.726338754644523769144108e-4, 0.477818471047398785369849e-6 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [60, 85].
/// </summary>
private static readonly double[] ErfImpMn = { -0.196457797609229579459841e-4, 0.157243887666800692441195e-5, 0.543902511192700878690335e-7, 0.317472492369117710852685e-9 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [60, 85].
/// </summary>
private static readonly double[] ErfImpMd = { 1, 0.052803989240957632204885, 0.000926876069151753290378112, 0.541011723226630257077328e-5, 0.535093845803642394908747e-15 };
/// <summary> Polynomial coefficients for a numerator in ErfImp
/// calculation for Erfc(x) in the interval [85, 110].
/// </summary>
private static readonly double[] ErfImpNn = { -0.789224703978722689089794e-5, 0.622088451660986955124162e-6, 0.145728445676882396797184e-7, 0.603715505542715364529243e-10 };
/// <summary> Polynomial coefficients for a denominator in ErfImp
/// calculation for Erfc(x) in the interval [85, 110].
/// </summary>
private static readonly double[] ErfImpNd = { 1, 0.0375328846356293715248719, 0.000467919535974625308126054, 0.193847039275845656900547e-5 };
/// <summary>
/// **************************************
/// COEFFICIENTS FOR METHOD ErfInvImp *
/// **************************************
/// </summary>
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0, 0.5].
/// </summary>
private static readonly double[] ErvInvImpAn = { -0.000508781949658280665617, -0.00836874819741736770379, 0.0334806625409744615033, -0.0126926147662974029034, -0.0365637971411762664006, 0.0219878681111168899165, 0.00822687874676915743155, -0.00538772965071242932965 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0, 0.5].
/// </summary>
private static readonly double[] ErvInvImpAd = { 1, -0.970005043303290640362, -1.56574558234175846809, 1.56221558398423026363, 0.662328840472002992063, -0.71228902341542847553, -0.0527396382340099713954, 0.0795283687341571680018, -0.00233393759374190016776, 0.000886216390456424707504 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.5, 0.75].
/// </summary>
private static readonly double[] ErvInvImpBn = { -0.202433508355938759655, 0.105264680699391713268, 8.37050328343119927838, 17.6447298408374015486, -18.8510648058714251895, -44.6382324441786960818, 17.445385985570866523, 21.1294655448340526258, -3.67192254707729348546 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.5, 0.75].
/// </summary>
private static readonly double[] ErvInvImpBd = { 1, 6.24264124854247537712, 3.9713437953343869095, -28.6608180499800029974, -20.1432634680485188801, 48.5609213108739935468, 10.8268667355460159008, -22.6436933413139721736, 1.72114765761200282724 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x less than 3.
/// </summary>
private static readonly double[] ErvInvImpCn = { -0.131102781679951906451, -0.163794047193317060787, 0.117030156341995252019, 0.387079738972604337464, 0.337785538912035898924, 0.142869534408157156766, 0.0290157910005329060432, 0.00214558995388805277169, -0.679465575181126350155e-6, 0.285225331782217055858e-7, -0.681149956853776992068e-9 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x less than 3.
/// </summary>
private static readonly double[] ErvInvImpCd = { 1, 3.46625407242567245975, 5.38168345707006855425, 4.77846592945843778382, 2.59301921623620271374, 0.848854343457902036425, 0.152264338295331783612, 0.01105924229346489121 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 3 and 6.
/// </summary>
private static readonly double[] ErvInvImpDn = { -0.0350353787183177984712, -0.00222426529213447927281, 0.0185573306514231072324, 0.00950804701325919603619, 0.00187123492819559223345, 0.000157544617424960554631, 0.460469890584317994083e-5, -0.230404776911882601748e-9, 0.266339227425782031962e-11 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 3 and 6.
/// </summary>
private static readonly double[] ErvInvImpDd = { 1, 1.3653349817554063097, 0.762059164553623404043, 0.220091105764131249824, 0.0341589143670947727934, 0.00263861676657015992959, 0.764675292302794483503e-4 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 6 and 18.
/// </summary>
private static readonly double[] ErvInvImpEn = { -0.0167431005076633737133, -0.00112951438745580278863, 0.00105628862152492910091, 0.000209386317487588078668, 0.149624783758342370182e-4, 0.449696789927706453732e-6, 0.462596163522878599135e-8, -0.281128735628831791805e-13, 0.99055709973310326855e-16 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 6 and 18.
/// </summary>
private static readonly double[] ErvInvImpEd = { 1, 0.591429344886417493481, 0.138151865749083321638, 0.0160746087093676504695, 0.000964011807005165528527, 0.275335474764726041141e-4, 0.282243172016108031869e-6 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 18 and 44.
/// </summary>
private static readonly double[] ErvInvImpFn = { -0.0024978212791898131227, -0.779190719229053954292e-5, 0.254723037413027451751e-4, 0.162397777342510920873e-5, 0.396341011304801168516e-7, 0.411632831190944208473e-9, 0.145596286718675035587e-11, -0.116765012397184275695e-17 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x between 18 and 44.
/// </summary>
private static readonly double[] ErvInvImpFd = { 1, 0.207123112214422517181, 0.0169410838120975906478, 0.000690538265622684595676, 0.145007359818232637924e-4, 0.144437756628144157666e-6, 0.509761276599778486139e-9 };
/// <summary> Polynomial coefficients for a numerator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x greater than 44.
/// </summary>
private static readonly double[] ErvInvImpGn = { -0.000539042911019078575891, -0.28398759004727721098e-6, 0.899465114892291446442e-6, 0.229345859265920864296e-7, 0.225561444863500149219e-9, 0.947846627503022684216e-12, 0.135880130108924861008e-14, -0.348890393399948882918e-21 };
/// <summary> Polynomial coefficients for a denominator of ErfInvImp
/// calculation for Erf^-1(z) in the interval [0.75, 1] with x greater than 44.
/// </summary>
private static readonly double[] ErvInvImpGd = { 1, 0.0845746234001899436914, 0.00282092984726264681981, 0.468292921940894236786e-4, 0.399968812193862100054e-6, 0.161809290887904476097e-8, 0.231558608310259605225e-11 };
/// <summary>Calculates the error function.</summary>
/// <param name="x">The value to evaluate.</param>
/// <returns>the error function evaluated at given value.</returns>
@ -199,10 +419,7 @@ namespace MathNet.Numerics
else
{
// Worst case absolute error found: 6.688618532e-21
double[] nc = { 0.00337916709551257388990745, -0.00073695653048167948530905, -0.374732337392919607868241, 0.0817442448733587196071743, -0.0421089319936548595203468, 0.0070165709512095756344528, -0.00495091255982435110337458, 0.000871646599037922480317225 };
double[] dc = { 1, -0.218088218087924645390535, 0.412542972725442099083918, -0.0841891147873106755410271, 0.0655338856400241519690695, -0.0120019604454941768171266, 0.00408165558926174048329689, -0.000615900721557769691924509 };
result = (z*1.125) + (z*Evaluate.Polynomial(z, nc)/Evaluate.Polynomial(z, dc));
result = (z*1.125) + (z*Evaluate.Polynomial(z, ErfImpAn)/Evaluate.Polynomial(z, ErfImpAd));
}
}
else if ((z < 110) || ((z < 110) && invert))
@ -213,105 +430,79 @@ namespace MathNet.Numerics
if (z < 0.75)
{
// Worst case absolute error found: 5.582813374e-21
double[] nc = { -0.0361790390718262471360258, 0.292251883444882683221149, 0.281447041797604512774415, 0.125610208862766947294894, 0.0274135028268930549240776, 0.00250839672168065762786937 };
double[] dc = { 1, 1.8545005897903486499845, 1.43575803037831418074962, 0.582827658753036572454135, 0.124810476932949746447682, 0.0113724176546353285778481 };
r = Evaluate.Polynomial(z - 0.5, nc)/Evaluate.Polynomial(z - 0.5, dc);
r = Evaluate.Polynomial(z - 0.5, ErfImpBn)/Evaluate.Polynomial(z - 0.5, ErfImpBd);
b = 0.3440242112F;
}
else if (z < 1.25)
{
// Worst case absolute error found: 4.01854729e-21
double[] nc = { -0.0397876892611136856954425, 0.153165212467878293257683, 0.191260295600936245503129, 0.10276327061989304213645, 0.029637090615738836726027, 0.0046093486780275489468812, 0.000307607820348680180548455 };
double[] dc = { 1, 1.95520072987627704987886, 1.64762317199384860109595, 0.768238607022126250082483, 0.209793185936509782784315, 0.0319569316899913392596356, 0.00213363160895785378615014 };
r = Evaluate.Polynomial(z - 0.75, nc)/Evaluate.Polynomial(z - 0.75, dc);
r = Evaluate.Polynomial(z - 0.75, ErfImpCn)/Evaluate.Polynomial(z - 0.75, ErfImpCd);
b = 0.419990927F;
}
else if (z < 2.25)
{
// Worst case absolute error found: 2.866005373e-21
double[] nc = { -0.0300838560557949717328341, 0.0538578829844454508530552, 0.0726211541651914182692959, 0.0367628469888049348429018, 0.00964629015572527529605267, 0.00133453480075291076745275, 0.778087599782504251917881e-4 };
double[] dc = { 1, 1.75967098147167528287343, 1.32883571437961120556307, 0.552528596508757581287907, 0.133793056941332861912279, 0.0179509645176280768640766, 0.00104712440019937356634038, -0.106640381820357337177643e-7 };
r = Evaluate.Polynomial(z - 1.25, nc)/Evaluate.Polynomial(z - 1.25, dc);
r = Evaluate.Polynomial(z - 1.25, ErfImpDn)/Evaluate.Polynomial(z - 1.25, ErfImpDd);
b = 0.4898625016F;
}
else if (z < 3.5)
{
// Worst case absolute error found: 1.045355789e-21
double[] nc = { -0.0117907570137227847827732, 0.014262132090538809896674, 0.0202234435902960820020765, 0.00930668299990432009042239, 0.00213357802422065994322516, 0.00025022987386460102395382, 0.120534912219588189822126e-4 };
double[] dc = { 1, 1.50376225203620482047419, 0.965397786204462896346934, 0.339265230476796681555511, 0.0689740649541569716897427, 0.00771060262491768307365526, 0.000371421101531069302990367 };
r = Evaluate.Polynomial(z - 2.25, nc)/Evaluate.Polynomial(z - 2.25, dc);
r = Evaluate.Polynomial(z - 2.25, ErfImpEn) /Evaluate.Polynomial(z - 2.25, ErfImpEd);
b = 0.5317370892F;
}
else if (z < 5.25)
{
// Worst case absolute error found: 8.300028706e-22
double[] nc = { -0.00546954795538729307482955, 0.00404190278731707110245394, 0.0054963369553161170521356, 0.00212616472603945399437862, 0.000394984014495083900689956, 0.365565477064442377259271e-4, 0.135485897109932323253786e-5 };
double[] dc = { 1, 1.21019697773630784832251, 0.620914668221143886601045, 0.173038430661142762569515, 0.0276550813773432047594539, 0.00240625974424309709745382, 0.891811817251336577241006e-4, -0.465528836283382684461025e-11 };
r = Evaluate.Polynomial(z - 3.5, nc)/Evaluate.Polynomial(z - 3.5, dc);
r = Evaluate.Polynomial(z - 3.5, ErfImpFn) /Evaluate.Polynomial(z - 3.5, ErfImpFd);
b = 0.5489973426F;
}
else if (z < 8)
{
// Worst case absolute error found: 1.700157534e-21
double[] nc = { -0.00270722535905778347999196, 0.0013187563425029400461378, 0.00119925933261002333923989, 0.00027849619811344664248235, 0.267822988218331849989363e-4, 0.923043672315028197865066e-6 };
double[] dc = { 1, 0.814632808543141591118279, 0.268901665856299542168425, 0.0449877216103041118694989, 0.00381759663320248459168994, 0.000131571897888596914350697, 0.404815359675764138445257e-11 };
r = Evaluate.Polynomial(z - 5.25, nc)/Evaluate.Polynomial(z - 5.25, dc);
r = Evaluate.Polynomial(z - 5.25, ErfImpGn) /Evaluate.Polynomial(z - 5.25, ErfImpGd);
b = 0.5571740866F;
}
else if (z < 11.5)
{
// Worst case absolute error found: 3.002278011e-22
double[] nc = { -0.00109946720691742196814323, 0.000406425442750422675169153, 0.000274499489416900707787024, 0.465293770646659383436343e-4, 0.320955425395767463401993e-5, 0.778286018145020892261936e-7 };
double[] dc = { 1, 0.588173710611846046373373, 0.139363331289409746077541, 0.0166329340417083678763028, 0.00100023921310234908642639, 0.24254837521587225125068e-4 };
r = Evaluate.Polynomial(z - 8, nc)/Evaluate.Polynomial(z - 8, dc);
r = Evaluate.Polynomial(z - 8, ErfImpHn) /Evaluate.Polynomial(z - 8, ErfImpHd);
b = 0.5609807968F;
}
else if (z < 17)
{
// Worst case absolute error found: 6.741114695e-21
double[] nc = { -0.00056907993601094962855594, 0.000169498540373762264416984, 0.518472354581100890120501e-4, 0.382819312231928859704678e-5, 0.824989931281894431781794e-7 };
double[] dc = { 1, 0.339637250051139347430323, 0.043472647870310663055044, 0.00248549335224637114641629, 0.535633305337152900549536e-4, -0.117490944405459578783846e-12 };
r = Evaluate.Polynomial(z - 11.5, nc)/Evaluate.Polynomial(z - 11.5, dc);
r = Evaluate.Polynomial(z - 11.5, ErfImpIn) /Evaluate.Polynomial(z - 11.5, ErfImpId);
b = 0.5626493692F;
}
else if (z < 24)
{
// Worst case absolute error found: 7.802346984e-22
double[] nc = { -0.000241313599483991337479091, 0.574224975202501512365975e-4, 0.115998962927383778460557e-4, 0.581762134402593739370875e-6, 0.853971555085673614607418e-8 };
double[] dc = { 1, 0.233044138299687841018015, 0.0204186940546440312625597, 0.000797185647564398289151125, 0.117019281670172327758019e-4 };
r = Evaluate.Polynomial(z - 17, nc)/Evaluate.Polynomial(z - 17, dc);
r = Evaluate.Polynomial(z - 17, ErfImpJn) /Evaluate.Polynomial(z - 17, ErfImpJd);
b = 0.5634598136F;
}
else if (z < 38)
{
// Worst case absolute error found: 2.414228989e-22
double[] nc = { -0.000146674699277760365803642, 0.162666552112280519955647e-4, 0.269116248509165239294897e-5, 0.979584479468091935086972e-7, 0.101994647625723465722285e-8 };
double[] dc = { 1, 0.165907812944847226546036, 0.0103361716191505884359634, 0.000286593026373868366935721, 0.298401570840900340874568e-5 };
r = Evaluate.Polynomial(z - 24, nc)/Evaluate.Polynomial(z - 24, dc);
r = Evaluate.Polynomial(z - 24, ErfImpKn) /Evaluate.Polynomial(z - 24, ErfImpKd);
b = 0.5638477802F;
}
else if (z < 60)
{
// Worst case absolute error found: 5.896543869e-24
double[] nc = { -0.583905797629771786720406e-4, 0.412510325105496173512992e-5, 0.431790922420250949096906e-6, 0.993365155590013193345569e-8, 0.653480510020104699270084e-10 };
double[] dc = { 1, 0.105077086072039915406159, 0.00414278428675475620830226, 0.726338754644523769144108e-4, 0.477818471047398785369849e-6 };
r = Evaluate.Polynomial(z - 38, nc)/Evaluate.Polynomial(z - 38, dc);
// Worst case absolute error found: 5.896543869e-24
r = Evaluate.Polynomial(z - 38, ErfImpLn) /Evaluate.Polynomial(z - 38, ErfImpLd);
b = 0.5640528202F;
}
else if (z < 85)
{
// Worst case absolute error found: 3.080612264e-21
double[] nc = { -0.196457797609229579459841e-4, 0.157243887666800692441195e-5, 0.543902511192700878690335e-7, 0.317472492369117710852685e-9 };
double[] dc = { 1, 0.052803989240957632204885, 0.000926876069151753290378112, 0.541011723226630257077328e-5, 0.535093845803642394908747e-15 };
r = Evaluate.Polynomial(z - 60, nc)/Evaluate.Polynomial(z - 60, dc);
r = Evaluate.Polynomial(z - 60, ErfImpMn) /Evaluate.Polynomial(z - 60, ErfImpMd);
b = 0.5641309023F;
}
else
{
// Worst case absolute error found: 8.094633491e-22
double[] nc = { -0.789224703978722689089794e-5, 0.622088451660986955124162e-6, 0.145728445676882396797184e-7, 0.603715505542715364529243e-10 };
double[] dc = { 1, 0.0375328846356293715248719, 0.000467919535974625308126054, 0.193847039275845656900547e-5 };
r = Evaluate.Polynomial(z - 85, nc)/Evaluate.Polynomial(z - 85, dc);
r = Evaluate.Polynomial(z - 85, ErfImpNn) /Evaluate.Polynomial(z - 85, ErfImpNd);
b = 0.5641584396F;
}
@ -398,10 +589,8 @@ namespace MathNet.Numerics
// long double: Max error found: 1.017064e-20
// Maximum Deviation Found (actual error term at infinite precision) 8.030e-21
const float y = 0.0891314744949340820313f;
double[] nc = { -0.000508781949658280665617, -0.00836874819741736770379, 0.0334806625409744615033, -0.0126926147662974029034, -0.0365637971411762664006, 0.0219878681111168899165, 0.00822687874676915743155, -0.00538772965071242932965 };
double[] dc = { 1, -0.970005043303290640362, -1.56574558234175846809, 1.56221558398423026363, 0.662328840472002992063, -0.71228902341542847553, -0.0527396382340099713954, 0.0795283687341571680018, -0.00233393759374190016776, 0.000886216390456424707504 };
double g = p*(p + 10);
double r = Evaluate.Polynomial(p, nc)/Evaluate.Polynomial(p, dc);
double r = Evaluate.Polynomial(p, ErvInvImpAn) /Evaluate.Polynomial(p, ErvInvImpAd);
result = (g*y) + (g*r);
}
else if (q >= 0.25)
@ -417,11 +606,9 @@ namespace MathNet.Numerics
// long double : Max error found: 6.084616e-20
// Maximum Deviation Found (error term) 4.811e-20
const float y = 2.249481201171875f;
double[] nc = { -0.202433508355938759655, 0.105264680699391713268, 8.37050328343119927838, 17.6447298408374015486, -18.8510648058714251895, -44.6382324441786960818, 17.445385985570866523, 21.1294655448340526258, -3.67192254707729348546 };
double[] dc = { 1, 6.24264124854247537712, 3.9713437953343869095, -28.6608180499800029974, -20.1432634680485188801, 48.5609213108739935468, 10.8268667355460159008, -22.6436933413139721736, 1.72114765761200282724 };
double g = Math.Sqrt(-2*Math.Log(q));
double xs = q - 0.25;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpBn) /Evaluate.Polynomial(xs, ErvInvImpBd);
result = g/(y + r);
}
else
@ -448,50 +635,40 @@ namespace MathNet.Numerics
{
// Max error found: 1.089051e-20
const float y = 0.807220458984375f;
double[] nc = { -0.131102781679951906451, -0.163794047193317060787, 0.117030156341995252019, 0.387079738972604337464, 0.337785538912035898924, 0.142869534408157156766, 0.0290157910005329060432, 0.00214558995388805277169, -0.679465575181126350155e-6, 0.285225331782217055858e-7, -0.681149956853776992068e-9 };
double[] dc = { 1, 3.46625407242567245975, 5.38168345707006855425, 4.77846592945843778382, 2.59301921623620271374, 0.848854343457902036425, 0.152264338295331783612, 0.01105924229346489121 };
double xs = x - 1.125;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpCn) /Evaluate.Polynomial(xs, ErvInvImpCd);
result = (y*x) + (r*x);
}
else if (x < 6)
{
// Max error found: 8.389174e-21
const float y = 0.93995571136474609375f;
double[] nc = { -0.0350353787183177984712, -0.00222426529213447927281, 0.0185573306514231072324, 0.00950804701325919603619, 0.00187123492819559223345, 0.000157544617424960554631, 0.460469890584317994083e-5, -0.230404776911882601748e-9, 0.266339227425782031962e-11 };
double[] dc = { 1, 1.3653349817554063097, 0.762059164553623404043, 0.220091105764131249824, 0.0341589143670947727934, 0.00263861676657015992959, 0.764675292302794483503e-4 };
double xs = x - 3;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpDn) /Evaluate.Polynomial(xs, ErvInvImpDd);
result = (y*x) + (r*x);
}
else if (x < 18)
{
// Max error found: 1.481312e-19
const float y = 0.98362827301025390625f;
double[] nc = { -0.0167431005076633737133, -0.00112951438745580278863, 0.00105628862152492910091, 0.000209386317487588078668, 0.149624783758342370182e-4, 0.449696789927706453732e-6, 0.462596163522878599135e-8, -0.281128735628831791805e-13, 0.99055709973310326855e-16 };
double[] dc = { 1, 0.591429344886417493481, 0.138151865749083321638, 0.0160746087093676504695, 0.000964011807005165528527, 0.275335474764726041141e-4, 0.282243172016108031869e-6 };
double xs = x - 6;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpEn) /Evaluate.Polynomial(xs, ErvInvImpEd);
result = (y*x) + (r*x);
}
else if (x < 44)
{
// Max error found: 5.697761e-20
const float y = 0.99714565277099609375f;
double[] nc = { -0.0024978212791898131227, -0.779190719229053954292e-5, 0.254723037413027451751e-4, 0.162397777342510920873e-5, 0.396341011304801168516e-7, 0.411632831190944208473e-9, 0.145596286718675035587e-11, -0.116765012397184275695e-17 };
double[] dc = { 1, 0.207123112214422517181, 0.0169410838120975906478, 0.000690538265622684595676, 0.145007359818232637924e-4, 0.144437756628144157666e-6, 0.509761276599778486139e-9 };
double xs = x - 18;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpFn) /Evaluate.Polynomial(xs, ErvInvImpFd);
result = (y*x) + (r*x);
}
else
{
// Max error found: 1.279746e-20
const float y = 0.99941349029541015625f;
double[] nc = { -0.000539042911019078575891, -0.28398759004727721098e-6, 0.899465114892291446442e-6, 0.229345859265920864296e-7, 0.225561444863500149219e-9, 0.947846627503022684216e-12, 0.135880130108924861008e-14, -0.348890393399948882918e-21 };
double[] dc = { 1, 0.0845746234001899436914, 0.00282092984726264681981, 0.468292921940894236786e-4, 0.399968812193862100054e-6, 0.161809290887904476097e-8, 0.231558608310259605225e-11 };
double xs = x - 44;
double r = Evaluate.Polynomial(xs, nc)/Evaluate.Polynomial(xs, dc);
double r = Evaluate.Polynomial(xs, ErvInvImpGn) /Evaluate.Polynomial(xs, ErvInvImpGd);
result = (y*x) + (r*x);
}
}

Write Preview
Loading…
Cancel
Save