Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Penetrating the Rectangular Potential Barrier, E < V0

  1. Jan 25, 2016 #1
    The following code compares my result to the known correct result. Can anyone find where my error is? If you run the octave code you will see that my answer has positive concavity whereas the correct one has negative concavity...I've checked this so many times and can't find my error. I wanted to see if I could correct it at this stage before trying to simplify to the known answer...thanks for any help.


    %I wrote this to check my work thus far on solving the GIANT algebra problem of penetrating the rectangular barrier, where E < V0
    % setting some constants

    hbar = 1.0545718*10^(-34);
    m = 9.10938356 * 10^(-31);
    a = 10^-20; %can be purely arbitrary

    V0 = 100; %can be purely arbitrary
    inc = 1;
    E = inc:inc:(V0-inc);

    %wavenumber formulas
    k1 = (sqrt(2*m*E))/hbar;
    k2 = (sqrt(2*m*(V0-E)))/hbar;

    %this section contains the solution
    TSoln = (1 + (((sinh(k2*a)).^2)./((4*E/V0).*(1-E/V0)))).^(-1)

    %this section contains my work thus far

    iVal = i*k2/k1;
    BLABLA1 = ( 1 - ((1-iVal).*(e.^(-2*k2*a) )./(1+iVal) ) ) + iVal.* (1 + ((1-iVal).*(e.^(-2*k2*a) )./(1+iVal)));
    BLABLA2 = ( 1 - ((1+iVal).*(e.^(2*k2*a) )./(1-iVal) ) ) - iVal.* (1 + ((1+iVal).*(e.^(2*k2*a) )./(1-iVal)));

    C = (e.^(-i*k1*a)).* ((e.^(-k2*a)./ ( BLABLA1 ) ) + (e.^(k2*a)./ ( BLABLA2 )));
    TChip = 4*(C.*conj(C))

    plot (E, [TSoln;TChip]);
  2. jcsd
  3. Jan 25, 2016 #2
    Why not numerically integrate the differential equation. If the sign of the second derivative is right, the concavity of the solution will also be right.
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook