на тему “программное обеспечение для моделирования элионных технологий” 5 (стр. 2 из 5)

В 1950х годах метод использовался для расчетов при разработке водородной бомбы. Основные заслуги в развитии метода в это время принадлежат сотрудникам лабораторий ВВС США и корпорации RAND.

В 1970х годах в новой области математики — теории вычислительной сложности было показано, что существует класс задач, сложность (количество вычислений, необходимых для получения точного ответа) которых растет с размерностью задачи экспоненциально. Иногда можно, пожертвовав точностью, найти алгоритм, сложность которого растет медленнее. Но для большинства задач этого нельзя сделать (например, задача определения объема выпуклого тела в n-мерном евклидовом пространстве) и метод Монте-Карло является единственной возможностью для получения достаточно точного ответа за приемлемое время. В настоящее время основные усилия исследователей направлены на создание эффективных Монте-Карло алгоритмов различных физических, химических и социальных процессов для параллельных вычислительных систем.

1.2 Основной принцип реализации метода Монте-Карло

Метода Монте-Карло - специальный метод изучения поведения заданной статистики при проведении многократных повторных выборок, существенно использующий вычислительные возможности современных компьютеров. При проведении анализа по методу Монте-Карло компьютер использует процедуру генерации псевдослучайных чисел для имитации данных из изучаемой генеральной совокупности. Процедура анализа по методу Монте-Карло состоит в построении выборки из генеральной совокупности в соответствии с указаниями пользователя, а затем производит следующие действия:

Для каждого повторения по методу Монте-Карло:

1. Имитирует случайную выборку из генеральной совокупности.

2. Проводит анализ выборки.

3. Сохраняет результаты.

После большого числа повторений, сохраненные результаты хорошо имитирует реальное распределение выборочной статистики. Метод Монте-Карло позволяет получить информацию о выборочном распределении в случаях, когда обычная теория выборочных распределений оказывается бессильной.

1.3 Генераторы случайных чисел

Существует много задач, в которых естественным образом присутствует случайный элемент. Например, многие типы сбоев в коммуникационных системах моделируются лучше всего как случайные эффекты. Есть и другие задачи, которые сами по себе не содержат элемента случайности, но где по тем или иным причинам полезно брать случайные выборки(например, полная перепись населения).

Еще одна задача, где случайный выбор оказывается удивительно полезным, - это грубая оценка объемов сложных областей в евклидовом пространстве более чем четырех измерений. Разумеется, сюда входит и приближенное вычисление и интегралов методом Монте-Карло.

Для любого из подобных приложений, связанных с ЭВМ, необходимо вырабатывать большие последовательности чисел, ведущих себя как случайные выборки из заданного распределения. Имеются три класса методов для выборки таких случайных чисел:

1. Можно использовать физический случайный процесс, например момент эмитирования электрона горячим катодом, интерпретируемый как фаза в неком временном цикле. Такие эффекты легко получить, но трудно добиться, чтобы связанные с ними результаты не содержали систематической ошибки, чему приходиться противодействовать другими методами. Кроме того, поскольку физический процесс случаен и невоспроизводим, машинную программу, основанную на нем, отлаживать очень нелегко.

2. Можно вычислять случайные числа отдельно от их использования и хранить их в виде файла на диске. Этот подход был широко распространен в не знавшую ЭВМ эпоху. Использование этого подхода с ЭВМ поднимает следующую проблему – проблему выбора случайной позиции в файле.

3. Наиболее распространенный метод заключается в применении алгоритма, вырабатывающего псевдослучайные числа, немедленно потребляемы программой.

Для данной программы был выбран алгоритм генерации псевдослучайных чисел, описанный в книге [5]. Достоинствами этого алгоритма является хорошая равномерность распределения случайных величин и скорость генерации. Для инициализации алгоритма используется значение, генерируемое стандартным генератором случайных величин.

1.4 Используемая физическая реализация метода Монте-Карло

Метод Монте-Карло (МК) широко используется для моделирования ионного легирования и физического распыления металлов, полупроводников и диэлектриков уже более двадцати лет [1]. Причем, в последнее время он получает все большее распространение в связи с гигантским скачком в производительности современных ПЭВМ, что отводит на второй план его основной недостаток − огромный объем вычислений.

Первой, по настоящему мощной программой реализации метода МК являлась программа TRIM 85 (Transport and Ranges of Ion in Matter, Fortran 77), написанная Бирзаком с соавторами [1]. В варианте, описанном в [1], мишень могла состоять из трех слоев, а каждый из них мог содержать до семи элементов. TRIM85 2D – базовая версия бескаскадного варианта программы стала основой всех последующих каскадных версий: TRIM88, TRSP2C, TRSPV1C, TRIM91, TRIM95, SRIM2001, 2003 (Stopping and Ranges of Ion in Matter), SRIM 2006, SRIM 2008.

