clear all; close all; clc; format long; %constantes me = 9.10955*10^(-31); eV = 1.6021895*10^(-19); hb = 1.0545896*10^(-34); mp = 1.672651*10^(-27); a = 1*10^(-10); function [p] = rho(E,m,hb) p = sqrt((2*m*E)/hb^2); endfunction function [k] = kha(E,m,hb,V) k = sqrt((2*m*(E-V))/hb^2); endfunction function [a] = alph(E,m,hb,V) a = sqrt((2*m*(V-E))/hb^2); endfunction %rho(3*eV,me,hb) %kha(3*eV,me,hb,2*eV) %alph(2*eV,me,hb,4*eV) function [r] = Ref(E,m,hb,V,a) if(E == V) p = rho(E,m,hb); r = p*p*a*a/4; r = r/(1+r); else p = rho(E,m,hb); k = kha(E,m,hb,V); r = (p^2-k^2)/(2*p*k); r = r*r; r = (r*(sin(k*a))^2)/(1+r*(sin(k*a))^2); endif endfunction function [t] = Tran(E,m,hb,V,a) if(E == V) p = rho(E,m,hb); t = p*p*a*a/4; t = 1/(1+t); else p = rho(E,m,hb); k = kha(E,m,hb,V); t = (p^2-k^2)/(2*p*k); t = t*t; t = 1/(1+t*(sin(k*a))^2); endif endfunction function [r] = Refc(E,m,hb,V,a) im = sqrt(-1); if(E == 0) r = 1; elseif(E == V) p = rho(E,m,hb); r = p*p*a*a/4; r = r/(1+r); else p = rho(E,m,hb); k = alph(E,m,hb,V); r = (p^2-(im*k)^2)/(2*p*(im*k)); r = r*r; r = (r*(im*sinh(k*a))^2)/(1+r*(im*sinh(k*a))^2); endif endfunction function [t] = Tranc(E,m,hb,V,a) im = sqrt(-1); if(E == 0) t = 0; elseif(E == V) p = rho(E,m,hb); t = p*p*a*a/4; t = 1/(1+t); else p = rho(E,m,hb); k = alph(E,m,hb,V); t = (p^2-(im*k)^2)/(2*p*(im*k)); t = t*t; t = 1/(1+t*(im*sinh(k*a))^2); endif endfunction N = 200; V0 = 2*eV; E = linspace(0*V0,10*V0,N); Rc = zeros(N); Tc = zeros(N); for i = 1:N if(E(i) < V0) Rc(i) = Refc(E(i),me,hb,V0,a); %Tc(i) = Tranc(E(i),me,hb,V0,a); else Rc(i) = Ref(E(i),me,hb,V0,a); %Tc(i) = Tran(E(i),me,hb,V0,a); endif Tc(i) = 1 - Rc(i); endfor Ex = E/eV; l1 = length(Ex); l2 = length(Rc); %coefficients R et T figure('name', 'Coefficients R et T, unité en eV:'); plot(Ex,Rc,Ex,Tc) %ouvre un fichier ou le créé fid = fopen('coeffs.txt','w'); fclose(fid); fid = fopen('coeffs.txt','at'); %écrit dans ce fichier fprintf(fid,'%s\t','x1'); fprintf(fid,'%s\t','y1'); fprintf(fid,'%s\t\n','y2'); for i = 1:N fprintf(fid,'%i\t',Ex(i)); fprintf(fid,'%i\t',Rc(i)); fprintf(fid,'%i\t\n',Tc(i)); endfor fclose(fid); % Fonctions de calcul des coefficients function [aa] = AA(E,m,hb,V,a) p = rho(E,m,hb); k = kha(E,m,hb,V); aa = i*exp(-i*p*a); aa = aa * (k*k - p*p)*sin(k*a); den = 2*p*k*cos(k*a) - i*(p*p + k*k)*sin(k*a); aa = aa/den; endfunction function [ba] = BA(E,m,hb,V,a,ind) p = rho(E,m,hb); k = kha(E,m,hb,V); ba = ind*p*(p+ind*k)*exp(-ind*i*a*(k+ind*p)/2); den = 2*p*k*cos(k*a)-i*(p*p+k*k)*sin(k*a); ba = ba/den; endfunction function [ca] = CA(E,m,hb,V,a) p = rho(E,m,hb); k = kha(E,m,hb,V); ca = 2*p*k*exp(-i*p*a); ca = ca * (2*p*k*cos(k*a) + i*(k*k+ p*p)*sin(k*a)); den = 4*p*p*k*k*cos(k*a)*cos(k*a)+(p*p+k*k)^2 * sin(k*a)*sin(k*a); ca = ca/den; endfunction function [] = Proba(E,m,hb,V,a,n,s) %définition des constantes p3 = rho(E,m,hb); k3 = kha(E,m,hb,V); %calcul des coefficients AAmp = AA(E,m,hb,V,a); BAp = BA(E,m,hb,V,a,1); BAm = BA(E,m,hb,V,a,-1); CAp = CA(E,m,hb,V,a); %Définition de l'intervalle X1x = linspace(-5*a,-a/2,n); X2x = linspace(-a/2,a/2,n); X3x = linspace(a/2,5*a,n); Xfx = [X1x X2x X3x]; for i = 1:n im = sqrt(-1); Y1y(i) = AAmp*exp(-1*im*p3*X1x(i)) + exp(im*p3*X1x(i)); Y1y(i) = Y1y(i)*conj(Y1y(i)); Yfy(i) = Y1y(i); endfor for i = 1:n im = sqrt(-1); Y2y(i) = BAm*exp(-1*im*k3*X2x(i)) + BAp*exp(im*k3*X2x(i)); Y2y(i) = Y2y(i)*conj(Y2y(i)); Yfy(i+n) = Y2y(i); endfor for i = 1:n Y3y(i) = CAp * conj(CAp); Yfy(i+2*n) = Y3y(i); endfor %norme de Phi au carré figure('name', s); plot(X1x,Y1y,"b",X2x,Y2y,"r",X3x,Y3y,"b"); %ouvre un fichier ou le créé %fid = fopen('test.txt','w'); %fclose(fid); %fid = fopen('test.txt','at'); %écrit dans ce fichier %fprintf(fid,'%s\t','x1'); %fprintf(fid,'%s\t','y1'); %fprintf(fid,'%s\t','x2'); %fprintf(fid,'%s\t','y2'); %fprintf(fid,'%s\t','x3'); %fprintf(fid,'%s\n','y3'); %for i = 1:n % fprintf(fid,'%i',X1x(i)); % fprintf(fid,'\t%i\t',Y1y(i)); % fprintf(fid,'%i',X2x(i)); % fprintf(fid,'\t%i\t',Y2y(i)); % fprintf(fid,'%i',X3x(i)); % fprintf(fid,'\t%i\n',Y3y(i)); %endfor %fclose(fid); endfunction %cas E=3eV Proba(3*eV,me,hb,2*eV,a,N,"E=3eV, V0 = 2eV"); %cas E=1eV Proba(1*eV,me,hb,2*eV,a,N,"E=1eV, V0 = 2eV"); Proba(1*eV,mp,hb,2*eV,a,N,"E=1eV, V0 = 2eV, proton"); %Proba Tunnel electron et 1eV: function [] = Probatunnel(E,m,hb,V,a,s) disp(strcat('Probabilite de passage par effet tunnel', s, ':')); PtunnelE = CA(E,m,hb,V,a); proba = PtunnelE * conj(PtunnelE); if(proba > 10^(-2)) prc = strcat(mat2str(round(proba*100000)/1000) , '%'); else prc = strcat(mat2str(proba*100) , '%'); endif disp(prc); endfunction %Probabilité de passage par effet tunnel pour l'électron et le proton Probatunnel(1*eV,me,hb,2*eV,a," pour l'electron"); Probatunnel(1*eV,mp,hb,2*eV,a," pour le proton");