LambertW.cin

// LambertW.cin is numerical implementation of function LambertW in C++.

//Copyleft 2012 by Dmitrii Kouznetsov.

z_type LambertWo(z_type z){ int n,m=48; z_type d=-z;
DB LambertWoCoe[66]={1.
,            1.00000000000000
,            1.50000000000000
,            2.66666666666667
,            5.20833333333333
,           10.80000000000000
,           23.34305555555556
,           52.01269841269841
,          118.62522321428571
,          275.57319223985894
,          649.78717234347448
,         1551.16051948051950
,         3741.44970295923895
,         9104.50024115801898
,        22324.30851270660423
,        55103.62197290384211
,       136808.86090394295752
,       341422.05066583841108
,       855992.96599660767242
,      2154990.20609108870849
,      5445552.92231446318328
,     13807330.00216663256288
,     35117044.98513923585415
,     89568002.56102797389030
,    229041684.61879497766495
,    587103504.11717987060547
,   1508256053.85779309272766
,   3882630161.29318952560425
,  10013943136.65483093261719
,  25873567362.65761184692383
,  66962097093.58074951171875
, 173571165959.91983032226562
, 450568046564.23547363281250
,1171223178256.48754882812500
,3048462517882.44677734375000
,7944240398206.81152343750000
,20726462893727.0546875000000
,54134309927745.6406250000000
,141537001323842.750000000000
,370420141783066.687500000000
,970343148988291.750000000000
,2544149870142830.00000000000
,6676186549945614.00000000000
,17533383884497606.0000000000
,46082787899842400.0000000000
,121208404038419968.000000000
,319031756001941888.000000000
,840290335422231168.000000000
,2214659548716029696.00000000
,5840571433838373888.00000000
,15412168778694733824.0000000
};
// #include "LambertWoCoe.inc"
z_type s=LambertWoCoe[m]*d; for(n=m-1;n>0;n--){ s+=LambertWoCoe[n]; s*=d;}
return z*(1.+s); }

z_type LambertWe(z_type z){ int n,m=100; z_type t=1./M_E+z; t*=2*M_E; t=sqrt(t); 
double LambertWeCoe[102]={ -1., 1.,
-0.3333333333333333333,    0.1527777777777777778,   -0.0796296296296296296, 0.0445023148148148148,
-0.0259847148736037625,    0.0156356325323339212,   -0.00961689202429943171, 0.00601454325295611786,
-0.00381129803489199923,   0.00244087799114398267,  -0.00157693034468678425,  0.00102626332050760715,
-0.000672061631156136204,  0.00044247306181462091,  -0.000292677224729627445, 0.000194387276054539318,
-0.000129574266852748819,  0.0000866503580520812717,-0.0000581136075044138168, 0.0000390766848674390516,
-0.0000263380647472310987, 0.0000177903458050795854,-0.0000120403527395599769, 8.16353198249661217e-6,
-5.54420320856735914e-6, 3.77109496110725341e-6,    -2.56870503905509544e-6, 1.7520067268263412e-6,
-1.19644530891572567e-6, 8.17994056528003472e-7,    -5.59855188137879574e-7, 3.83566385149181379e-7,
-2.63037861927186308e-7, 1.80544727751016444e-7,    -1.24027544004224703e-7, 8.52703516168582825e-8,
-5.86686308977225918e-8, 4.03947301280155627e-8,    -2.78316209626026845e-8, 1.91881310685738554e-8,
-1.32371244669424267e-8, 9.13710905680806547e-9,    -6.3105381623506781e-9, 4.3606925472094428e-9,
-3.01485057758130614e-9, 2.08539250292067727e-9,    -1.44315359778020971e-9, 9.9915217179721807e-10,
-6.92049499350375257e-10, 4.79536537580206161e-10,  -3.32413217764048006e-10, 2.30515597891712908e-10,
-1.59912196300685515e-10, 1.10972718669265744e-10,  -7.70368786717788538e-11, 5.34962929680470565e-11,
-3.71609045080338343e-11, 2.58215148167552937e-11,  -1.79475645333454047e-11, 1.24782428944225017e-11,
-8.67803499992131116e-12, 6.03678282054492658e-12,  -4.20051193774766886e-12, 2.92353237458081712e-12,
-2.03525702096338716e-12, 1.41720636451186404e-12,  -9.87066370244645293e-13, 6.87632082200663056e-13,
-4.79136913828980365e-13, 3.33929038076548621e-13,  -2.327754844257851e-13, 1.62295447584966799e-13,
-1.13177248655861478e-13, 7.89393171526768803e-14,  -5.50689596325770156e-14, 3.84235573773504502e-14,
-2.68141146603903671e-14, 1.87155512160818899e-14,  -1.30651139440572381e-14, 9.12207036595587938e-15,
-6.37003114897626912e-15, 4.44893399727651941e-15,  -3.10767106994306571e-15, 2.17108719363861928e-15,
-1.51698449707275936e-15, 1.06009611885002865e-15,  -7.40914663180577404e-16, 5.17903258763594061e-16,
-3.62063975903811412e-16, 2.53149433731430854e-16,  -1.77020013639351241e-16, 1.23799924052967824e-16,
-8.65904219941911638e-17, 6.05716926896423085e-17,  -4.23758948025618507e-17, 2.96494297432517426e-17,
-2.07472769988816835e-17, 1.45195179941344395e-17,  -1.01622223292926686e-17, 7.1132762840052161e-18};
// #include "LambertWeCoe.cin"
z_type s=LambertWeCoe[m]*t; for(n=m-1; n>0; n--) { s+=LambertWeCoe[n]; s*=t;} return -1.+s;}

z_type LambertW(z_type z){ 
if( abs(z)<.2 ) return LambertWo(z);
if( abs(1./M_E+z)<.2 ) return LambertWe(z);
return Tania(log(z)-1.);} // Except the negative part of the real axis, Tania does the LambertW well.

//End of the numerical implementation of LambertW. // // //