View Full Version : Gamma function
danbaron
31-07-2011, 09:58
The factorial function (!) only works for integers, it is discrete.
The Gamma function works for integers (n),
gamma(n) = (n-1)!
but also for values between the integers, it is continuous,
gamma(n + 0.32457), etc.
Here is a C program which has an approximation of the Gamma function, from,
http://en.wikipedia.org/wiki/Gamma_function (the equation just above the heading, "History")
If you run the program, you see that,
gamma(2) = 1!
gamma(2) < gamma(2.5) < gamma(3)
gamma(3) = 2!
gamma(3) < gamma(3.5) < gamma(4)
gamma(4) = 3!
gamma(4) < gamma(4.5) < gamma(5)
gamma(5) = 4!
etc.
#include <stdio.h>
#include <math.h>
double gamma(double x)
// rough approximation of Gamma function
{
static double pi;
pi = 4 * atan(1);
static double e;
e = exp(1);
int i;
double sum = 0;
static double c[5];
c[0] = 1;
c[1] = 1.0/12;
c[2] = 1.0/288;
c[3] = -139.0/51840;
c[4] = -571.0/2488320;
for(i=0;i<5;i++) sum += c[i]/pow(x, i);
return sum *= pow(x,x-0.5) * pow(e,-x) * sqrt(2*pi);
}
int main()
{
char c;
double i;
for(i = 0.5; i<=40; i += 0.5) printf("%04.1f %060.10f\n", i, gamma(i));
c = getchar();
return 0;
}
danbaron
01-08-2011, 06:26
I guess it's pretty stupid to write the gamma function, when C already includes it (tgamma).
Remember that, for integers (n),
gamma(n) = (n-1)!.
(Below, you can see that as n grows, so does the error.)
(Also notice that, the output simulates the appearance of rain falling from thunderstorm clouds, when viewed from a distance. (I don't know if that is a mathematical property.))
' code ------------------------------------------------------------------------------------------------
#include <stdio.h>
#include <math.h>
int main()
{
char c;
double i;
for(i = 0.5; i<=40; i += 0.5) printf(".1f 0.10f\n", i, tgamma(i));
c = getchar();
return 0;
}
' output ----------------------------------------------------------------------------------------------
00.5 0000000000000000000000000000000000000000000000001.7724538509
01.0 0000000000000000000000000000000000000000000000001.0000000000
01.5 0000000000000000000000000000000000000000000000000.8862269255
02.0 0000000000000000000000000000000000000000000000001.0000000000
02.5 0000000000000000000000000000000000000000000000001.3293403882
03.0 0000000000000000000000000000000000000000000000002.0000000000
03.5 0000000000000000000000000000000000000000000000003.3233509704
04.0 0000000000000000000000000000000000000000000000006.0000000000
04.5 0000000000000000000000000000000000000000000000011.6317283966
05.0 0000000000000000000000000000000000000000000000024.0000000000
05.5 0000000000000000000000000000000000000000000000052.3427777846
06.0 0000000000000000000000000000000000000000000000120.0000000000
06.5 0000000000000000000000000000000000000000000000287.8852778150
07.0 0000000000000000000000000000000000000000000000720.0000000000
07.5 0000000000000000000000000000000000000000000001871.2543057978
08.0 0000000000000000000000000000000000000000000005040.0000000000
08.5 0000000000000000000000000000000000000000000014034.4072934834
09.0 0000000000000000000000000000000000000000000040320.0000000000
09.5 0000000000000000000000000000000000000000000119292.4619946090
10.0 0000000000000000000000000000000000000000000362880.0000000001
10.5 0000000000000000000000000000000000000000001133278.3889487856
11.0 0000000000000000000000000000000000000000003628800.0000000009
11.5 0000000000000000000000000000000000000000011899423.0839622486
12.0 0000000000000000000000000000000000000000039916800.0000000075
12.5 0000000000000000000000000000000000000000136843365.4655658574
13.0 0000000000000000000000000000000000000000479001600.0000001077
13.5 0000000000000000000000000000000000000001710542068.3195733000
14.0 0000000000000000000000000000000000000006227020800.0000007451
14.5 0000000000000000000000000000000000000023092317922.3142378032
15.0 0000000000000000000000000000000000000087178291200.0000104308
15.5 0000000000000000000000000000000000000334838609873.5564574599
16.0 0000000000000000000000000000000000001307674368000.0003223500
16.5 0000000000000000000000000000000000005189998453040.1263327700
17.0 0000000000000000000000000000000000020922789888000.0051576600
17.5 0000000000000000000000000000000000085634974475162.0572060300
18.0 0000000000000000000000000000000000355687428096000.0585764600
18.5 0000000000000000000000000000000001498612053315336.0709547900
19.0 0000000000000000000000000000000006402373705728001.1475086200
19.5 0000000000000000000000000000000027724322986333718.2905000000
20.0 0000000000000000000000000000000121645100408832033.2812000000
20.5 0000000000000000000000000000000540624298233507433.9061000000
21.0 0000000000000000000000000000002432902008176640607.4166000000
21.5 0000000000000000000000000000011082798113786906003.9520000000
22.0 0000000000000000000000000000051090942171709448099.1363000000
22.5 0000000000000000000000000000238280159446418474544.0000000000
23.0 0000000000000000000000000001124000727777607971802.0000000000
23.5 0000000000000000000000000005361303587544413749128.0000000000
24.0 0000000000000000000000000025852016738884984515607.0000000000
24.5 0000000000000000000000000125990634307293761521577.0000000000
25.0 0000000000000000000000000620448401733239552410000.0000000000
25.5 0000000000000000000000003086770540528697529220000.0000000000
26.0 0000000000000000000000015511210043330988264640000.0000000000
26.5 0000000000000000000000078712648783481796272090000.0000000000
27.0 0000000000000000000000403291461126605793833730000.0000000000
27.5 0000000000000000000002085885192762267589569090000.0000000000
28.0 0000000000000000000010888869450418354972400000000.0000000000
28.5 0000000000000000000057361842800962331239100000000.0000000000
29.0 0000000000000000000304888344611713895574200000000.0000000000
29.5 0000000000000000001634812519827426644042100000000.0000000000
30.0 0000000000000000008841761993739703670144000000000.0000000000
30.5 0000000000000000048226969334909066557884200000000.0000000000
31.0 0000000000000000265252859812191200035000000000000.0000000000
31.5 0000000000000001470922564714727341197000000000000.0000000000
32.0 0000000000000008222838654177925782278000000000000.0000000000
32.5 0000000000000046334060788513915613293000000000000.0000000000
33.0 0000000000000263130836933693625032901000000000000.0000000000
33.5 0000000000001505856975626701932920000000000000000.0000000000
34.0 0000000000008683317618811891588840000000000000000.0000000000
34.5 0000000000050446208683494507567950000000000000000.0000000000
35.0 0000000000295232799039604142308230000000000000000.0000000000
35.5 0000000001740394199580549821257590000000000000000.0000000000
36.0 0000000010333147966386148254900000000000000000000.0000000000
36.5 0000000061783994085110651212700000000000000000000.0000000000
37.0 0000000371993326789906364865600000000000000000000.0000000000
37.5 0000002255115784106522798538200000000000000000000.0000000000
38.0 0000013763753091226549819111800000000000000000000.0000000000
38.5 0000084566841903994709253311100000000000000000000.0000000000
39.0 0000523022617466599549516000000000000000000000000.0000000000
39.5 0003255823413303802954033000000000000000000000000.0000000000
40.0 0020397882081197472289204000000000000000000000000.0000000000
actually it's interesting to find new ways to approximate the gamma function, when I get around to it I will post some implementations. :)
danbaron
04-08-2011, 06:03
http://en.wikipedia.org/wiki/Gamma_function
http://en.wikipedia.org/wiki/Euler%E2%80%93Mascheroni_constant
(http://en.wikipedia.org/wiki/Euler%E2%80%93Mascheroni_constant)
I can say this.
Before stumbling onto the realization that C now includes the Gamma function, I tried to implement it using the infinite product definition of Euler and Weierstrass, using the Euler-Mascheroni constant.
In C, I had real problems getting the constant to converge, and with the accuracy of the partial product.
My experience really made me wonder how the intrinsic C Gamma function is so accurate.
the Lanczos approximation is a popular method but there are other interesting formulas, my idea is to round up the newer methods and see how they perform.
here'a Lanczos implementation
Uses "console"
Function gamma(ByVal y As Ext) As Ext
Dim As Ext Pi = 3.1415926535897932385
Dim As Ext sq2pi = 2.50662827463100050241577 'sqrt(2Pi)
Dim As Ext g = 607/128 ' best resu'ts when 4<=g<=5
Dim As Ext t, w, gam, x = y-1
Dim As Integer i, cg = 14
Dim c(15) As Ext
'Lanczos approximation for the complex plane
'calculated using vpa digits(256)
'the best set of coeffs was selected from
'a solution space of g=0 to 32 with 1 to 32 terms
'these coeffs really give superb performance
'of 15 sig. digits for 0<=real(z)<=171
'coeffs should sum to about g*g/2+23/24
'http://www.numericana.com/answer/info/godfrey.htm
c( 1) = 0.99999999999999709182 'thiBasic arrays start at 1 ?
c( 2) = 57.156235665862923517
c( 3) = -59.597960355475491248
c( 4) = 14.136097974741747174
c( 5) = -0.49191381609762019978
c( 6) = 0.33994649984811888699/10000
c( 7) = 0.46523628927048575665/10000
c( 8) = -0.98374475304879564677/10000
c( 9) = 0.15808870322491248884/1000
c(10) = -0.21026444172410488319/1000
c(11) = 0.21743961811521264320/1000
c(12) = -0.16431810653676389022/1000
c(13) = 0.84418223983852743293/10000
c(14) = -0.26190838401581408670/10000
c(15) = 0.36899182659531622704/100000
If ( x < 0 ) Then
x=-x
If Frac(x)=0 Then
gam=10^4932
Else
t = c(1)
For i=1 To cg
t = t + c(i+1)/(x+i)
Next
w = x + g + 0.5
gam=w^(x+0.5)*Exp(-w)*sq2pi*t
gam = Pi*x/(gam*Sin(Pi*x))
End If
Else
t = c(1)
For i=1 To cg
t = t + c(i+1)/(x+i)
Next
w = x + g + 0.5
gam=w^(x+0.5)*Exp(-w)*sq2pi*t
End If
Function = gam
End Function
Dim As Ext y, x = 1
While x>0
Console_Write "x "
x=Console_ReadLine()
If x=0 Then Exit While
y = gamma (x)
Console_WriteLine "Gamma "<rim$(Str$(x))&" = "&Format$(y,17)
Wend
Console_WriteLine "All done. Press any key to finish"
Console_WaitKey
danbaron
06-08-2011, 06:08
Wow, it works good, Johan.
I didn't know anything about the Lanczos implementation.
I had heard of him before.
http://www.gap-system.org/~history/Biographies/Lanczos.html (http://www.gap-system.org/%7Ehistory/Biographies/Lanczos.html)
100! = 93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
From your implementation,
Gamma 101 = 9.332621544394413E+157.
So, the first 15 digits are correct.
(I'll look at the implementation of the Gamma function, in "Numerical Recipes in C", now.)
Dan
danbaron
06-08-2011, 06:58
For 100!, the function from, "Numerical Recipes in C", is only correct for 9 digits, so the Lanczos implementation beats it.
' code --------------------------------------------------------------------------------------
#include <stdio.h>
#include <math.h>
double gammaln(double xx)
// From "Numerical Recipes in C", 2nd ed., p.214.
// Returns ln(Gamma(xx)), for xx > 0.
{
int j;
double x,y,tmp,ser;
static double cof[6];
cof[0] = 76.18009172947146;
cof[1] = -86.50532032941677;
cof[2] = 24.01409824083091;
cof[3] = -1.231739572450155;
cof[4] = 0.1208650973866179e-2;
cof[5] = -0.5395239384953e-5;
y=x=xx;
tmp=x+5.5;
tmp -= (x+0.5)*log(tmp);
ser=1.000000000190015;
for(j=0;j<=5;j++) ser += cof[j]/++y;
return -tmp+log(2.5066282746310005*ser/x);
}
int main()
{
char c;
double i;
i = 101;
printf("%05.1f %25.20e\n", i, exp(gammaln(i)));
c = getchar();
return 0;
}
' output ------------------------------------------------------------------------------------
101.0 9.33262154537299415097e+157
danbaron
06-08-2011, 07:24
On the other hand, "Super-Mathematica" will give the exact result for Gamma(1001) (maybe it calculates 1000!).
http://www.wolframalpha.com/input/?i=Gamma[1001]
Result:
402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
Check result, from Racket.
Welcome to DrRacket, version 5.1 [3m].
Language: racket.
> (fact 1000)
402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
>
danbaron
06-08-2011, 07:39
For Gamma(1001 + 1/1000000),
http://www.wolframalpha.com/input/?i=Gamma[1001.000001]
this is all I could get from "Super-Mathematica".
4.0239003988057584997902173976163592011969374179001511239151002642392471030979284415068711496796131858912114129894951019307516229727006305338006965443711395900922510707344833143982161291576196145549524252116586498402388174031102033111615654658920197684590636911143385941614831004326546290077112736324575126945175550119098139642713208122214518177314109490782649292284077778713228049034887764970619206691682137254901666889850542038221053151508464693805173317329282324843144943943599436514279514538607552755619053756475634096657267703608254014755211290420692809243460322686928748132961796450943572995373926927227559737480726681980983218643224702815506486387652208117712019519946798982740183845603220534925666831744612502114908382710346769016257450862275257405860597569351553374047531992774021108553153132441252879047714632798432033625895109963398433028869133545373648578084478742621362263016745157981846648908850249353140681346060239970206896405377407332686209471584545211937434617377658459807959756083415727145958537520435373456645804904729989329964598572790374526220294313751473179336564051622005259346833431900944629570839237229481513254432209634074253387839700361771632471703417078441912584157269029731250759827936047129505797040335623306987799596245979371220717534340058326409719830787501939605783323783893179283319665132503930707413454484359916812366340570334427980972071481375201194203796140336701067494124741436407947480809292277983651740542286852851769396310225165686126061205580214600926014228591882855121520450907878434749341566595088164908352656635137221802845387984029365296643760097823289428298591031940194745031433675760122709191816281215841944892670582736814231986110261890538652395205193148570945253573720796741616911796790452097412733314334762831996654388583958842108704860923901748565161209141596165235353480433962920389547091119284361184464231936462484553514420363877766052084394699696456390950010286890909749030319272522477519054298031886134749961363250122571370745518417715230196582522727560938697643549511273538629747079520406413534104121981790847012314953152367591394890854899650620559354790016382878566029849228881525569110021829561650569024779437891761222842924279083000460739467674535783721585102675448762616885235579170233920356944428551215896798274345929092227557674918376075932337861270314616782583093758957562631739284704672146519207109744977914313663591015850751595443459206763794303939770317473799857171590222662301548722568155715675746100133944439923277239113204532536037954523472602398569689151629072755297077462874071783240172319640630147138674416129455126541858280762275787308341634326062407520397616731640922794849171635359556882358169551264574095302587695705568308648876430372219809533500063843080428965101561654802023673373158391310104669887045308869837699763815934065102121965101508796092056855661749281613039002338281146882969916854109020038123529481436584940630288994678359267874045978361539072677171117303849094719953611327095943318932698734386917593119024538031854596120942550157270174457336010032617790660722732015345076544097871929876417549941543216004206034777057170912796699287536720531688865073593195894207758319489176583240052580932558196302682862695669489650826090073859613164147950379619784492273985221444482402596660180723828594850186630616315085145571389650601927654671155352412964610239592652562298130511431240532574515358034040225214190832152608962048874360481181163117706057960456785858703174930508635662497576032855328818899495348962384480754728207869918433136816825060953289916424884475186026464979296732458581572805992968768306412122770943206407590057950197962062011156536819841696174893179842028679782235569275374134554732020123806556821231345176... × 10^2567
(And, my experience with "Super-Mathematica" is that it usually gives the answer for almost any input, NOW.)
(The answers for Gamma(1001) and Gamma(1001.000001) diverge after only 4 digits. It seems strange. I guess it must be correct.)
For 100!, the function from, "Numerical Recipes in C", is only correct for 9 digits, so the Lanczos implementation beats it.
Dan, both are Lanczos approximations, the one I posted uses more terms, but there's something wrong with thinBasic,
it's supposed to use extended precision in it's calculations but
Val("1e+4932")
gives
1.7014118346046923E+38
here's the Stieltjes' approximation.
uses "console"
Function gamma ( x As Ext ) As Ext
'Stieltjes' approximation to the Gamma function
'x!=exp(log(2*Pi)/2-x+(x+1/2)*log(x)+cf(x))
'cf(x) = c1
' ---------------------
' c2
' x + ----------------
' c3
' x + ------------
' c4
' x + --------
' x + '
' '
' '
Dim As Ext cf, logpi2half, xt, gam, pp
Dim As Integer i
logpi2half = 0.91893853320467274178032973640562
' COEFFICIENTS FOR CONTINUED FRACTION
Dim As Ext C(14)
c( 1) = 0.0833333333333333333333333333333333
c( 2) = 0.0333333333333333333333333333333333
c( 3) = 0.252380952380952380952380952380952
c( 4) = 0.525606469002695417789757412398922
c( 5) = 1.01152306812684171174737212473062
c( 6) = 1.51747364915328739842849151949548
c( 7) = 2.26948897420495996090915067220988
c( 8) = 3.00991738325939817007314073420772
c( 9) = 4.02688719234390122616887595318144
c(10) = 5.00276808075403005168850241227619
c(11) = 6.28391137081578218007266315495165
c(12) = 7.49591912238403392975235470826751
c(13) = 9.04066023436772669953113936026048
c(14) = 10.4893036545094822771883713045926
' CONTINUED FRACTION
xt=x+8
CF = C(9)+Xt
For i = 8 To 2 Step -1
cf=xt+c(i)/cf
Next
cf=c(1)/cf
' POCHHAMMER POLYNOMIAL
PP = X * (X+1) * (X+2) * (X+3) * (X+4) * (X+5) * (X+6) * (X+7) * (X+8)
gam = Exp(logpi2half-xt+(xt+.5)*Log(xt)+cf)
Return gam/pp
End Function
Dim As Ext y, x = 1
While x>0
Console_Write "x "
x=Console_ReadLine()
If x=0 Then Exit While
y = gamma (x)
Console_WriteLine "Gamma "<rim$(Str$(x))&" = "&Format$(y,17)
Wend
Console_WriteLine "All done. Press any key to finish"
Console_WaitKey
danbaron
07-08-2011, 06:18
x 101
Gamma 101 = 9.3326215443944156E+157
x 1001
Gamma 1001 = 4.023872600770938E+2567
x
Welcome to DrRacket, version 5.1 [3m].
Language: racket.
> (fact 100)
93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
> (fact 1000)
402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
>
(Your script needs a, 'uses "console"'.)
Wow, that is really good, Johan.
You know a lot of math.
I just was looking at continued fractions, for the first time the other night.
http://en.wikipedia.org/wiki/Continued_fraction
Considering Stieltjes, all I know is that integrals are named after him.
http://en.wikipedia.org/wiki/Riemann–Stieltjes_integral
(http://en.wikipedia.org/wiki/Riemann%E2%80%93Stieltjes_integral)http://en.wikipedia.org/wiki/Lebesgue–Stieltjes_integration (http://en.wikipedia.org/wiki/Lebesgue%E2%80%93Stieltjes_integration)
I have a book, by Daniel Stroock, "A Concise Introduction to the Theory of Integration", 2nd ed., Birkhauser. But, it seems to be beyond my current capacity/motivation.
Dan
-----------------------------------------------------------------
I agree with you about Val("1e+4932").
But, try this.
Uses "console"
Dim e As Ext
e = Val("9.9999999999999999999999999999999999e+4931")
PrintL e
Console_ReadLine
Petr Schreiber
07-08-2011, 09:00
Hi Jack,
I can confirm that, interesting. It seems there is some kind of overflow(?) taking place.
VAL("1e4931") is still okay.
Anyway - if you need to specify 1e4932, you can do it like:
Dim myNumber As Ext
myNumber = 10 ^ 4932
... and the assigned value will be correct.
Petr
P.S. One wish - if you find some odd behavior, please do not hesitate to report the problem via Support in the top forum bar. It is easier to track and fix problems then. Thanks!
Dan, both are Lanczos approximations, the one I posted uses more terms, but there's something wrong with thinBasic,
it's supposed to use extended precision in it's calculations but
Val("1e+4932")
gives
Dan, added the missing uses "console" guess I was in a hurry.
P.S. One wish - if you find some odd behavior, please do not hesitate to report the problem via Support in the top forum bar. It is easier to track and fix problems then. Thanks!
agree Petr. :)