clear; clc;
ln = @(x) log(x);
scale = 1e-10; hbar = 1.055*1e-34; me = 9.11*1e-31; eV2J = 1.6*1e-19; J2eV = eV2J^(-1);
ma = 0.067; mb = 0.15;
b = 17; a = 45; c = 100;
dz = 1;
eV = 0;
deV = eV/(b+a+b);
rangeU = [
0
1 - deV*b/2
0 - deV*(b + a/2)
1 - deV*(b + a + b/2)
0 - eV
];
rangeM = [ ma, mb, ma, mb, ma ];
U( 1 : c ) = rangeU(1); m( 1 : c ) = rangeM(1)*me;
U( c+1 : c+b ) = rangeU(2); m( c+1 : c+b ) = rangeM(2)*me;
U( c+b+1 : c+b+a ) = rangeU(3); m( c+b+1 : c+b+a ) = rangeM(3)*me;
U( c+b+a+1 : c+b+a+b ) = rangeU(4); m( c+b+a+1 : c+b+a+b ) = rangeM(4)*me;
U( c+b+a+b+1 : c+b+a+b+c ) = rangeU(5); m( c+b+a+b+1 : c+b+a+b+c ) = rangeM(5)*me;
L = c + b + a + b + c;
z = 0:dz:L;
syms Zi Ui mi Ez t r;
ki = (2*mi*(Ez-Ui)*eV2J).^(0.5)./hbar;
Mi = [
exp( 1i*ki*Zi ), exp(-1i*ki*Zi )
ki./mi.*exp( 1i*ki*Zi ), -ki./mi.*exp(-1i*ki*Zi )
];
invMi = inv(Mi);
productM = [1 0; 0 1];
for index = length(z)-1:-1:2
localProduct = subs(invMi, {mi, Ui, Zi}, {m(index), U(index), z(index)*scale})*subs(Mi, {mi, Ui, Zi}, {m(index-1), U(index-1), z(index)*scale});
productM = productM*localProduct;
end
outWave = [t, 0]'; inWave = [1, r]';
[t, r] = solve( outWave == productM*inWave, t, r);
outAmp = limit(t, Ez, 0)
outAmp = double(outAmp);
outAmp =
limit((((exp(-(conj(((2315911741694958125423557058355*Ez)/118571099379011784113736688648896417641748464297615937576404566024103044751294464)^(1/2))*683982766270....
Error using symengine
DOUBLE cannot convert the input expression into a double array.
Error in sym/double (line 616)
Xstr = mupadmex('symobj::double', S.s, 0);
Error in test (line 54)
outAmp = double(outAmp);
