/*

(d4) This program computes the transition curve through


the origin (N = 0) in Mathieu's equation using a perturbation


method.  To call it, type


                     MATHIEU0()

*/


mathieu0():=(input(),setup1(),setup2(),
	 for i thru m do (step1(),step2a(),step3a()),output())$
input():=m:read("ENTER DEGREE OF TRUNCATION")$
setup1():=(w:0,for i thru m do w:w+k[i]*e^i,x:1,
       for i thru m do x:x+y[i](t)*e^i)$
setup2():=(temp1:diff(x,t,2)+x*(w+e*cos(t)),temp1:taylor(temp1,e,0,m),
       for i thru m do eq[i]:coeff(temp1,e,i))$
step1():=temp1
      :expand(trigreduce(expand(ev(eq[i],makelist([e[j],f[j]],j,1,i-1),
				   diff))))$
step2a():=(temp1:ode2(temp1,y[i](t),t),temp1:ev(temp1,%k1:0,%k2:0))$
step3a():=(f[i]:solve(coeff(expand(rhs(temp1)),t,2),k[i]),e[i]:ev(temp1,f[i]))$
output():=(print("delta=",ev(w,makelist([f[j]],j,1,m))),print(" "))$