function y=rtwell(ERG);
%transmission and reflection for a resonant tunneling barrier
global ENG k1a k2a exp_k1a exp_k2a  k1b k2b exp_k1b exp_k2b V_0 X RR
global AA BB


RR=AA./BB;
k1b =sqrt(ERG);
k2b= sqrt(ERG-V_0);
exp_k1b=exp(i*k1b);
exp_k2b=exp(i*k2b);
k1a =RR .*sqrt(ERG);
k2a= RR .*sqrt(ERG-V_0);
exp_k1a=exp(i*k1a);
exp_k2a=exp(i*k2a);
M= [ -exp_k1a,    1./exp_k2a,   exp_k2a,    0, 0, 0, 0, 0;
k1a .*exp_k1a,  k2a ./exp_k2a,  -k2a.*exp_k2a, 0, 0, 0, 0, 0;
0, -1 ./exp_k2b,  -exp_k2b,    1 ./exp_k1b,   exp_k1b, 0, 0, 0 ;
0, -k2b ./exp_k2b,  k2b .*exp_k2b, k1b ./exp_k1b, -k1b.* exp_k1b, 0, 0, 0 ;
0,0,0,  -exp_k1b, -1 ./exp_k1b, exp_k2b, 1 ./exp_k2b, 0;
0,0,0,  -k1b .*exp_k1b, k1b ./exp_k1b, k2b .*exp_k2b, -k2b ./exp_k2b, 0;
0,0,0,0,0,   -exp_k2a,      -1 ./exp_k2a, exp_k1a;
0,0,0,0,0, -k2a.*exp_k2a, k2a./exp_k2a, k1a.*exp_k1a ];
A=[1 ./exp_k1a, k1a./exp_k1a, 0, 0,0,0,0,0].';
y=inv(M)*A;