Моделирование СТМ-профилограмм

Федеральное государственное бюджетное образовательное учреждение
высшего профессионального образования
«Ижевский государственный технический университет»
им. М.Т.Калашникова
Факультет «Информатика и вычислительная техника»
Кафедра «Вычислительная техника»

Пояснительная записка
к курсовой работе по дисциплине
«Моделирование»
на тему «Моделирование СТМ-профилограмм»
Выполнил: студент гр.7-78-3

Принял: д.т.н., профессор
Шелковников Е.Ю.



Ижевск 2013 
Оглавление
1. Задание на курсовую работу 3
2. Введение 4
3. Принцип работы сканирующего туннельного микроскопа 6
4. Метод статических испытаний (Монте-Карло) для вычисления интегралов 9
5. Тестирование генератора случайных чисел 15
6. Листинг программы тестирования генератора случайных чисел 18
7. Формирование туннельного тока между конусообразным острием и подложкой 20
8. Методика построения СТМ-профилограмм нанообъектов для изучения пространственного разрешения туннельного ТПАК 25
9. Листинг программы и получение СТМ-профилограммы 30
10. Заключение 42
11. Список литературы 43

 
Задание на курсовую работу
Вариант №13
1. Выполнить исследование разрешающей способности сканирующего туннельного микроскопа (СТМ) при сканировании исследуемой поверхности острием иглы конусообразной формы. Построить графики распределения туннельного тока вдоль подложки при различных параметрах туннельного промежутка и острия иглы (туннельный зазор Z0 = 5Å, 7Å, 10Å, 15Å; туннельное напряжение UT = 0,01В, 0,1В; локальная работа электронов φ0 = 4.5 эВ; уровень Ферми Ef = 5.71 эВ; угол конуса острия =10, 30, 60, 90).
2. Провести моделирование работы отрицательной обратной связи туннельного микроскопа и построить СТМ-профилограммы для следующей поверхности:

Рисунок 1 – Исходная поверхность
 
Введение
Данная курсовая работа посвящена изучению работы сканирующего туннельного микроскопа (СТМ).
СТМ- программно-аппаратный комплекс, программно-управляемый от ЭВМ через ЦАП и АЦП и состоящий из датчика (бесконтактного считывающего элемента на основе зондирующей иглы, разделенной от исследуемой поверхности туннельным зазором), системы обратной связи (ОС) и пьезопреобразователей с высоковольтными усилителями по осям X,Y,Z для получения измерительной информации о нанорельефе поверхности, локальной работе выхода электронов и спектре электронных состояний образца с атомным разрешением.
Туннельный программно-аппаратный комплекс (ТПАК) открывает новые возможности для изучения наночастиц вследствие своего высокого пространственного разрешения, возможности работы в различных средах (вакууме, воздухе, жидкостях) и применения для решения задач нанотехнологии и неразрушающего контроля. Необходимо отметить, что в наночастицах резко проявляются размерные эффекты физико-химических свойств, которые могут быть использованы в практике различных технологических процессов и при создании материалов с новыми необычными для макроскопических тел свойствами. Новые свойства наночастиц проявляются, когда их размеры становятся соизмеримыми с характерным корреляционным масштабом того или иного физического явления. При этом для каждого типа свойств имеется своя граница «ультрамалости» наночастицы, что позволяет целенаправленно изменять их физико-химические и физико-механические свойства. Таким образом, неразрушающий контроль параметров наночастиц (в первую очередь, размеров) с помощью ТПАК является актуальной проблемой, от успешного решения которой зависит повышение эффективности разработок новых перспективных наноматериалов.
Следует отметить, что ТПАК позволяет производить неразрушающий контроль параметров поверхности с уникальным пространственным разрешением: до нескольких сотых ангстрема по нормали к поверхности и единицы ангстрем – вдоль поверхности. Атомное разрешение ТПАК определяется, в первую очередь, остротой зондирующей иглы и туннельным током с ее острия к поверхности в форме электронного конусного луча.
В ходе выполнения курсовой работы необходимо сгенерировать случайную последовательность чисел и проверить ее на равномерность распределения. Проверка делается с помощью специальных тестов. Случайные последовательности вырабатываются до тех пор, пока одна из них не будет удовлетворять условиям всех тестов. Полученная последовательность равномерно распределенных случайных чисел (РРСЧ) необходима для того, чтобы вычислять определенные интегралы методом Монте-Карло. Интегралы, в свою очередь, необходимы для того, чтобы находить токи, идущие от некоторой поверхности к точке на острие иглы СТМ, плотности этих токов и т.д. Используя вышесказанное, необходимо построить профилограммы движения острия точечной иглы относительно заданной поверхности, а также выполнить исследование разрешающей способности СТМ при сканировании поверхности острием иглы конусообразной формы. 
Принцип работы сканирующего туннельного микроскопа
В основе работы микроскопа лежит явление квантовомеханического туннелирования электрона сквозь барьер, высота которого превышает энергию электрона. Если приблизить два проводника (электрода), разделённые диэлектриком (вакуум, газ или жидкость) на достаточно малое расстояние (∼ 4-20Å) и создать между ними смещающее напряжение UT, то в цепи потечёт туннельный ток IT, пропорциональный приложенному напряжению и экспоненциально зависящий от расстояния между электродами. [1]
Для реализации туннельного микроскопа проводящий образец используется в качестве основного электрода по отношению к металлическому острию, служащему сканирующим электродом. Сканирующий электрод установлен над основным на расстоянии нескольких ангстрем, представляющим собой туннельный барьер. Металлическое остриё служит в качестве источника электронов или, так называемого, холодного катода. Туннельный эффект применяется только для высвобождения электронов из металла контакта в диэлектрик. Под влиянием сильного электрического поля электроны высвобождаются из эмиттирующего острия и ускоряются для формирования изображения. Каждый из электродов может двигаться по трём осям X, Y, Z с помощью трёх пьезопреобразователей (ПП), обеспечивающих грубое перемещение образца и точное перемещение острия. Электронный блок управления подсоединен к электродам (т.е. к образцу и к острию) и к ПП, обеспечивающим малые перемещения электродов, которые также должны быть определёнными и воспроизводимыми. Точное положение острия в трёх координатах становится известным на основании значений напряжений, обеспечивающих соответствующий шаг ПП. Различные конструкции ПП позволяют перемещать остриё зонда на расстояния от долей ангстрема до нескольких десятков микрометров.
Структурная схема СТМ представлена на рис.2, где 1 – образец; 2 – примесные атомы; И – игла; ППX,Y,Z – пьезопреобразователь; СОС – система обратной связи; ЦАП и АЦП – цифро-аналоговые и аналого-цифровые преобразователи; ВУX,Y,Z – высоковольтные усилители.

Рисунок 2–Структурная схема СТМ
Туннельный ток измеряется и поддерживается постоянным с помощью точной настройки ПП по оси Z посредством системы ООС. Строчная развёртка осуществляется пьезоприводом по оси Х, кадровая – по оси Y. Поскольку положение пьезопривода пропорционально пьезонапряжению, напряжения трёх пьезоэлектрических двигателей определяют положение острия на каждой оси. В целом координаты являются декартовыми с тремя ортогональными осями. Анализ информации осуществляется в виде анализа трёхмерного изображения поверхности образца. [3]
СОС функционирует таким образом, чтобы при постоянном туннельном напряжении поддерживать действительный (измеренный) туннельный ток равным заданному. Исследованиями установлено, что при поддерживании системой СОС заданного значения туннельного тока с точностью не хуже 2% величина туннельного зазора остаётся также постоянной с ошибкой, не превышающей 0,1Å. Если производить сканирование остриём поверхности образца, подавая на пьезоприводы X,Y развёртывающие напряжения подобно телевизионной развёртке, то система СОС будет стремиться поддерживать заданным туннельный ток, а значит и постоянный зазор игла-остриё, изменяя напряжение на Z-пьезоприводе в соответствии с рельефом образца.
Таким образом, напряжение на Z-пьезоприводе острия будет линейно отражать текущую высоту рельефа поверхности образца. Для управления процессом сканирования, задания необходимых туннельных тока и напряжения, запоминания и хранения измерительной информации, а также её преобразования, обработки и отображения в требуемом виде в составе СТМ необходима ЭВМ с ЦАП, АЦП и устройствами отображения.
Первичное СТМ-изображение обычно формируется в виде карты распределения измеряемого параметра (высоты рельефа, работы выхода электрона и т.д.) по плоскости сканирования. Строго говоря, регистрируется не рельеф, а электронное состояние поверхности, но в большинстве случаев на однородных по составу образцах эти зависимости подобны. Для корректного использования таких изображений требуется тщательная калибровка пьезоприводов на объектах с известной топографией.
 
Метод статических испытаний (Монте-Карло) для вычисления интегралов
В СТМ реализуется косвенный метод измерений. В результате косвенных измерений измеряются значения физических величин, функционально связанных с другими, являющимися конечной целью измерения. Для численных исследований параметров СТМ вместо физических экспериментальных данных можно воспользоваться экспериментальными данными, полученными в опытах над моделью объекта, т.е. с помощью вычислительного эксперимента. [6]
Задачи, связанные с вычислением многократных интегралов для исследования параметров СТМ, могут быть эффективно решены методом статистических испытаний, так как в этом случае к точности результатов расчёта не предъявляются очень жёсткие требования. Рассмотрим простейший пример для случая одного измерения.
Пусть ξ – непрерывная случайная величина, принимающая свои значения xiв некоторой области Ω на оси OX. Закон распределения задан плотностью вероятностей f(x) в Ω. Рассмотрим задачу об определении вероятности попадания случайной величины ξ в интервал ω с фиксированными границами a и b, содержащийся в Ω. Если обозначить искомую вероятность P(a ≤ ξ < b) = p, то она выражается в виде интеграла:
(4.1)
Этот интеграл можно вычислить по методу статистических испытаний (эксперимент с получением случайных значений случайной величины ξ). Если появившееся при данном испытании значение xi находится внутри интервала ω, данное испытание будем считать удачным. После проведения N испытаний подсчитаем число m удачных испытаний и вычислим частоту p попадания случайной величины ξ в интервал ω:
(4.2)
Располагая частотой p, мы можем приближённо оценить искомую вероятность p на основании закона больших чисел.
Для этого воспользуемся теоремой Бернулли: если событие А имеет вероятность p и если m – число наступления событий А при N независимых испытаниях, то каково бы ни было постоянное ε>0.
(4.3)
При достаточно большом числе испытаний в качестве оценки для интеграла можно взять частоту m/N, т.е.
(4.4)
Моделирование эксперимента включает:
из совокупности случайных чисел с законом распределения f(x) извлекается число xi;
случайное число xi сравнивается с границами a и b интервала ω. Результаты сравнения отмечаются специальным признаком β, равным единице, если выполнено неравенство:
a ≤ xi<b (4.5)
и равным нулю в противном случае;
полученная величина β прибавляется к содержимому «счётчику числа удачных испытаний»;
к содержимому «счётчика количества испытаний» прибавляется единица.
После проведения N испытаний определяется приближённое значение искомой вероятности:
(4.6)
Описанная процедура не требует запоминания всех случайных чисел, извлекаемых в процессе счёта. По ходу вычисления запоминаются только число испытаний N и число удачных испытаний m.
Для того, чтобы метод статистических испытаний можно было считать практически приемлемым, необходимо оценить точность равенства (4.1) и на этом основании определить число испытаний N для вычисления интеграла с достаточной точностью. Представление о точности можно получить, рассматривая p, как случайную величину. Она имеет математическое ожидание:
(4.7)
и дисперсию:
(4.8)
Поэтому средняя квадратичная ошибка равенства (4.4) будет равна:
(4.9)
Видно, что максимум достигается при р = 0,5.
Обсудим вопрос о точности метода более подробно. Равенство (4.4) имеет точность ε с надежностью α, если для неравенства справедливо соотношение:
(4.10)
Свяжем величины ε и α с числом испытаний N. Первую ориентировку в этом вопросе можно получить из неравенства Чебышева, справедливого для любой случайной величины. В наших обозначениях оно имеет вид:
(4.11)

Сопоставляя выражения (4.10) и (4.11), можно принять:
(4.12)
Если подставить вместо выражение из (4.9), получим:
(4.13)
Отсюда:
(4.14)
Формула (4.14), полученная на основании неравенства Чебышева, даёт сильно завышенное значение N. Более точную оценку для N можно получить в том случае, если использовать закон распределения случайной величины p ̅из равенства (4.10). Величина имеет при асимптотически нормальное распределение. На этом основании при больших N равенство (4.10) можно записать в виде:
(4.15)
гдеta– величина критического интервала, которая выбирается из таблиц нормального распределения по заданной надежности α. Сравнивая это соотношение с (27) мы видим, что или
(4.16)
Отсюда можно определить значение N:
(4.17)
Для того, чтобы наглядно представить порядок величины N, обеспечивающей заданную точность, положим, как и прежде, α=0,95 и вычислим N для различных значений p и ε. Результаты расчёта даны в табл.1. При увеличении требований к точности расчёта (уменьшение ε) существенно увеличивается необходимое число испытаний N. Это обстоятельство является одним из серьёзных ограничений применимости метода статистических испытаний. Метод целесообразно применять для решения тех задач, для которых требования к точности не являются слишком жёсткими.
Таблица 1 - Значения N при α=0,95 и различных p и ε
P 
0.05 0.01 0.005 0.001
0.1 0.9 140 3600 14000 360000
0.2 0.8 250 6200 25000 620000
0.3 0.7 330 8400 33000 840000
0.4 0.6 380 9400 38000 940000
0.5 390 9800 39000 980000

Рассмотрим простой вариант метода Монте-Карло для вычисления одно-мерных интегралов. Выберем N случайных чисел 0≤Ui≤1 и осуществим пере-ход к отрезку [a,b] по формуле ui = a+Ui(b-a).
Рассмотрим величину:
(4.18)
которая представляет собой среднее арифметическое значений функции f, умноженное на b-a. Такой способ вычисления интеграла можно интерпретировать как статистический вариант метода прямоугольников, когда в качестве узла берется случайное число, равномерно распределенное на интервале интегрирования [a,b]. Вследствие случайности узла погрешность интеграла также будет носить случайный характер. Погрешность будет уменьшаться с ростом числа испытаний N по формуле ε ≈ N-1/2.
Для кратных интегралов можно записать:
(4.19)
ГдеVc - k-мерный объем области интегрирования.
 
