Solving the Schrodinger equation using sage (Well)

savecancel

Hello, @Grim! I am not quite sure about Schrödinger's Equation, but I know how to obtain $k=\frac{n\pi}{a}$. I believe this will result a helpful version of your code:

k, a = var('k, a', domain='positive')
psi = function('psi')(x)
eqs = diff(psi, x, 2) + k^(2)*psi == 0
EQS = desolve(eqs, psi, ivar=x)
A, B = var('A, B')
EQS = EQS.subs(_K1=B, _K2=A)
show(LatexExpr(r'\psi(x) = '), EQS)
eq1 = EQS.subs(x=0) == 0  # this is psi(0) = 0
eq2 = EQS.subs(x=a) == 0  # this is psi(a) = 0
show(eq1)
show(eq2)

The first line defines the symbolic variables k and a to be positive, so you can avoid using the assume command. The second line defines psi as the unknown function; the third line defines the ordinary differential equation; and the fourth one solves the differential equation.

Now, things get a little complicated. If you write print(EQS), you will notice that there are a couple of new variables called _K1 and _K2. They are integration variables (byproducts of the resolution of the ODE). However, they are not symbolic variables that you have defined, so you won't be able to use them outside of EQS (for example, to use the solve command or any other). So, I replaced them with A and B in the sixth line. The seventh line (the first show command) is there only to check that I have obtained the same definition of $\psi$ as you in your code (it will show you $\psi(x)=A\cos(kx)+B\sin(kx)$). Then, eq1 and eq2 represent your boundary conditions $\psi(0)=0$ and $\psi(a)=0$, respectively.

At this point, there are many ways in which you can obtain the values of $A$, $B$ and $k$. The first method is manually solving equations eq1 and eq2, which are shown by the last two lines of my code. You will see that

$$A=0$$ and $$A\cos(ak)+B\sin(ak)=0$$

If you replace the first equation in the second, you obtain $B\sin(ax)=0$. Now, since $B>0$ (according to your code), you must have that $\sin(ak)=0$, which implies that $k=\frac{n\pi}{a}$. Finally, in this case, $B$ can be any positive number (regardless of its value, your ODE and boundary conditions are satisfied.)

The second way to solve this is using Sage itself. You have two equations, namely eq1 and eq2, so you will need two variables. You can write

sol = solve([eq1, eq2], A, k)
show(sol)

just after the previous code. You will get something like

$$\left[\left[A = 0, k = \frac{\pi + 2\pi z_{1879}}{a}\right], \left[A = 0, k = \frac{2 \pi z_{1900}}{a}\right]\right]$$

The variables z1879 and z1900 are similar to the _K1 and _K2 that I mentioned above. More precisely, they are integer parameters; just mentally replace them with $m$ and $n$. Then you have two possible solutions. After further analysis, you will notice that these two solutions can be written more simply as $A=0$ and $k=\frac{n\pi}{a}$. Once again, $B$ has no restrictions, besides being positive.

You can also use the Sympy package to solve these equations:

sol = solve([eq1, eq2], A, k, algorithm='sympy')
show(sol)

You will get

$$\left[\left\{A : 0, k : 0\right\}, \left\{A : 0, k : \frac{\pi}{a}\right\}\right]$$

These are just two particular solutions that you can easily generalize using reduction formulas and properties of the trigonometric functions.

What happens if you use solve([eq1, eq2], A, B) or solve([eq1, eq2], B, k)? You will get $[[A=0, B=0]]$ and $[]$, respectively. The first is a correct solution to both equations, but not one that you are looking for. The second is Sage telling you that could not find a solution. If you use SymPy, you will get $\left[\left\{B: -\frac{A}{\tan(ak)}\right\}\right]$, which is a way of admission of defeat, but it is more useful, because it hints you that it is Sage's (actually SymPy's) fault, while the $[]$ before could lead you to think that there is no solution.

I would like to point out a few interesting teachings derived from your question:

1. Do not depend on any software to solve a problem; it's just a tool.
2. SageMath is excellent in solving problems, but it is not perfect. Try different approaches (e.g., as we did when we used SymPy), especially if the problem is complex.
3. You can use a software to help you, but the more of the resolution you can do manually (using ingenuity), the more sure you will be of your result. In the end, if you can do it manually, why bother using Sage? In the end, there are things that are easy for a human mind, but very difficult for a software.