5. Моделирование движения иглы туннельного микроскопа
РАСЧЕТНЫЕ ФОРМУЛЫ
Туннельный ток, возникающий между иглой и поверхностью металла, вычисляется по формуле:
где Z(x,y) – туннельный зазор, вычисляется по формуле:
(x0,0,z0) – координата описывающая положение иглы. Игла движется вдоль оси OX.
(x,y,z) – точка лежащая на исследуемой поверхности металла.
Формула для вычисления функции j(Z(x,y)) (далее j(Z)):
где U = 0.01 В, уровень Ферми Ef = 5.71 эл/В.
где k – диэлектрическая постоянная (k=1), φ – локальный выход электронов (φ = 4.5).
;
ЛИСТИНГ ПРОГРАММЫ
Program Modeling_of_STM;
Uses graph,crt;
var
gd,gm,z0,StepGraph:integer;
const
cx=70;
cy=100;
Ef=5.71; { Константы }
f0=4.5;
U=0.1;
eU=0.1;
var
x0,y0,z1,z2,zr,zn,Step,p:real;
Ir,I1,I2,Iet:real;
{************************************************************************}
x,y,z:real;
int,in1,in2,in3,in4,in5,in6,in7:real;
{************************************************************************}
function J_current(z,f:extended):extended;
var a,b:extended;
begin
b:=-1.025*z*sqrt(f); { Вычисляем плотность тока }
if b<-500 then b:=-500; { Проверяем допутимые диапозоны }
if b>-10 then b:=-10;
a:=1620*U*Ef*exp(b);
J_current:=a;
end;
{************************************************************************}
function Fi(z:extended):extended;
var a,b,f1,f2,s1,s2,f:extended;
begin
s1:=3/f0; { вычисление потенциала }
a:=23/(z*3*f0+10-2*eU*z);
if a>=1 then a:=0.999;
if a<0.01 then a:=0.01;
s2:=z*(1-a)+s1;
if s2<s1 then s2:=s1+1;
f:=f0-eU/(2*z)*(s1+s2);
f:=f-2.86/(s2-s1);
if f<0 then f:=0.001;
f1:=s2*(z-s1);
f2:=s1*(z-s2);
if f1<0 then f1:=0.01;
if f2<0 then f2:=0.01;
f:=f*ln(f1/f2);
if f<0 then f:=0;
Fi:=f;
end;
{***************************************************************************}
function J_computing(z:real):extended;
var a,f:extended;
begin
f:=Fi(z); { Потенциал }
a:=J_current(z,f); { Плотность тока }
J_computing:=a;
end;
{***************************************************************************}
function Integral_compute(x1,y1,x2,y2,k:integer;popravka:integer;zX:real):extended;
const
q=2000;
var i,j,dx,dy:integer;
buf_integ,buf:real;
begin
dx:=x2-x1;
dy:=y2-y1;
buf_integ:=0;
i:=q;
j:=1;
while i>j do begin
x:=x1+dx*(random(65535)/65535);
y:=y1+dy*(random(65535)/65535);
case k of
1:z:=sqrt(sqr(x-x0)+sqr(y-y0)+sqr(zX-popravka));
2:z:=sqrt(sqr(x-x0)+sqr(y-y0)+sqr(zX-(y-popravka)));
3:if y0>=40 then z:=31
else z:=sqrt(sqr(x-x0)+sqr(popravka-y0)+sqr(-y+zX));
4:z:=sqrt(sqr(x-x0)+sqr(y-y0)+sqr(zX-(popravka-y)));
end;
if z>30 then begin {30-destination}
dec(i);
buf:=0;
end
else buf:=J_computing(z);
buf_integ:=buf_integ+buf;
inc(j);
end;
if i>q/2 then Integral_compute:=buf_integ*dx*dy/i else Integral_compute:=0;
end;
{***************************************************************************}
function Integral(zX:real):extended; {суммарный интеграл}
begin
in1:=Integral_compute(0,0,10,10,1,0,zX);
in2:=Integral_compute(0,10,10,20,2,10,zX);
in3:=Integral_compute(0,20,10,24,1,10,zX);
in4:=Integral_compute(0,24,10,34,4,34,zX);
in5:=Integral_compute(0,34,10,40,1,0,zX);
in6:=Integral_compute(0,0,10,16,3,40,zX);
in7:=Integral_compute(0,40,10,50,1,16,zX);
int:=2*(in1+in2+in3+in4+in5+in6+in7);
Integral:=int;
end;
{***************************************************************************}
procedure put_surface; {рисуем поверхность}
begin
setcolor(15);
line(0+cx,300+cy,100+cx,300+cy);
line(100+cx,300+cy,200+cx,200+cy);
line(200+cx,200+cy,240+cx,200+cy);
line(240+cx,200+cy,340+cx,300+cy);
line(340+cx,300+cy,400+cx,300+cy);
line(400+cx,300+cy,400+cx,140+cy);
line(400+cx,140+cy,501+cx,140+cy);
line(69,99,571,99);
line(69,99,69,419);
SetFillStyle(4,LightBlue);
floodfill(320+cx+1,300+cy+1,15);
setColor(7);
line(69,99,571,99);
line(69,99,69,419);
SetColor(lightred);
end;
{***************************************************************************}
Procedure GraphInterface;
begin
gd:=detect;
Initgraph(gd,gm,'c:\bp\bgi');
end;{GraphInterface}
{***************************************************************************}
Procedure GraphWindow(Px1,Py1,Px2,Py2:integer; WindowName:string{;bmp:pointer});
Label l1,l2;
begin
SetViewPort(Px1,Py1,Px2,Py2,true{false});
SetFillStyle(1,LightGray);
FloodFill(1,1,15);
SetColor(White);
Line(0,0,Px2-Px1,0);
Line(0,0,0,Py2-Py1);
SetColor(DarkGray);
Line(0,Py2-Py1,Px2-Px1,Py2-Py1);
Line(Px2-Px1,0,Px2-Px1,Py2-Py1);
SetFillStyle(1,Blue);
Bar(2,2,Px2-Px1-2,15);
SetColor(White);
SetTextStyle(2,0,4);
OutTextXY(6,2,WindowName);
SetLineStyle(0,0,1);
SetColor(DarkGray);
Line(2,15+5,Px2-Px1-2,15+5);
Line(Px2-Px1-4-10,13,Px2-Px1-4,13);
Line(Px2-Px1-4,3,Px2-Px1-4,3);
SetColor(White);
Line(2,15+6,Px2-Px1-2,15+6);
SetFillStyle(1,LightGray);
Bar(Px2-Px1-4-9,4,Px2-Px1-4,13);
Line(Px2-Px1-4-10,4,Px2-Px1-4,4);
Line(Px2-Px1-4-10,4,Px2-Px1-4-10,13);
SetColor(Black);
SetLineStyle(0,0,1);
Line(Px2-Px1-4-8,6,Px2-Px1-6,12);
Line(Px2-Px1-4-7,6,Px2-Px1-5,12);
Line(Px2-Px1-4-8,12,Px2-Px1-6,6);
Line(Px2-Px1-4-7,12,Px2-Px1-5,6);
end;{GraphWindow}
{***************************************************************************}
Procedure GraphRead(Left_Margin,Upper_Margin:integer;var varRead:integer);
var
ch:char;
str1:string;
int1,i,i1,LM:integer;
function deci(i:integer):longint;
var i1,di:integer;
begin
di:=1;
for i1:=1 to i do
di:=di*10;
di:=di div 10;
deci:=di;
end;
label l1;
begin
varRead:=0;
str1:='';
LM:=Left_Margin;
repeat
l1: ch:=readkey;
if (ch>='0')and(ch<='9')or(ord(ch)=13)then
if (ord(ch)<>13) then
begin
str1:=str1+ch;
outTextXY(LM,Upper_Margin,ch);
LM:=LM+TextWidth(ch);
end
else
begin
if str1=''then
begin
goto l1;
end;
end;
until ord(ch)=13;
for i:=1 to length(str1) do
begin
ch:=str1[i];
val(ch,int1,i1);
varRead:=varRead+int1*deci(i);
end;
end;
{***************************************************************************}
Procedure _Enter_Window(Left_Margin,Upper_Margin,Right_Margin:integer;
CommText:string;var int_var:integer);
begin
SetFillStyle(1,15);
bar(Left_Margin,Upper_Margin,Right_Margin,Upper_Margin+18);
SetFillStyle(1,8);
bar(Left_Margin+1,Upper_Margin+1,Right_Margin,Upper_Margin+18);
SetFillStyle(1,7);
bar(Left_Margin+1,Upper_Margin+1,Right_Margin-1,Upper_Margin+18-1);
SetColor(15);
OutTextXY(Left_Margin+1+5,Upper_Margin+3,CommText);
SetColor(15);
SetFillStyle(1,15);
bar(Left_Margin+TextWidth(CommText)+10,Upper_Margin+2,Right_Margin-2,Upper_Margin+18-2);
SetFillStyle(1,8);
bar(Left_Margin+TextWidth(CommText)+10,Upper_Margin+2,Right_Margin-3,Upper_Margin+18-3);
SetFillStyle(1,7);
bar(Left_Margin+TextWidth(CommText)+11,Upper_Margin+3,Right_Margin-3,Upper_Margin+18-3);
GraphRead(Left_Margin+TextWidth(CommText)+11,Upper_Margin+3,int_var);
end; {_Enter_Window}
{***************************************************************************}
Begin
GraphInterface;
GraphWindow(1,1,640,480,'The Monte-Carlo Method');
SetFillStyle(1,7);
SetFillStyle(1,15);
bar(69,99,572,420);
SetFillStyle(1,0);
bar(68,98,571,419);
SetFillStyle(1,7);
bar(69,99,571,419);
SetFillStyle(1,0);
bar(70,100,570,418);
_Enter_Window(40,480-30,140,'Enter Z0:',z0);
zr:=z0;
_Enter_Window(260,480-30,360,'Enter step:',StepGraph);
step:=stepGraph;
randomize;
Put_surface;
y0:=0;x0:=0;zn:=zr;
putpixel(cx,300+cy-round(zr*10),14);
Iet:=integral(zr);
while y0<=(50-Step) do
begin
if keypressed then break;
y0:=y0+Step; {продвигаемся по поверхности}
Ir:=integral(zr);
if abs(Ir-Iet)>=(Iet/100) then
begin
p:=0.11;
repeat
p:=p+0.11;
if Ir>Iet then begin z1:=zr;z2:=zr+p;end else begin z2:=zr;z1:=zr-p;end;
until (integral(z2)<=Iet)and(integral(z1)>=Iet);
while abs(Ir-Iet)>=(Iet/100) do
begin
zr:=(z2+z1)/2; Ir:=integral(zr); {подбор zr}
if Ir<Iet then z2:=zr else z1:=zr;
if keypressed then break;
end;
end;
line(round((y0-Step)*10)+cx,300+cy-round(zn*10),round(y0*10)+cx,300+cy-round(zr*10));
zn:=zr;
end;
sound(2000);
delay(100);
nosound;
readln;
end.
РЕЗУЛЬТАТЫ РАБОТЫ:
Z0=3Å; Шаг=1Å;
