Subject: Si, Ci, Ei, also more pattern matching examples.
From: Richard Fateman
Date: Thu, 11 Dec 2008 17:10:13 -0800
Here's a set of rules for integrating some stuff with Si, Ci, and Ei.
Is there another place where they are defined etc?
/* Integration of some forms including cos(x)/x^n and sin(x)/x^n in
terms of Ci(x)=integrate(cos(x)/x,x) and Si(x)=integrate(sin(x)/x,x)*/
(matchdeclare(w,true,nn, nonnegintegerp),
tellsimp('integrate(cos(w)/w^nn,w), cosint[nn](w)),
tellsimp('integrate(sin(w)/w^nn,w), sinint[nn](w)),
/*cosint[n](w) is integral of cos(w)/w^n */
cosint[n](w):= block([k:1-n], (cos(w)*w^k +sinint[-k](w))/k),
sinint[n](w):= block([k:1-n], (sin(w)*w^k -cosint[-k](w))/k),
cosint[0](w):= sin(w),
sinint[0](w):= -cos(w),
cosint[1](w):= Ci(w),
sinint[1](w):= Si(w),
gradef(Ci(w), cos(w)/w),
gradef(Si(w), sin(w)/w))
/* you may be surprised [I was] that this will find closed forms for
integrate (a*cos(x)*x^4+b*sin(x))/x^2,x);
integrate(cos(sqrt(x))/x^2,x);
integrate(cos(x+1)/(x+1)^2,x);
integrate(cos(x^4)/x,x);
integrate(x*cos(sqrt(x^2+2))/(x^2+2), x);
integrate(cos(x)*cos(sin(x))/sin(x), x);
integrate(-sin(x)*cos(cos(x))*cos(sin(cos(x)))/sin(cos(x)),x);
integrate(Ci(x)*cos(x)/x,x);
integrate(Ci(x^2)*cos(x^2)/x,x);
*/
/* more rules, for x^n*Ci(x), x^n*Si(x), n is pos or neg but not -1. */
(matchdeclare (notm1, lambda([r],is (r # -1))),
tellsimp('integrate(Ci(w)*w^notm1,w),
block([k:1+notm1],(Ci(w)*w^k -integrate(w^(k-1)*cos(w),w))/k)),
(tellsimp('integrate(Si(w)*w^notm1,w),
block([k:1+notm1],(Si(w)*w^k -integrate(w^(k-1)*sin(w),w))/k))))
/*exponential integral, Ei(x)= integral(exp(x)/x,x)
matches don't work so easily, since the matcher has to look inside
both %e^x and x^n to see which is which, and back up. Here's a
way around. Another hazard is that 'integrate(z,w) simplifies to z*w.
Is that a bug? To get around this we introduce mintegrate(). */
(matchdeclare([z,w,v],true,nnn, negintegerp),
negintegerp(n):=is (n<0) and integerp(n),
tellsimpafter(mintegrate(z,v), intwrtf(v,z)),
defmatch(exp_mat_,exp(w)*w^nnn,w), /* if we know w, we can find exp(w),
w^n */
intwrtf(v,z):=block([ans:exp_mat_(z,v)],print(ans),
if (ans= false) then integrate(z,v) else
ev(expint[-nnn](v),ans) ),
/*expint[n](w) is integral of exp(w)*w^n */
expint[n](w):= block([k:1-n], (exp(w)*w^k -expint[-k](w))/k),
expint[1](w):=Ei(w),
expint[0](w):=exp(w),
gradef(Ei(w), exp(w)/w))