Решение задачи о туннелировании

Наконец, следуя нашему принципу обучения посредством решения одних и тех же задач разными способами, повторим в Mathematica пример расчета прохождения частицы через одномерный прямоугольный потенциальный барьер [2, с. 101–102], разобранный в Теме 1. Программа имеет следующий вид (вместо E используем V).

u0=40;m=1/2;h=1;

(V> 0) ™True

(V < u0) ™True

Обратим внимание на формат оператора, присваивающего свойства переменным – аналог оператора assume в Maple. Далее описываем величину ε и уравнения:

epsilon=(2*m/h^2)*V;

eq1=f1''[x]+(2 m/h^2) V f1[x]Š0

eq2=f2''[x]+(2 m/h^2) (EE-u0) f2[x]Š0

eq3=f3''[x]+(2 m/h^2) EE f3[x]Š0

Очень важное наблюдение: в Mathematica есть три вида «знаков равенства»: обычный знак равенства =, двоеточие со знаком равенства:= и двойной знак равенства ==. Первый является оператором присваивания для символов (переменных), второй – оператором присваивания для функций, а третий является просто знаком равенства в соотношениях и уравнениях (т.е. когда он непосредственно ничего не присваивает). В Maple имеется лишь два таких знака, причем:= – это оператор присваивания и для символов, и для уравнений а = – это знак равенства в уравнениях.

r=DSolve[{eq1,eq2,eq3,f1[1]Šf2[1],f1'[1]Šf2'[1],f2[0]Šf3[0],f2'[0]Šf3'[0],f1[0]+I*f1[Pi/2/Sqrt[epsilon]]==2,f3[0]-I*f3[Pi/2/Sqrt[epsilon]]==0},{f1[x],f2[x],f3[x]},x]

Нельзя не обратить внимания, насколько компактнее, по сравнению с Maple, система Mathematica выдала решение! А ведь это то же самое решение, в чем мы сейчас убедимся. Это является проявлением уже упоминавшейся несколько большей оптимальности расчетов в Mathematica.

Выполняем операторы присваивания:

f1[x]=f1[x]/.r;f2[x]=f2[x]/.r;f3[x]=f3[x]/.r;

(структура операторов получилась несколько странная с точки зрения логики, но в каждой системе иногда надо подстраиваться под общее строение и синтаксис языка).

Вычисляем коэффициент прохождения (напомним, имея в виду (1.10), можно получить амплитуду прошедшей сквозь барьер волны S, просто подставив в и далее вычисляя коэффициент прохождения как квадрат модуля S, умножая S на его комплексно сопряженное (оператор conjugate):

S=f3[x]/.x®0

T=N[S Conjugate[S]];

Проверим, получили ли мы тот же результат, что и раньше, для чего вычислим коэффициент прохождения для энергии V=38:

T/.V®38

С помощью оператора Plot получаем хорошо знакомый график для коэффициента прохождения в зависимости от энергии

Plot[T/.V®Ex,{Ex,0,4*u0}]

Рис. 3.4.Изменение коэффициента прохождения в зависимости от энергии для потенциального барьера, показанного на рисунке 1.4, полученное в Mathematica

Наконец, проверим еще раз соотношение (1.11), для чего вычислим коэффициент отражения (обозначим, как и выше, через ). Функция ft представляет собой сумму коэффициентов прохождения и отражения. Так же как это было сделано в MATLAB, для ознакомления с циклами в Mathematica, произведем проверку не путем построения графика, а составим простейший цикл, выводя значения ft для пяти различных значений аргумента E:

r=(f1[x]-1)/.x®0; r2=N[r Conjugate[r]];

ft[q_]:=N[(T+r2)/.V®q];

nn=5;For[k=1,k<=nn,Print[ft[k*4*u0/nn]];k++]

{1. +0. ä}

.

Дата добавления: 2015-07-20