Numerically solve recurrence relation

  • Context: Mathematica 
  • Thread starter Thread starter Poirot
  • Start date Start date
  • Tags Tags
    Recurrence Relation
Join the discussion
Registration is free. Ask a follow-up in this thread, or start your own.
5 replies · 2K views
Poirot
Messages
87
Reaction score
2
I am trying to numerically solve the recurrence relation for:
$$ T[n]=\frac{sm}{2\hbar^2i k_n}(T[n-1]+T[n+1]) \\
k_n=(\frac{2m(E-\hbar \omega n)}{\hbar^2})^{1/2}$$
My code is:
Code:
RecurrenceTable[{Tn[
    n] == (s*
      m/(2*\[HBar]^2*
         I*(2*m*(En - \[HBar]*\[Omega]*n)/\[HBar]^2)^(1/2)))*(Tn[
       n - 1] + Tn[n + 1]), Tn[-2] == 1, Tn[2] == 1}, Tn, {n, -2, 2}]
Taking all the constants to be 1, except E which I have plugged to be zero, to ensure non negative root.

I think I need my seed values to be $$T_{-2}$$ and $$T_{2} (=1)$$ so it's at the end of the boundary but I keep getting the error that the recursion depth of 1024 has been exceeded and I've tried "Block" and set it as high as possible without crashing the kernel to no avail. If I plug $$T_{-2}=T_{-1}=1$$ I actually get an answer but if I use interpolating to see what this looks like it's blowing up exponentially (I'm taking the Abs of it as complex). I was told to look for convergence but to me, this seems like divergence, which may be due to not being able to set the end points.

How can I actually plug end points into get a reasonable answer?

Thanks in advance, and apologies for my lack of coding ability, I'm very new to Mathematica.
 
Physics news on Phys.org
I'm not sure what you mean by "solve" here.

But this looks something like the wave function of a harmonic oscillator. If that's what it is, think what you know about such a system. Why do you say
E which I have plugged to be zero, to ensure non negative root.
when you have
2m(E−ℏωn)ℏ2
? It looks as though you should be getting the square root of a negative number now unless n is always negative. (Are negative values of n allowed?) I'd multiply by the kn to get them out of the denominator. (This is a clue.)
 
John Park said:
I'm not sure what you mean by "solve" here.

But this looks something like the wave function of a harmonic oscillator. If that's what it is, think what you know about such a system. Why do you say when you have ? It looks as though you should be getting the square root of a negative number now unless n is always negative. (Are negative values of n allowed?) I'd multiply by the kn to get them out of the denominator. (This is a clue.)
Sorry I made a mistake, I meant to say En is set to 10 to ensure non negative for now.

This system is actually for an oscillating delta potential $$V(x,t) = s cos(\omega t) \delta(x)$$ and the T I'm looking at here is for the transmission coefficient.
n runs from -∞ to ∞ (as a consequence of a Fourier transform and Floquet theorem) but I'm taking a finite range for now, I think to see what it tends to.

I tried moving the kn factor over to the other side and evaluating but get the same error.

Thanks for your help, by the way, greatly appreciated!
 
Sorry I made a mistake, I meant to say En is set to 10 to ensure non negative for now.

I don't think I was very helpful, but what struck me was the term involving (E − ℏ ω n ). In the simple-harmonic-osciilator problem, a very similar term in a recurrence relation vanishes if E has the appropriate value; all higher terms also vanish. This automatically gives the energy levels of the SHO, but it is also crucial to the convergence of the wave function. Unless the series terminates at some n, the wave function does increase exponentially at large distances.

Edit: This factor (power series) in the wave function increases exponentially; the wave function itself "merely" becomes un-normalisable.

But do you see parallels with your problem?

Also, if your (E − ℏ ω n ) becomes negative, the whole factor involving kn becomes real, if s is real: what is the significance of that? (Not a trick question: I don't know.)
 
Last edited:
Thank you for your replies, I found another way of doing this by using the "Do" function and "NSolve", which allowed me to define values outside of the interval I was looking at (which also means I can force it to tail off for large n/-n to make it physically realisable.