Все перечисленные программы, как аналитические теории и большинство других программ расчетов методом Монте-Карло, базируется на парных столкновениях ионов с первоначально покоящимися атомами мишени. На протяжении последних лет был значительно усовершенствован основной физический подход к моделированию с помощью программы TRIM. Главной отличительной чертой программ SRIM(TRIM) является использование для расчета угла рассеяния в СЦМ так называемой “магической” формулы.

Широкое использование TRIM(SRIM)-алгоритма в задачах моделирования движения ускоренных ионов в твердом теле во многом обусловлено простотой и вместе с тем достаточной корректностью описания процесса упругого рассеяния ионов на атомах мишени. Основная проблема рассмотрения процесса упругого рассеяния в рамках метода МК состоит в многократном расчете угла рассеяния q налетающей частицы в СЦМ. Так при моделировании 10⁴ ионных траекторий в режиме учета движения атомов отдачи (10³…10⁴) приходится рассчитывать до 10¹⁰ таких интегралов, что для ПЭВМ конца прошлого столетия было невозможно. В большинстве известных программ данный вопрос решается путем использования аналитических формул, основанных на ряде грубых приближений, часто неприемлемых, как в случае тяжелых ионов и низких энергий. Эти недостатки свойственны наиболее известным ранним версиям программам реализации метода МК – PIBER (Program for Ion Beam Exposure of Resist) и MARLOWE. Программа TRIM отличается от всех остальных применением несложного аналитического выражения для расчета q, которое часто называют “магической” формулой.

Для расчета Q во всех программах TRIM(SRIM) используется соотношение:

cos(Q/2) = (r + p + d)/(r + r0),

(1)

где r = r₁ + r₂; d = d₁ + d₂; r₁, r₂ – радиусы кривизны траекторий; d₁, d₂ – небольшие поправочные параметры; p – прицельный параметр; r₀ – расстояние наибольшего сближения. Вывод данной формулы вытекает непосредственно из геометрии рассеяния (рис.1.)

Расстояние наибольшего сближения определяется по формуле (2).

1 – [ V(r₀)/E_отн)] – ( b²/r²₀) =0 (2)

Рисунок 1 – Геометрия рассеяния в системе центра масс

Из этой формулы методом Ньютона за два-три шага можно определить r₀ c точностью порядка 0,1%. Радиус кривизны r в СЦМ определяется исходя из определения центробежной силы f_c .

Было бы удобнее выразить E_r в величинах Z₁·Z₂·e²/a и длину в единицах а. Таким образом, вводится понятие приведенной энергии:

e = a·E_r/(Z₁·Z₂·e²), (3)

где а – длина экранирования.

В результате всего перечисленного, формула (1) принимает вид:

cos(Q/2) = (B + R_c + D)/(R₀ + R_c), (4)

где B = p/a; R_c = r/a; R₀ = r₀/a; D = d/a – величины, выраженные в единицах длины экранирования.

При высоких энергиях атомные столкновения могут быть описаны при помощи неэкранированного потенциала Кулона, то есть резерфордовского рассеяния. Таким образом, формула для Q должна асимптотически приближаться к резерфордовскому результату с увеличением приведенной энергии e. Для достижения этого записываем формулу для Q при резерфордовском рассеянии, а затем при помощи подгоночных параметров переходим к случаю с меньшей энергией:

D = A·(R₀ – B)/(1 + G),

где A = 2·a·e·B^b; G = g·[(1 + A²)^1/2 – A]^-1; a = 1 + C₁·e^-1/2; b = (C₂ + e^1/2)/(C₃ + e^1/2); g = (C₄ + e)/(C₅ + e). (5)

Величины С₁…С₅ – подгоночные параметры, определяемые потенциалом ион-атомного взаимодействия. Значения C_i рассчитаны для трех наиболее распространенных потенциалов: универсального, Kr-C и мольеровского и приведены в табл.1.1.

Отметим, что в базовой версии TRIM-85 [1] допускается возможность выбора потенциала ион-атомного взаимодействия. В дальнейших версиях авторы применяют только универсальный потенциал, что, видимо, не совсем правильно.

Все вышесказанное не используется для больших значений e. При e > 10 для увеличения скорости вычислений используется кулоновский потенциал, так как все заметные отклонения и передача энергии происходят при достаточно близких столкновениях. Тогда: