# Solution of equation with power series (perturbation)

So I need to use Mathematica to find the solution of $y=x- \epsilon \sin(2y)$ as a power series in terms of $\epsilon$. I'd assume I'd need to create an equation $f=x-y- \epsilon \sin(2y)$, then express $y(\epsilon)=\sum_n a_n \epsilon^n$, then input into series, but I can't seem to get it to work. Some help would be appreciated.

-
Hello Lenny! Welcome to the Mathematica StackExchange community. Have you tried to do anything by yourself? I mean, is there any piece of code which you've tried to use? –  Rod Jun 13 '13 at 10:13
I have. What I have tried so far: eq[eps_,y_]:=x-eps Sin[2y]; y[eps_]:=Sum[Subscript[a,n] eps^n,{n,0,20}]; Series[eq[eps,y[eps]],{eps,0,2}]; –  Lenny Jun 13 '13 at 10:16
Is eps any Mathematica function? Or is it just an argument of a function? –  Rod Jun 13 '13 at 10:19
eps is just epsilon –  Lenny Jun 13 '13 at 10:20
Hint: let $x=3\eta/2-\xi, y=\eta/2$, then your equation will reduce to Kepler's Equation. –  Silvia Jun 13 '13 at 11:03

Let define the equation to solve $f=x-y\epsilon\sin(2 x)\equiv 0$ and series expansion of $y$ in powers of $\varepsilon$.

f = x - y - \[Epsilon] Sin[2 y];
ord = 3;
y = x + Sum[a[n] \[Epsilon]^n, {n, ord}] + O[\[Epsilon]]^(ord + 1);


Then expand the equation and solve for any value of $\varepsilon$ parameter

eqs = Normal[Series[f, {\[Epsilon], 0, ord}]];
SolveAlways[eqs == 0, \[Epsilon]]
(* => {{a[3] -> 2 (a[1] - a[1]^3 + a[1] Cos[4 x]), a[2] -> -2 a[1] Cos[2 x], Sin[2 x] -> -a[1]}} *)


If we are interested in the expansion coefficients $a_n$ in terms of $x$ then following code will refine the solution

sol = First@Solve[Equal@@@First@SolveAlways[eqs == 0, \[Epsilon]], a/@Range[ord]]
(* => {a[1] -> -Sin[2 x], a[2] -> 2 Cos[2 x] Sin[2 x],  a[3] -> 2 Sin[2 x] (-1 - Cos[4 x] + Sin[2 x]^2)} *)


One should also do check

f/.sol//Simplify
(* => O[\[Epsilon]]^4 *)

-

First, one can express $\varepsilon$ as a function of $y$:

Clear[eps]
eps[y_] = ε /. Solve[x == y - ε Sin[2 y], ε][[1]]


$$-(x-y) \csc (2 y).$$

Then, expand into series around $y=x$, since it is a solution for $\varepsilon=0$. Applying InverseSeries after that gives:

res = InverseSeries[Series[eps[y], {y, x, 3}]] /. y -> ε // FullSimplify


$$x+\varepsilon \sin (2 x)+\varepsilon ^2 \sin (4 x)+\varepsilon ^3 \sin (2 x) (3 \cos (4 x)+1)+O\left(\varepsilon ^4\right).$$

Check:

Series[y - ε Sin[2 y] /. y -> res, {ε, 0, 3}] // FullSimplify


$$x+O\left(\varepsilon ^4\right).$$

-