2
$\begingroup$

I'm using the following code for calculation of reflection and transmission coefficients through square well potential:

minE = 5;
maxE = 200;
nE = 200;
energies = Subdivide[minE, maxE, nE - 1];
reflections = ConstantArray[0, nE];
transmissions = ConstantArray[0, nE];

k = Sqrt[2 energies];
Do[
 solution = 
  NDSolve[{psi''[x] == -2 (energies[[i]] - v[x]) psi[x], 
    psi[xMax] == 1, psi'[xMax] == I*k[[i]]}, psi, {x, -xMax, xMax}];
 max = FindMaximum[
    Abs[psi[x]] /. solution, {x, -.8 xMax, -.2 xMax}][[1]];
 min = FindMinimum[
    Abs[psi[x]] /. solution, {x, -.8 xMax, -.2 xMax}][[1]];
 reflections[[i]] = ((max - min)/(max + min))^2;,
 {i, nE}]

reflFunc = ListInterpolation[reflections, {{minE, maxE}}];
transFunc = ListInterpolation[1 - reflections, {{minE, maxE}}];
Plot[{reflFunc[x], transFunc[x]}, {x, minE, maxE}]

Now, this is kind of problem where I loop over some variable (energy) with manually defined accuracy (by nE and Subdivide). There are two questions:

  1. Is there a way to do it with functions instead of loop and lists (how to directly create function I finally get through ListInterpolation) and how to control accuracy/ "energy step"/nE in that case

  2. If there is another more elegant way, I'd be most obliged to see it

$\endgroup$
  • $\begingroup$ As far as I know, for 1D Schrödinger square potential problem, one can find the analytical expression of reflection/transmission as a function of energy. Wolfram can facilitate the derivation work. $\endgroup$ – Αλέξανδρος Ζεγγ Jul 23 at 0:25
0
$\begingroup$

After some reading I'm down to following code, it somehow gets everything into function, but it is surely not most elegant:

ffunc[energy_] := (
  k = Sqrt[2 energy];
  solution = 
   NDSolve[{psi''[x] == -2 (energy - v[x]) psi[x], psi[xMax] == 1, 
     psi'[xMax] == I*k}, psi, {x, -xMax, xMax}];
  max = FindMaximum[
     Abs[psi[x]] /. solution, {x, -.8 xMax, -.2 xMax}][[1]];
  min = FindMinimum[
     Abs[psi[x]] /. solution, {x, -.8 xMax, -.2 xMax}][[1]];
  ((max - min)/(max + min))^2)

v[x_] := If[x > 0 && x < 1, 50, 0];
xMax = 10;
minE = 5;
maxE = 200;
nE = 200;
energies = Subdivide[minE, maxE, nE - 1];
reflections = ffunc /@ energies;
reflFunc = ListInterpolation[reflections, {{minE, maxE}}];
transFunc = ListInterpolation[1 - reflections, {{minE, maxE}}];
ListPlot[reflections]

Anyway, this example evaluates in 9.36s which is a bit longer than original code at 9.02s. Plot with MaxRecursion->0 takes 10.86

I am aware this change does nothing functional, but at least body of Do command is in function now.

| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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