1/* Filename mathieu0.mac 2 3 *************************************************************** 4 * * 5 * <package name> * 6 * <functionality description> * 7 * * 8 * from: Computer Algebra in Applied Math. * 9 * by Rand (Pitman,1984) * 10 * Programmed by Richard Rand * 11 * These files are released to the public domain * 12 * * 13 *************************************************************** 14*/ /* 15(d4) This program computes the transition curve through 16 17 18the origin (N = 0) in Mathieu's equation using a perturbation 19 20 21method. To call it, type 22 23 24 MATHIEU0() 25 26*/ 27 28mathieu0():=(input(),setup1(),setup2(), 29 for i thru m do (step1(),step2a(),step3a()),output())$ 30input():=m:read("ENTER DEGREE OF TRUNCATION")$ 31setup1():=(w:0,for i thru m do w:w+k[i]*e^i,x:1, 32 for i thru m do x:x+y[i](t)*e^i)$ 33setup2():=(temp1:diff(x,t,2)+x*(w+e*cos(t)),temp1:taylor(temp1,e,0,m), 34 for i thru m do eq[i]:coeff(temp1,e,i))$ 35step1():=temp1 36 :expand(trigreduce(expand(ev(eq[i],makelist([e[j],f[j]],j,1,i-1), 37 diff))))$ 38step2a():=(temp1:ode2(temp1,y[i](t),t),temp1:ev(temp1,%k1:0,%k2:0))$ 39step3a():=(f[i]:solve(coeff(expand(rhs(temp1)),t,2),k[i]),e[i]:ev(temp1,f[i]))$ 40output():=(print("delta=",ev(w,makelist([f[j]],j,1,m))),print(" "))$ 41