# Solving coupled eigenvalue differential equations

I am trying to solve an equation of the form as follows

$\left(\begin{array}{cc} -\frac{\hbar^{2}}{2m}\frac{\delta^{2}}{\delta z^{2}}+\sin^{2}\left(z\right) & z\\ z & \frac{\hbar^{2}}{2m}\frac{\delta^{2}}{\delta z^{2}}-\sin^{2}\left(z\right) \end{array}\right)\left(\begin{array}{c} u_{n}(z)\\ v_{n}(z) \end{array}\right)==\epsilon_{n}\left(\begin{array}{c} u_{n}(z)\\ v_{n}(z) \end{array}\right)$

Assuming required no. of boundary conditions but with no knowledge of the eigenvalues, please tell how to solve such a system of equations in Mathematica. Here $n$ is an index for the eigenfunction and the corresponding eigenvalue. Take $z$ as real and assume that $z$ runs on a finite lattice of fixed size of length $L=Na$ where $N$ is the number of sites on the lattice and a is the lattice constant. If you take $a==1$, take $N$ of the order of $10000$. There is a specified value of the eigenfunctions at the boundary of the lattice which you can take zero for the simplest case. Please multiply everything by units such that dimensions are consistent- Here the derivative term is the kinetic energy term while sin2z term is the potential term.

-
Is this a Mathematica related question? If so do you have any Mathematica code. It would save the time to type the equation.... –  user21 Jun 20 '13 at 7:03
I want to know how to solve this in Mathematica. There are algorithms to solve the coupled eigenvalue differential equations but I have not been able to find a Mathematica command(s) that directly utilizes them for solving such a set of equations. And that or a combination of commands to help solve this, is what I am looking for. –  cleanplay Jun 20 '13 at 10:57
What kind of solution is expected? A closed form may be not attainable. Where seems to be no a ready procedure to find eigenvalues numerically, so probably an implementation of an algorithm you're mentioning would be needed. –  Andrew Jun 20 '13 at 11:41
You're units are dimensionally inconsistent. The derivative term has units of $s^{-1}$, but $\sin^2(z)$ is unitless. Also, what are the boundary conditions? They effect the eigenvalues of the system. Also, as this appears to be a Hamiltonian, is $z$ real? If not, shouldn't the matrix be Hermitian? –  rcollyer Jun 20 '13 at 13:24
@rcollyer I somtimes find in my own research that $2=1$ :D –  sebhofer Jun 20 '13 at 15:26

Sorry, I formulate my comments here, because there is not enough space in the comment section.

1) The dimension question cannot be ignored. The original Schrödinger equation should look like

$-\frac{h^2}{2m} \left(\frac{d}{\text{dx}}\right)^2 u(x)+A \sin ^2\left(\frac{2 \pi x}{L}\right) u(x)+B x v(x)=a u(x)$

Letting $x \to q z$ and $u \to p U, v \to p V$ the equation becomes

$-\frac{h^2 p}{2m q^2} \left(\frac{d}{\text{dz}}\right)^2 U(z)+A p \sin ^2\left(\frac{2 \pi q z}{L}\right) U(z)+\left(B p q z\right) V(z)=\left(a p\right)U(z)$

Now we can set two factors equal to 1 letting $q =\left (\frac{h^2}{2 m B} \right)^\frac {1}{3}$ and $p =\left (\frac{2m}{h^2 B^2} \right)^\frac {1}{3}$ resulting in

$- \left(\frac{d}{\text{dz}}\right)^2 U(z)+S \sin ^2\left(R z\right) U(z)+V(z)=T U(z)$

an equation with two parameters S and R and eigenvalue T. The second equation has a similar form.

2) The system seems to describe the motion of a charged particle in a periodic potential in a constant electric field. (I would be helpful if you describe the phsical problem). The equations are hence a mixture between Mathieu- and Airy differential equations.

3) Both are tricky when it comes to Eigenvalues. Airy is easier, it corresponds to the motion of a charged particle in a constant electric field; the eigenvalues are related to the zeroes of the Airy-function. Mathieu (coming from the Sin^2-terms) is much more complicated. It leads, roughly speaking, to energy "bands". You might wish to consult a book on quantum mechanics.

4) It's always a good idea to simplify the equations. Here it means neglecting first the z-Term -> Mathieu, and then the sin^2-term -> Airy, and study the two areas separately.

5) Once you have clearly formulated your physical problem you can go and solve the system in MMA numerically. I have done it with simple boundary conditions like {u,v}[0]==0 and {u',v'}[0] = {1,1} and varying a and b. No problem but of no use without a proper question.

Best regards, Wolfgang

-