Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I am trying to solve the 1-D Schrödinger equation with the Crank-Nicholson method. I use the atomic units:

  • time < 10^-14 => 413 a.u.t.
  • length 30 * 10^-9 m => 567 Bohr radii

a = 50; 
σ = a/10;
p[x_, 0] := Exp[-(x^2)/(2 σ^2)]; (* граничные условия *)
p[\[PlusMinus]a, t_] := 0; (* н.у. *)
NN = 100;
h = 1/NN;
repl = x -> ((i*6*a)/NN - 3*a); 

eq = D[Subscript[p, i][t], t] == I/2(Subscript[p, i - 1][t] + 2 Subscript[p, i][t] - 
  Subscript[p, i + 1][t])/h^2(*+D[u[x],x]*Subscript[p,i][t]*)/. repl; 

Table[eq, {i, 0, NN}];
boundary = {Subscript[p, 0][t] == 0, Subscript[p, NN][t] == 0}; 
p0[x_] := Exp[-(x^2)/(2 σ^2)]; 
Cauchy = Table[Subscript[p, i][0] == p0[x] /. repl, {i, 1, NN - 1}]; (*задача Коши*)

eqns = Join[Table[eq, {i, 1, NN - 1}], boundary, Cauchy];
sol = NDSolve[N[eqns], 
N[Table[Subscript[p, i][t], {i, 0, NN}]], {t, 0, 413}];

and get this error:

NDSolve::mconly: For the method IDA, only machine real code is available. Unable to continue with complex values or beyond floating-point exceptions. >>

I want to get something like this:

sol = NDSolve[{I D[u[t, x], t] == (-1/2)D[u[t, x], {x, 2}], 
     u[0., x] == Exp[-(x^2./(2*σ^2))], u[t, a] == 0, 
     u[t, -a] == 0}, u, {t, 0, 4130}, {x, -a, a}, (*MaxStepSize->0.01,*)
     AccuracyGoal -> 3, PrecisionGoal -> 3];

Animate[Plot[Evaluate[Abs[u[t, x] /. First[sol]]^2], {x, -a, a}, 
     PlotRange -> {0, 1}], {t, 0, 413}]

What's wrong with my code and how to fix it?

share|improve this question
    
Back in the late 90s I recall Ron Knapp writing an example for CN for Schrodinger. Please search the Wolfram library and I think you should be able to locate it. –  Mike Honeychurch Dec 13 '13 at 2:18
    
@MikeHoneychurch I tried to search but can't find Ron Knapp's example, could you search if you have time? I'm very interested in this. But I find this example by Terry Robb, which solves the TDSE using a external Fortran subroutine SCHROED. And SCHROED uses a Crank-Nicholson scheme. –  xslittlegrass Mar 6 '14 at 4:51

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.