0
$\begingroup$

I have a code given by

V =((1/2 - (3 r)/2 + r^2) (r^8 (-1 + (3 r)/2) + 
4 (2 - (9 r)/2 + 2 r^2) + 2 r^4 (3 - 9 r + 6 r^2)))/(r^6 (2 + 
r^4)^2);

X= r + 5/2 ArcTan[4 (-(3/4) + r)] + 3/4 Log[1/2 - (3 r)/2 + r^2];

zz = X /. r -> 3/2

fV[z0_?NumericQ] :=With[{z = SetPrecision[z0 + zz, 100 + 1]},If[Abs[z] 
<= 35,Re[V /. FindRoot[X == z, {r, 10000001/10000000}, 
  MaxIterations -> 10000, WorkingPrecision -> 100]], 0]]

My problem: consider this formula*

Subscript[\[Psi], m, n, p] = 
Subscript[\[Psi], m - 1, n - 1, p] + Subscript[\[Psi], m, n, p - 1] -
Subscript[\[Psi], m - 1, n - 1, p - 1] - 
h^2/8 (fV[(p - (n - 1))/2 h] Subscript[\[Psi], m - 1, n - 1, p] + 
 fV[((p - 1) - n)/2 h] Subscript[\[Psi], m, n, p - 1])

with the following initial conditions

Subscript[\[Psi], 0, 0, 0] = 1;
Subscript[\[Psi], 1, 1, 0] = 1;
Subscript[\[Psi], 2, 2, 0] = 1;
Subscript[\[Psi], 3, 3, 0] = 1;
Subscript[\[Psi], 4, 4, 0] = 1;
Subscript[\[Psi], 5, 5, 0] = 1;
Subscript[\[Psi], 6, 6, 0] = 1;
Subscript[\[Psi], 7, 7, 0] = 1;
Subscript[\[Psi], 8, 8, 0] = 1;
Subscript[\[Psi], 9, 9, 0] = 1;
Subscript[\[Psi], 10, 10, 0] = 1;
Subscript[\[Psi], 0, 0, p] = Exp[-0.25 (p h)^2];
h = .5;

I need a table given by

Table[{t, Subscript[\[Psi], t/h, t/h, t/h]}, {t, 0, 5, h}]

as the result of calculations. How can I tell the Mathematica to calculate Subscript[\[Psi], m, n, p] from formula* and provide me the mentioned table? I will be thankful if someone helps.

$\endgroup$

1 Answer 1

3
$\begingroup$

You could define a function instead of using Subscript (this is in general a good idea) and employ memoization (in order to preventing that values have to be computed multiple times:

ClearAll[fV];
fV[z0_?NumericQ] := fV[z0] = With[{z = SetPrecision[z0 + zz, 100 + 1]}, 
   If[Abs[z] <= 35, 
    Re[V /. FindRoot[X == z, {r, 10000001/10000000}, 
       MaxIterations -> 10000, WorkingPrecision -> 100]], 0]];

h = 1/2;
ClearAll[Ψ];
Ψ[0, 0, 0] = 1;
Ψ[1, 1, 0] = 1;
Ψ[2, 2, 0] = 1;
Ψ[3, 3, 0] = 1;
Ψ[4, 4, 0] = 1;
Ψ[5, 5, 0] = 1;
Ψ[6, 6, 0] = 1;
Ψ[7, 7, 0] = 1;
Ψ[8, 8, 0] = 1;
Ψ[9, 9, 0] = 1;
Ψ[10, 10, 0] = 1;
Ψ[0, 0, p_] := Ψ[0, 0, p] = Exp[- 1/4 (p h)^2];
Ψ[m_, n_, p_] := Ψ[m, n, p] = Ψ[m - 1, n - 1, p] + Ψ[m, n,  - 1] - Ψ[m - 1, n - 1, p - 1] - h^2/8 (fV[(p - (n - 1))/2 h] Ψ[m - 1, n - 1, p] + fV[((p - 1) - n)/2 h] Ψ[m, n, p - 1])
Table[{t, Ψ[t/h, t/h, t/h]}, {t, 0, 5, h}]

{{0, 1}, {1/2, 0.936178}, {1, 0.766748}, {3/2, 0.54557}, {2, 0.330802}, {5/2, 0.161054}, {3, 0.0478391}, {7/2, -0.0171044}, {4, -0.0494697}, {9/2, -0.0633942}, {5, -0.068163}}

$\endgroup$
3
  • $\begingroup$ Thank you very much. When I employ it for higher values of t, say t=21, with h=1/10 and \[CapitalPsi][m_, n_, 0] := \[CapitalPsi][m, n, 0] = 1;, I cannot make its table and it gives some FindRoot errors. Could you please tell me why? $\endgroup$
    – Mehrab
    Feb 22, 2019 at 13:56
  • $\begingroup$ Consider h=1/50 $\endgroup$
    – Mehrab
    Feb 22, 2019 at 14:30
  • $\begingroup$ I obtain plenty of recursion limit errors. Probably your initial conditions are not sufficient; the recursive evaluation of Ψ does not end up in one of the cases Ψ[0, 0, 0], Ψ[1, 1, 0],...., or Ψ[0, 0, p]. To Ψ calls itseld recursively over and over again infinitely often. Of course, a computer cannot do that so Mathematica breaks this vicious circle at some point and throws an error. $\endgroup$ Feb 22, 2019 at 15:19

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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