[q-lang-users] taylor series stuff
Brought to you by:
agraef
Menu
▾
▴
[q-lang-users] taylor series stuff
|
From: <ed...@ri...> - 2007-06-29 19:44:53
|
<p> <br />Hello all,</p><p>This is something I've been playing with he= re and there. My intention is to find a Taylor's series for certain functio= ns (mostly statistical ones, probability density functions) integrate them = term by term to get an approximation of the cumulative density functions, a= nd then maybe, add some equations to convert the integrated polynomials to = Horner's scheme for quick approximations.</p><p>I'm after comments, critici= sm, and suggestions.</p><p>Eddie </p><p>source: </p><p>/* Derivat= ives and Taylor's series stuff by Eddie Rucker <br />Fri Jun 29 14:53:09 CD= T 2007<br /><br />Use at your own risk!<br /><br />All derivatives are with= respect to (WRT) some symbol.<br />Also, expressions must be in parenthesi= s, i.e. "dd (Expression) symbol."<br />dd_lambda will return the = derivative of a lambda.<br />dd_lambda_nth will return the nth derivative o= f a lambda.<br />taylor F A N will return a lambda that is a series about p= oint A with N terms.<br />Make sure you use analytic functions with the tay= lors series expansion or it <br />will blow up on you!<br />I'm sure there = are some warts and please overlook the hackery of the algebra<br />simplifi= cations.<br /><br /> Examples:<br /> <br /> &n= bsp; =3D=3D> dd (2*X^3+3*X*Y^2+4) X<br /> 6*X^2+3*Y^2<br />&= nbsp; <br /> =3D=3D> dd (2*X^3+3*X*Y^2+4) Y<br /> = ; 6*X*Y<br /> <br /> =3D=3D> dd (exp (2*X)= ) X<br /> 2*exp (2*X)<br /><br /> =3D=3D> dd (co= s (2*X^2+X*Y)^2) X<br /> -2*(4*X+Y)*sin (2*X^2+X*Y)*cos (2*X^2+= X*Y)<br /><br /> =3D=3D> dd (cos (2*X^2+X*Y)^2) Y<br /> = ; -2*X*sin (2*X^2+X*Y)*cos (2*X^2+X*Y)<br /><br /> =3D=3D= > dd (ln (2*X)) X<br /> 1/X<br /><br /> =3D=3D&g= t; dd_nth (X^(1%2)) X 3 // take the third derivative<br /> = ; 3%8*X^(-5%2)<br /><br /> =3D=3D> dd_nth (ln (2*X)) X= 2<br /> -X^(-2)<br /><br /> =3D=3D> dd (2^(2*X)= ) X<br /> 1.38629436111989*2^(2*X)<br /><br />  = ; =3D=3D> dd_lambda (\X . sin X) 3.141597<br /> -0.999= 999999990554<br /><br /> =3D=3D> taylor (\X . sin X) 0 7 // = taylor series of sin about the point 0, terms up to X^7<br /> &n= bsp; \X1 . X1 + -0.166666666666667*X1^3+0.00833333333333333*X1^5 + -0.00019= 8412698412698*X1^7<br /><br /> =3D=3D> var f =3D taylor (\X = . sin X) 0 7 // assign the series to a function<br /> =3D=3D>= ; f 2<br /> 0.907936507936508<br /> =3D=3D>= ; sin 2<br /> 0.909297426825682<br /> =3D=3D&= gt; var f =3D taylor (\X . sin X) 2 7 // better approximation about 2<br />= =3D=3D> f 2<br /> 0.909297426825682<br />= <br /> =3D=3D> taylor (\X . exp X) 0 7<br /> &nbs= p; \X1 . 1.0+X1+0.5*X1^2+0.166666666666667*X1^3+0.0416666666666667*X1^4+0.0= 0833333333333333*X1^5+0.00138888888888889*X1^6+0.000198412698412698*X1^7<br= /><br /> =3D=3D> var g =3D taylor (\X . ln X) 2 7 // be car= eful not to make a series about 0 but 2 is okay<br /> =3D=3D>= ; g 2 <br /> 0.693147180559945<br /><br /> = =3D=3D> ln 2<br /> 0.693147180559945<br /><br /> = =3D=3D> taylor (\X . exp (X^2)) 0 9<br /> \X1 .= 1.0+X1^2+0.5*X1^4+0.166666666666667*X1^6+0.0416666666666667*X1^8<br />*/<b= r /><br />// Some of these simplification rules are in the Q in a nutshell<= br /><br />(X^N)^M =3D X^(N*M); // added<br /= >X^N*X^M =3D X^(N+M); // added<br />Y^1 =3D Y= ; &n= bsp; // added<br />Y^1.0 =3D Y; &n= bsp; // added<br />Y^0 =3D _FAIL_= if (isnum Y) and then (Y=3D0); /= / added<br /> =3D 1 otherwi= se;<br />Y^0.0 =3D _FAIL_ if (isnum Y) and then (Y=3D0.0);  = ; // this isn't the proper way to deal with floats!<br /> = =3D 1 otherwise;<br />X*(Y= *Z) =3D (X*Y)*Z;<br />X+(Y+Z) =3D (X+Y)+Z;<br />X*Y+X*Z =3D X*(Y+Z); &= nbsp; // changed direction<br />Y*X+Z*X =3D (Y+Z)*X; = ; // changed direction<br />0*X =3D 0;<br />0.0*X = =3D 0.0; = // this isn't the proper way to deal with floats!<br />X*0 =3D 0;<br />X*0.= 0 =3D 0.0;  = ; // this isn't the proper way to deal with floats!<br />1*X =3D X;<br />1.= 0*X =3D X;<br />X*1 =3D X;<br />X*1.0 =3D X;<br />0+X =3D X;<br />0.0+X =3D= X;<br />X+0 =3D X;<br />X+0.0 =3D X;<br />0-X =3D -X;<br />0.0-X =3D X;<br= />X-0 =3D X;<br />X-0.0 =3D X;<br />(-1)*X =3D -X; = // added and maybe questionable<= br />(-1.0)*X =3D -X;<br />X*(-1) =3D -X; &nbs= p; // same here<br />X*(-1.0) =3D -X;<br />X*= A =3D A*X if isnum A; // added<br />-(-X) =3D X; &nb= sp; // added<br= />X*-Y =3D -X*Y; &nbs= p; // added<br />X/(X*Y) =3D 1/Y; = // added<br />X/(Y*X) =3D 1/Y; &n= bsp; // added<br />N*X+X =3D (N+1)*X;&n= bsp; // added<br />X+N*X =3D (N+1)*X;&n= bsp; // added<br /><br />// some derivative functio= ns -- taken from my calculus textbook.<br /><br />dd (exp X) WRT =3D (dd X = WRT) * exp X;<br /><br />dd (ln X) WRT =3D (dd X WRT) * 1/X;<br /><br />dd = (sin X) WRT =3D (dd X WRT) * cos X;<br /><br />dd (cos X) WRT =3D - (dd X W= RT) * sin X;<br /><br />dd (tan X) WRT =3D (dd X WRT) / ((cos X)^2);<br /><= br />dd (sinh X) WRT =3D (dd X WRT) / (cosh X);<br /><br />dd (cosh X) WRT = =3D (dd X WRT) / (sinh X);<br /><br />dd (tanh X) WRT =3D (dd X WRT) / ((co= sh X)^2);<br /><br />dd (asin X) WRT =3D (dd X WRT) / sqrt (1 - X^2);<br />= <br />dd (acos X) WRT =3D - (dd X WRT) / sqrt (1 - X^2);<br /><br />dd (ata= n X) WRT =3D (dd X WRT) / (1 + X^2);<br /><br />dd (asinh X) WRT =3D (dd X = WRT) / sqrt (X^2 + 1);<br /><br />dd (acosh X) WRT =3D (dd X WRT) / sqrt (X= ^2 - 1);<br /><br />dd (atanh X) WRT =3D (dd X WRT) / (1 - X^2);<br /><br /= >dd (U*V) WRT =3D (dd U WRT)*V + (dd V WRT)*U; // product rule<br /><br />d= d (U/V) WRT =3D dd (U*V^(-1)) WRT; // product rule<br /><br />dd (U+V) WRT = =3D (dd U WRT) + (dd V WRT); // addition rule<br /><br />dd (U-V) WRT =3D (= dd U WRT) - (dd V WRT); // addition rule<br /><br />dd (U^V) WRT =3D (dd V = WRT) * ln U * U^V if isnum U;<br /> &nb= sp; =3D (dd U WRT) * V * U^(V-1); // power and chain rule<br /><br />dd (-U= ) WRT =3D - (dd U WRT); // derivative function is odd<br /><br />dd U WRT = =3D 1 if U =3D=3D WRT; // constant or symbol<br /> = =3D 0 otherwise;<br /><br />// higher order derivatives<br />dd_nth U WRT N= =3D dd U WRT &n= bsp; if N =3D 1; <br />&nbs= p; =3D dd (dd_nth U WRT (N= -1)) WRT otherwise;<br /><br />dd_lambda (\X. Expr) =3D \X. `(dd (Expr) X);= <br /><br />dd_lambda_nth (\X. Expr) N =3D \X. `(dd_nth (Expr) X N);<br /><= br />taylor F A N =3D \X . `(taylor_ F A 0 N 1 0);<br />taylor_ F A I N Fac= t Sum<br /> =3D (Sum + (F A)/Fac= t * (X-A)^I) if I =3D N;<br /> = =3D taylor_ (DD) A (I+1) N (Fact*(I+1)) (Sum + (F A)/Fact * (X-A)^I)<br />&= nbsp; where DD =3D dd_lambda F;<br /><= br /> =0D </p><BR>= |