Тестирование генератора случайных чисел
Модели с вероятностными процессами обычно строят на основе генераторов РРСЧ, которые определяют качество систем моделирования. Однако полученные программным путем числа являются псевдослучайными, т.к. задав конкретное число можно повторить всю последовательность. Поэтому самой важной характеристикой качества генерации СЧ является интервал апериодичности – длина последовательности неповторяющихся случайных чисел. Чем этот интервал больше, тем выше качество генератора. Необходима тщательная проверка генераторов с помощью специальных количественных критериев – тестов. Каждый тест применяется к последовательности статистических испытаний <Un>=U0, U1… 0, которые должны быть распределены равномерно между 0 и 1. [7]
Некоторые тесты предназначены для проверки целочисленных последовательностей, в них проверке подвергается вспомогательная последовательность <Yn>=Y0, Y1…, которая определяется при помощи Yn=целая часть(d*Un), входящие в эту последовательность целые числа распределены равномерно между 0 и (d-1), где d – любое целое число.
Приведем некоторые виды тестов РРСЧ.
Тест проверки частот.
Пусть мы имеем последовательность целых случайных <Yn> чисел в промежутке [0,d-1]. Для каждого целого числа от 0 до d-1 подсчитывается их количество в последовательности <Yn>. Подсчет ведут в массиве Yd размерностью [0, d-1]. Затем к массиву Yd применяют критерий хи-квадрат со степенями свободы v = d-1 и вероятностью p=1/d.
Тест проверки комбинаций (покер-тест).
В этом тесте величина d берется небольшой и поэтому диапазон значений элементов последовательности <Yn>=d*<Un> так же невелик. В классическом "покер-тесте" рассматривается N групп из 5 следующих друг за другом целых чисел. Выделяется 7 типов комбинаций отличающихся различным содержанием цифр:
1) abcde; 2) aabcd; 3) aabbcd; 4) aaabc; 5) aaabb; 6) aaaab; 7) aaaaa.
Затем применяют критерий хи-квадарат.
Тест серий.
В тесте серий проверяется равномерность и независимость пар следующих друг за другом случайных чисел. Для этого подсчитывается сколько раз встретилась каждая пара (Yrj,Yrj+1)=(q,r) при 0≤j<n. Величины q и r могут принимать любые значения от 0 до d. Затем применяется критерий хи-квадрат с числом степеней свободы v=d2, и с вероятностью p=1/d2.
Тест проверки интервалов.
В этом тесте проверяется длина интервалов между появлениями значений Uj принадлежащих некоторому заданному отрезку. Если α и β это два действительных числа, причем 0≤α<β≤1, то подсчитываются длины последовательностей Uj, Uj+1, … , Uj+r в которых только Uj+r лежат между α и β. Такая последовательность из (r+1) чисел определяет интервал r. Пусть в массив count[0..t-1] содержит число интервалов с длинами 0,1,…,(t-1) и count[t] содержит число интервалов длиной >(t-1). Затем к массиву count применяется критерий хи-квадрат со степенями свободы v=t и вероятностями p0=p, p1=p(1-p) p2=p(1-p)2,… pt=p(1-p)t.
Тест коллекционера.
Пусть нам дана последовательность целых случайных чисел:
Y0,Y1,… , где 0≤Yj<d.
Далее определяются длины сегментов Yj+1,Yj+2,…,Yj+r, необходимых для того, чтобы собрать "полный набор" целых чисел от 0 до d-1. Затем длины сегментов обрабатываются с помощью критерия хи-квадрат с v=t-d+1 степенями свободы, после этого соответствующие вероятности равны:
(30)
где d£r<t. В фигурных скобках указаны числа Стирлинга II-го рода. Их можно вычислить по формуле:
(31)
Тест проверки перестановок.
Разделим исходную последовательность случайных чисел на n групп по t элементам в каждой: (Ujt, Ujt+1, … , Ujt+t-1), 0≤j<n.
В каждой группе возможно всего t! вариантов относительного расположения чисел. Далее подсчитывается сколько раз встречается каждое конкретное относительное расположение. После этого применяется критерий хи-квадрат с числом степеней свободы v=t! и p=1/t!.
Для тестирования последовательностей случайных чисел на равномерность будем использовать тесты 1-4, а именно тест проверки частот, покер-тест, тест серий и тест проверки интервалов. Генерацию и проверку последовательности чисел будем проводить до тех пор, пока очередная случайная последовательность не будет удовлетворять всем тестам. В результате, удачная последовательность чисел будет записана во внешний файл.
 
