next up previous
Next: C Progam for calculating Up: Mathieu's Equation, solution, and Previous: Hill's Method solution

Sträng's recursion formula for $ \Delta (0)$

First we note that by the symmetry of $ \Delta (0)$ , $ \gamma_{-n} = \gamma_{n}$ . Following Sträng, we define

$\displaystyle A_i = \left(\begin{array}{ccccccccc}
&&&&&&0&\gamma_{2i}&1\end{array}\right)$ (2.10)

We have $ \Delta_i=det(A_i) \ni \Delta(0)=\lim_{i\rightarrow \infty}\Delta_i.$ We can decompose $ A_i$ in terms of $ A_{i-1}$ ,

$\displaystyle A_i = \left(\begin{array}{ccccc}
\end{array}\right)$ (2.11)

A Laplace decomposition yields

$\displaystyle det(A_i)=\left\vert\begin{array}{cccc}
\end{array}\right\vert$ (2.12)

Here $ rA_{i-1}$ represents $ A_{i-1}$ with its left most column chopped off. Again, following Sträng we define $ lA$ as the matrix A with its rightmost column removed, $ uA$ the matrix A with its lowest row removed, $ lA$ the matrix A with its upper most row removed. Ultimately, $ uldr(A_{i-1})=A_{i-2}$ , and given the symmetry involved, det$ (rd(A_{i-1}))=$det$ (ul(A_{i-1}))$ . Following this procedure we find

$\displaystyle \Delta_i=\Delta_{i-1}-2\gamma_{2i} \gamma_{2(i-1)}$det$\displaystyle (rd(A_{2(i-1)}))+(\gamma_{2i} \gamma_{2(i-1)})^2 \Delta_{i-2}$ (2.13)

We also note, similarly using Lapalcian decomposition,

$\displaystyle \Omega_i=$det$\displaystyle (ul(A_i))=$det$\displaystyle (rd(A_i))  \rightarrow  \Omega_i=det(A_{i-1})-\gamma_{2i} \gamma_{2(i-1)} \Omega_{i-2}$

so that

$\displaystyle \frac{\Delta_{i-1}-\Omega_i}{\gamma_{2i} \gamma_{2(i-1)}}=\Omega_{i-1}=$det$\displaystyle (rd(A_{i-1}))$


$\displaystyle \Delta_i = \Delta_{i-1} + 2(\Omega_i-\Delta_{i-1})+(\gamma_{2i} \gamma_{2(i-1)})^2\Delta_{i-2}$

$\displaystyle \frac{\Delta_i + \Delta_{i-1} -(\gamma_{2i} \gamma_{2(i-1)})^2\De...
... \Delta_{i-2} -(\gamma_{2(i-1)} \gamma_{i-2})^2\Delta_{i-3})}{2} = \Omega_{i-1}$

Plugging this into Equation 2.13,

$\displaystyle \Delta_i = (1-\gamma_{2i} \gamma_{2(i-1)})\Delta_{i-1}+\left( (\g...
...i-2} +\gamma_{2i} \gamma_{2(i-1)}(\gamma_{2(i-1)}\gamma_{2(i-2)})^2\Delta_{i-3}$

Define $ \alpha_{2i} = \gamma_{2i} \gamma_{2(i-1)}$ and $ 1-\alpha_{2i}=\beta_{2i}$ and find,

$\displaystyle \Delta_i = \beta_{2i}\Delta_{i-1}-\alpha_{2i} \beta_{2i} \Delta_{i-2} + \alpha_{2i}\alpha_{2(i-1)}^2\Delta_{i-3}$ (2.14)

We can recursively solve for $ \Delta(0)=\lim_{i\rightarrow \infty}\Delta_i$ to as much accuracy as nessessary, though the program presented below found convergence to a fair tolerance quite quickly. We first must ``seed'' the recurssion with the first three $ \Delta_i$ . This can be done by hand, though we have deferred to the kindess of our computer algebraic program Maple instead.

Maple finds,

dc := -2*e2^2*e0*e4^2*e6+e2^2*e4^2-2*e4^2*e2*e0*e6^2+2*e2*e4^2*e6

A:= matrix([[1,e4,0,0,0],[e2,1,e2,0,0],[0,e0,1,e0,0],
da := 1-2*e2*e4-2*e2*e0+2*e2^2*e0*e4+e2^2*e4^2

db := 1-2*e2*e0

Our program seeks to find all stable values of $ \mu $ , i.e. those that satisfy Equation 2.9 as real values (i.e. all iso-$ \mu $ for which $ \mu $ is exclusively imaginary.

Our code finds all such iso-$ \mu $ by looping through the a and q axis. If our $ \mu $ formula returns ``nan'' which is the C language's way of saying not a real number, then we set the value of $ \mu $ to zero, though of course it is only the imaginary part of $ \mu $ which is actually zero. We perform a contour plot on our data output and find the elegant avian like image of the stability region of Mathieu's equation (Figure 2.1).

Figure 2.1: Stability diagram for Mathieu's general equation

For the quadropole field, the rf linear Paul trap, we have the following stability regime (Figure 2.2). The original stability diagram is simply reflected about the x-axis as $ a \rightarrow -a$ between the two.

Figure 2.2: Stability diagram for linear rf Paul trap; stable regions are those in which the two stability diagrams intersect.

Figure 2.3: Stability diagram for linear rf Paul trap, closer view of lowest region of mutal stability (intersection).

For the chamber rf Paul trap, we recall that for the z direction we must allow for $ (a,q) \rightarrow (2a,2q)$ (Equation 1.3), and so the lowest region of stability (and the largest region at that) has a slightly different look (Figure 2.4).

Figure 2.4: Stability diagram for rf Paul trap, lowest region of stability (intersection of isobars).

next up previous
Next: C Progam for calculating Up: Mathieu's Equation, solution, and Previous: Hill's Method solution
tim jones 2008-07-07