12
$\begingroup$

Is there a function in Mathematica that can be used to find the perturbation solution of an equation like $x^2 − 1 = \epsilon \,x$, $x − 2 = \epsilon \cosh(x)$ or $x^2 − 1 = \epsilon\, e^x$?

$\endgroup$
  • $\begingroup$ What would you expect the output to look like? $\endgroup$ – bill s Dec 6 '16 at 3:09
  • $\begingroup$ @bills, Something like $x=a_0+a_1\epsilon +a_2\epsilon^{2}+...$ $\endgroup$ – MrDi Dec 6 '16 at 3:10
  • 8
    $\begingroup$ There's a tutorial about this exact thing. $\endgroup$ – rcollyer Dec 6 '16 at 3:36
13
$\begingroup$

Decide up to which power you would like to expand:

pow = 4;

Let's do one of the equations you mentioned as an example (bring all terms to one side and save as a single expression, eq in this case):

eq = x - 2 - e Cosh[x];

Write an ansatz for the x solution with unknown coefficients:

ansatz = Sum[a[i] e^i, {i, 0, pow}];

Substitute the ansatz into your expression and do a series expansion:

expand = Series[eq /. x -> ansatz, {e, 0, pow}];

Solve the constraints of overall factors vanishing:

vars = Table[a[i], {i, 0, pow}];
sols = Solve[expand == 0, vars][[1]] // Simplify // Quiet;

Finally, insert the solution into your ansatz to obtain the result:

res = ansatz /. sols

2 + e Cosh[2] + e^2 Cosh[2] Sqrt[-1 + Cosh[2]^2] + 1/3 e^4 Cosh[2] Sqrt[-1 + Cosh[2]^2] (-3 + 8 Cosh[2]^2) + e^3 (-Cosh[2] + (3 Cosh[2]^3)/2)

Don't forget to test numerically, whether your result is actually correct at the end:

enum = 10^-10;
xnum = x /.FindRoot[(eq /. e -> enum) == 0, {x, 2}, WorkingPrecision -> 60];
(res /. e -> enum) - xnum

-3.627605747396*10^-47

This shows that expanding to fourth order with an e= 10^-10 indeed consistently matches the result up to about 10^-50 accuracy, so the expansion was correct. Rinse and repeat for the other examples.

PS: In cases where your equation admits several solutions you might have to be a little bit more careful, but the principle still stays the same.

$\endgroup$
6
$\begingroup$

For these simple examples you could use InverseSeries. For example:

InverseSeries[Series[(x^2 - 1)/x, {x, 1, 10}], e]
InverseSeries[Series[(x - 2)/Cosh[x], {x, 2, 10}], e]
InverseSeries[Series[(x^2 - 1)/E^x, {x, 1, 10}], e]

You need to solve for e, and then do a series around the value of x when e is zero.

$\endgroup$
3
$\begingroup$

The new in M12 function AsymptoticSolve can be used to find these perturbation expansions:

AsymptoticSolve[x^2 - 1 == ϵ x, {x, 1}, {ϵ, 0, 7}]

{{x -> 1 + ϵ/2 + ϵ^2/8 - ϵ^4/128 + ϵ^6/ 1024}}

AsymptoticSolve[x - 2 == ϵ Cosh[x], {x, 2}, {ϵ, 0, 7}]

{{x -> 2 + ϵ Cosh[2] + ϵ^2 Cosh[2] Sinh[2] + 1/2 ϵ^3 (Cosh[2]^3 + 2 Cosh[2] Sinh[2]^2) + 1/3 ϵ^4 (5 Cosh[2]^3 Sinh[2] + 3 Cosh[2] Sinh[2]^3) + 1/24 ϵ^5 (13 Cosh[2]^5 + 88 Cosh[2]^3 Sinh[2]^2 + 24 Cosh[2] Sinh[2]^4) + 1/15 ϵ^6 (47 Cosh[2]^5 Sinh[2] + 100 Cosh[2]^3 Sinh[2]^3 + 15 Cosh[2] Sinh[2]^5) + 1/720 ϵ^7 (541 Cosh[2]^7 + 7746 Cosh[2]^5 Sinh[2]^2 + 7800 Cosh[2]^3 Sinh[2]^4 + 720 Cosh[2] Sinh[2]^6)}}

AsymptoticSolve[x^2 - 1 == ϵ Exp[x], {x, 1}, {ϵ, 0, 7}]
AsymptoticSolve[x^2 - 1 == ϵ Exp[x], {x, -1}, {ϵ, 0, 7}]

{{x -> 1 + (E ϵ)/2 + (E^2 ϵ^2)/8 + (E^3 ϵ^3)/16 + ( 13 E^4 ϵ^4)/384 + (E^5 ϵ^5)/48 + (69 E^6 ϵ^6)/ 5120 + (841 E^7 ϵ^7)/92160}}

{{x -> -1 - ϵ/(2 E) + (3 ϵ^2)/(8 E^2) - (7 ϵ^3)/( 16 E^3) + (235 ϵ^4)/(384 E^4) - (121 ϵ^5)/(128 E^5) + ( 7959 ϵ^6)/(5120 E^6) - (245953 ϵ^7)/(92160 E^7)}}

$\endgroup$
2
$\begingroup$

According to standard perturbation theory for static Hamiltonians of this type, the ground state of the system may be approximated by:

\begin{align} |0\rangle_V &= |0\rangle_H - \epsilon \sum_{\alpha\neq 0} \frac{U_{\alpha 0}}{\alpha} |\alpha\rangle_H \\ \notag &\simeq N^{-\frac12} \left ( |0\rangle_H + \beta_1 |1\rangle_H + \beta_2 |2\rangle_H \right)\,, \end{align} where $N=1+\beta_1^2+\beta_2^2$, $\beta_k = -\epsilon U_{\alpha 0}/\alpha$ and $U_{\alpha 0}={}_H\langle \alpha|U| 0 \rangle_H$, $\alpha=1,2$.

Hope this helps to formulate your problem.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.