Листинг программы тестирования генератора случайных чисел

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
int main(int argc, const char * argv[]){
int result = 0;
int n=5000,a[n],i,b[n/10],t=0,r=0;
setlocale(0,"");
while (result == 0)
{ result = 1;
float aa[7]={2.088,3.325,5.899,8.343,11.39,16.92,21.67};
srand(time(NULL));
for(i=0;i<n;i++)
a[i]=rand()%10;
while(t<10){
for(i=0;i<n/2;i++){
if(a[i]==t)r++;
}b[t]=r;r=0;
t++;}
int np;float V=0,s=0;
np=n/20;
for(i=0;i<10;i++){
V=((b[i]-np)*(b[i]-np));
V=V/np;
s=s+V; }
r=0;
printf("проверкачастотV=%2.2f",s);
printf("
");
for(i=0;i<7;i++)
if (aa[i]<s)r++;
if(r==0)
printf("плохойр-т 0<s<1 ");
if(r==7)
printf("плохойр-т 99<s<100");
if(r==6)
printf("подозр-йр-т 95<s<99");
if(r==1)
printf("подозр-йр-т 1<s<5");
if(r==2)
printf("почтиподозр-йр-т 5<s<25 ");
if(r==5)
printf("почтиподозр-йр-т 75<s<95");
if(r==3)
printf("хорошийр-т 25<s<50 ");
if(r==4)
printf("хорошийр-т 50<s<75 ");
if ((r<3)||(r>4))
result *= 0;
s=0;


///#2
printf("
");
float aa1[7]={0.8721,1.635,3.455,5.348,7.841,12.59,16.81};
int z[5]={0,0,0,0,0},c[5],j,k=0,x=1,d,q;
for(i=0;i<n/50;i++){
for(j=0;j<5;j++)
c[j]=a[5*i+j];d=5;
while(k<5){ j=k+1;
while(j<d){
if(c[k]==c[j]){ x++;
for(q=j;q<d-1;q++) c[q]=c[q+1];
d--;
j--;}
j++;}
k++; }
k=0;
if(x==1)z[0]=z[0]+1;
if(x==2)z[1]=z[1]+1;
if(x==3)z[2]=z[2]+1;
if(x==4)z[3]=z[3]+1;
if(x==5)z[4]=z[4]+1;
x=1; }
float zz[5]={3,5,1.5,0.12,0.001};
for(i=0;i<4;i++){np=zz[i]*10;
V=((z[i]-np)*(z[i]-np));
V=V/np;
s=s+V; }
printf("
");
printf("покертест V=%2.2f",s);
printf("
");r=0;
for(i=0;i<7;i++){
if (aa1[i]<s)r++; }
if(r==0)
printf("плохойр-т 0<s<1 ");
if(r>=7)
printf("плохойр-т 99<s<100");
if(r==6)
printf("подозр-йр-т 95<s<99");
if(r==1)
printf("подозр-йр-т 1<s<5");
if(r==2)
printf("почтиподозр-йр-т 5<s<25 ");
if(r==5)
printf("почтиподозр-йр-т 75<s<95");
if(r==3)
printf("хорошийр-т 25<s<50 ");
if(r==4)
printf("хорошийр-т 50<s<75 ");
if ((r<3)||(r>4))
result *= 0;
s=0;
//№3
d=10;int g=0,l=0,bb[d*d];j=0;
for(i=0;i<d*d;i++){
for(k=0;k<n/2;k++){
if((a[2*k]==g) && (a[2*k+1]==j))
l++;}
bb[i]=l;l=0;j++;
if(j==d){g++;j=0;}
} printf("
");
float npp;
npp=25;
for(i=0;i<d*d;i++){
V=((bb[i]-npp)*(bb[i]-npp));
V=V/npp;
s=s+V;}
printf("
");
printf("тестсерий V=%2.2f",s);
printf("
");r=0;
float aa2[7]={69.23,77.046,89.185,99.341,110.189,125.458,136.971};
for(i=0;i<7;i++){
if (aa2[i]<s)r++; }
if(r==0)
printf("плохойр-т 0<s<1 ");
if(r>=7)
printf("плохойр-т 99<s<100");
if(r==6)
printf("подозр-йр-т 95<s<99");
if(r==1)
printf("подозр-йр-т 1<s<5");
if(r==2)
printf("почтиподозр-йр-т 5<s<25 ");
if(r==5)
printf("почтиподозр-йр-т 75<s<95");
if(r==3)
printf("хорошийр-т 25<s<50 ");
if(r==4)
printf("хорошийр-т 50<s<75 ");
if ((r<3)||(r>4))
result *= 0;
s=0;

//#4
t=14;
int nn=400,kk[t+1];
for(i=0;i<=t;i++)
kk[i]=0;
i=0;j=0;r=0;
while(j<nn){
if(a[i]!=5)
r++;
if ((a[i]==5)||(r==t)){
kk[r]=kk[r]+1;
j++;
r=0; }
i++; }
printf("
");
float p=0.1,pp,p1;
for(i=0;i<t;i++){
p1=powf(1-p,i);
pp=p*p1;
npp=nn*pp;
V=((kk[i]-npp)*(kk[i]-npp));
V=V/(npp);
s=s+V; }
p1=powf(1-p,15);
npp=nn*p1;
V=((kk[t]-npp)*(kk[t]-npp));
V=V/(npp);
s=s+V;
printf("
");
printf("тестинтерваловV=%2.2f
",s);
float aa3[7]={5.229,7.261,11.04,14.34,18.25,25,30.58};
for(i=0;i<7;i++){
if (aa3[i]<s)r++; }
if(r==0)
printf("плохойр-т 0<s<1 ");
if(r>=7)
printf("плохойр-т 99<s<100");
if(r==6)
printf("подозр-йр-т 95<s<99");
if(r==1)
printf("подозр-йр-т 1<s<5");
if(r==2)
printf("почтиподозр-йр-т 5<s<25 ");
if(r==5)
printf("почтиподозр-йр-т 75<s<95");
if(r==3)
printf("хорошийр-т 25<s<50 ");
if(r==4)
printf("хорошийр-т 50<s<75 ");
if ((r<3)||(r>4))
result *= 0;}
char *FN = "good_test.txt";
FILE *file;
file = fopen(FN, "w+");
int ii = 0;
while(ii < n){
fprintf(file, "%d
", a[ii]);
ii++;}
fclose;
system("pause");
return 0;

 
Формирование туннельного тока между конусообразным острием и подложкой
ПР СТМ определяется размерами области пространства, в котором в основном протекают токи между остриём и подложкой. С целью упрощения задачи определения туннельного тока в СТМ рассмотрим её в приближении условно изолированных атомов, когда на поверхности острия (или подложки) расположены эмитирующие точки – независимые точечные источники тока с плотностью J, экспоненциально зависящей от расстояния между рассматриваемыми точками острия и подложки и отражающей величины как нормальной к эмитирующей поверхности, так и тангенциальной компонент начальных скоростей туннелирующих электронов. При этом размер источника тока на поверхности сравним с диаметром атома, а общий ток IT, текущий через туннельный промежуток, является суммой всех элементарных токов J. Следует отметить, что ВКБ-приближение, справедливое для поля плоского конденсатора, т.е. для перпендикулярной к подложке линии плотности тока, дополняется в данной модели другими линиями плотности тока с целью учёта ухудшения ПР за счёт тангенциальной составляющей начальных скоростей туннелирующих электронов.[2]
Проведём анализ туннельного тока для острия иглы в форме кругового конуса, т.к. такая форма достаточно точно аппроксимирует рабочие микровыступы острия при ориентации (111) монокристалла, определяющие атомное разрешение СТМ.
Рассмотрим формирование туннельного тока с поверхности подложки ко всей поверхности острия иглы. Выберем на плоскости т.А с координатами (х; 0) и т.В на поверхности конуса с радиусом r и углом  между радиусом r и фронтальной плоскостью (рисунок 2a), где  – полуугол раствора конуса острия; Z0 – расстояние между остриём и подложкой; Xp – радиус окружности, полученной пересечением конической поверхности (порождаемой движением отрезка образующей конуса острия между его вершиной и подложкой) с поверхностью подложки; Hк – высота т.В, b – расстояние между т.А и проекцией т.В; Хm, Rm – соответственно, предельные значения х и r. Для определения тока , протекающего с подложки на всю поверхность конуса, выберем в окрестности т.А элементарную площадку dxdk = xddxна поверхности подложки и в окрестности т.В элементарную площадку dl1dl2 = (r/sin)ddr на боковой поверхности конуса (рисунок 2б). Для удобства расчётов представим ток в виде двух составляющих: тока , протекающего с поверхности подложки, ограниченной окружностью с радиусом Хp, и тока , протекающего с остальной части подложки. Плотность тока ( ) на поверхности подложки для первой (второй) составляющей численно равна току, протекающему через единичную площадку в окрестностях т.А с координатой xXp(xXp) ко всей поверхности острия:
; , (7.1)
где p = arcsin(Xp/x)+/2 – предельное значение угла для x>Xp. Из рис.3.3б определим расстояние Zпк между точками А и В:
(7.2)
С учётом (8.1) запишем выражение для общего тока :
(7.3)
Рассмотрим формирование туннельного тока с боковой поверхности конуса ко всей поверхности подложки. Выберем на поверхности конуса т.А с радиусом r и на поверхности подложки т.В с координатами (x; 0) согласно рисунку 3а, где  – полуугол раствора конуса острия; Z0 – расстояние между остриём и подложкой; Xp – радиус окружности, полученной пересечением конической поверхности (порождаемой движением отрезка образующей конуса острия между его т.А и подложкой) с поверхностью подложки; Hк – высота т.А; Хm, Rm – соответственно, предельные значения х и r. Для определения тока , протекающего из т.А на всю поверхность подложки выберем в окрестности т.В элементарную площадку dxdk = xddx(рисунок 3б). Плотность тока на поверхности острия численно равна току, протекающему через единичную площадку в окрестностях т.А ко всей поверхности подложки. Представим суммарный ток из т.А в виде двух составляющих, физической границей которых является коническая поверхность, порождаемая вращением отрезка образующей конуса острия между т. А и подложкой вокруг нормали из т.А к плоскости подложки. Плотности тока для обеих составляющих запишутся следующим образом:
; (7.4)
Из рис.3.4 определим расстояние Zкп между точками А и В:
(7.5)
Плотность тока Jкп между точками А и В определяется из (7.4). Для определения всего тока с конуса на подложку выберем в окрестности т.А элементарную площадку на боковой поверхности конуса dl1dl2 = (r/sin)ddr(рисунок 3а). С учетом (7.3) запишем выражение для общего тока :
(7.6)
 

Рисунок 3 -Формирование туннельного тока с плоской подложки на конусообразное остриё

Рисунок 4 - Формирование туннельного тока с поверхности конусообразного острия к плоской подложке
 
Методика построения СТМ-профилограмм нанообъектов для изучения пространственного разрешения туннельного ТПАК
Недостатком, ограничивающим потенциальные возможности ТПАК, является то, что электроны туннелируют к поверхности расходящимся пучком, а это приводит к значительным искажениям СТМ-изображения и ухудшению ПР ТПАК. Для исследования возможностей улучшения ПР, как одной из важнейших метрологических характеристик нанопрофилометра, необходима информация о механизме формирования СТМ-изображения.[5]
На основании проведённых численных исследований электронного пучка острия иглы создана методика построения СТМ-изображения при движении острия иглы известных формы и размеров над исследуемой поверхностью, суть которой заключается в следующем. Выбираются геометрические параметры и материал иглы, геометрические параметры и материал образца, рабочие параметры сканирования. Программно моделируется процесс работы ТПАК, при котором туннельные токи рассчитываются с учётом геометрии этих поверхности, а также с использованием обобщенной формулы Симмонса для потенциального барьера произвольной формы и метода статистических испытаний.
Алгоритм формирования СТМ-изображения с учетом полученных выражений для туннельного тока представлен на рис.5. Перед началом сканирования устанавливаются необходимые параметры: величина работы выхода электронов материалов иглы и образца, туннельный зазорZ0, угол конуса α. Затем вычисляется эталонный туннельный токIТэт, который виртуальная СОС ТПАК поддерживает постоянным.. Вычисление IТэт проводится для заданного туннельного зазора Z0 и плоской поверхности. Компьютерное сканирование начинается с установки ЗО в начальную точку (x0,y0) растра исходной поверхности, СТМ-изображение которой необходимо получить. Далее находится плотность тока в каждой точке растра (xi,yi), определяется «видимость» ЗО в данной точке и вычисляется туннельный ток IТ(xi,yi) = JТ(xi,yi) •S(xi,yi) •sinγ, если точка «видима» (где S(xi,yi) = (xi+1 – xi-1)(yi+1 – yi-1)/4 – площадь поверхности в окрестностях данной точки (xi,yi);γ – угол между нормалью в точке(xi,yi) и отрезком прямой, являющимся минимальным расстоянием между точкой(xi,yi) и ЗО). Затем вычисляется суммарный ток со сканируемой поверхности на ЗО путем суммирования всех найденных IT(xi,yi) и сравнивается с эталонным током IТэт. Если токи IT(xi,yi) и IТэт не равны, то выполняется коррекция туннельного зазора следующим образом. Если IT(xi,yi)>IТэт, то туннельный зазор Z0 увеличивается (при IT(xi,yi)<IТэт, Z0 уменьшается) на определенную величину и повторно вычисляется IT(xi,yi). Далее с использованием линейной интерполяции уточняется значение туннельного зазора Z0 и снова вычисляется IT(xi,yi). Уточнение выполняется до тех пор, пока IT(xi,yi) не станет равен IТэт с заданной погрешностью δT. После этого в формируемый массив высот СТМ-изображения заносится сумма текущей высоты исходной поверхности и текущего значения туннельного зазора Z0. Затем центр ЗО перемещается в следующую точку сканирования (пока не будет отсканирована вся поверхность) и процедура «виртуальной» стабилизации туннельного тока повторяется. Таким образом, в результате компьютерного сканирования исходной поверхности формируется новый массив высот исследуемой поверхности.
Иглу нужно поддерживать на такой высоте, чтобы туннельный ток оставался постоянным. Для этого нужно находить такую новую высоту, чтобы выполнялось равенство:
 
Рисунок 5 -Алгоритм компьютерного сканирования поверхности с вычислением туннельного тока 
(8.1)
Это уравнение проще всего решить методом секущих или модифицированным методом Ньютона последовательными итерациями по формуле:
(8.2)
Как видно, что перед тем, как приступить к последовательным итерациям, нам нужно знать два предыдущих приближения решения уравнения.
Пусть, и – первое и второе приближение соответственно. Принимая во внимание итерационную формулу метода секущих можно записать:
(8.3)
где i=2,3,4,….
Выше было сказано, что для метода секущих требуется знать два приближения z1 и z2. При этом очень важно, чтобы эти числа были близко расположены к решению уравнения (8.1), в противном случае алгоритм будет работать неустойчиво, возможно, бесконечное зацикливание, или будет найдена совсем другая высота z, которая тоже является решением уравнения (8.1). Чтобы избежать этих неприятностей, нужно выбирать малый шаг движения иглы в направлении оси OX, чтобы z1 близко находилось к решению уравнения (8.1). z2 следует выбирать, как
z2 = z1 + δ ,еслиI1>Iэт,
z2 = z1 – δ, еслиI1<Iэт (8.4)
 
У этого метода имеются следующие недостатки:
малый шаг вдоль оси OX означает большое количество вычислений, а следовательно проигрыш в скорости вычислений;
требование малого шага – является необходимым, но недостаточным условием сходимости метода секущих т.к. по мере продвижения иглы вдоль оси OX всеже возможны большие скачки значений туннельного тока из-за неравномерности поверхности. Это может привести к затягиванию итерационного процесса, что неизбежно влечет потерю производительности.
 
Листинг программы и получение СТМ-профилограммы

unit CTM;
interface

uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, StdCtrls, ExtCtrls, ComCtrls;

const N = 10000;
U = 0.01;
Ef = 5.71;
k = 1;
Np = 385;
fi0 = 4.5;
type
Trndm = array[1..2]of real;
Tprfgrm = array[0..Np]of TPoint;
TForm1 = class(TForm)
Button1: TButton;
Image1: TImage;
Label1: TLabel;
Memo1: TMemo;
Label2: TLabel;
Label3: TLabel;
Label4: TLabel;
ComboBox1: TComboBox;
procedure Button1Click(Sender: TObject);
procedure FormCreate(Sender: TObject);
private
procedure ReadFromFile;
function jz(x:real; y:real; z:real):real;
function Integral(xn,xm,yn,ym,z:real):real;
function GetJ(xp:real; zp:real):real;
function GetPrfgrm: Tprfgrm;
procedure InitPov;
public
{ Public declarations }
end;

var
form1: TForm1;
z0 : real;
Z: array[1..2]of real;
J: array[1..2]of real;
pov: array[1..12]of TPoint;
prfgrm: Tprfgrm;
M:array[1..n]of Trndm;
koef: integer;
xOy: TPoint;

implementation
{$R *.dfm}

function TForm1.jz(x:real; y:real; z:real):real;
var S1,S2,fi,Zxy:real;
begin
Zxy:=sqrt(sqr(z)+sqr(x)+sqr(y));
S1:=3/(k*fi0);
S2:=Zxy*(1-23/(3*fi0*k*Zxy+10-2*U*k*Zxy))+S1;
fi:=fi0-U*(S1+S2)/(2*Zxy)-2.86/(k*(S2-S1))*ln(S2*(Zxy-S1)/(S1*(Zxy-S2)));
result:=1620*U*Ef*exp((-1.025)*Zxy*sqrt(fi));
end;
function TForm1.Integral(xn,xm,yn,ym,z:real):real;
var a:array[1..2]of real;
b:array[1..2]of real;
sum:real;
V:real;
i,j,k:integer;
x:array[1..2]of real;
p:boolean;
begin
b[1]:=xm; a[1]:=xn;
b[2]:=ym; a[2]:=yn;
V:=1; sum:=0;
k:=0;
for i:=1 to 2 do
V:=V*(b[i]-a[i]);
for i:=1 to N do
begin
for j:=1 to 2 do
x[j]:=a[j]+(b[j]-a[j])*M[i,j];
p:=true;
for j:=1 to 2 do
if(x[j]<a[j])and(x[j]>b[j]) then
p:=false;
if p then
sum:=sum+jz(x[1],x[2],z);
end;
sum:=sum*V/N;
result:=sum;
end;

function TForm1.GetJ(xp:real; zp:real):real;
const Max = 10;
var xn,xm,yn,ym,z,J:real;
vidim: array[1..11]of boolean;
i:integer;
begin
yn:=-Max; ym:=Max;
J:=0;
for i:=1 to 11 do begin
vidim[i]:=false;
if ((abs(pov[i].x-xp)<=Max)and((abs(pov[i].y)-zp)<=Max))or
((abs(pov[i+1].x-xp)<=Max)and((abs(pov[i+1].y)-zp)<=Max))then
vidim[i]:=true;
end;
if vidim[1] then
if zp>0 then
begin
xn:=pov[1].x-xp; if xn>Max then xn:=Max;
xm:=pov[2].x-xp; if xm>Max then xm:=Max;
z:=zp;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[2] then
if zp>xp-pov[2].x then
begin
z:=(zp-(xp-pov[2].x))*cos(PI/4);
xn:=(pov[2].x-(xp+z*cos(PI/4)))/cos(PI/4);
xm:=(pov[3].x-(xp+z*cos(PI/4)))/cos(PI/4);
if xn>Max then xn:=Max;
if xm>Max then xm:=Max;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[3] then
if zp>pov[3].y then
begin
xn:=pov[3].x-xp; if xn>Max then xn:=Max;
xm:=pov[4].x-xp; if xm>Max then xm:=Max;
z:=zp-pov[3].y;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[4] then
if zp>xp-pov[4].x then
begin
z:=(zp-(xp-pov[4].x))*cos(PI/4);
xn:=(pov[4].x-(xp+z*cos(PI/4)))/cos(PI/4);
xm:=(pov[5].x-(xp+z*cos(PI/4)))/cos(PI/4);
if xn>Max then xn:=Max;
if xm>Max then xm:=Max;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[5] then
if zp>pov[5].y then
begin
xn:=pov[5].x-xp; if xn>Max then xn:=Max;
xm:=pov[6].x-xp; if xm>Max then xm:=Max;
z:=zp-pov[5].y;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[6] then
if xp<pov[6].x then
begin
xn:=pov[6].y-zp; if xn>Max then xn:=Max;
xm:=pov[7].y-zp; if xm>Max then xm:=Max;
z:=pov[6].x-xp;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[7] then
if zp>pov[7].y then
begin
xn:=pov[7].x-xp; if xn>Max then xn:=Max;
xm:=pov[8].x-xp; if xm>Max then xm:=Max;
z:=zp-pov[7].y;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[8] then
if xp<pov[8].x then
begin
xn:=pov[8].y-zp; if xn>Max then xn:=Max;
xm:=pov[9].y-zp; if xm>Max then xm:=Max;
z:=pov[8].x-xp;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[9] then
if zp>pov[9].y then
begin
xn:=pov[9].x-xp; if xn>Max then xn:=Max;
xm:=pov[10].x-xp; if xm>Max then xm:=Max;
z:=zp-pov[9].y;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[10] then
if xp<pov[10].x then
begin
xn:=pov[10].y-zp; if xn>Max then xn:=Max;
xm:=pov[11].y-zp; if xm>Max then xm:=Max;
z:=pov[10].x-xp;
J:=J+Integral(xn,xm,yn,ym,z);
end;
if vidim[11] then
if zp>pov[11].y then
begin
xn:=pov[11].x-xp; if xn>Max then xn:=Max;
xm:=pov[12].x-xp; if xm>Max then xm:=Max;
z:=zp-pov[11].y;
J:=J+Integral(xn,xm,yn,ym,z);
end;
result:=J;
end;

function TForm1.GetPrfgrm: Tprfgrm;
const e = 0.05;
step = 0.1;
type Point=record
x:real;
y:real;
end;
var J0:real;
z:array[1..3]of real;
J:array[1..2]of real;
prfgrmr:array[0..Np]of Point;
prfgrm: Tprfgrm;
i,k:integer;
process:real;
begin
prfgrmr[0].x:=0; prfgrmr[0].y:=z0;
J0:=0;
for i:=1 to 10 do
J0:=J0+GetJ(0,z0);
J0:=J0/10;
for i:=1 to Np do
begin
prfgrmr[i].x:=prfgrmr[i-1].x+0.1;
z[3]:=prfgrmr[i-1].y;
J[2]:=GetJ(prfgrmr[i].x,z[3]);
if (J[2]<(J0-e*J0))or(J[2]>(J0+e*J0))then
begin
z[2]:=z[3]+step*((J[2]-J0)/abs(J[2]-J0)); //step снужнымзнаком
J[1]:=GetJ(prfgrmr[i].x,z[2]);
if (J[1]<(J0-e*J0))or(J[1]>(J0+e*J0))then
begin
k:=0;
repeat
z[1]:=z[2]-(J[1]-J0)*(z[2]-z[3])/(J[1]-J[2]);
J[2]:=J[1];
J[1]:=GetJ(prfgrmr[i].x,z[1]);
z[3]:=z[2];
z[2]:=z[1];
inc(k);
if k>10000 then break;
until (J[1]>=(J0-e*J0))and(J[1]<=(J0+e*J0));
prfgrmr[i].y:=z[1];
if k>10000 then prfgrmr[i].y:=prfgrmr[i-1].y;
end
else prfgrmr[i].y:=z[2];
end
else prfgrmr[i].y:=z[3];
process:=100*i/Np;
form1.memo1.Lines.Add(floattostr(trunc(process))+ %);
end;

for i:=0 to Np do
begin
prfgrm[i].x:=round(xOy.x+koef*prfgrmr[i].x);
prfgrm[i].y:=round(xOy.y-koef*prfgrmr[i].y);
end; {}
Result:=prfgrm;
end;

procedure TForm1.InitPov;
begin
pov[1].x:=-10; pov[1].y:=0; pov[2].x:=10; pov[2].y:=0;
pov[3].x:=14; pov[3].y:=4; pov[4].x:=17; pov[4].y:=4;
pov[5].x:=20; pov[5].y:=7; pov[6].x:=25; pov[6].y:=7;
pov[7].x:=25; pov[7].y:=0; pov[8].x:=27; pov[8].y:=0;
pov[9].x:=27; pov[9].y:=4; pov[10].x:=29; pov[10].y:=4;
pov[11].x:=29; pov[11].y:=8; pov[12].x:=39; pov[12].y:=8;
form1.Memo1.Text:=;
end;

procedure TForm1.Button1Click(Sender: TObject);
var i:integer;
SaveF: TextFile;
s: string;
begin
z0:=StrToFloat(combobox1.Text);
if z0<5 then z0:=5;

//---------------Профилограмма---------------
With Image1.Canvas do
begin
Pen.Color := RGB(0,0,0);
Pen.Width := 3;
Pen.Style := psSolid;
prfgrm:=GetPrfgrm;
Polyline(prfgrm);
end;
end;

procedure TForm1.FormCreate(Sender: TObject);
var i,j:integer;
begin
InitPov;
for i:=1 to N do
for j:= 1 to 2 do M[i,j]:=random; {}

xOy.x:=15; xOy.Y:=Image1.Height-15;
koef:=trunc((Image1.Width-xOy.x)/pov[12].x-5);

//---------------Оси---------------
With Image1.Canvas do
begin
Pen.Color := RGB(0,0,0);
Pen.Width := 1;
MoveTo(xOy.x,0);
LineTo(xOy.x-3,15);LineTo(xOy.x+3,15);LineTo(xOy.x,0);
LineTo(xOy.X,xOy.y);
LineTo(Image1.Width-150,xOy.y);
LineTo(Image1.Width-165,xOy.Y-3);
LineTo(Image1.Width-165,xOy.Y+3);
LineTo(Image1.Width-150,xOy.y);
//---------------Подписьосей---------------
for i:=1 to 40 do begin
MoveTo(xOy.x+i*koef,xOy.y-2);LineTo(xOy.x+i*koef,xOy.y+2);
TextOut( xOy.x+i*koef-4, xOy.y+2, IntToStr(i));
end;
for i:=1 to 35 do begin
MoveTo(xOy.x-2,xOy.y-i*koef);LineTo(xOy.x+2,xOy.y-i*koef);
TextOut( xOy.x-15, xOy.y-i*koef-4, IntToStr(i));
end;
//-------------------------------------------

//---------------Сетка---------------
Pen.Color := RGB(150,150,150);
Pen.Style := psDot;
for i:=1 to 40 do begin
MoveTo(xOy.x+i*koef,xOy.y);
LineTo(xOy.x+i*koef,0);
end;
for i:=1 to 35 do begin
MoveTo(xOy.x,xOy.y-i*koef);
LineTo(Image1.Width,xOy.y-i*koef);
end;

//---------------Поверхность---------------
Pen.Color := RGB(250,0,0);
Pen.Width := 2;
Pen.Style := psSolid;
MoveTo(xOy.x+koef*pov[1].x,xOy.y-koef*pov[1].y);
for i := 2 to 12 do
begin
LineTo(xOy.x+koef*pov[i].x,xOy.y-koef*pov[i].y);
end;
end;
end;

procedure TForm1.ReadFromFile;
var
ch: real;
i,j: Integer;
f: File of real;
begin
AssignFile(f, goodrnd.mas);
for i := 1 to n do
for j:= 1 to 2 do
M[i,j] := 0;
Reset(f);
for i:=1 to N do
for j:= 1 to 2 do
begin
read(f,ch);
M[i,j]:=ch;
end;
CloseFile(f);
end;
end.


unit CTM;
interface

uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, StdCtrls, ExtCtrls, ComCtrls, Math;

Const N = 100; Ef = 5.71; k = 1; Np = 150;
Rm = 200; tp = 1; fi0 = 4.5; R = 5;
type
Tprfgrm = array[-Np..Np]of TPoint;
Trndm = array[1..2]of real;
TForm1 = class(TForm)
Button1: TButton;
Image1: TImage;
Label1: TLabel;
Label4: TLabel;
Label5: TLabel;
Button2: TButton;
Memo1: TMemo;
Label2: TLabel;
CombZ0: TComboBox;
CombU: TComboBox;
CombAlf: TComboBox;
Label3: TLabel;
procedure Button1Click(Sender: TObject);
procedure FormCreate(Sender: TObject);
procedure Button2Click(Sender: TObject);
procedure Clear;
private
function GetZ(rr:real; x:real; bet:real; tp:integer):real;
function GetJz(rr:real; x:real; bet:real; tp:integer):real;
function Integral1(betn,betm:real; rr,xx:real; tp:integer):real;
function Integral2(rrn,rrm,xx:real;tp:integer;p:boolean):real;
function GetJ(xx:real):real;
function GetPrfgrm: Tprfgrm;
public
{ Public declarations }
end;

var
Form1: TForm1;
Z0, alf, U, Xp, vmt :real;
M:array[1..N]of Trndm;
prfgrm: Tprfgrm;
koefx, koefy: integer;
xOy: TPoint;

implementation
{$R *.dfm}

function TForm1.GetZ(rr, x, bet: real; tp:integer): real;
var Z:real;
begin
If tp=0 then Z:=sqrt(sqr(Z0+R-sqrt(abs(R*R-rr*rr)))+x*x+rr*rr-2*x*rr*cos(bet))
else if tp=1 then Z:=sqrt(sqr(Z0+R*(1-1/sin(alf))+rr/tan(alf))+x*x+rr*rr-2*x*rr*cos(bet));
Result:=Z;
end;

function TForm1.GetJz(rr:real; x:real; bet:real; tp:integer):real;
var S1,S2,fi,Zxy:real;
begin
Zxy:=GetZ(rr,x,bet,tp);
S1:=3/(k*fi0); S2:=Zxy*(1-23/(3*fi0*k*Zxy+10-2*U*k*Zxy))+S1;
fi:=fi0-U*(S1+S2)/(2*Zxy)-2.86/(k*(S2-S1))*ln(S2*(Zxy-S1)/(S1*(Zxy-S2)));
result:=1620*U*Ef*exp((-1.025)*Zxy*sqrt(abs(fi)));
end;

function TForm1.Integral1(betn,betm:real; rr,xx:real; tp:integer):real;
var sum:real;
V,bet:real;
i:integer;
begin
sum:=0; V:=betm-betn;
for i:=1 to N do
begin
bet:=betn+(betm-betn)*M[i,1];
if(bet>=betn)and(bet<=betm) then
if tp=0 then sum:=sum+GetJz(rr,xx,bet,tp)*R*rr/sqrt(abs(sqr(R)-sqr(rr)))
else if tp=1 then sum:=sum+GetJz(rr,xx,bet,tp)*rr/sin(alf);
end;
sum:=sum*V/N;
result:=sum;
end;

function TForm1.Integral2(rrn,rrm,xx:real;tp:integer;p:boolean):real;
var i:integer;
sum, V, bp, Xd, rr: real;
begin
sum:=0; V:=rrm-rrn;
for i:=1 to N do
begin
rr:=rrn+(rrm-rrn)*M[i,1];
if(rr>=rrn)and(rr<=rrm) then
begin
if p then
begin
if xx=0 then bp:=0
else if tp=1 then bp:=arccos(Xp/xx)
else if tp=0 then
begin
Xd:=(R*(R/sqrt(sqr(R)-sqr(rr))-1)-Z0)*sqrt(sqr(R)-sqr(rr))/rr;
if Xd>0 then bp:=arccos(Xd/xx)
else bp:=arcsin(-Xd/xx)+PI/2;
end;
end
else bp:=PI;
sum:=sum+Integral1(-bp,bp,rr,xx,tp);
end;
end;
sum:=sum*V/N;
result:=sum;
end;

function TForm1.GetJ(xx:real):real;
const Max = 5;
var i:integer;
a1, a2, J, Rp1, Rp2, Rps: real;
begin
J:=0;
if xx=R then
begin
a1:=0; a2:=0
end
else
begin
a1:=(sqr(xx)*(R+Z0)-R*xx*sqrt(sqr(xx)+2*R*Z0+sqr(Z0)))/(sqr(xx)-sqr(R));
a2:=(sqr(xx)*(R+Z0)+R*xx*sqrt(sqr(xx)+2*R*Z0+sqr(Z0)))/(sqr(R)-sqr(xx));
end;
Rp1:=R*sqrt(1-sqr(R/(R+Z0-a1)));
if Rp1>Rm then Rp1:=Rm;
Rp2:=R*sqrt(1-sqr(R/(R+Z0+a2)));
if Rp2>Rm then Rp2:=Rm;
Rps:=R*cos(alf);
if Rps>Rm then Rps:=Rm;
if vrnt=1 then
begin
if xx<Xp thenJ:=J+Integral2(0,Rp1,xx,0,False)+Integral2(Rp1,Rp2,xx,0,True)
elseJ:=J+Integral2(0,Rp1,xx,0,False)+Integral2(Rp1,Rps,xx,0,True)+Integral2(Rps,Rm,xx,1,True);
end
else if vrnt=2 then
begin
if xx<Xp then J:=J+Integral2(0,Rps,xx,0,False)+Integral2(Rps,Rm,xx,1,False)
elseJ:=J+Integral2(0,Rp1,xx,0,False)+Integral2(Rp1,Rps,xx,0,True)+Integral2(Rps,Rm,xx,1,True);
end;
result:=2*PI*J;
end;

function TForm1.GetPrfgrm: Tprfgrm;
const e = 0.05; step = 0.1;
type Point=record
x, y:real;
y:real;
end;
varprfgrmr:array[-Np..Np]of Point;
prfgrm: Tprfgrm;
i:integer;
J0, max,process:real;
begin
max:=0;
for i:=0 to Np do
begin
prfgrmr[i].x:=0.1*i; prfgrmr[-i].x:=-prfgrmr[i].x; prfgrmr[i].y:=GetJ(prfgrmr[i].x);
if prfgrmr[i].y>max then max:=prfgrmr[i].y;
prfgrmr[-i].y:=prfgrmr[i].y;
form1.memo1.Lines.add(floattostr(trunc(100*i/Np))+ %);
end;
for i:=-Np to Np do
begin
prfgrm[i].x:=round(xOy.x+koefx*prfgrmr[i].x); prfgrm[i].y:=round(xOy.y-koefy*prfgrmr[i].y/max);
process:=50*i/Np+50; form1.Memo1.Lines.Add(floattostr(trunc(process))+ %);
end;
Result:=prfgrm;
end;


procedure TForm1.Button1Click(Sender: TObject);
var i:integer;
begin
z0:=StrToFloat(CombZ0.text);
if combU.ItemIndex=0 then U:=0.1;
if combU.Itemindex=1 then U:=0.01;
alf:=StrToFloat(CombAlf.Text)*PI/180; alf:=alf/2;

Xp:=abs(R+Z0-R/sin(alf))*tan(alf);
if (R+Z0)<(R/sin(alf)) then vrnt:=1
else vrnt:=2;
With Image1.Canvas do
begin
Pen.Color := RGB(0,0,250); //Профилограмма
Pen.Width := 2; prfgrm:=GetPrfgrm; Polyline(prfgrm);
end;
end;

procedure TForm1.FormCreate(Sender: TObject);
var i,j:integer;
begin
for i:=1 to N do
for j:= 1 to 2 do
M[i,j]:=random;
xOy.x:=round(Image1.Width/2); xOy.Y:=Image1.Height-15;
koefy:=trunc((xOy.y-20)); koefx:=trunc((Image1.Width-40)/round(Np/5));
Clear; form1.Memo1.Text:=;
end;

procedure TForm1.Button2Click(Sender: TObject);
begin
Clear; form1.memo1.Text:=;
end;

procedure TForm1.Clear;
var i,j:integer;
begin
With Image1.Canvas do
begin
Brush.Color := RGB(255,255,255); Pen.Color := RGB(255,255,255);
Rectangle(0, 0, Image1.Width, Image1.Height);
Pen.Color := RGB(0,0,0); //Оси
Pen.Width := 1; Pen.Style := psSolid;
MoveTo(xOy.x,0);
LineTo(xOy.x-3,15); LineTo(xOy.x+3,15); LineTo(xOy.x,0);
LineTo(xOy.x,xOy.y); MoveTo(Image1.Width,xOy.y);
LineTo(Image1.Width-15,xOy.Y-3); LineTo(Image1.Width-15,xOy.Y+3); LineTo(Image1.Width,xOy.y); LineTo(0,xOy.y);
Pen.Color := RGB(150,150,150); //Подписи осей и сетка
Pen.Style := psDot;
for i:=-round(Np/10) to round(Np/10) do
begin
if i<>0 then
begin
MoveTo(xOy.x+i*koefx,xOy.y); LineTo(xOy.x+i*koefx,0);
end;
TextOut( xOy.x+i*koefx-4, xOy.y+2, FloatToStr(i));
end;
for i:=1 to 10 do
begin
MoveTo(0,round(xOy.y-i*0.1*koefy)); LineTo(Image1.Width,round(xOy.y-i*0.1*koefy));
TextOut( xOy.x-15, round(xOy.y-i*0.1*koefy)-4, FloatToStr(0.1*i));
end;
TextOut( xOy.x+6, 4, J/Jmax);
end;
end;
END.




Построим с помощью разработанной программы профилограммы:

Рисунок6 – Профилограммы, построенныев разработанной программе
Внизу построена заданная поверхность, далее снизу вверх соответственно изображены профилограммы для Z0 = 5; Z0 = 7; Z0 = 12; Z0 = 15.Все размеры на осях указаны в ангстремах. 
Построим с помощью программы графики отношения очередной плотности тока к максимальной плотности тока :

Рисунок 7 – Семейство графиков описывающих разрешающую способность СТМ при Z0 = 5; Z0 = 7; Z0 = 10; Z0 = 15 (снизу вверх соответственно)
По оси абсцисс расположено смещение на подложке относительно оси иглы (в ангстремах), по оси ординат – относительные единицы.

Рисунок 8 – Семейство графиков описывающих разрешающую способность СТМ при угле конуса острияα = 10; α = 30;α = 60; α = 90 (снизу вверх соответственно) 
Заключение
В данной курсовой работе были разработаны:
программа для тестирования равномерности распределения последовательности случайных чисел;
программа для построения профилограмм, описывающих движение точечного острия сканирующего туннельного микроскопа над заданной поверхностью;
программа для построения графиков с целью исследования разрешающей способности поверхности конусообразного сканирующего туннельного микроскопа.
Результаты работы программы для тестирования равномерности распределения последовательности случайных чисел использовались в программах 2 и 3 для высчитывания определенных интегралов.
Вторая программа имитирована работу СТМ. Поверхность, по которой «проходила» вершина иглы сканирующего туннельного микроскопа получилась гладкой, без резких перепадов, что говорит о близости моделирования работы реального туннельного микроскопа.
Третья программа позволяет оценить зависимость разрешающей способности СТМ от параметровострия, например, при увеличении радиуса закругления острия и при прочих равных параметрах наблюдается существенное снижение разрешающей способности.
Графики зависимостей плотностей токов от отдаленности от острия конусообразной иглы, остаются примерно одинаковыми при разных входных параметрах (угла вершины конуса, высоты над плоской поверхностью и т.д.). Причем, плотность тока достигает своего максимума на вершине конусообразной иглы. На высоте порядка 4 ангстрем ток практически перестает существовать.

Список литературы
Эдельман В.С. Сканирующая туннельная микроскопия // ПТЭ, 1989.– №5.– С.25-49.
Кук И., Силверман П. Растровая туннельная микроскопия // ПНИ, 1989.– №2.– С.3-22.
Елинсон М.И., Васильев Г.Ф. Автоэлектронная эмиссия.- М.: Физмат. ГИЗ, 1958.– 272с.
Симмонс Дж.Г. Силы изображения в туннельных переходах металл-окисел-металл // Туннельные явления в твёрдых телах / Под ред. Э.Бурштейна и С.Лундквиста.– М.: Мир. 1973.– С.131-142.
Бусленко Н.П., Шрейдер Ю.А. Метод статистических испытаний.– М.: ГИФМЛ, 1961.– 226с.
Шелковников Е.Ю. Исследование параметров туннельного микроскопа с применением метода Монте-Карло // Информационно-измерительные системы на базе наукоёмких технологий. - Ижевск, 1998. - С.95-102.
Кнут Д. Искусство программирования для ЭВМ.– М.: Мир, 1977.–Т.2.– 724с.
Шелковников Е.Ю.Методические указания к выполнению курсового проектастудентов по дисциплине «Моделирование»